diff --git a/README.md b/README.md
index 45911f7e5..8c7e85214 100644
--- a/README.md
+++ b/README.md
@@ -44,6 +44,7 @@ As a tour guide that uses this binder, you can watch the `bioptim` workshop that
 
 A GUI is available to run all the current examples. To run it you can use the following command, from the root folder of the project:
 ```bash
+conda install -c conda-forge pyqt pyqtgraph
 python -m bioptim.examples
 ```
 Please refer to section [Examples](#examples) for more information on how to run the examples.
diff --git a/bioptim/examples/__main__.py b/bioptim/examples/__main__.py
index 995472e3b..824d592eb 100644
--- a/bioptim/examples/__main__.py
+++ b/bioptim/examples/__main__.py
@@ -23,9 +23,12 @@
 # Avoid clash with module name
 examples_ = OrderedDict(
     [
-        ("acados", OrderedDict([("Static arm", "static_arm.py"), ("Cube", "cube.py"), ("Pendulum", "pendulum.py")])),
         (
-            "fatigue",
+            "toy_examples/acados",
+            OrderedDict([("Static arm", "static_arm.py"), ("Cube", "cube.py"), ("Pendulum", "pendulum.py")]),
+        ),
+        (
+            "toy_examples/fatigue",
             OrderedDict(
                 [
                     ("Pendulum with fatigue", "pendulum_with_fatigue.py"),
@@ -37,32 +40,39 @@
             "getting_started",
             OrderedDict(
                 [
-                    ("Custom Bounds", "custom_bounds.py"),
-                    ("Custom constraint", "custom_constraint.py"),
-                    ("Custom initial guess", "custom_initial_guess.py"),
                     ("Custom objectives", "custom_objectives.py"),
                     ("Custom parameters", "custom_parameters.py"),
                     ("Custom phase transitions", "custom_phase_transitions.py"),
                     ("Custom plotting", "custom_plotting.py"),
-                    ("Example cyclic movement", "example_cyclic_movement.py"),
                     ("Example external forces", "example_external_forces.py"),
                     ("Example inequality constraint", "example_inequality_constraint.py"),
-                    ("Example mapping", "example_mapping.py"),
                     ("Example multiphase", "example_multiphase.py"),
-                    ("Example optimal time", "example_optimal_time.py"),
                     ("Example simulation", "example_simulation.py"),
-                    ("Pendulum", "pendulum.py"),
+                    ("Example cyclic movement", "example_cyclic_movement.py"),
+                    ("Basic OCP", "basic_ocp.py"),
+                    ("Example constraint weight", "custom_constraint_weights.py"),
+                    ("How to plot", "how_to_plot.py"),
+                ]
+            ),
+        ),
+        (
+            "toy_examples/feature_examples",
+            OrderedDict(
+                [
+                    ("Custom Bounds", "custom_bounds.py"),
+                    ("Custom constraint", "custom_constraint.py"),
+                    ("Custom initial guess", "custom_initial_guess.py"),
                 ]
             ),
         ),
         (
-            "moving_horizon_estimation",
+            "toy_examples/moving_horizon_estimation",
             OrderedDict(
                 [("Cyclic nmpc", "cyclic_nmpc.py"), ("Mhe", "mhe.py"), ("Multi cyclic nmpc", "multi_cyclic_nmpc.py")]
             ),
         ),
         (
-            "muscle_driven_ocp",
+            "toy_examples/muscle_driven_ocp",
             OrderedDict(
                 [
                     ("Muscle activations tracker", "muscle_activations_tracker.py"),
@@ -73,7 +83,7 @@
             ),
         ),
         (
-            "muscle_driven_with_contact",
+            "toy_examples/muscle_driven_with_contact",
             OrderedDict(
                 [
                     ("Contact forces inequality constraint muscle", "contact_forces_inequality_constraint_muscle.py"),
@@ -86,7 +96,7 @@
             ),
         ),
         (
-            "optimal_time_ocp",
+            "toy_examples/optimal_time_ocp",
             OrderedDict(
                 [
                     ("Multiphase time constraint", "multiphase_time_constraint.py"),
@@ -96,7 +106,7 @@
             ),
         ),
         (
-            "symmetrical_torque_driven_ocp",
+            "toy_examples/symmetrical_torque_driven_ocp",
             OrderedDict(
                 [
                     ("Symmetry by constraint", "symmetry_by_constraint.py"),
@@ -105,7 +115,7 @@
             ),
         ),
         (
-            "torque_driven_ocp",
+            "toy_examples/torque_driven_ocp",
             OrderedDict(
                 [
                     ("Maximize predicted height center_of_mass", "maximize_predicted_height_CoM.py"),
@@ -126,7 +136,7 @@
             ),
         ),
         (
-            "track",
+            "toy_examples/tracking",
             OrderedDict(
                 [
                     ("Track marker on segment", "track_marker_on_segment.py"),
diff --git a/bioptim/examples/getting_started/example_multistart.py b/bioptim/examples/getting_started/example_multistart.py
index 7de6125f7..bd21d98e7 100644
--- a/bioptim/examples/getting_started/example_multistart.py
+++ b/bioptim/examples/getting_started/example_multistart.py
@@ -23,6 +23,7 @@
     PhaseDynamics,
     SolutionMerge,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
@@ -196,7 +197,8 @@ def prepare_multi_start(
 def main():
     # --- Prepare the multi-start and run it --- #
 
-    bio_model_path = ["models/pendulum.bioMod"]
+    example_path = ExampleUtils.folder
+    bio_model_path = [example_path + "/models/pendulum.bioMod"]
     final_time = [1]
     n_shooting = [30, 40, 50]
     seed = [0, 1, 2, 3]
@@ -208,7 +210,7 @@ def main():
         "seed": seed,
     }
 
-    save_folder = "./temporary_results"
+    save_folder = example_path + "/./temporary_results"
     multi_start = prepare_multi_start(
         combinatorial_parameters=combinatorial_parameters,
         save_folder=save_folder,
diff --git a/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_pendulum.py b/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_pendulum.py
index 2564a215c..4b4f882f7 100644
--- a/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_pendulum.py
+++ b/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_pendulum.py
@@ -14,6 +14,7 @@
     VariationalTorqueBiorbdModel,
     VariationalOptimalControlProgram,
 )
+from bioptim.examples.utils import ExampleUtils
 import numpy as np
 
 
@@ -92,7 +93,8 @@ def main():
     n_shooting = 100
 
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp(bio_model_path="models/pendulum.bioMod", final_time=1, n_shooting=n_shooting)
+    biorbd_model_path = ExampleUtils.folder + "/models/pendulum.bioMod"
+    ocp = prepare_ocp(bio_model_path=biorbd_model_path, final_time=1, n_shooting=n_shooting)
 
     # --- Print ocp structure --- #
     ocp.print(to_console=False, to_graph=False)
diff --git a/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_with_holonomic_constraints_pendulum.py b/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_with_holonomic_constraints_pendulum.py
index 4a7b83b1e..6d0728151 100644
--- a/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_with_holonomic_constraints_pendulum.py
+++ b/bioptim/examples/toy_examples/discrete_mechanics_and_optimal_control/example_variational_integrator_with_holonomic_constraints_pendulum.py
@@ -17,6 +17,7 @@
     VariationalTorqueBiorbdModel,
     VariationalOptimalControlProgram,
 )
+from bioptim.examples.utils import ExampleUtils
 import numpy as np
 
 
@@ -111,7 +112,8 @@ def main():
     n_shooting = 100
 
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp(bio_model_path="models/pendulum_holonomic.bioMod", final_time=1, n_shooting=n_shooting)
+    biorbd_model_path = ExampleUtils.folder + "/models/pendulum_holonomic.bioMod"
+    ocp = prepare_ocp(bio_model_path=biorbd_model_path, final_time=1, n_shooting=n_shooting)
 
     # --- Print ocp structure --- #
     ocp.print(to_console=False, to_graph=False)
diff --git a/bioptim/examples/toy_examples/feature_examples/custom_constraint_weights.py b/bioptim/examples/toy_examples/feature_examples/custom_constraint_weights.py
index 3917c093e..ac0b60a42 100644
--- a/bioptim/examples/toy_examples/feature_examples/custom_constraint_weights.py
+++ b/bioptim/examples/toy_examples/feature_examples/custom_constraint_weights.py
@@ -32,6 +32,7 @@
     PhaseDynamics,
     ConstraintWeight,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def custom_weight(node: int, n_nodes: int) -> float:
@@ -172,6 +173,7 @@ def main():
     """
 
     nodes_to_test = [Node.START, Node.INTERMEDIATES, Node.ALL_SHOOTING]
+    biorbd_model_path = ExampleUtils.folder + "/models/cube.bioMod"
 
     for interpolation_type in InterpolationType:
         for node in nodes_to_test:
@@ -183,7 +185,11 @@ def main():
 
             print(f"Solving problem using {interpolation_type} weight applied at {node} nodes.")
             ocp = prepare_ocp(
-                "models/cube.bioMod", n_shooting=30, final_time=2, interpolation_type=interpolation_type, node=node
+                biorbd_model_path=biorbd_model_path,
+                n_shooting=30,
+                final_time=2,
+                interpolation_type=interpolation_type,
+                node=node,
             )
             sol = ocp.solve()
             print("\n")
diff --git a/bioptim/examples/toy_examples/feature_examples/example_multinode_objective.py b/bioptim/examples/toy_examples/feature_examples/example_multinode_objective.py
index 13ff7e7ad..4def3fab3 100644
--- a/bioptim/examples/toy_examples/feature_examples/example_multinode_objective.py
+++ b/bioptim/examples/toy_examples/feature_examples/example_multinode_objective.py
@@ -19,6 +19,7 @@
     MultinodeObjectiveList,
     CostType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def multinode_min_controls(controllers: list[PenaltyController]):
@@ -125,7 +126,8 @@ def main():
 
     # --- Prepare the ocp --- #
     n_shooting = 30
-    ocp = prepare_ocp(biorbd_model_path="models/pendulum.bioMod", final_time=1, n_shooting=n_shooting)
+    biorbd_model_path = ExampleUtils.folder + "/models/pendulum.bioMod"
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path, final_time=1, n_shooting=n_shooting)
     ocp.add_plot_penalty(CostType.ALL)
 
     # --- Solve the ocp --- #
diff --git a/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums.py b/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums.py
index e62689941..daad7e3d6 100644
--- a/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums.py
+++ b/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums.py
@@ -25,6 +25,7 @@
     SolutionMerge,
     OdeSolver,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def compute_all_states(sol, bio_model: HolonomicTorqueBiorbdModel):
@@ -197,8 +198,8 @@ def main():
     Runs the optimization and animates it
     """
 
-    model_path = "models/two_pendulums.bioMod"
-    ocp, bio_model = prepare_ocp(biorbd_model_path=model_path)
+    biorbd_model_path = ExampleUtils.folder + "/models/two_pendulums.bioMod"
+    ocp, bio_model = prepare_ocp(biorbd_model_path=biorbd_model_path)
 
     # --- Solve the program --- #
     sol = ocp.solve(Solver.IPOPT())
@@ -211,7 +212,7 @@ def main():
     if viewer == "bioviz":
         import bioviz
 
-        viz = bioviz.Viz(model_path)
+        viz = bioviz.Viz(biorbd_model_path)
         viz.load_movement(q)
         viz.exec()
 
@@ -219,7 +220,7 @@ def main():
         import pyorerun
 
         viz = pyorerun.PhaseRerun(t_span=np.concatenate(sol.decision_time()).squeeze())
-        viz.add_animated_model(pyorerun.BiorbdModel(model_path), q=q)
+        viz.add_animated_model(pyorerun.BiorbdModel(biorbd_model_path), q=q)
 
         viz.rerun("double_pendulum")
 
diff --git a/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums_algebraic.py b/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums_algebraic.py
index cb80d5c0e..28b3ff856 100644
--- a/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums_algebraic.py
+++ b/bioptim/examples/toy_examples/holonomic_constraints/two_pendulums_algebraic.py
@@ -28,6 +28,7 @@
     CostType,
     OdeSolver,
 )
+from bioptim.examples.utils import ExampleUtils
 from .custom_dynamics import (
     ModifiedHolonomicTorqueBiorbdModel,
     constraint_holonomic,
@@ -232,8 +233,8 @@ def main():
     Runs the optimization and animates it
     """
 
-    model_path = "models/two_pendulums.bioMod"
-    ocp, bio_model = prepare_ocp(biorbd_model_path=model_path)
+    biorbd_model_path = ExampleUtils.folder + "/models/two_pendulums.bioMod"
+    ocp, bio_model = prepare_ocp(biorbd_model_path=biorbd_model_path)
     ocp.add_plot_penalty(CostType.ALL)
 
     # --- Solve the program --- #
@@ -248,14 +249,14 @@ def main():
     if viewer == "bioviz":
         import bioviz
 
-        viz = bioviz.Viz(model_path)
+        viz = bioviz.Viz(biorbd_model_path)
         viz.load_movement(q)
         viz.exec()
 
     if viewer == "pyorerun":
         from pyorerun import BiorbdModel as PyorerunBiorbdModel, PhaseRerun
 
-        pyomodel = PyorerunBiorbdModel(model_path)
+        pyomodel = PyorerunBiorbdModel(biorbd_model_path)
         viz = PhaseRerun(t_span=np.concatenate(sol.decision_time()).squeeze())
         viz.add_animated_model(pyomodel, q=q)
 
diff --git a/bioptim/examples/toy_examples/inverse_optimal_control/double_pendulum_torque_driven_IOCP.py b/bioptim/examples/toy_examples/inverse_optimal_control/double_pendulum_torque_driven_IOCP.py
index 05066d725..c3a7c90aa 100755
--- a/bioptim/examples/toy_examples/inverse_optimal_control/double_pendulum_torque_driven_IOCP.py
+++ b/bioptim/examples/toy_examples/inverse_optimal_control/double_pendulum_torque_driven_IOCP.py
@@ -23,12 +23,13 @@
     PhaseDynamics,
     SolutionMerge,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
     weights,
     coefficients,
-    biorbd_model_path="models/double_pendulum.bioMod",
+    biorbd_model_path,
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     n_threads: int = 4,
     expand_dynamics: bool = True,
@@ -119,7 +120,8 @@ def fitness(self, weights):
         """
         global i_inverse
         i_inverse += 1
-        ocp = prepare_ocp(weights, self.coefficients)
+        biorbd_model_path = ExampleUtils.folder + "/models/double_pendulums.bioMod"
+        ocp = prepare_ocp(weights, self.coefficients, biorbd_model_path=biorbd_model_path)
         sol = ocp.solve(self.solver)
         print(
             f"+++++++++++++++++++++++++++ Optimized the {i_inverse}th ocp in the inverse algo +++++++++++++++++++++++++++"
@@ -148,7 +150,8 @@ def main():
 
     # Generate data using OCP
     weights_to_track = [0.4, 0.3, 0.3]
-    ocp_to_track = prepare_ocp(weights=weights_to_track, coefficients=[1, 1, 1])
+    biorbd_model_path = ExampleUtils.folder + "/models/double_pendulums.bioMod"
+    ocp_to_track = prepare_ocp(weights=weights_to_track, coefficients=[1, 1, 1], biorbd_model_path=biorbd_model_path)
     ocp_to_track.add_plot_penalty(CostType.ALL)
     solver = Solver.IPOPT()
     # solver.set_linear_solver("ma57")  # Much faster, but necessite libhsl installed
@@ -163,7 +166,7 @@ def main():
     for i in range(len(weights_to_track)):
         weights_pareto = [0, 0, 0]
         weights_pareto[i] = 1
-        ocp_pareto = prepare_ocp(weights=weights_pareto, coefficients=[1, 1, 1])
+        ocp_pareto = prepare_ocp(weights=weights_pareto, coefficients=[1, 1, 1], biorbd_model_path=biorbd_model_path)
         sol_pareto = ocp_pareto.solve(solver)
         # sol_pareto.animate()
         coefficients.append(sol_pareto.cost)
@@ -215,7 +218,7 @@ def main():
     # Compare the kinematics
     import biorbd
 
-    ocp_final = prepare_ocp(weights=pop_weights, coefficients=coefficients)
+    ocp_final = prepare_ocp(weights=pop_weights, coefficients=coefficients, biorbd_model_path=biorbd_model_path)
     sol_final = ocp_final.solve(solver)
     states = sol_final.decision_states(to_merge=SolutionMerge.NODES)
     controls = sol_final.decision_controls(to_merge=SolutionMerge.NODES)
diff --git a/bioptim/examples/toy_examples/moving_horizon_estimation/cyclic_nmpc.py b/bioptim/examples/toy_examples/moving_horizon_estimation/cyclic_nmpc.py
index 9fbf198d1..fab123b9d 100644
--- a/bioptim/examples/toy_examples/moving_horizon_estimation/cyclic_nmpc.py
+++ b/bioptim/examples/toy_examples/moving_horizon_estimation/cyclic_nmpc.py
@@ -22,6 +22,7 @@
     Solution,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 class MyCyclicNMPC(CyclicNonlinearModelPredictiveControl):
@@ -76,14 +77,14 @@ def prepare_nmpc(
 
 
 def main():
-    model_path = "models/arm2.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/arm2.bioMod"
     torque_max = 50
 
     cycle_duration = 1
     cycle_len = 20
     n_cycles = 3
 
-    nmpc = prepare_nmpc(model_path, cycle_len=cycle_len, cycle_duration=cycle_duration, max_torque=torque_max)
+    nmpc = prepare_nmpc(biorbd_model_path, cycle_len=cycle_len, cycle_duration=cycle_duration, max_torque=torque_max)
 
     def update_functions(_nmpc: CyclicNonlinearModelPredictiveControl, cycle_idx: int, _sol: Solution):
         return cycle_idx < n_cycles  # True if there are still some cycle to perform
diff --git a/bioptim/examples/toy_examples/moving_horizon_estimation/mhe.py b/bioptim/examples/toy_examples/moving_horizon_estimation/mhe.py
index 6e74ad6c3..6a803d399 100644
--- a/bioptim/examples/toy_examples/moving_horizon_estimation/mhe.py
+++ b/bioptim/examples/toy_examples/moving_horizon_estimation/mhe.py
@@ -34,6 +34,7 @@
     PhaseDynamics,
     SolutionMerge,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def states_to_markers(bio_model, states):
@@ -182,7 +183,7 @@ def get_solver_options(solver):
 
 
 def main():
-    biorbd_model_path = "models/cart_pendulum.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/cart_pendulum.bioMod"
     bio_model = TorqueBiorbdModel(biorbd_model_path)
 
     solver = Solver.IPOPT()  # or Solver.ACADOS()  # If ACADOS is used, it must be manually installed
diff --git a/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc.py b/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc.py
index a4eb53254..829174971 100644
--- a/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc.py
+++ b/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc.py
@@ -25,6 +25,7 @@
     SolutionMerge,
     OdeSolver,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 class MyCyclicNMPC(MultiCyclicNonlinearModelPredictiveControl):
@@ -66,7 +67,7 @@ def prepare_nmpc(
     u_bounds = BoundsList()
     u_bounds["tau"] = [-max_torque] * model.nb_q, [max_torque] * model.nb_q
 
-    new_objectives = Objective(ObjectiveFcn.Lagrange.MINIMIZE_STATE, key="q")
+    new_objectives = Objective(ObjectiveFcn.Lagrange.MINIMIZE_CONTROL, key="tau")
 
     # Rotate the wheel and force the marker of the hand to follow the marker on the wheel
     wheel_target = np.linspace(-2 * np.pi * n_cycles_simultaneous, 0, cycle_len * n_cycles_simultaneous + 1)[
@@ -97,7 +98,7 @@ def prepare_nmpc(
 
 
 def main():
-    model_path = "models/arm2.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/arm2.bioMod"
     torque_max = 50
 
     cycle_duration = 1
@@ -107,7 +108,7 @@ def main():
     n_cycles = 4
 
     nmpc = prepare_nmpc(
-        model_path,
+        biorbd_model_path,
         cycle_len=cycle_len,
         cycle_duration=cycle_duration,
         n_cycles_to_advance=n_cycles_to_advance,
diff --git a/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc_with_parameters.py b/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc_with_parameters.py
index 94ba99e55..878f67b7a 100644
--- a/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc_with_parameters.py
+++ b/bioptim/examples/toy_examples/moving_horizon_estimation/multi_cyclic_nmpc_with_parameters.py
@@ -39,6 +39,7 @@
     VariableScaling,
     ContactType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 class MyCyclicNMPC(MultiCyclicNonlinearModelPredictiveControl):
@@ -212,7 +213,7 @@ def prepare_nmpc(
 
 
 def main():
-    model_path = "models/arm2.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/arm2.bioMod"
     torque_max = 50
 
     cycle_duration = 1
@@ -222,7 +223,7 @@ def main():
     n_cycles = 4
 
     nmpc = prepare_nmpc(
-        model_path,
+        biorbd_model_path,
         cycle_len=cycle_len,
         cycle_duration=cycle_duration,
         n_cycles_to_advance=n_cycles_to_advance,
diff --git a/bioptim/examples/toy_examples/muscle_driven_ocp/muscle_excitations_tracker.py b/bioptim/examples/toy_examples/muscle_driven_ocp/muscle_excitations_tracker.py
index 81bf4dedd..4da438b0d 100644
--- a/bioptim/examples/toy_examples/muscle_driven_ocp/muscle_excitations_tracker.py
+++ b/bioptim/examples/toy_examples/muscle_driven_ocp/muscle_excitations_tracker.py
@@ -34,6 +34,7 @@
     DynamicsOptions,
 )
 from bioptim.optimization.optimization_variable import OptimizationVariableContainer
+from bioptim.examples.utils import ExampleUtils
 
 
 def generate_data(
@@ -331,9 +332,8 @@ def main():
 
     # Define the problem
     use_residual_torque = True
-    bio_model = MusclesWithExcitationsBiorbdModel(
-        "models/arm26_muscle_driven_ocp.bioMod", with_residual_torque=use_residual_torque
-    )
+    biorbd_model_path = ExampleUtils.folder + "/models/arm26_muscle_driven_ocp.bioMod"
+    bio_model = MusclesWithExcitationsBiorbdModel(biorbd_model_path, with_residual_torque=use_residual_torque)
     final_time = 0.5
     n_shooting_points = 30
     phase_dynamics = PhaseDynamics.SHARED_DURING_THE_PHASE
@@ -345,9 +345,7 @@ def main():
 
     # Track these data
     # To allow for non free variable, the model must be reloaded
-    bio_model = MusclesWithExcitationsBiorbdModel(
-        "models/arm26_muscle_driven_ocp.bioMod", with_residual_torque=use_residual_torque
-    )
+    bio_model = MusclesWithExcitationsBiorbdModel(biorbd_model_path, with_residual_torque=use_residual_torque)
     ocp = prepare_ocp(
         bio_model,
         final_time,
diff --git a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle.py b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle.py
index cbd8307c8..84addaa80 100644
--- a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle.py
+++ b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle.py
@@ -24,6 +24,7 @@
     SolutionMerge,
     ContactType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(biorbd_model_path, phase_time, n_shooting, min_bound, max_bound, expand_dynamics=True):
@@ -112,7 +113,7 @@ def prepare_ocp(biorbd_model_path, phase_time, n_shooting, min_bound, max_bound,
 
 
 def main():
-    biorbd_model_path = "models/2segments_4dof_2contacts_1muscle.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_4dof_2contacts_1muscle.bioMod"
     t = 0.3
     ns = 10
     dt = t / ns
@@ -128,7 +129,7 @@ def main():
     sol = ocp.solve(Solver.IPOPT(show_online_optim=platform.system() == "Linux"))
 
     nlp = ocp.nlp[0]
-    nlp.model = BiorbdModel(biorbd_model_path)
+    nlp.model = MusclesBiorbdModel(biorbd_model_path)
 
     states = sol.decision_states(to_merge=SolutionMerge.NODES)
     controls = sol.decision_controls(to_merge=SolutionMerge.NODES)
diff --git a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle_excitations.py b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle_excitations.py
index 99a2792d7..48d9c51b1 100644
--- a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle_excitations.py
+++ b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inequality_constraint_muscle_excitations.py
@@ -26,6 +26,7 @@
     ContactType,
     MusclesWithExcitationsBiorbdModel,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(biorbd_model_path, phase_time, n_shooting, min_bound, ode_solver=OdeSolver.RK4(), expand_dynamics=True):
@@ -117,7 +118,7 @@ def prepare_ocp(biorbd_model_path, phase_time, n_shooting, min_bound, ode_solver
 
 
 def main():
-    biorbd_model_path = "models/2segments_4dof_2contacts_1muscle.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_4dof_2contacts_1muscle.bioMod"
     t = 0.3
     ns = 10
     dt = t / ns
diff --git a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_constraint_muscle.py b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_constraint_muscle.py
index 659d50256..c3f2a0b84 100644
--- a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_constraint_muscle.py
+++ b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_constraint_muscle.py
@@ -33,6 +33,7 @@
     SolutionIntegrator,
     ContactType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def contact_velocity_all_points(controller):
@@ -188,7 +189,7 @@ def prepare_ocp(
 
 
 def main():
-    biorbd_model_path = "models/2segments_4dof_2contacts_1muscle.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_4dof_2contacts_1muscle.bioMod"
     t = 0.1
     ns = 100
 
diff --git a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_soft_contacts_muscle.py b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_soft_contacts_muscle.py
index 8c8fcc905..eee4f3c5c 100644
--- a/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_soft_contacts_muscle.py
+++ b/bioptim/examples/toy_examples/muscle_driven_with_contact/contact_forces_inverse_dynamics_soft_contacts_muscle.py
@@ -26,6 +26,7 @@
     ContactType,
     DefectType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(biorbd_model_path, phase_time, n_shooting, expand_dynamics=True):
@@ -118,7 +119,7 @@ def main():
     # This example does not converge, but it is a good example of how to set up the problem
     # And the dynamics seems fine (inf_pr = 5.66e-07) when restoration failed
 
-    biorbd_model_path = "models/2segments_4dof_2soft_contacts_1muscle.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_4dof_2soft_contacts_1muscle.bioMod"
     t = 1
     ns = 100
     ocp = prepare_ocp(
diff --git a/bioptim/examples/toy_examples/muscle_driven_with_contact/muscle_activations_contacts_tracker.py b/bioptim/examples/toy_examples/muscle_driven_with_contact/muscle_activations_contacts_tracker.py
index 4ba4121a6..c5bc1121d 100644
--- a/bioptim/examples/toy_examples/muscle_driven_with_contact/muscle_activations_contacts_tracker.py
+++ b/bioptim/examples/toy_examples/muscle_driven_with_contact/muscle_activations_contacts_tracker.py
@@ -23,6 +23,8 @@
     Node,
     ContactType,
 )
+from bioptim.examples.utils import ExampleUtils
+
 
 # Load track_segment_on_rt
 spec = importlib.util.spec_from_file_location(
@@ -99,13 +101,13 @@ def prepare_ocp(
 
 def main():
     # Define the problem
-    model_path = "models/2segments_4dof_2contacts_1muscle.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_4dof_2contacts_1muscle.bioMod"
     final_time = 0.7
     ns = 20
 
     # Generate data using another optimization that will be feedback in as tracking data
     ocp_to_track = data_to_track.prepare_ocp(
-        biorbd_model_path=model_path,
+        biorbd_model_path=biorbd_model_path,
         phase_time=final_time,
         n_shooting=ns,
         min_bound=50,
@@ -128,7 +130,7 @@ def main():
 
     # Track these data
     ocp = prepare_ocp(
-        biorbd_model_path=model_path,
+        biorbd_model_path=biorbd_model_path,
         phase_time=final_time,
         n_shooting=ns,
         muscle_activations_ref=muscle_activations_ref,
diff --git a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_muscle_driven.py b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_muscle_driven.py
index 7692d0ef3..90cc36393 100644
--- a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_muscle_driven.py
+++ b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_muscle_driven.py
@@ -47,8 +47,9 @@
     ContactType,
     BiMapping,
 )
-from bioptim.examples.stochastic_optimal_control.arm_reaching_torque_driven_implicit import ExampleType
-from bioptim.examples.stochastic_optimal_control.models.leuven_arm_model import LeuvenArmModel
+from bioptim.examples.utils import ExampleUtils
+from bioptim.examples.toy_examples.stochastic_optimal_control.arm_reaching_torque_driven_implicit import ExampleType
+from bioptim.examples.toy_examples.stochastic_optimal_control.models.leuven_arm_model import LeuvenArmModel
 
 
 def sensory_reference(
@@ -638,7 +639,7 @@ def main():
     plot_sol_flag = False
     vizualise_sol_flag = False
 
-    biorbd_model_path = "models/LeuvenArmModel.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/LeuvenArmModel.bioMod"
 
     hand_initial_position = np.array([0.0, 0.2742])  # Directly from Tom's version
     hand_final_position = np.array([9.359873986980460e-12, 0.527332023564034])  # Directly from Tom's version
diff --git a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_collocations.py b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_collocations.py
index a97ba3c31..b18608299 100644
--- a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_collocations.py
+++ b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_collocations.py
@@ -28,6 +28,7 @@
     SolutionMerge,
     DynamicsOptions,
 )
+from bioptim.examples.utils import ExampleUtils
 
 from bioptim.examples.toy_examples.stochastic_optimal_control.arm_reaching_torque_driven_implicit import ExampleType
 
@@ -272,7 +273,7 @@ def main():
     use_sx = True
     vizualize_sol_flag = True
 
-    biorbd_model_path = "models/LeuvenArmModel.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/LeuvenArmModel.bioMod"
 
     hand_final_position = np.array([9.359873986980460e-12, 0.527332023564034])  # Directly from Tom's version
 
diff --git a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_explicit.py b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_explicit.py
index 7045763b1..007f78e30 100644
--- a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_explicit.py
+++ b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_explicit.py
@@ -41,8 +41,9 @@
     PhaseDynamics,
     BiMapping,
 )
-from bioptim.examples.stochastic_optimal_control.arm_reaching_torque_driven_implicit import ExampleType
-from bioptim.examples.stochastic_optimal_control.common import (
+from bioptim.examples.utils import ExampleUtils
+from bioptim.examples.toy_examples.stochastic_optimal_control.arm_reaching_torque_driven_implicit import ExampleType
+from bioptim.examples.toy_examples.stochastic_optimal_control.common import (
     dynamics_torque_driven_with_feedbacks,
 )
 
@@ -594,7 +595,7 @@ def main():
     use_sx = False
     vizualize_sol_flag = True
 
-    biorbd_model_path = "models/LeuvenArmModel.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/LeuvenArmModel.bioMod"
 
     hand_final_position = np.array([9.359873986980460e-12, 0.527332023564034])  # Directly from Tom's version
 
diff --git a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_implicit.py b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_implicit.py
index cf42c3e65..bf0e6251d 100644
--- a/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_implicit.py
+++ b/bioptim/examples/toy_examples/stochastic_optimal_control/arm_reaching_torque_driven_implicit.py
@@ -36,6 +36,7 @@
     ControlType,
     VariableScalingList,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 class ExampleType(Enum):
@@ -343,7 +344,7 @@ def main():
     with_cholesky = True
     with_scaling = True
 
-    biorbd_model_path = "models/LeuvenArmModel.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/LeuvenArmModel.bioMod"
 
     hand_final_position = np.array([9.359873986980460e-12, 0.527332023564034])  # Directly from Tom's version
 
diff --git a/bioptim/examples/toy_examples/stochastic_optimal_control/rockit_matrix_lyapunov.py b/bioptim/examples/toy_examples/stochastic_optimal_control/rockit_matrix_lyapunov.py
index 39c0025c5..cd0fdf75a 100644
--- a/bioptim/examples/toy_examples/stochastic_optimal_control/rockit_matrix_lyapunov.py
+++ b/bioptim/examples/toy_examples/stochastic_optimal_control/rockit_matrix_lyapunov.py
@@ -27,12 +27,12 @@
     PhaseDynamics,
     SolutionMerge,
 )
-from bioptim.examples.stochastic_optimal_control.models.rockit_model import (
+from bioptim.examples.toy_examples.stochastic_optimal_control.models.rockit_model import (
     RockitModel,
     RockitDynamicsOCP,
     RockitDynamicsSOCP,
 )
-from bioptim.examples.stochastic_optimal_control.common import (
+from bioptim.examples.toy_examples.stochastic_optimal_control.common import (
     test_matrix_semi_definite_positiveness,
     test_eigen_values,
     reshape_to_matrix,
diff --git a/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_constraint.py b/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_constraint.py
index 6b71ec9aa..4711ed396 100644
--- a/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_constraint.py
+++ b/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_constraint.py
@@ -30,10 +30,11 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
-    biorbd_model_path: str = "models/cubeSym.bioMod",
+    biorbd_model_path,
     ode_solver: OdeSolverBase = OdeSolver.RK4(),
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     expand_dynamics: bool = True,
@@ -114,7 +115,8 @@ def main():
     Solves an ocp where the symmetry is enforced by constraints, and animates it
     """
 
-    ocp = prepare_ocp()
+    biorbd_model_path = ExampleUtils.folder + "/models/cubeSym.bioMod"
+    ocp = prepare_ocp(biorbd_model_path)
 
     # Objective and constraints plots
     ocp.add_plot_penalty()
diff --git a/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_mapping.py b/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_mapping.py
index 8bb8e1ad8..297044c8e 100644
--- a/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_mapping.py
+++ b/bioptim/examples/toy_examples/symmetrical_torque_driven_ocp/symmetry_by_mapping.py
@@ -44,10 +44,11 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
-    biorbd_model_path: str = "models/cubeSym.bioMod",
+    biorbd_model_path,
     ode_solver: OdeSolverBase = OdeSolver.RK4(),
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     expand_dynamics: bool = True,
@@ -151,7 +152,8 @@ def main():
     Solves an ocp where the symmetry must be respected, and animates it
     """
 
-    ocp = prepare_ocp()
+    biorbd_model_path = ExampleUtils.folder + "/models/cubeSym.bioMod"
+    ocp = prepare_ocp(biorbd_model_path)
 
     # --- Solve the program --- #
     sol = ocp.solve(Solver.IPOPT(show_online_optim=platform.system() == "Linux"))
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/example_multi_biorbd_model.py b/bioptim/examples/toy_examples/torque_driven_ocp/example_multi_biorbd_model.py
index 9638d5b58..fe46b867d 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/example_multi_biorbd_model.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/example_multi_biorbd_model.py
@@ -15,11 +15,12 @@
     PhaseDynamics,
     SolutionMerge,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
-    biorbd_model_path: str = "models/triple_pendulum.bioMod",
-    biorbd_model_path_modified_inertia: str = "models/triple_pendulum_modified_inertia.bioMod",
+    biorbd_model_path,
+    biorbd_model_path_modified_inertia,
     n_shooting: int = 40,
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     expand_dynamics: bool = True,
@@ -75,8 +76,14 @@ def prepare_ocp(
 
 
 def main():
+    example_folder = ExampleUtils.folder
+    biorbd_model_path = example_folder + "/models/triple_pendulum.bioMod"
+    biorbd_model_path_modified_inertia = example_folder + "/models/triple_pendulum_modified_inertia.bioMod"
+
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp()
+    ocp = prepare_ocp(
+        biorbd_model_path=biorbd_model_path, biorbd_model_path_modified_inertia=biorbd_model_path_modified_inertia
+    )
     ocp.add_plot_penalty()
 
     # --- Solve the program --- #
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/example_pendulum_time_dependent.py b/bioptim/examples/toy_examples/torque_driven_ocp/example_pendulum_time_dependent.py
index 08403806e..f38a1ba4e 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/example_pendulum_time_dependent.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/example_pendulum_time_dependent.py
@@ -33,6 +33,7 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 class TimeDependentModel(TorqueBiorbdModel):
@@ -174,8 +175,10 @@ def main():
     If pendulum is run as a script, it will perform the optimization and animates it
     """
 
+    biorbd_model_path = ExampleUtils.folder + "/models/pendulum.bioMod"
+
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp(biorbd_model_path="models/pendulum.bioMod", final_time=1, n_shooting=30)
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path, final_time=1, n_shooting=30)
 
     # Custom plots
     ocp.add_plot_penalty(CostType.ALL)
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/example_rigid_contact.py b/bioptim/examples/toy_examples/torque_driven_ocp/example_rigid_contact.py
index 82d3aaa53..f3b8b81c2 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/example_rigid_contact.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/example_rigid_contact.py
@@ -34,6 +34,7 @@
     PenaltyController,
     ExternalForceSetVariables,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def custom_com_over_contact(controller: PenaltyController) -> MX:
@@ -147,7 +148,7 @@ def main():
     """
     Defines a multiphase ocp and animate the results
     """
-    biorbd_model_path = "../torque_driven_ocp/models/3segments_4dof_1contact.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/3segments_4dof_1contact.bioMod"
     n_shooting = 30
     final_time = 1
     defects_type = DefectType.QDDOT_EQUALS_FORWARD_DYNAMICS
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/example_soft_contact.py b/bioptim/examples/toy_examples/torque_driven_ocp/example_soft_contact.py
index 0721c9120..72572e753 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/example_soft_contact.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/example_soft_contact.py
@@ -27,6 +27,7 @@
     SolutionMerge,
     ContactType,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_single_shooting(
@@ -194,7 +195,7 @@ def main():
     """
     Defines a multiphase ocp and animate the results
     """
-    biorbd_model_path = "../torque_driven_ocp/models/soft_contact_sphere.bioMod"
+    biorbd_model_path = ExampleUtils.folder + "/models/soft_contact_sphere.bioMod"
     ode_solver = OdeSolver.RK8()
 
     # Prepare OCP to reach the second marker
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/minmax_torque_by_extra_parameter_multiphase.py b/bioptim/examples/toy_examples/torque_driven_ocp/minmax_torque_by_extra_parameter_multiphase.py
index afc197bf0..06f211781 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/minmax_torque_by_extra_parameter_multiphase.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/minmax_torque_by_extra_parameter_multiphase.py
@@ -29,6 +29,7 @@
     PenaltyController,
     ParameterObjectiveList,
 )
+from bioptim.examples.utils import ExampleUtils
 from matplotlib import pyplot as plt
 
 
@@ -51,10 +52,10 @@ def my_parameter_function(bio_model: BioModel, value: MX):
 
 
 def prepare_ocp(
+    biorbd_model_path: str,
     parameter_option: int = 0,
-    bio_model_path: str = "models/double_pendulum.bioMod",
 ) -> OptimalControlProgram:
-    bio_model = (TorqueBiorbdModel(bio_model_path), TorqueBiorbdModel(bio_model_path))
+    bio_model = (TorqueBiorbdModel(biorbd_model_path), TorqueBiorbdModel(biorbd_model_path))
 
     # Problem parameters
     n_shooting = (30, 30)
@@ -69,7 +70,7 @@ def prepare_ocp(
     constraints = ConstraintList()
 
     # Define the parameter to optimize
-    parameters = ParameterList()
+    parameters = ParameterList(use_sx=False)
     parameter_init = InitialGuessList()
     parameter_bounds = BoundsList()
     parameter_objectives = ParameterObjectiveList()
@@ -184,6 +185,8 @@ def prepare_ocp(
 
 def main():
     # --- Prepare the ocp --- #
+    biorbd_model_path = ExampleUtils.folder + "/models/double_pendulum.bioMod"
+
     fig, axs = plt.subplots(1, 3)
     axs[0].set_title("Joint coordinates")
     axs[0].set_ylabel("q [°]")
@@ -199,7 +202,7 @@ def main():
         axs[ax].grid(True)
 
     for i, linestyle in enumerate(linestyles):
-        ocp = prepare_ocp(parameter_option=i)
+        ocp = prepare_ocp(biorbd_model_path=biorbd_model_path, parameter_option=i)
 
         # --- Solve the ocp --- #
         sol = ocp.solve()
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/ocp_mass_with_ligament.py b/bioptim/examples/toy_examples/torque_driven_ocp/ocp_mass_with_ligament.py
index 5e17a6edd..cc3f371fc 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/ocp_mass_with_ligament.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/ocp_mass_with_ligament.py
@@ -14,6 +14,7 @@
     TorqueBiorbdModel,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
@@ -99,8 +100,8 @@ def prepare_ocp(
 
 
 def main():
-    model_path = "./models/mass_point_with_ligament.bioMod"
-    ocp = prepare_ocp(biorbd_model_path=model_path)
+    biorbd_model_path = ExampleUtils.folder + "/models/mass_point_with_ligament.bioMod"
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path)
 
     # --- Solve the program --- #
     sol = ocp.solve()
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/phase_transition_uneven_variable_number_by_bounds.py b/bioptim/examples/toy_examples/torque_driven_ocp/phase_transition_uneven_variable_number_by_bounds.py
index 33905766d..628992e99 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/phase_transition_uneven_variable_number_by_bounds.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/phase_transition_uneven_variable_number_by_bounds.py
@@ -15,18 +15,19 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
-    biorbd_model_path_with_translations: str = "models/double_pendulum_with_translations.bioMod",
+    biorbd_model_path: str,
     n_shooting: tuple = (40, 40),
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     expand_dynamics: bool = True,
 ) -> OptimalControlProgram:
 
     bio_model = (
-        TorqueBiorbdModel(biorbd_model_path_with_translations),
-        TorqueBiorbdModel(biorbd_model_path_with_translations),
+        TorqueBiorbdModel(biorbd_model_path),
+        TorqueBiorbdModel(biorbd_model_path),
     )
 
     # Problem parameters
@@ -117,7 +118,8 @@ def prepare_ocp(
 
 def main():
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp()
+    biorbd_model_path = ExampleUtils.folder + "/models/double_pendulum_with_translations.bioMod"
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path)
 
     # --- Solve the program --- #
     sol = ocp.solve(Solver.IPOPT(show_online_optim=False, show_options=dict(show_bounds=True)))
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/spring_load.py b/bioptim/examples/toy_examples/torque_driven_ocp/spring_load.py
index 52516ebe5..2174d9f96 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/spring_load.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/spring_load.py
@@ -25,6 +25,7 @@
     PhaseDynamics,
     SolutionMerge,
 )
+from bioptim.examples.utils import ExampleUtils
 
 # scenarios are based on a Mayer term (at Tf)
 # 0: maximize upward speed - expected kinematics: negative torque to get as low as possible and release
@@ -166,7 +167,7 @@ def dynamics(
 
 
 def prepare_ocp(
-    biorbd_model_path: str = "models/mass_point.bioMod",
+    biorbd_model_path: str,
     phase_dynamics: PhaseDynamics = PhaseDynamics.SHARED_DURING_THE_PHASE,
     expand_dynamics: bool = True,
     phase_time: float = 0.5,
@@ -237,11 +238,15 @@ def prepare_ocp(
 def main():
     phase_time = 0.5
     n_shooting = 30
+    biorbd_model_path = ExampleUtils.folder + "/models/mass_point.bioMod"
+
     fig, axs = plt.subplots(1, 3)
 
     for scenario in range(8):  # in [1]: #
         print(scenarios[scenario]["label"])
-        ocp = prepare_ocp(phase_time=phase_time, n_shooting=n_shooting, scenario=scenario)
+        ocp = prepare_ocp(
+            biorbd_model_path=biorbd_model_path, phase_time=phase_time, n_shooting=n_shooting, scenario=scenario
+        )
 
         ocp.print(to_console=True, to_graph=False)
 
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/torque_activation_driven.py b/bioptim/examples/toy_examples/torque_driven_ocp/torque_activation_driven.py
index d0b0f0233..3601babfc 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/torque_activation_driven.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/torque_activation_driven.py
@@ -17,6 +17,7 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
@@ -112,9 +113,10 @@ def main():
     """
     Prepares and solves an ocp with torque actuators, the animates it
     """
+    biorbd_model_path = ExampleUtils.folder + "/models/2segments_2dof_2contacts.bioMod"
 
     ocp = prepare_ocp(
-        biorbd_model_path="models/2segments_2dof_2contacts.bioMod",
+        biorbd_model_path=biorbd_model_path,
         n_shooting=30,
         final_time=2,
     )
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/torque_driven_free_floating_base.py b/bioptim/examples/toy_examples/torque_driven_ocp/torque_driven_free_floating_base.py
index 8b5cc2973..409d1d157 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/torque_driven_free_floating_base.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/torque_driven_free_floating_base.py
@@ -19,6 +19,7 @@
     Solver,
     Node,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(biorbd_model_path: str):
@@ -164,7 +165,6 @@ def prepare_ocp(biorbd_model_path: str):
         bio_model=bio_model,
         n_shooting=n_shooting,
         phase_time=phase_time,
-        dynamics=dynamics,
         x_bounds=x_bounds,
         u_bounds=u_bounds,
         x_init=x_initial_guesses,
@@ -176,7 +176,8 @@ def prepare_ocp(biorbd_model_path: str):
 
 if __name__ == "__main__":
     # --- Prepare the ocp --- #
-    ocp = prepare_ocp(biorbd_model_path="models/trunk_and_2arm.bioMod")
+    biorbd_model_path = ExampleUtils.folder + "/models/trunk_and_2arm.bioMod"
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path)
 
     # --- Solve the ocp --- #
     solver = Solver.IPOPT(show_online_optim=platform.system() == "Linux")
diff --git a/bioptim/examples/toy_examples/torque_driven_ocp/track_markers_with_torque_actuators.py b/bioptim/examples/toy_examples/torque_driven_ocp/track_markers_with_torque_actuators.py
index 409c51bb5..1caefbba5 100644
--- a/bioptim/examples/toy_examples/torque_driven_ocp/track_markers_with_torque_actuators.py
+++ b/bioptim/examples/toy_examples/torque_driven_ocp/track_markers_with_torque_actuators.py
@@ -26,6 +26,7 @@
     Solver,
     PhaseDynamics,
 )
+from bioptim.examples.utils import ExampleUtils
 
 
 def prepare_ocp(
@@ -129,8 +130,9 @@ def main():
     """
     Prepares and solves an ocp with torque actuators, the animates it
     """
+    biorbd_model_path = ExampleUtils.folder + "/models/cube_with_actuators.bioMod"
 
-    ocp = prepare_ocp("models/cube.bioMod", n_shooting=30, final_time=2, actuator_type=2)
+    ocp = prepare_ocp(biorbd_model_path=biorbd_model_path, n_shooting=30, final_time=2, actuator_type=2)
 
     # --- Solve the program --- #
     sol = ocp.solve(Solver.IPOPT(show_online_optim=platform.system() == "Linux"))
diff --git a/bioptim/optimization/vector_utils.py b/bioptim/optimization/vector_utils.py
index e37cfad25..112de5738 100644
--- a/bioptim/optimization/vector_utils.py
+++ b/bioptim/optimization/vector_utils.py
@@ -106,7 +106,7 @@ def _compute_value_for_node(
     """
     collapsed_values = np.ndarray((variable_container.shape, 1))
 
-    real_keys = [key for key in defined_values.keys() if key is not "None"]
+    real_keys = [key for key in defined_values.keys() if key != "None"]
     for key in real_keys:
 
         if defined_values[key].type == InterpolationType.ALL_POINTS:
diff --git a/tests/shard1/test_prepare_all_examples.py b/tests/shard1/test_prepare_all_examples.py
index 4704ec875..d431ba71c 100644
--- a/tests/shard1/test_prepare_all_examples.py
+++ b/tests/shard1/test_prepare_all_examples.py
@@ -572,8 +572,7 @@ def test__torque_driven_ocp__phase_transition_uneven_variable_number_by_bounds()
     bioptim_folder = TestUtils.bioptim_folder()
 
     ocp_module.prepare_ocp(
-        biorbd_model_path_with_translations=bioptim_folder
-        + "/examples/models/double_pendulum_with_translations.bioMod",
+        biorbd_model_path=bioptim_folder + "/examples/models/double_pendulum_with_translations.bioMod",
         n_shooting=(5, 5),
         phase_dynamics=PhaseDynamics.SHARED_DURING_THE_PHASE,
         expand_dynamics=False,
diff --git a/tests/shard2/test_global_nmpc_final.py b/tests/shard2/test_global_nmpc_final.py
index 49ed69d59..80a445d1f 100644
--- a/tests/shard2/test_global_nmpc_final.py
+++ b/tests/shard2/test_global_nmpc_final.py
@@ -25,7 +25,6 @@
     ],
 )
 def test_multi_cyclic_nmpc_get_final(phase_dynamics, ode_solver):
-
     if platform.system() == "Windows" or platform.system() == "Darwin":
         pytest.skip("This test is skipped on Windows and macOS because sensitive")
 
@@ -71,12 +70,12 @@ def update_functions(_nmpc, cycle_idx, _sol):
         npt.assert_almost_equal(q[:, -1], np.array([1.37753244e-40, 1.04359174e00, 1.03625065e00]))
 
         # initial and final velocities
-        npt.assert_almost_equal(qdot[:, 0], np.array((6.28293718, 2.5617072, -0.00942694)))
-        npt.assert_almost_equal(qdot[:, -1], np.array([6.28293718, 2.41433059, -0.59773899]), decimal=5)
+        npt.assert_almost_equal(qdot[:, 0], np.array([6.28318531, 2.57685385, -0.06542029]))
+        npt.assert_almost_equal(qdot[:, -1], np.array([6.28318531, 2.40501008, -0.63136038]), decimal=5)
 
         # initial and final controls
-        npt.assert_almost_equal(tau[:, 0], np.array((0.00992505, 4.88488618, 2.4400698)))
-        npt.assert_almost_equal(tau[:, -1], np.array([-0.00992505, 5.19414727, 2.34022319]), decimal=4)
+        npt.assert_almost_equal(tau[:, 0], np.array([1.34399907e-11, 4.92501368e00, 2.53191286e00]))
+        npt.assert_almost_equal(tau[:, -1], np.array([2.82405932e-12, 4.96328211e00, 2.20504597e00]), decimal=4)
 
         # check time
         n_steps = nmpc.nlp[0].dynamics_type.ode_solver.n_integration_steps
@@ -122,11 +121,11 @@ def update_functions(_nmpc, cycle_idx, _sol):
             npt.assert_almost_equal(q[:, 0], np.array((-12.56637061, 1.04359174, 1.03625065)))
             npt.assert_almost_equal(q[:, -1], np.array([0.0, 1.04359174, 1.03625065]))
             # initial and final velocities
-            npt.assert_almost_equal(qdot[:, 0], np.array([6.28810582, 2.55280178, 0.02627301]), decimal=5)
-            npt.assert_almost_equal(qdot[:, -1], np.array([6.28810582, 2.42201137, -0.57799103]), decimal=5)
+            npt.assert_almost_equal(qdot[:, 0], np.array([6.28318531, 2.5769604, -0.06617214]), decimal=5)
+            npt.assert_almost_equal(qdot[:, -1], np.array([6.28318531, 2.40569762, -0.63227535]), decimal=5)
             # initial and final controls
-            npt.assert_almost_equal(tau[:, 0], np.array([-0.19682043, 4.84862251, 2.37825343]), decimal=4)
-            npt.assert_almost_equal(tau[:, -1], np.array([0.19682043, 5.35831142, 2.43094827]), decimal=4)
+            npt.assert_almost_equal(tau[:, 0], np.array([-1.57225142e-10, 4.92697679e00, 2.53359257e00]), decimal=4)
+            npt.assert_almost_equal(tau[:, -1], np.array([1.75585447e-10, 4.96765087e00, 2.20617819e00]), decimal=4)
 
         # check time
         n_steps = nmpc.nlp[0].dynamics_type.ode_solver.polynomial_degree
@@ -203,9 +202,9 @@ def update_functions(_nmpc, cycle_idx, _sol):
     # check some result of the third structure
     assert len(sol[2]) == 3
 
-    npt.assert_almost_equal(sol[2][0].cost.toarray().squeeze(), 0.0002)
-    npt.assert_almost_equal(sol[2][1].cost.toarray().squeeze(), 0.0002)
-    npt.assert_almost_equal(sol[2][2].cost.toarray().squeeze(), 0.0002)
+    npt.assert_almost_equal(sol[2][0].cost.toarray().squeeze(), 0.0)
+    npt.assert_almost_equal(sol[2][1].cost.toarray().squeeze(), 0.0)
+    npt.assert_almost_equal(sol[2][2].cost.toarray().squeeze(), 0.0)
 
 
 @pytest.mark.parametrize("phase_dynamics", [PhaseDynamics.SHARED_DURING_THE_PHASE])
diff --git a/tests/shard3/test_global_torque_driven_ocp.py b/tests/shard3/test_global_torque_driven_ocp.py
index 1f8dd6cda..1d35df53e 100644
--- a/tests/shard3/test_global_torque_driven_ocp.py
+++ b/tests/shard3/test_global_torque_driven_ocp.py
@@ -461,7 +461,7 @@ def test_phase_transition_uneven_variable_number_by_bounds(phase_dynamics):
     biorbd_model_path_with_translations = bioptim_folder + "/examples/models/double_pendulum_with_translations.bioMod"
 
     ocp = ocp_module.prepare_ocp(
-        biorbd_model_path_with_translations=biorbd_model_path_with_translations,
+        biorbd_model_path=biorbd_model_path_with_translations,
         n_shooting=(10, 10),
         phase_dynamics=phase_dynamics,
         expand_dynamics=True,