csu-hmc / opty

A library for using direct collocation in the optimization of dynamic systems.
http://opty.readthedocs.io
Other
86 stars 20 forks source link

Problem with kinetic energy as objective function #157

Closed Peter230655 closed 1 month ago

Peter230655 commented 1 month ago

When I use $\int\text{kinetic energy} \ dt$ as objective function, I get this error message:; 24-05-13 error messagw with kinetic energy.txt

When I use another objective function, say, $\int F^2 \ dt$, where F is the driving force, no error happens. I ran across this same problem with several opty simulations.

NOTE: to get obj and obj_grad I always use a function, which I call Timo's function, (I got it from @tjstienstra from one of his programs) so conceivably the problem could be there?

Below is some code, the shortest simulation I have: A pendulum is rotably attached to a car. This car may be moved horizontally (along the x axis) driven by a force F. The task is to move the car in such a way, that the pendulum becomes upright.

In line 206 of this code ther is a variable kinetik. If kinetik = True, the kinetic energy will be used as objective function (and will cause the error) If kinetik = False, the other objective function is used and all seems fine.

Any help is greatly appreciated!

#A pendulum is attach 
to the center of cart which rolls in the X direction.\
#A force acts on the cart in X direction, the gravity acts in the negative Y direction.\
#The goal is to get the pendulum in the upright position by suitably moving the cart to the left and to the right.\
#This is 99.9% a copy of the example of the revered pendulum from the docs of opty and of Timo's simulation sent to me.
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
import sympy.physics.mechanics as me
from collections import OrderedDict
import time
import numpy as np
import sympy as sm
from opty.direct_collocation import Problem
from opty.utils import building_docs
import matplotlib.pyplot as plt
import matplotlib
import matplotlib.animation as animation
from IPython.display import HTML
matplotlib.rcParams['animation.embed_limit'] = 2**128
from matplotlib import patches
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# function to create the obj and obj_grad
# I call it Timo's function because I found it in his program
def create_objective_function(
    objective, state_symbols, input_symbols, unknown_symbols,
    N, node_time_interval=1.0, integration_method="backward euler"):
    """Creates function to evaluate the objective and objective gradient.

    Parameters
    ----------
    objective : sympy.Expr
        The objective function to be minimized, which is a function of the
        states and inputs. The objective function can contain non-nested
        indefinite integrals of time, e.g. ``Integral(f(t)**2, t)``.
    state_symbols : iterable of symbols
        The state variables.
    input_symbols : iterable of symbols
        The input variables.
    unknown_symbols : iterable of symbols
        The unknown parameters.
    N : int
        Number of collocation nodes, i.e. the number of time steps.
    node_time_interval : float
        The value of the time interval. The default is 1.0, as this term only
        appears in the objective function as a scaling factor.
    integration_method : str, optional
        The method used to integrate the system. The default is "backward
        euler".

    """
    def lambdify_function(expr, multiplication_array, take_sum):
        if take_sum:
            def integration_function(x):
                return node_time_interval * np.sum(x * multiplication_array)
        else:
            def integration_function(x):
                return node_time_interval * x * multiplication_array
        return sm.lambdify(
            (states, inputs, params), expr,
            modules=[{int_placeholder.name: integration_function}, "numpy"],
            cse=True)

    def parse_expr(expr, in_integral=False):
        if not expr.args:
            return expr
        if isinstance(expr, sm.Integral):
            if in_integral:
                raise NotImplementedError("Nested integrals are not supported.")
            if expr.limits != ((me.dynamicsymbols._t,),):
                raise NotImplementedError(
                    "Only indefinite integrals of time are supported.")
            return int_placeholder(parse_expr(expr.function, True))
        return expr.func(*(parse_expr(arg) for arg in expr.args))

    # Parse function arguments
    states = sm.ImmutableMatrix(state_symbols)
    inputs = sm.ImmutableMatrix(sort_sympy(input_symbols))
    params = sm.ImmutableMatrix(sort_sympy(unknown_symbols))
    if states.shape[1] > 1 or inputs.shape[1] > 1 or params.shape[1] > 1:
        raise ValueError(
            'The state, input, and unknown symbols must be column matrices.')
    n, q = states.shape[0], inputs.shape[0]
    i_idx, r_idx = n * N, (n + q) * N

    # Compute analytical gradient of the objective function
    objective_grad = sm.ImmutableMatrix([objective]).jacobian(
        states[:] + inputs[:] + params[:])

    # Replace the integral with a custom function
    int_placeholder = sm.Function("_IntegralFunction")
    objective = parse_expr(objective)
    objective_grad = tuple(parse_expr(objective_grad))

    # Replace zeros with an array of zeros, otherwise lambdify will return a
    # scalar zero instead of an array of zeros.
    objective_grad = tuple(
        np.zeros(N) if grad == 0 else grad for grad in objective_grad[:n + q]
        ) + tuple(objective_grad[n + q:])

    # Define evaluation functions based on the integration method
    if integration_method == "backward euler":
        obj_expr_eval = lambdify_function(
            objective, np.hstack((0, np.ones(N - 1))), True)
        obj_grad_time_expr_eval = lambdify_function(
            objective_grad[:n + q], np.hstack((0, np.ones(N - 1))), False)
        obj_grad_param_expr_eval = lambdify_function(
            objective_grad[n + q:], np.hstack((0, np.ones(N - 1))), True)
        def obj(free):
            states = free[:i_idx].reshape((n, N))
            inputs = free[i_idx:r_idx].reshape((q, N))
            return obj_expr_eval(states, inputs, free[r_idx:])

        def obj_grad(free):
            states = free[:i_idx].reshape((n, N))
            inputs = free[i_idx:r_idx].reshape((q, N))
            return np.hstack((
                *obj_grad_time_expr_eval(states, inputs, free[r_idx:]),
                obj_grad_param_expr_eval(states, inputs, free[r_idx:])
             ))

    elif integration_method == "midpoint":
        obj_expr_eval = lambdify_function(
            objective, np.ones(N - 1), True)
        obj_grad_time_expr_eval = lambdify_function(
            objective_grad[:n + q], np.hstack((0.5, np.ones(N - 2), 0.5)),
            False)
        obj_grad_param_expr_eval = lambdify_function(
            objective_grad[n + q:], np.ones(N - 1), True)
        def obj(free):
            states = free[:i_idx].reshape((n, N))
            states_mid = 0.5 * (states[:, :-1] + states[:, 1:])
            inputs = free[i_idx:r_idx].reshape((q, N))
            inputs_mid = 0.5 * (inputs[:, :-1] + inputs[:, 1:])
            return obj_expr_eval(states_mid, inputs_mid, free[r_idx:])

        def obj_grad(free):
            states = free[:i_idx].reshape((n, N))
            states_mid = 0.5 * (states[:, :-1] + states[:, 1:])
            inputs = free[i_idx:r_idx].reshape((q, N))
            inputs_mid = 0.5 * (inputs[:, :-1] + inputs[:, 1:])
            return np.hstack((
                *obj_grad_time_expr_eval(states, inputs, free[r_idx:]),
                obj_grad_param_expr_eval(states_mid, inputs_mid, free[r_idx:])
            ))

    else:
        raise NotImplementedError(
            f"Integration method '{integration_method}' is not implemented.")

    return obj, obj_grad

def sort_sympy(seq):
    """Returns a sorted list of the symbols."""
    seq = list(seq)
    try:  # symbols
        seq.sort(key=lambda x: x.name)
    except AttributeError:  # functions
        seq.sort(key=lambda x: x.__class__.__name__)
    return seq
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# set up Kane's equations
start = time.time()

N, A      = sm.symbols('N A', cls = me.ReferenceFrame)
t         = me.dynamicsymbols._t
O, P1, P2 = sm.symbols('O P1 P2', cls = me.Point)                   # O fixed point in N, P1 center of mass of the cart, P2 center of mass of the pendulum 
O.set_vel(N, 0)                                                     # Velocity of the fixed point is zero   
q1, q2, u1, u2, F = me.dynamicsymbols('q1 q2 u1 u2 F')              # Generalized coordinates of the cart and pendulum; F is the force applied to the cart
l, m1, m2, g, iZZ = sm.symbols('l, m1, m2, g, iZZ')                 # Length of the pendulum, masses of the pendulum and cart, gravitational constant, and moment of inertia of the pendulum

A.orient_axis(N, q2, N.z)
A.set_ang_vel(N, u2 * N.z)

P1.set_pos(O, q1 * N.x)
P1.set_vel(N, u1 * N.x)
P2.set_pos(P1, l * A.x)
P2.v2pt_theory(P1, N, A)

P1a = me.Particle('P1a', P1, m1)                                    # Particle for the cart

I = me.inertia(A, 0, 0, iZZ)
P2a = me.RigidBody('P2a', P2, A, m2, (I, P2))                       # Rigid body for the pendulum
BODY = [P1a, P2a]

FL = [(P1, F * N.x - m1*g*N.y), (P2, -m2*g*N.y)]                    # Force applied to the cart and gravitational force on the pendulum
kd = sm.Matrix([q1.diff() - u1, q2.diff() - u2])                    # Kinematic differential equations

q_ind = [q1, q2]
u_ind = [u1, u2]

KM = me.KanesMethod(N, q_ind=q_ind, u_ind=u_ind, kd_eqs=kd)
(fr, frstar) = KM.kanes_equations(BODY, FL)
EOM = kd.col_join(fr + frstar)
EOM
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# preparation for opty
state_symbols = tuple((*q_ind, *u_ind))
constant_symbols = (l, m1, m2, g, iZZ)
specified_symbols = (F,)
unknown_symbols = []
methode = "backward euler"

#=====================================================================================================================
#=====================================================================================================================
kinetik = False     # if True the kinetic energy is the objective function, else the force**2 is the objective function
#=====================================================================================================================
#=====================================================================================================================

target_angle = np.pi /2.0                       # The desired angle for the pendulum to swing up to.
duration = 15.0                                 # The duration of the swing up maneuver.
num_nodes = 250                                 # The number of nodes to use in the direct collocation problem.
save_animation = False

interval_value = duration / (num_nodes - 1)

# Specify the known system parameters.
par_map = OrderedDict()
par_map[l]   = 2.0
par_map[m1]  = 1.0
par_map[m2]  = 1.0
par_map[g]   = 9.81
par_map[iZZ] = 2.0

kin_energie = sum([body.kinetic_energy(N) for body in BODY])
print('kinetic energy FS', kin_energie.free_symbols)
print('kinetic energy DS', me.find_dynamicsymbols(kin_energie))
print('kinetic energy', kin_energie)
if kinetik == False:
    objective = sm.Integral(F**2, t)
else:
    objective = sm.Integral(kin_energie, t)

# Specify the objective function and it's gradient.
obj, obj_grad = create_objective_function(
    objective, 
    state_symbols, 
    specified_symbols, 
    unknown_symbols,
    num_nodes, 
    interval_value, 
    methode)

# This portion is copied from Timo's code
t0, tf = 0.0, duration

# Opty want instance constraints, basically things like `q2(2.0) - 3.141` to specify that q2 should be 3.141 when t=2.0
# I personally like to rewrite it to initial state and final state constraints and then translate them into the set of equations opty wants.
initial_state_constraints = {q1: 0.0, u1: 0.0, q2: -np.pi/2., u2: 0.0}
final_state_constraints = {q1: 5., q2: np.pi/2., u2: 0.0}
instance_constraints = (
) + tuple(
    xi.subs({t: t0}) - xi_val for xi, xi_val in initial_state_constraints.items()
) + tuple(
    xi.subs({t: tf}) - xi_val for xi, xi_val in final_state_constraints.items()
)
print('instance constraints:', instance_constraints)
bounds = {F: (-15., 15.)}
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# solve the optimization problem
prob = Problem(obj, obj_grad, EOM, state_symbols, num_nodes, interval_value,
               known_parameter_map=par_map,
               instance_constraints=instance_constraints,
               bounds=bounds)

# Use a random positive initial guess.
initial_guess = np.random.randn(prob.num_free)

# Find the optimal solution.
solution, info = prob.solve(initial_guess)
print('message from optimizer:', info['status_msg'])

#prob.plot_constraint_violations(solution)
prob.plot_objective_value()
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#Plot some results
fig, ax = plt.subplots(5, 1, figsize=(10, 25), sharex=True)
times = np.linspace(0.0, duration, num=num_nodes)

for i, j in enumerate(state_symbols + specified_symbols):
    ax[i].plot(times, solution[i * num_nodes:(i + 1) * num_nodes])
    ax[i].grid()
    ax[i].set_ylabel(str(j))
ax[-1].set_xlabel('time')
ax[0]. set_title('Generalized coordinates and speeds')
ax[-1]. set_title('Control force');
moorepants commented 1 month ago

Error message from the text file:

ypeError                                 Traceback (most recent call last)
Cell In[3], line 274
    271 initial_guess = np.random.randn(prob.num_free)
    273 # Find the optimal solution.
--> 274 solution, info = prob.solve(initial_guess)
    275 print('message from optimizer:', info['status_msg'])
    277 #prob.plot_constraint_violations(solution)

File c:\Users\Peter\anaconda3\envs\sympy-dev\Lib\site-packages\cyipopt\cython\ipopt_wrapper.pyx:658, in ipopt_wrapper.Problem.solve()

File c:\Users\Peter\anaconda3\envs\sympy-dev\Lib\site-packages\cyipopt\cython\ipopt_wrapper.pyx:871, in ipopt_wrapper.objective_cb()

File c:\Users\Peter\anaconda3\envs\sympy-dev\Lib\site-packages\sympy\core\expr.py:340, in Expr.__float__(self)
    338 if result.is_number and result.as_real_imag()[1]:
    339     raise TypeError("Cannot convert complex to float")
--> 340 raise TypeError("Cannot convert expression to float")

TypeError: Cannot convert expression to float
moorepants commented 1 month ago

The program runs fine for me, but on Linux.

moorepants commented 1 month ago

What happens when you run:

obj(initial_guess)

or

obj_grad(initial_guess)
Peter230655 commented 1 month ago

The program runs fine for me, but on Linux.

You mean the program that gives me the errors on Windows runs fine under Linux? (Similar what happened with sm.pi, I recall)

moorepants commented 1 month ago

Yes.

Peter230655 commented 1 month ago

What happens when you run:

obj(initial_guess)

or

obj_grad(initial_guess)

Frankly, I am not sure what you mean. Should I use opty.utils.create_objective_function?

moorepants commented 1 month ago

obj, and obj_grad are two functions you define in your program. I am suggesting to execute those functions and see if they raise an error.

Peter230655 commented 1 month ago

obj, and obj_grad are two functions you define in your program. I am suggesting to execute those functions and see if they raise an error.

My objective function is a sympy function, call it obj. So, I would form its gradient, call in obj_grad, then lambdify both functions?

moorepants commented 1 month ago

That is not what this line says:

obj, obj_grad = create_objective_function(
    objective, 
    state_symbols, 
    specified_symbols, 
    unknown_symbols,
    num_nodes, 
    interval_value, 
    methode)
Peter230655 commented 1 month ago

That is not what this line says:

obj, obj_grad = create_objective_function(
    objective, 
    state_symbols, 
    specified_symbols, 
    unknown_symbols,
    num_nodes, 
    interval_value, 
    methode)

this is what I call Timo's function. It is at the very top of my program, right after importing sympy, numpy, etc. It is this function that I use on all my opty simulations. I got it from Timo, probably it was in a simulation he sent we when I started to pay around with opty.

moorepants commented 1 month ago

I know what the function is. You may be misunderstanding what it does.

Peter230655 commented 1 month ago

I thought it returns obj and obj_grad, as needed in Problem? As it worked before, I never gave any thought what exactly it does.

Peter230655 commented 1 month ago

I wanted to see, if the create_objective_function from opty may work. So I wrote:

from opty.utils import parse_free from opty.utils import create_objective_function

Both are is section 1.3.3 of the documentation. While parse free works fine, with the second import line I get this error: _ImportError: cannot import name 'create_objectivefunction' from 'opty.utils' (c:\Users\Peter\anaconda3\envs\sympy-dev\Lib\site-packages\opty\utils.py)

What am I doing wrong? Could this be related to the fact that I use Timo's symy_dev environment? I do not think so, because I get the same error with my 'normal' environment. In both environments I have opty version 1.2.0 Thanks!

moorepants commented 1 month ago

You are attempting to import features from the development version of opty when you do not have it installed.

Peter230655 commented 1 month ago

You are attempting to import features from the development version of opty when you do not have it installed.

Clear. So, I guess what I call Timo's function will be similar to the utils.create_objective_function once the development version becomes the released version?

moorepants commented 1 month ago

The functions are likely identical. You can see the one in opty here: https://github.com/csu-hmc/opty/blob/master/opty/utils.py#L255

Peter230655 commented 1 month ago

obj, and obj_grad are two functions you define in your program. I am suggesting to execute those functions and see if they raise an error.

Now it finally clicked, and I ran obj(initial_guess) obj_grad(initial gues)

I case of the objective = kinetic energy, this sympy function depends on the constant_symbols, which I set in the dictionary par_map. When I print obj(initial_guess) something is printed, which contains the constant symbols. When I looked at Timo's function, it does not seem to ask for par_map, so how would it know them? obj_grad(initial_guess) throws an error: name l (this is one key of par_map, one of the constant_symbols) is not defined.

My other objective does not depend on the constant_symbols, nor did any of the other objectives in the other simulations I tested.

moorepants commented 1 month ago

It is helpful if you post what you input and what the output (or errors) are here.

Peter230655 commented 1 month ago

It is helpful if you post what you input and what the output (or errors) are here.

This is the input ( it is from the simulation I posted here earlier)

par_map = OrderedDict()
par_map[l]   = 2.0
par_map[m1]  = 1.0
par_map[m2]  = 1.0
par_map[g]   = 9.81
par_map[iZZ] = 2.0

kin_energie = sum([body.kinetic_energy(N) for body in BODY])
print('kinetic energy FS', kin_energie.free_symbols)
print('kinetic energy DS', me.find_dynamicsymbols(kin_energie))
print('kinetic energy', kin_energie)
if kinetik == False:
    objective = sm.Integral(F**2, t)
else:
    objective = sm.Integral(kin_energie, t)

# Specify the objective function and it's gradient.
obj, obj_grad = create_objective_function(
    objective, 
    state_symbols, 
    specified_symbols, 
    unknown_symbols,
    num_nodes, 
    interval_value, 
    methode)

initial_guess = np.random.randn(prob.num_free)
print('obj',obj(initial_guess))
print('obj_grad',obj_grad(initial_guess))

This is the output of the print, when kinetik = True, that is the kinetic energy is the objective. You can see, that sm.symbols, the constant_symbols like iZZ, m1, m2 are there.

obj 7.08275041081572*iZZ + 6.61823431244662*m1 + 0.0301204819277108*m2*(1.37096022470083e-9*l**2 + 4.23801010364116e-5*l + 0.708304255011999) + 0.0301204819277108*m2*(1.65933266568474e-8*l**2 - 0.000290092923805568*l + 1.281868764177) + 0.0301204819277108*m2*(3.15186980395656e-6*l**2 + 0.000379117416554737*l + 0.384058211065885) + 0.0301204819277108*m2*(2.39440095180923e-5*l**2 - 0.00072290019469842*l + 0.00560355131506347) + 0.0301204819277108*m2*(3.79817241492966e-5*l**2 + 0.0114667063491362*l + 1.29894922177963) + 0.0301204819277108*m2*(7.71710409400169e-5*l**2 + 0.000271800159458009*l + 0.0283893070100778) + 0.0301204819277108*m2*(0.000253613146936964*l**2 + 0.0437547616458357*l + 6.02462382723217) + 0.0301204819277108*m2*(0.000270038290596149*l**2 - 0.00764989721371629*l + 0.110206165640946) + 0.0301204819277108*m2*(0.000305428881611059*l**2 - 0.00994036640582995*l + 0.350285025378736) + 0.0301204819277108*m2*(0.000330008119004537*l**2 - 0.0061580190972643*l + 1.65482175209466) + 0.0301204819277108*m2*(0.000659344752412932*l**2 - 0.0276473025422879*l + 0.462369874459908) + 0.0301204819277108*m2*(0.000659510181127682*l**2 - 0.0708791300907772*l + 1.92991202347619) + 0.0301204819277108*m2*(0.000711479000513146*l**2 - 0.00386471780836697*l + 0.235948335669646) + 0.0301204819277108*m2*(0.000712915843477491*l**2 - 0.0470846451491102*l + 1.03002095367928) + 0.0301204819277108*m2*(0.000745420394753938*l**2 + 0.048957635886274*l + 0.808768230029053) + 0.0301204819277108*m2*(0.00215159121345072*l**2 - 0.0465807967800967*l + 2.07434702507136) + 0.0301204819277108*m2*(0.00247619623342746*l**2 + 0.00244740216013135*l + 0.00108551241658223) + 0.0301204819277108*m2*(0.00261726156621969*l**2 - 0.0410009015196672*l + 0.941984154846162) + 0.0301204819277108*m2*(0.00287881841871155*l**2 + 0.0584739601271916*l + 0.344809543961222) + 0.0301204819277108*m2*(0.00303940214763494*l**2 + 0.0617146109056688*l + 3.23004756287086) + 0.0301204819277108*m2*(0.00318397683828226*l**2 - 0.00974611255623183*l + 1.12452955009995) + 0.0301204819277108*m2*(0.00390282810079114*l**2 - 0.0102746255327791*l + 0.0450734675202643) + 0.0301204819277108*m2*(0.00412803841575544*l**2 + 0.112912742451682*l + 1.19898365802908) + 0.0301204819277108*m2*(0.0061414002773282*l**2 - 0.00130700666915758*l + 0.12413845466657) + 0.0301204819277108*m2*(0.00649393180506241*l**2 - 0.0558621141157683*l + 0.527224983681373) + 0.0301204819277108*m2*(0.00658136139389158*l**2 - 0.00116994153173267*l + 0.000147717929751299) + 0.0301204819277108*m2*(0.00844285749677711*l**2 + 0.0123888151876511*l + 0.00454501173576538) + 0.0301204819277108*m2*(0.00917787705508385*l**2 - 0.031424256983172*l + 0.128959925426516) + 0.0301204819277108*m2*(0.0101135675015842*l**2 - 0.220378839181579*l + 1.20227129990794) + 0.0301204819277108*m2*(0.0102348134392688*l**2 + 0.0103649267779385*l + 0.018122630804559) + 0.0301204819277108*m2*(0.0104115572843454*l**2 - 0.0720495751704606*l + 0.632189997823937) + 0.0301204819277108*m2*(0.0106553337253021*l**2 - 0.0283721130376657*l + 0.0222219285890145) + 0.0301204819277108*m2*(0.0109112278909702*l**2 - 0.0695602227769736*l + 0.738511822856318) + 0.0301204819277108*m2*(0.0123649552817868*l**2 - 0.0191206983017037*l + 0.0238709774827859) + 0.0301204819277108*m2*(0.0127683252682951*l**2 - 0.106462633691716*l + 1.29904779375213) + 0.0301204819277108*m2*(0.0149646516814596*l**2 + 0.0756322252370686*l + 0.958945935292591) + 0.0301204819277108*m2*(0.0162190507391938*l**2 - 0.128810681012421*l + 1.34359357567967) + 0.0301204819277108*m2*(0.0193916987820967*l**2 - 0.0466150954242617*l + 0.685045986238947) + 0.0301204819277108*m2*(0.0222362351656691*l**2 + 0.226506836053932*l + 0.642196534390226) + 0.0301204819277108*m2*(0.0224342027643725*l**2 - 0.380859443854557*l + 4.68868726519033) + 0.0301204819277108*m2*(0.0225355797551932*l**2 - 0.125949269667664*l + 0.92422588937501) + 0.0301204819277108*m2*(0.022668022725448*l**2 - 0.539302056892285*l + 7.04359905518912) + 0.0301204819277108*m2*(0.0228690727268977*l**2 - 0.00172431863312592*l + 6.70134917887404e-5) + 0.0301204819277108*m2*(0.0230353003595574*l**2 + 0.0436371227030878*l + 0.926254569247128) + 0.0301204819277108*m2*(0.0245809357245742*l**2 + 0.0616791919351306*l + 0.380346784017252) + 0.0301204819277108*m2*(0.026141303821014*l**2 + 0.380455122731175*l + 1.41456562036976) + 0.0301204819277108*m2*(0.026970224296859*l**2 + 0.0421741052949763*l + 0.0164897207779384) + 0.0301204819277108*m2*(0.0281686934701908*l**2 + 0.0403380044665294*l + 0.0345883419398707) + 0.0301204819277108*m2*(0.0299114345409024*l**2 - 0.00190048272658487*l + 0.021119789257449) + 0.0301204819277108*m2*(0.0318422461588351*l**2 + 0.0197344407613327*l + 0.0178879746941327) + 0.0301204819277108*m2*(0.0325664056332679*l**2 - 0.33172879369327*l + 1.459008308485) + 0.0301204819277108*m2*(0.0363183372835989*l**2 + 0.402398986640922*l + 2.18225816551529) + 0.0301204819277108*m2*(0.036972019455468*l**2 + 0.035894348709945*l + 0.562895163197006) + 0.0301204819277108*m2*(0.0447747300068689*l**2 - 0.0982096937284161*l + 0.252887626243105) + 0.0301204819277108*m2*(0.0511363750233591*l**2 + 0.00323740348727349*l + 0.33924397092965) + 0.0301204819277108*m2*(0.052118423105867*l**2 - 0.0147362381428838*l + 0.00334608299152447) + 0.0301204819277108*m2*(0.0533598933991968*l**2 + 0.00662485209811963*l + 0.000232014167853789) + 0.0301204819277108*m2*(0.0553528944114719*l**2 - 0.657736079470405*l + 2.67303354483774) + 0.0301204819277108*m2*(0.0608388363632886*l**2 + 0.125761029005016*l + 0.100305535923108) + 0.0301204819277108*m2*(0.0608915567094739*l**2 - 0.191445614864388*l + 0.170330221168256) + 0.0301204819277108*m2*(0.0625735860436686*l**2 - 0.119014275348076*l + 0.861692467168893) + 0.0301204819277108*m2*(0.0626251465505362*l**2 - 0.067415581404847*l + 0.0723960557571266) + 0.0301204819277108*m2*(0.0643495498689684*l**2 + 0.251233663150443*l + 0.443209004604748) + 0.0301204819277108*m2*(0.0737767215575459*l**2 + 1.02005251887947*l + 3.54174457994057) + 0.0301204819277108*m2*(0.0737882045151557*l**2 - 0.103689411772581*l + 0.0513621602129999) + 0.0301204819277108*m2*(0.0754923473766136*l**2 - 0.665896041461111*l + 1.50452298238393) + 0.0301204819277108*m2*(0.0770620720534784*l**2 - 0.130536635340903*l + 0.108753222149204) + 0.0301204819277108*m2*(0.0773528426254158*l**2 - 0.561848486543818*l + 1.2940039641335) + 0.0301204819277108*m2*(0.0886583274138098*l**2 + 0.741425911861996*l + 1.93244889510292) + 0.0301204819277108*m2*(0.0942233608310071*l**2 - 0.138216814298681*l + 2.31919268332503) + 0.0301204819277108*m2*(0.100553144196389*l**2 - 0.63636938700589*l + 2.21085322903382) + 0.0301204819277108*m2*(0.101200733326992*l**2 + 0.635224279064149*l + 0.998035396135556) + 0.0301204819277108*m2*(0.119726739413426*l**2 + 1.43991076205258*l + 4.41820320093888) + 0.0301204819277108*m2*(0.124735936269166*l**2 - 0.272047351360634*l + 0.506378179076027) + 0.0301204819277108*m2*(0.124814015854552*l**2 + 0.0406734817175383*l + 0.0114877527966646) + 0.0301204819277108*m2*(0.125437457547957*l**2 + 0.295486196428763*l + 0.243286992192846) + 0.0301204819277108*m2*(0.125800352089473*l**2 - 0.620899832848617*l + 0.776493067439683) + 0.0301204819277108*m2*(0.12764298083469*l**2 + 0.14481071821561*l + 0.155008374681286) + 0.0301204819277108*m2*(0.128278888701751*l**2 - 1.29611838090318*l + 4.09835452995269) + 0.0301204819277108*m2*(0.129689693782827*l**2 + 0.505177455439501*l + 0.947212667494443) + 0.0301204819277108*m2*(0.131569587219619*l**2 - 0.241166363239798*l + 0.115433739338364) + 0.0301204819277108*m2*(0.132603154930277*l**2 - 0.021952114309859*l + 0.00287618631643805) + 0.0301204819277108*m2*(0.133812957271161*l**2 + 0.0614748211368591*l + 0.0150751262286752) + 0.0301204819277108*m2*(0.134534385268784*l**2 + 0.00797551619337831*l + 0.000239473328306061) + 0.0301204819277108*m2*(0.141965269445816*l**2 + 0.436258970752153*l + 0.382797061356946) + 0.0301204819277108*m2*(0.144954043619145*l**2 - 0.962245666133865*l + 3.47430835118538) + 0.0301204819277108*m2*(0.147693846979137*l**2 - 0.228110558767207*l + 0.378442731076156) + 0.0301204819277108*m2*(0.149995424910081*l**2 - 0.0824793667612512*l + 0.212172041895239) + 0.0301204819277108*m2*(0.151745615881445*l**2 + 0.00817210930667576*l + 1.02404515591469) + 0.0301204819277108*m2*(0.156262747594926*l**2 + 1.70592568394975*l + 5.19383478399391) + 0.0301204819277108*m2*(0.159021992089636*l**2 - 0.170249921805809*l + 0.208499579390378) + 0.0301204819277108*m2*(0.174610246859131*l**2 + 0.746268035928062*l + 2.39425054574447) + 0.0301204819277108*m2*(0.182507908485376*l**2 - 0.280195983548149*l + 0.138452472764186) + 0.0301204819277108*m2*(0.187498202246873*l**2 - 0.998596232134618*l + 1.70316595806976) + 0.0301204819277108*m2*(0.191256943190781*l**2 - 0.887052689870352*l + 1.03848246847259) + 0.0301204819277108*m2*(0.197815396940617*l**2 - 0.067789286281133*l + 0.00744365279870785) + 0.0301204819277108*m2*(0.198208165404656*l**2 - 0.285752117215271*l + 0.106046259050734) + 0.0301204819277108*m2*(0.203608211867754*l**2 - 0.177415015062798*l + 0.068450587595239) + 0.0301204819277108*m2*(0.20959567069045*l**2 - 0.194030168557906*l + 0.442124158803954) + 0.0301204819277108*m2*(0.226573175883361*l**2 - 0.244429642849245*l + 0.0778177044302685) + 0.0301204819277108*m2*(0.230750384158729*l**2 - 0.20983843466048*l + 0.0494129216563038) + 0.0301204819277108*m2*(0.245385849796199*l**2 + 0.0308102150681052*l + 0.520927207900351) + 0.0301204819277108*m2*(0.275481262192461*l**2 - 0.61226152011676*l + 0.422476204661143) + 0.0301204819277108*m2*(0.283208529848993*l**2 - 0.60728838634765*l + 0.592838175497362) + 0.0301204819277108*m2*(0.283384143804916*l**2 - 0.928355896391765*l + 3.57196353273703) + 0.0301204819277108*m2*(0.289598184914034*l**2 + 0.13554234179672*l + 1.36715949670619) + 0.0301204819277108*m2*(0.300531015458334*l**2 + 0.196023357979333*l + 0.0509318143777542) + 0.0301204819277108*m2*(0.325988130289768*l**2 - 0.0329236614194518*l + 0.0556207867692121) + 0.0301204819277108*m2*(0.329731350806745*l**2 - 0.501918190730618*l + 1.06982077589328) + 0.0301204819277108*m2*(0.343136174880681*l**2 + 0.0848752107916888*l + 1.434761710674) + 0.0301204819277108*m2*(0.348886890970498*l**2 - 0.00205644363841442*l + 0.000469648062573292) + 0.0301204819277108*m2*(0.356633369636093*l**2 + 0.159739541812927*l + 0.472228701693582) + 0.0301204819277108*m2*(0.360360218916456*l**2 - 0.0282015312579339*l + 0.68514682435757) + 0.0301204819277108*m2*(0.367549738349028*l**2 + 0.475874822118955*l + 0.954974750560148) + 0.0301204819277108*m2*(0.379245398008886*l**2 + 0.0064749883439895*l + 0.000141047680387917) + 0.0301204819277108*m2*(0.390079017287836*l**2 + 0.1820422087055*l + 0.0371362203599001) + 0.0301204819277108*m2*(0.394797866985169*l**2 - 0.244909725748765*l + 0.52674911617785) + 0.0301204819277108*m2*(0.394806135617756*l**2 - 0.523267655616229*l + 0.235160903624199) + 0.0301204819277108*m2*(0.398458058550227*l**2 - 0.632795644930882*l + 0.537230677512729) + 0.0301204819277108*m2*(0.41399949111788*l**2 - 0.729755852146709*l + 0.774279683218096) + 0.0301204819277108*m2*(0.415361481850441*l**2 - 0.0224061740164679*l + 0.0525697770207983) + 0.0301204819277108*m2*(0.436672541905971*l**2 + 0.172082837839294*l + 0.0697779027605327) + 0.0301204819277108*m2*(0.447138262692456*l**2 - 0.0181823905427556*l + 0.592640059743925) + 0.0301204819277108*m2*(0.448305528275037*l**2 + 0.238888851200015*l + 0.0531127977531318) + 0.0301204819277108*m2*(0.452642346796372*l**2 - 0.142916680342148*l + 0.0700032355036384) + 0.0301204819277108*m2*(0.45770071276561*l**2 - 0.496167261811478*l + 0.265119459018913) + 0.0301204819277108*m2*(0.457996678833234*l**2 - 0.312027501854529*l + 1.4865815029476) + 0.0301204819277108*m2*(0.458344527584803*l**2 - 0.588175807357479*l + 0.203624451407457) + 0.0301204819277108*m2*(0.474736379069997*l**2 - 0.122882464247322*l + 0.0610154017022102) + 0.0301204819277108*m2*(0.481785338657585*l**2 - 0.0345444243682733*l + 1.10201455974279) + 0.0301204819277108*m2*(0.481792745437333*l**2 + 0.354800945954665*l + 0.368644756890051) + 0.0301204819277108*m2*(0.484199150627016*l**2 - 0.579717963951604*l + 0.267577334666263) + 0.0301204819277108*m2*(0.484854988038492*l**2 + 0.1169197660718*l + 0.010086585171041) + 0.0301204819277108*m2*(0.491634340705604*l**2 + 0.0755502676208614*l + 0.645114222471071) + 0.0301204819277108*m2*(0.492822305128626*l**2 + 1.47386277691497*l + 3.58148106853646) + 0.0301204819277108*m2*(0.508909242581344*l**2 - 0.550041329324119*l + 0.149399587879751) + 0.0301204819277108*m2*(0.518913964522314*l**2 + 0.267918458556401*l + 0.868452449284215) + 0.0301204819277108*m2*(0.519411504260071*l**2 + 0.824127702047716*l + 0.8090532760367) + 0.0301204819277108*m2*(0.55057781202221*l**2 + 0.0108987371196025*l + 0.0113045918093708) + 0.0301204819277108*m2*(0.556338707548029*l**2 + 0.346361936882871*l + 0.0614493031607012) + 0.0301204819277108*m2*(0.561380275692216*l**2 + 0.568700301174827*l + 3.24923706032237) + 0.0301204819277108*m2*(0.561978117804303*l**2 + 0.561083850480972*l + 1.66052100404395) + 0.0301204819277108*m2*(0.573514960560652*l**2 + 0.867411862851288*l + 0.407736828124429) + 0.0301204819277108*m2*(0.588924985687355*l**2 + 0.902927057495917*l + 2.59852012436397) + 0.0301204819277108*m2*(0.598959775900144*l**2 - 1.1105736279876*l + 0.792343408416711) + 0.0301204819277108*m2*(0.603464741993695*l**2 + 0.163352354976369*l + 0.78691113680164) + 0.0301204819277108*m2*(0.624691583777247*l**2 + 0.488156474230515*l + 0.50155520912995) + 0.0301204819277108*m2*(0.635309976422854*l**2 + 1.0509862937621*l + 2.77510147397331) + 0.0301204819277108*m2*(0.63649268991924*l**2 - 1.01740124320764*l + 0.585720915905752) + 0.0301204819277108*m2*(0.636588128422039*l**2 + 0.0205183169204711*l + 0.437801096053157) + 0.0301204819277108*m2*(0.660551578199027*l**2 - 0.565498441390817*l + 0.662004046158495) + 0.0301204819277108*m2*(0.671909300261682*l**2 + 1.63316485167686*l + 1.55507058839825) + 0.0301204819277108*m2*(0.678676054467725*l**2 + 0.32170500471052*l + 0.102634579356275) + 0.0301204819277108*m2*(0.698854659746744*l**2 - 0.865364136867816*l + 0.270527919194651) + 0.0301204819277108*m2*(0.716979733732061*l**2 - 1.3186101759248*l + 0.996567507880505) + 0.0301204819277108*m2*(0.736379596003541*l**2 + 1.69750005695366*l + 1.46676536251534) + 0.0301204819277108*m2*(0.741071028481318*l**2 + 1.47205960569168*l + 0.734146210233477) + 0.0301204819277108*m2*(0.745918367018672*l**2 - 0.809032012122499*l + 0.636921978217837) + 0.0301204819277108*m2*(0.7510899675911*l**2 + 0.184795261894334*l + 0.137842712018193) + 0.0301204819277108*m2*(0.756934246538041*l**2 - 0.466420118157217*l + 0.0726876947738683) + 0.0301204819277108*m2*(0.773130615273342*l**2 - 1.92883718842278*l + 4.15295741459772) + 0.0301204819277108*m2*(0.789898575761953*l**2 + 0.943804479832147*l + 0.487372234736823) + 0.0301204819277108*m2*(0.809115848348826*l**2 + 0.52953344345489*l + 0.0998691606136252) + 0.0301204819277108*m2*(0.80922505030744*l**2 + 0.327356329029082*l + 0.143457716667336) + 0.0301204819277108*m2*(0.810837281201528*l**2 + 0.308741806226561*l + 0.211018710624879) + 0.0301204819277108*m2*(0.811004929038135*l**2 - 0.400114241057312*l + 0.263035496771458) + 0.0301204819277108*m2*(0.812582137603469*l**2 - 0.113396837998203*l + 0.40249958435924) + 0.0301204819277108*m2*(0.812660082204572*l**2 + 0.156649188258613*l + 0.0229814916984944) + 0.0301204819277108*m2*(0.825211878115664*l**2 + 0.00653548749045388*l + 2.38712106279792) + 0.0301204819277108*m2*(0.846222620286103*l**2 - 0.464519150872453*l + 0.725363110336508) + 0.0301204819277108*m2*(0.854208582833716*l**2 - 0.00795803443553801*l + 6.27478072793531e-5) + 0.0301204819277108*m2*(0.873703005685823*l**2 - 0.111488227249747*l + 0.038850131450676) + 0.0301204819277108*m2*(0.893644944545664*l**2 + 0.25319746907418*l + 0.0205609281485535) + 0.0301204819277108*m2*(0.895552745033009*l**2 + 2.02111610794695*l + 1.407413292872) + 0.0301204819277108*m2*(0.931200882237973*l**2 + 0.319069794462965*l + 1.65347092170158) + 0.0301204819277108*m2*(0.936712381503346*l**2 - 0.0645829825900696*l + 0.0938427265079706) + 0.0301204819277108*m2*(0.96824675873476*l**2 - 0.407673011625742*l + 0.16613395884864) + 0.0301204819277108*m2*(0.969508486504755*l**2 - 1.21094829281329*l + 0.380163337020004) + 0.0301204819277108*m2*(0.976973066107582*l**2 - 0.199857303314436*l + 0.114053232012741) + 0.0301204819277108*m2*(0.985503500129417*l**2 + 0.46345311134669*l + 1.00909103058072) + 0.0301204819277108*m2*(0.994637607625096*l**2 - 0.136733209670741*l + 0.00537688390983427) + 0.0301204819277108*m2*(1.02753678469026*l**2 - 2.01414424732987*l + 1.05531484024757) + 0.0301204819277108*m2*(1.0695005343037*l**2 - 1.34599499149329*l + 2.19501504844494) + 0.0301204819277108*m2*(1.07586197585779*l**2 + 1.78943674020001*l + 2.03800427120714) + 0.0301204819277108*m2*(1.09146784412901*l**2 + 1.40712500973707*l + 0.478539757856129) + 0.0301204819277108*m2*(1.10242571599159*l**2 + 0.086480271647605*l + 1.15802733517656) + 0.0301204819277108*m2*(1.10480830076656*l**2 - 2.43720148420066*l + 3.89801336146915) + 0.0301204819277108*m2*(1.13090890788169*l**2 + 0.142557013794384*l + 0.635568014220426) + 0.0301204819277108*m2*(1.13369387277927*l**2 - 0.17236873830086*l + 0.24936615302305) + 0.0301204819277108*m2*(1.19150199766931*l**2 + 0.195355715949399*l + 0.212421378471926) + 0.0301204819277108*m2*(1.20250526777958*l**2 - 2.10721044349849*l + 1.19404497876124) + 0.0301204819277108*m2*(1.21649599085607*l**2 + 0.00524979657676285*l + 0.0105628927811466) + 0.0301204819277108*m2*(1.2258053664495*l**2 + 2.22671097443634*l + 1.12754005785476) + 0.0301204819277108*m2*(1.28268973317703*l**2 - 0.60531334116776*l + 0.0714970709786448) + 0.0301204819277108*m2*(1.29337060818131*l**2 - 3.00182017798056*l + 5.14874742676025) + 0.0301204819277108*m2*(1.29730907274323*l**2 + 0.548742611907275*l + 0.226501971180922) + 0.0301204819277108*m2*(1.3089397913749*l**2 + 3.44588386830323*l + 2.56765432755135) + 0.0301204819277108*m2*(1.31408004638314*l**2 + 0.0207197216550482*l + 0.0242947488760451) + 0.0301204819277108*m2*(1.47082352167101*l**2 + 1.66675232290335*l + 0.534070822025748) + 0.0301204819277108*m2*(1.51006618053085*l**2 - 0.29805118150224*l + 0.139352414948119) + 0.0301204819277108*m2*(1.55389589060909*l**2 + 0.304684308326072*l + 0.0769958305424621) + 0.0301204819277108*m2*(1.65043467406627*l**2 + 1.85423418459618*l + 1.55068873592805) + 0.0301204819277108*m2*(1.65156154676568*l**2 + 1.02344045010824*l + 0.15855156204782) + 0.0301204819277108*m2*(1.71561576052715*l**2 - 2.25791680478442*l + 0.800315823218847) + 0.0301204819277108*m2*(1.71593542399204*l**2 + 0.531982499013133*l + 0.19506544163939) + 0.0301204819277108*m2*(1.72165025701379*l**2 - 0.0997502862994697*l + 0.0545627749919387) + 0.0301204819277108*m2*(1.75727355518775*l**2 - 0.973996222474864*l + 0.484093149961652) + 0.0301204819277108*m2*(1.76033541571913*l**2 - 1.30128870926804*l + 0.36212995841669) + 0.0301204819277108*m2*(1.85117954841089*l**2 - 6.28717045170271*l + 6.12906050374592) + 0.0301204819277108*m2*(1.90760768935142*l**2 + 1.03684610077284*l + 2.13541386100234) + 0.0301204819277108*m2*(1.95600165123866*l**2 - 0.0674229913150292*l + 0.00509827604304091) + 0.0301204819277108*m2*(1.96847862743878*l**2 - 1.54801022586405*l + 0.901391349802683) + 0.0301204819277108*m2*(2.0017197138234*l**2 + 1.2731443869249*l + 0.390452069318419) + 0.0301204819277108*m2*(2.03974359787578*l**2 + 3.02894606656761*l + 1.12507157421903) + 0.0301204819277108*m2*(2.05373897227361*l**2 + 0.90209249829607*l + 0.485090186567332) + 0.0301204819277108*m2*(2.06264904206331*l**2 - 0.17832338010376*l + 0.0220151648789403) + 0.0301204819277108*m2*(2.08365215901323*l**2 + 0.393384806177849*l + 0.0199217030990186) + 0.0301204819277108*m2*(2.12078182425756*l**2 + 3.37835791242028*l + 2.40675561054744) + 0.0301204819277108*m2*(2.15284006433945*l**2 - 0.219204643567509*l + 0.411391374752928) + 0.0301204819277108*m2*(2.16742421585616*l**2 - 2.55967494569968*l + 0.763208554961176) + 0.0301204819277108*m2*(2.26274754918832*l**2 + 0.764695448822515*l + 0.2586193462406) + 0.0301204819277108*m2*(2.30312619732623*l**2 - 4.04763285955175*l + 2.93268637703749) + 0.0301204819277108*m2*(2.37831051733834*l**2 - 1.17671872913411*l + 1.38717222517796) + 0.0301204819277108*m2*(2.51107660509336*l**2 - 0.152533416757132*l + 0.0230231059423946) + 0.0301204819277108*m2*(2.55462769462537*l**2 + 0.896365107770482*l + 0.0932929819988125) + 0.0301204819277108*m2*(2.65607997360035*l**2 - 1.15642583308946*l + 0.386081701259664) + 0.0301204819277108*m2*(2.96053535999947*l**2 + 0.208517414462969*l + 0.636703996861026) + 0.0301204819277108*m2*(2.98169407002421*l**2 - 3.15795593788294*l + 2.48482943853216) + 0.0301204819277108*m2*(3.07852368653619*l**2 + 1.83266332867313*l + 0.543042963306605) + 0.0301204819277108*m2*(3.12569730440801*l**2 + 2.71856001990821*l + 0.72427667079162) + 0.0301204819277108*m2*(3.21560850228717*l**2 + 1.28848333429209*l + 0.299254343267957) + 0.0301204819277108*m2*(3.21776375965121*l**2 + 0.815570068213069*l + 0.0832971488629137) + 0.0301204819277108*m2*(3.5355883648463*l**2 - 0.103650468477893*l + 0.0015715266411021) + 0.0301204819277108*m2*(3.55189301668637*l**2 - 3.1146886235898*l + 0.684201591311954) + 0.0301204819277108*m2*(3.60698942976626*l**2 - 0.237936321189115*l + 0.051086751273537) + 0.0301204819277108*m2*(3.62281823216439*l**2 - 0.0170816858633584*l + 4.30361381620563e-5) + 0.0301204819277108*m2*(3.78026232298605*l**2 - 0.69788199012343*l + 0.146357735004042) + 0.0301204819277108*m2*(4.14317719897579*l**2 - 2.49516378768941*l + 1.42152645458621) + 0.0301204819277108*m2*(4.29529443602412*l**2 - 1.82883529067181*l + 0.219492730287726) + 0.0301204819277108*m2*(4.37421223424354*l**2 + 3.75004941431772*l + 0.805502338307209) + 0.0301204819277108*m2*(4.42743210018134*l**2 - 0.543511815892691*l + 1.43466760352947) + 0.0301204819277108*m2*(4.68238712738811*l**2 - 0.0130038715994134*l + 0.26303254384268) + 0.0301204819277108*m2*(4.81303645913592*l**2 - 0.672259950773602*l + 0.476972842962155) + 0.0301204819277108*m2*(4.92234108059468*l**2 + 4.97311533678152*l + 1.29741381682889) + 0.0301204819277108*m2*(5.28208704255195*l**2 + 0.0291614169608416*l + 0.00061660149702531) + 0.0301204819277108*m2*(5.32716287350701*l**2 + 0.666035709258146*l + 0.957557277947415) + 0.0301204819277108*m2*(5.93594311193423*l**2 + 0.71281938816025*l + 0.190624993819451) + 0.0301204819277108*m2*(9.58865525428093*l**2 - 1.21063338613784*l + 0.13750935230048) + 0.0301204819277108*m2*(10.8829519541931*l**2 + 5.3317727316327*l + 1.72335209820097)

NameError Traceback (most recent call last) Cell In[5], line 10 8 initial_guess = np.random.randn(prob.num_free) 9 print('obj',obj(initial_guess)) ---> 10 print('obj_grad',obj_grad(initial_guess)) 11 # Find the optimal solution. 12 solution, info = prob.solve(initial_guess)

Cell In[2], line 94 91 states = free[:i_idx].reshape((n, N)) 92 inputs = free[i_idx:r_idx].reshape((q, N)) 93 return np.hstack(( ---> 94 *obj_grad_time_expr_eval(states, inputs, free[r_idx:]), 95 obj_grad_param_expr_eval(states, inputs, free[r_idx:]) 96 ))

File :9, in _lambdifygenerated(_Dummy_89, _Dummy_90, _Dummy_91) 7 x2 = _Dummy_86 8 x3 = 2x1 ----> 9 x4 = lsin(x2) 10 x5 = (1/2)m2 11 return (numpy.array((0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)), _IntegralFunction(-lm2x0x1cos(x2)), _IntegralFunction(m1x0 + x5(2x0 - x3x4)), _IntegralFunction(iZZx1 + x5*(l*2x3 - 2x0x4)), numpy.array((0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)))

NameError: name 'l' is not defined

moorepants commented 1 month ago

So it seems that the issue is that create_objective_function() does not support symbolic objectives that include any known parameters or known specified trajectories (as you surmised).

tjstienstra commented 1 month ago

Refer to the docstring. They should be substituted beforehand.

moorepants commented 1 month ago

The program runs fine for me, but on Linux.

This statement now becomes and oddity.

moorepants commented 1 month ago

Ah, I just copied and pasted the above code and it ran, but kinetic was set to False.

Peter230655 commented 1 month ago

Refer to the docstring. They should be substituted beforehand.

Fair enough! I should have read this!! :-( :-(

Peter230655 commented 1 month ago

Can I close the issue?

tjstienstra commented 1 month ago

Fair enough! I should have read this!! :-( :-(

😂 Familiar problem

Peter230655 commented 1 month ago

Fair enough! I should have read this!! :-( :-(

😂 Familiar problem

with me for sure! Just info: I did the substitution and all worked fine!