Skip to content

Commit

Permalink
Merge pull request #58 from Kev1CO/refactoring_ocp_creation
Browse files Browse the repository at this point in the history
Structure refactor
  • Loading branch information
Kev1CO authored Apr 23, 2024
2 parents b51d8bd + 5c7b763 commit 85fe804
Show file tree
Hide file tree
Showing 78 changed files with 3,116 additions and 2,959 deletions.
2 changes: 0 additions & 2 deletions .github/workflows/run_tests_win.yml
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,6 @@ jobs:

- name: Print mamba info
run: |
mamba config --show
mamba info
mamba list
Expand Down Expand Up @@ -90,7 +89,6 @@ jobs:

- name: Print mamba info
run: |
mamba config --show
mamba info
- name: Install bioptim on Windows
Expand Down
9 changes: 8 additions & 1 deletion README.md
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,14 @@ You can also compute the models form initial value problem.
For that, use the IvpFes class to build the computed problem.

```python
ocp = IvpFes(model=DingModelFrequency(), ...)
from cocofest import IvpFes, DingModelFrequencyWithFatigue

fes_parameters = {"model": DingModelFrequencyWithFatigue(), "n_stim": 10}
ivp_parameters = {"n_shooting": 20, "final_time": 1}

ivp = IvpFes(fes_parameters, ivp_parameters)

result, time = ivp.integrate()
```

## Summation truncation
Expand Down
8 changes: 4 additions & 4 deletions cocofest/custom_constraints.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
"""
This custom constraint are for the functional electrical stimulation frequency and intensity.
This class regroups all the custom constraints that are used in the optimization problem.
"""

from casadi import MX, SX
Expand All @@ -13,7 +13,7 @@ def pulse_time_apparition_as_phase(controller: PenaltyController) -> MX | SX:
return controller.time.cx - controller.parameters["pulse_apparition_time"].cx[controller.phase_idx]

@staticmethod
def pulse_time_apparition_bimapping(controller: PenaltyController) -> MX | SX:
def equal_to_first_pulse_interval_time(controller: PenaltyController) -> MX | SX:
if controller.ocp.n_phases <= 1:
RuntimeError("There is only one phase, the bimapping constraint is not possible")

Expand All @@ -25,7 +25,7 @@ def pulse_time_apparition_bimapping(controller: PenaltyController) -> MX | SX:
return first_phase_tf - current_phase_tf

@staticmethod
def pulse_duration_bimapping(controller: PenaltyController) -> MX | SX:
def equal_to_first_pulse_duration(controller: PenaltyController) -> MX | SX:
if controller.ocp.n_phases <= 1:
RuntimeError("There is only one phase, the bimapping constraint is not possible")
return (
Expand All @@ -34,7 +34,7 @@ def pulse_duration_bimapping(controller: PenaltyController) -> MX | SX:
)

@staticmethod
def pulse_intensity_bimapping(controller: PenaltyController) -> MX | SX:
def equal_to_first_pulse_intensity(controller: PenaltyController) -> MX | SX:
if controller.ocp.n_phases <= 1:
RuntimeError("There is only one phase, the bimapping constraint is not possible")
return (
Expand Down
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
Original file line number Diff line number Diff line change
Expand Up @@ -27,10 +27,8 @@
n_stim=n_stim,
n_shooting=10,
final_time=1,
time_min=0.01,
time_max=0.1,
time_bimapping=True,
custom_objective=objective_functions,
pulse_event={"min": 0.01, "max": 0.1, "bimapping": True},
objective={"custom": objective_functions},
with_residual_torque=True,
)

Expand Down
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
"""
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Bakir's 2022 intensity work.
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Hmed's 2018 intensity work.
This ocp was build to maintain an elbow angle of 90 degrees.
The stimulation frequency will be optimized between 1 and 10 Hz as well as the pulse intensity between minimal
sensitivity threshold and 130mA to satisfy the maintained elbow. No residual torque is allowed.
Expand Down Expand Up @@ -48,15 +48,14 @@
n_stim=n_stim,
n_shooting=5,
final_time=1,
time_min=0.05,
time_max=1,
time_bimapping=True,
pulse_intensity_min=minimum_pulse_intensity,
pulse_intensity_max=130,
pulse_intensity_bimapping=False,
pulse_event={"min": 0.05, "max": 1, "bimapping": True},
pulse_intensity={
"min": minimum_pulse_intensity,
"max": 130,
"bimapping": False,
},
with_residual_torque=True,
custom_objective=objective_functions,
q_tracking=track_q,
objective={"custom": objective_functions, "q_tracking": track_q},
use_sx=True,
)

Expand Down
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
"""
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Bakir's 2022 intensity work.
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Hmed's 2018 intensity work.
This ocp was build to maintain an elbow angle of 90 degrees.
The stimulation frequency will be optimized between 1 and 10 Hz as well as the pulse intensity between minimal
sensitivity threshold and 130mA to satisfy the maintained elbow. No residual torque is allowed.
Expand Down Expand Up @@ -55,14 +55,14 @@
n_stim=n_stim,
n_shooting=10,
final_time=1,
time_min=0.1,
time_max=1,
time_bimapping=True,
pulse_intensity_min=minimum_pulse_intensity,
pulse_intensity_max=130,
pulse_intensity_bimapping=False,
pulse_event={"min": 0.1, "max": 1, "bimapping": True},
pulse_intensity={
"min": minimum_pulse_intensity,
"max": 130,
"bimapping": False,
},
with_residual_torque=False,
custom_objective=objective_functions,
objective={"custom": objective_functions},
)

sol = ocp.solve(Solver.IPOPT(_max_iter=1000))
Expand Down
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
"""
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Bakir's 2022 intensity work.
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Hmed's 2018 intensity work.
This ocp was build to maintain an elbow angle of 90 degrees.
The stimulation frequency will be optimized between 1 and 10 Hz as well as the pulse intensity between minimal
sensitivity threshold and 130mA to satisfy the maintained elbow. No residual torque is allowed.
Expand Down Expand Up @@ -61,14 +61,14 @@
n_stim=n_stim,
n_shooting=10,
final_time=1,
time_min=0.1,
time_max=1,
time_bimapping=True,
pulse_intensity_min=minimum_pulse_intensity,
pulse_intensity_max=130,
pulse_intensity_bimapping=False,
pulse_event={"min": 0.1, "max": 1, "bimapping": True},
pulse_intensity={
"min": minimum_pulse_intensity,
"max": 130,
"bimapping": False,
},
with_residual_torque=False,
custom_objective=objective_functions,
objective={"custom": objective_functions},
)

sol = ocp.solve(Solver.IPOPT(_max_iter=1000))
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -47,13 +47,11 @@
n_stim=n_stim,
n_shooting=n_shooting,
final_time=1,
time_min=0.01,
time_max=1,
time_bimapping=False,
pulse_event={"min": 0.01, "max": 1, "bimapping": False},
with_residual_torque=True,
custom_objective=objective_functions,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=False,
objective={"custom": objective_functions},
activate_force_length_relationship=True,
activate_force_velocity_relationship=False,
minimize_muscle_fatigue=True,
)

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -43,13 +43,15 @@
n_stim=n_stim,
n_shooting=10,
final_time=1,
pulse_duration_min=minimum_pulse_duration,
pulse_duration_max=0.0006,
pulse_duration_bimapping=False,
pulse_duration={
"min": minimum_pulse_duration,
"max": 0.0006,
"bimapping": False,
},
with_residual_torque=False,
custom_objective=objective_functions,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=False,
objective={"custom": objective_functions},
activate_force_length_relationship=True,
activate_force_velocity_relationship=False,
minimize_muscle_fatigue=True,
)

Expand Down
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
"""
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Bakir's 2022 intensity work.
This example will do a 10 stimulation example with Ding's 2003 frequency model associated to Hmed's 2018 intensity work.
Those ocp were build to move the elbow from 0 to 90 degrees angle.
The stimulation frequency will be set to 10Hz and pulse intensity will be optimized to satisfy the motion and to minimize the overall muscle fatigue.
Intensity can be optimized from sensitivity threshold to 130mA. No residual torque is allowed.
Expand Down Expand Up @@ -45,13 +45,15 @@
n_stim=n_stim,
n_shooting=10,
final_time=1,
pulse_intensity_min=minimum_pulse_intensity,
pulse_intensity_max=130,
pulse_intensity_bimapping=False,
pulse_intensity={
"min": minimum_pulse_intensity,
"max": 130,
"bimapping": False,
},
with_residual_torque=False,
custom_objective=objective_functions,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=False,
objective={"custom": objective_functions},
activate_force_length_relationship=True,
activate_force_velocity_relationship=False,
minimize_muscle_fatigue=True,
)

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@

sol_list = []
sol_time = []
muscle_force_length_relationship = [False, True]
activate_force_length_relationship = [False, True]

for i in range(2):
ocp = OcpFesMsk.prepare_ocp(
Expand All @@ -28,10 +28,10 @@
n_stim=10,
n_shooting=10,
final_time=1,
pulse_duration=0.00025,
pulse_duration={"fixed": 0.00025},
with_residual_torque=False,
muscle_force_length_relationship=muscle_force_length_relationship[i],
muscle_force_velocity_relationship=muscle_force_length_relationship[i],
activate_force_length_relationship=activate_force_length_relationship[i],
activate_force_velocity_relationship=activate_force_length_relationship[i],
use_sx=False,
)
sol = ocp.solve(Solver.IPOPT(_max_iter=1000))
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@
ConstraintList,
Node,
Solver,
SolutionMerge,
)

from cocofest import DingModelFrequencyWithFatigue, OcpFesMsk
Expand Down Expand Up @@ -88,23 +89,22 @@
n_stim=n_stim,
n_shooting=n_shooting,
final_time=1,
time_min=0.01,
time_max=0.1,
time_bimapping=False,
pulse_event={"min": 0.01, "max": 0.1, "bimapping": False},
with_residual_torque=False,
custom_constraint=constraint,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=True,
activate_force_length_relationship=True,
activate_force_velocity_relationship=True,
minimize_muscle_fatigue=True if pickle_file_list[i] == "minimize_muscle_fatigue.pkl" else False,
minimize_muscle_force=True if pickle_file_list[i] == "minimize_muscle_force.pkl" else False,
use_sx=False,
)

sol = ocp.solve(Solver.IPOPT(_max_iter=10000)).merge_phases()
time = sol.time
states = sol.states
controls = sol.controls
parameters = sol.parameters
sol = ocp.solve(Solver.IPOPT(_max_iter=10000))

time = sol.decision_time(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
states = sol.decision_states(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
controls = sol.decision_controls(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
parameters = sol.decision_parameters()

dictionary = {
"time": time,
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@
ConstraintList,
Node,
Solver,
SolutionMerge,
)

from cocofest import DingModelIntensityFrequencyWithFatigue, OcpFesMsk
Expand Down Expand Up @@ -91,23 +92,26 @@
n_stim=n_stim,
n_shooting=n_shooting,
final_time=1,
pulse_intensity_min=minimum_pulse_intensity,
pulse_intensity_max=80,
pulse_intensity_bimapping=False,
pulse_intensity={
"min": minimum_pulse_intensity,
"max": 80,
"bimapping": False,
},
with_residual_torque=False,
custom_constraint=constraint,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=True,
activate_force_length_relationship=True,
activate_force_velocity_relationship=True,
minimize_muscle_fatigue=True if pickle_file_list[i] == "minimize_muscle_fatigue.pkl" else False,
minimize_muscle_force=True if pickle_file_list[i] == "minimize_muscle_force.pkl" else False,
use_sx=False,
)

sol = ocp.solve(Solver.IPOPT(_max_iter=10000)).merge_phases()
time = sol.time
states = sol.states
controls = sol.controls
parameters = sol.parameters
sol = ocp.solve(Solver.IPOPT(_max_iter=10000))

time = sol.decision_time(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
states = sol.decision_states(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
controls = sol.decision_controls(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
parameters = sol.decision_parameters()

dictionary = {
"time": time,
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@
ConstraintList,
Node,
Solver,
SolutionMerge,
)

from cocofest import DingModelPulseDurationFrequencyWithFatigue, OcpFesMsk
Expand Down Expand Up @@ -76,7 +77,6 @@
)

for i in range(len(pickle_file_list)):

ocp = OcpFesMsk.prepare_ocp(
biorbd_model_path="../../msk_models/arm26.bioMod",
bound_type="start_end",
Expand All @@ -85,23 +85,26 @@
n_stim=n_stim,
n_shooting=n_shooting,
final_time=1.5,
pulse_duration_min=minimum_pulse_duration,
pulse_duration_max=0.0006,
pulse_duration_bimapping=False,
pulse_duration={
"min": minimum_pulse_duration,
"max": 0.0006,
"bimapping": False,
},
with_residual_torque=False,
custom_constraint=constraint,
muscle_force_length_relationship=True,
muscle_force_velocity_relationship=True,
activate_force_length_relationship=True,
activate_force_velocity_relationship=True,
minimize_muscle_fatigue=True if pickle_file_list[i] == "minimize_muscle_fatigue.pkl" else False,
minimize_muscle_force=True if pickle_file_list[i] == "minimize_muscle_force.pkl" else False,
use_sx=False,
)

sol = ocp.solve(Solver.IPOPT(_max_iter=10000)).merge_phases()
time = sol.time
states = sol.states
controls = sol.controls
parameters = sol.parameters
sol = ocp.solve(Solver.IPOPT(_max_iter=10000))

time = sol.decision_time(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
states = sol.decision_states(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
controls = sol.decision_controls(to_merge=[SolutionMerge.PHASES, SolutionMerge.NODES])
parameters = sol.decision_parameters()

dictionary = {
"time": time,
Expand Down
Loading

0 comments on commit 85fe804

Please sign in to comment.