Navigation

The optimization study introduced in the preceding lesson used AnyBody’s builtin
facilities for optimizing. Sometimes that is not enough, either because the
objective functions depends on data that is not easily included in AnyBody, or
because a different solver is needed.

In this tutorial we use an external optimizer together with AnyBody. The example
is based on the model from the previous lesson but uses an
optimizer from the Scipy python library. The same
could also have been achived with other optimization frameworks (like NLopt, or languages (like MatLab).

We also need one additional Python library (AnyPyTools) which will make it
easier to work with AnyBody from Python. AnyPyTools can be easily
installed from the command prompt. Open the Anaconda command prompt and type:

importmathimportscipyfromanypytoolsimportAnyPyProcessfromanypytools.macro_commandsimportLoad,OperationRun,Dump,SetValuedefrun_model(saddle_height,saddle_pos,silent=False):"""Run the AnyBody model and return the metabolism results"""macro=[Load("BikeModel2D.main.any"),SetValue("Main.BikeParameters.SaddleHeight",saddle_height),SetValue("Main.BikeParameters.SaddlePos",saddle_pos),OperationRun("Main.Study.InverseDynamics"),Dump("Main.Study.Output.Pmet"),Dump("Main.Study.Output.Abscissa.t"),]app=AnyPyProcess(silent=silent)results=app.start_macro(macro)returnresults[0]defobjfun(designvars):"""Calculate the objective function value"""saddle_height=designvars[0]saddle_pos=designvars[1]result=run_model(saddle_height,saddle_pos,silent=True)if"ERROR"inresult:raiseValueError("Failed to run model")pmet=scipy.integrate.trapz(result["Pmet"],result["Abscissa.t"])returnfloat(pmet)defseat_distance_constraint(designvars):"""Compute contraint value which must be larger than zero"""returnmath.sqrt(designvars[0]**2+designvars[1]**2)-0.66constaints={"type":"ineq","fun":seat_distance_constraint}bounds=[(0.61,0.69),(-0.22,-0.05)]initial_guess=(0.68,-0.15)solution=scipy.optimize.minimize(objfun,initial_guess,constraints=constaints,bounds=bounds,method="SLSQP")print(solution)

A copy of the file can be downloaded here.
For now you can place the optimize.py next to your main file BikeModel2D.main.any.
If you didn’t complete the model from lesson 2, you can download the
finshed model here.

For the external optimizers to work, we need a way to run AnyBody models from
Python and record the results of the simulations, so we need to create a
function to do this. There are more information on how to do this in the
documentaion for AnyPyTools. So here we will
just show how the code looks and not discuss the details.

8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23

defrun_model(saddle_height,saddle_pos,silent=False):"""Run the AnyBody model and return the metabolism results"""macro=[Load("BikeModel2D.main.any"),SetValue("Main.BikeParameters.SaddleHeight",saddle_height),SetValue("Main.BikeParameters.SaddlePos",saddle_pos),OperationRun("Main.Study.InverseDynamics"),Dump("Main.Study.Output.Pmet"),Dump("Main.Study.Output.Abscissa.t"),]app=AnyPyProcess(silent=silent)results=app.start_macro(macro)returnresults[0]result=run_model(0.66,-0.16)print(result["Main.Study.Output.Pmet"])

The function run_model() takes saddle_height and saddle_pos as input
and return the Pmet metabolism output from AnyBody.

If you use an interactive Python environment (like IPython) you could try calling the function directly to to
test it and investigate the results:

The next step is to define the objective function. The objective function should
take a list of design values as input and return the objective function value.
In Lesson 2 the objective function was the time integral of the
metabolism variable. So we will do the same here with Scipy’s
scipy.integrate.trapz(): function.

23
24
25
26
27
28
29
30
31
32
33

defobjfun(x):saddle_height=x[0]saddle_pos=x[1]result=run_model(saddle_height,saddle_pos,silent=True)if"ERROR"inresult:raiseValueError("Failed to run model")# Integrate Pmetpmet=scipy.integrate.trapz(result["Pmet"],result["Abscissa.t"])returnfloat(pmet)

Note

The function also checks the results for errors reported
by AnyBody and raises a ValueError exception if that happens.
There could be ways of handle error without failing but it is important to
handle model failures, otherwise they may go unnoticed or mess with the
optimization results.

Again, we can run this function interactively to test it:

In [9]: pmet=objfun([0.66,-0.16])In [10]: print(pmet)505.329399532772

Now we get the time integral of the Pmet variable as as single value,
and we are now ready to define the optimization process.

defseat_distance_constraint(x):""" Compute contraint value which must be larger than zero"""return(math.sqrt(x[0]**2+x[1]**2)-0.66)constaints={"type":"ineq","fun":seat_distance_constraint}bounds=[(0.61,0.69),(-0.22,-0.05)]initial_guess=(0.68,-0.15)

And there we have it!
We can now take advantage of the many different algorithms and settings available for scipy.optimize.minimize().
We could also use a different package or customize our own algorithm, constraints etc.
The possibilities are practically endless.

For more information regarding the AnyPyTools python package follow this link.
You can also check out this webcast.