Acceleration Qontrol#

This example demonstrates basic acceleration control.

The following example solves a qp problem expressed at the joint acceleration level such that:

\[\begin{split}\begin{array}{ccc}\boldsymbol{\ddot{q}}^{opt} = & \underset{\boldsymbol{\ddot{q}}}{\mathrm{argmin}} & ||J(\boldsymbol{q})\boldsymbol{\ddot{q}} + \dot{J}(\boldsymbol{q},\boldsymbol{\dot{q}})\boldsymbol{\dot{q}} - \boldsymbol{\dot{v}}^{target} || + \omega || \boldsymbol{\ddot{q}} ||^2\\& \textrm{s.t.} & \boldsymbol{\dot{q}^{min}} \leq \boldsymbol{\dot{q}}(\boldsymbol{\ddot{q}}) \leq \boldsymbol{\dot{q}^{max}}. \\ & & \boldsymbol{q}^{min} \leq \boldsymbol{q}(\boldsymbol{\ddot{q}}) \leq \boldsymbol{q}^{max} \end{array}\end{split}\]

The robot main tasks consists in following a simple trajectory defined in Cartesian space. The mujoco library is used to simulate the robot behaviour.

Full Code:

You can find the source code of this example here.

  1// Copyright 2021 DeepMind Technologies Limited
  2//
  3// Licensed under the Apache License, Version 2.0 (the "License");
  4// you may not use this file except in compliance with the License.
  5// You may obtain a copy of the License at
  6//
  7//     http://www.apache.org/licenses/LICENSE-2.0
  8//
  9// Unless required by applicable law or agreed to in writing, software
 10// distributed under the License is distributed on an "AS IS" BASIS,
 11// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 12// See the License for the specific language governing permissions and
 13// limitations under the License.
 14
 15#include "mujoco/mujoco_sim.h"
 16#include "Qontrol/Qontrol.hpp"
 17#include "trajectory_generation/trajectory_generation.h"
 18
 19
 20using namespace Qontrol;
 21
 22class MujocoQontrol : public MujocoSim 
 23{
 24public:
 25//------------------------------------------- simulation -------------------------------------------
 26std::shared_ptr<Qontrol::Model::RobotModel> model;
 27std::shared_ptr<Qontrol::JointAccelerationProblem> acceleration_problem;
 28std::shared_ptr<Qontrol::Task::CartesianAcceleration<ControlOutput::JointAcceleration>> main_task;
 29
 30Qontrol::RobotState robot_state;
 31
 32TrajectoryGeneration* traj;
 33std::string resource_path;
 34
 35void initController() override
 36{
 37  model =
 38      Model::RobotModel::loadModelFromFile(resource_path+"robot.urdf");
 39  
 40  acceleration_problem = std::make_shared<Qontrol::JointAccelerationProblem>(model);  
 41  main_task = acceleration_problem->task_set->add<Task::CartesianAcceleration>("MainTask"); 
 42  auto regularisation_task = acceleration_problem->task_set->add<Task::JointAcceleration>("RegularisationTask",1e-5); 
 43
 44  auto joint_configuration_constraint = acceleration_problem->constraint_set->add<Constraint::JointConfiguration>("JointConfigurationConstraint");
 45  auto joint_velocity_constraint = acceleration_problem->constraint_set->add<Constraint::JointVelocity>("JointVelocityConstraint");
 46  
 47  mju_copy(d->qpos, m->key_qpos, m->nu);
 48  robot_state.joint_position.resize(model->getNrOfDegreesOfFreedom());
 49  robot_state.joint_velocity.resize(model->getNrOfDegreesOfFreedom());
 50
 51  traj = new TrajectoryGeneration(resource_path + "trajectory.csv",
 52                                  m->opt.timestep);
 53}
 54
 55void updateController() override
 56{
 57
 58  for (int i=0; i<model->getNrOfDegreesOfFreedom() ; ++i)
 59  {
 60    robot_state.joint_position[i] = d->qpos[i];
 61    robot_state.joint_velocity[i] = d->qvel[i];
 62  }
 63  model->setRobotState(robot_state);
 64
 65  traj->update();
 66  pinocchio::SE3 traj_pose(traj->pose.matrix());
 67  
 68  pinocchio::SE3 current_pose(model->getFramePose(model->getTipFrameName()).matrix());
 69  const pinocchio::SE3 tipMdes = current_pose.actInv(traj_pose);
 70  auto err = pinocchio::log6(tipMdes).toVector();
 71  Eigen::Matrix<double, 6, 1> p_gains;
 72  p_gains << 1000, 1000, 1000, 1000, 1000, 1000;
 73  Eigen::Matrix<double, 6, 1> d_gains = 2.0 * p_gains.cwiseSqrt();
 74  Eigen::Matrix<double, 6, 1> xdd_star =
 75    p_gains.cwiseProduct(err) +
 76    d_gains.cwiseProduct(traj->velocity - 
 77    model->getFrameVelocity(model->getTipFrameName())) + 
 78    traj->acceleration;
 79  
 80  main_task->setTargetAcceleration(xdd_star);
 81  
 82  acceleration_problem->update(m->opt.timestep);
 83
 84  if (acceleration_problem->solutionFound())
 85  {
 86    sendJointVelocity(acceleration_problem->getJointVelocityCommand());
 87  }
 88}
 89 
 90};
 91
 92int main(int argc, const char** argv) {
 93  MujocoQontrol mujoco_qontrol;
 94  Qontrol::Log::Logger::parseArgv(argc, argv);
 95
 96  mjvCamera cam;
 97  mjv_defaultCamera(&cam);
 98
 99  mjvOption opt;
100  mjv_defaultOption(&opt);
101
102  mjvPerturb pert;
103  mjv_defaultPerturb(&pert);
104
105  // simulate object encapsulates the UI
106  auto sim = std::make_unique<mj::Simulate>(
107      std::make_unique<mj::GlfwAdapter>(),
108      &cam, &opt, &pert, /* is_passive = */ false
109  );
110
111  std::string robot = argv[1];
112  std::string mujoco_scene = "./resources/"+robot+"/scene.xml";
113  mujoco_qontrol.resource_path = "./resources/"+robot+"/";
114
115  // start physics thread
116  std::thread physicsthreadhandle( &MujocoQontrol::PhysicsThread, mujoco_qontrol, sim.get(), mujoco_scene.c_str());
117
118  // start simulation UI loop (blocking call)
119  sim->RenderLoop();
120  physicsthreadhandle.join();
121
122  return 0;
123}

Explanation of the code

Declaration

First we declare all the objects that will be used to define our problem.

std::shared_ptr<Qontrol::Model::RobotModel> model;

We use pinocchio for our model library.

std::shared_ptr<Qontrol::JointAccelerationProblem> acceleration_problem;

The output of our qp controller is at the acceleration level.

std::shared_ptr<Qontrol::Task::CartesianAcceleration<ControlOutput::JointAcceleration>> main_task;

The main task is expressed as a Cartesian acceleration task.

Initialization

void initController() override
{
  model =
      Model::RobotModel::loadModelFromFile(resource_path+"robot.urdf");

During initialization we instantiate the model with the robot urdf.

std::shared_ptr<Qontrol::JointAccelerationProblem> acceleration_problem;

We initialize the problem by giving it the model. By default, the qpmad library is used.

  auto regularisation_task = acceleration_problem->task_set->add<Task::JointAcceleration>("RegularisationTask",1e-5); 

We then fill the task set of acceleration_problem with the main task and the regularisation task. Each tasks is given a name and a relative weight \(\omega\). This weight can be modified at any time. In this example, the regularisation tasks is defined as a joint acceleration task. Its means that this task will minimize the overall robot joint acceleration.

  auto joint_velocity_constraint = acceleration_problem->constraint_set->add<Constraint::JointVelocity>("JointVelocityConstraint");

We then fill the constraint set of acceleration_problem with the two pre-implemented constraints. Each constraint is given a name. These constraints will automatically be updated during the update of Qontrol.

  robot_state.joint_velocity.resize(model->getNrOfDegreesOfFreedom());

We also create the robot state and resize it.

#include "trajectory_generation/trajectory_generation.h"

We create a simple trajectory that has been precalculated and store in a csv file. This trajectory start at the robot current Cartesian pose and does a translation of (-0.1, -0.1, -0.1) m.

Update

void updateController() override
{

  for (int i=0; i<model->getNrOfDegreesOfFreedom() ; ++i)
  {
    robot_state.joint_position[i] = d->qpos[i];
    robot_state.joint_velocity[i] = d->qvel[i];
  }
  model->setRobotState(robot_state);

  traj->update();

The update function is called every milliseconds. At the beginning of each update we fill the new robot state according to the simulated robot. We also update the trajectory so that it gives the next Cartesian pose to reach in 1 ms.

  Eigen::Matrix<double, 6, 1> xdd_star =

We then compute the desired Cartesian acceleration using a simple derivate proportionnal controller. Pinocchio is used to compute the error between the desired Cartesian pose and the current Cartesian pose. This is done by the log6 function. The p_gains are the proportionnal gains of the controller and the d_gains are the derivate gains of the controller.

  main_task->setTargetAcceleration(xdd_star);

The desired Cartesian acceleration is then fed to the main task.

  
  acceleration_problem->update(m->opt.timestep);

  if (acceleration_problem->solutionFound())
  {
    sendJointVelocity(acceleration_problem->getJointVelocityCommand());
  }
}
 
};

Once we updated the necassary tasks and constraints we can update the whole problem. If a solution to the problem exist we can then get it and send it to the simulated robot.

Main function

int main(int argc, const char** argv) {
  MujocoQontrol mujoco_qontrol;
  Qontrol::Log::Logger::parseArgv(argc, argv);

  mjvCamera cam;
  mjv_defaultCamera(&cam);

  mjvOption opt;
  mjv_defaultOption(&opt);

  mjvPerturb pert;
  mjv_defaultPerturb(&pert);

  // simulate object encapsulates the UI
  auto sim = std::make_unique<mj::Simulate>(
      std::make_unique<mj::GlfwAdapter>(),
      &cam, &opt, &pert, /* is_passive = */ false
  );

  std::string robot = argv[1];
  std::string mujoco_scene = "./resources/"+robot+"/scene.xml";
  mujoco_qontrol.resource_path = "./resources/"+robot+"/";

  // start physics thread
  std::thread physicsthreadhandle( &MujocoQontrol::PhysicsThread, mujoco_qontrol, sim.get(), mujoco_scene.c_str());

  // start simulation UI loop (blocking call)
  sim->RenderLoop();
  physicsthreadhandle.join();

  return 0;

The main function fetches the robot name given in argv and starts the Mujoco simulation.

Full Code:

You can find the source code of this example here.

  1#!/usr/bin/env python3
  2"""
  3Non-Interactive Acceleration Control with Trajectory Following
  4
  5This example demonstrates:
  6- Loading a robot model in both MuJoCo and Qontrol
  7- Joint velocity control via joint acceleration QP
  8- Following a pre-defined Cartesian trajectory from CSV file
  9- Cartesian acceleration tracking with joint velocity commands
 10
 11Mirrors the C++ example: accelerationQontrol.cpp
 12
 13Usage:
 14    python acceleration_qontrol.py <robot_name>
 15    
 16Example:
 17    python acceleration_qontrol.py panda
 18    python acceleration_qontrol.py panda --verbose
 19"""
 20
 21import numpy as np
 22import qontrol
 23import pinocchio as pin
 24import argparse
 25import os
 26import time
 27from mujoco_helper import MujocoSimulator, compute_log6_error
 28from trajectory_generation import TrajectoryGeneration
 29
 30
 31class QontrolAccelerationController:
 32    """Acceleration controller with trajectory following using Qontrol"""
 33    
 34    def __init__(self, robot_name: str, resource_path: str, mujoco_sim: MujocoSimulator, verbose: bool = False):
 35        """
 36        Initialize the Qontrol acceleration controller
 37        
 38        Args:
 39            robot_name: Name of the robot (e.g., 'panda')
 40            resource_path: Path to robot resources folder
 41            mujoco_sim: MuJoCo simulator instance
 42            verbose: Print detailed diagnostics each iteration (default: False)
 43        """
 44        self.robot_name = robot_name
 45        self.resource_path = resource_path
 46        self.sim = mujoco_sim
 47        
 48        # =================================================================
 49        # Qontrol Setup - Robot Model and QP Problem
 50        # =================================================================
 51        urdf_path = os.path.join(resource_path, "robot.urdf")
 52        if not os.path.exists(urdf_path):
 53            raise FileNotFoundError(f"URDF file not found: {urdf_path}")
 54        
 55        print(f"Loading Qontrol model from: {urdf_path}")
 56        self.qontrol_model = qontrol.RobotModel.load_from_file(urdf_path)
 57        self.ndof = self.qontrol_model.get_nr_of_degrees_of_freedom()
 58        print(f"Robot has {self.ndof} degrees of freedom")
 59        
 60        # Create QP solver
 61        self.solver = qontrol.create_qpmad_solver()
 62        
 63        # Create acceleration-level problem
 64        self.acceleration_problem = qontrol.JointAccelerationProblem(self.qontrol_model, self.solver)
 65        
 66        # Add Cartesian acceleration task for end-effector tracking
 67        self.main_task = self.acceleration_problem.task_set.add_cartesian_acceleration("MainTask", 1.0)
 68        
 69        # Add joint acceleration regularization (minimize joint accelerations)
 70        self.regularization_task = self.acceleration_problem.task_set.add_joint_acceleration("RegularizationTask", 1e-5)
 71        
 72        # Add constraints
 73        self.joint_config_constraint = self.acceleration_problem.constraint_set.add_joint_configuration("JointConfigurationConstraint")
 74        self.joint_vel_constraint = self.acceleration_problem.constraint_set.add_joint_velocity("JointVelocityConstraint")
 75        
 76        # Store verbose flag
 77        self.verbose = verbose
 78
 79        # Robot state object
 80        self.robot_state = qontrol.RobotState()
 81        self.robot_state.resize(self.ndof)
 82        
 83        # Get tip frame name
 84        self.tip_frame = self.qontrol_model.get_tip_frame_name()
 85        print(f"Controlling frame: {self.tip_frame}")
 86        
 87        # Get initial joint state and update model
 88        qpos, qvel = self.sim.get_joint_state()
 89        self.robot_state.joint_position = qpos
 90        self.robot_state.joint_velocity = qvel
 91        self.qontrol_model.set_robot_state(self.robot_state)
 92        
 93        print("dt:", self.sim.dt)
 94       
 95        # =================================================================
 96        # Trajectory Setup
 97        # =================================================================
 98        trajectory_path = os.path.join(resource_path, "trajectory.csv")
 99        if not os.path.exists(trajectory_path):
100            raise FileNotFoundError(f"Trajectory file not found: {trajectory_path}")
101        
102        print(f"Loading trajectory from: {trajectory_path}")
103        self.trajectory = TrajectoryGeneration(trajectory_path, self.sim.dt)
104        
105        # =================================================================
106        # Control Parameters
107        # =================================================================
108        # PD gains for Cartesian space
109        self.p_gains = np.array([1000.0, 1000.0, 1000.0, 1000.0, 1000.0, 1000.0])
110        self.d_gains = 2.0 * np.sqrt(self.p_gains)  # Critical damping
111        
112        # Statistics tracking
113        self.update_times = []
114        self.iteration_count = 0
115        self.last_stats_time = time.perf_counter()
116    
117    def update(self):
118        """Main control update: compute and apply joint velocity commands"""
119        start_time = time.perf_counter()
120        
121        # 1. Update Qontrol model with current state from simulation
122        qpos, qvel = self.sim.get_joint_state()
123        self.robot_state.joint_position = qpos
124        self.robot_state.joint_velocity = qvel
125        self.qontrol_model.set_robot_state(self.robot_state)
126        
127        # 2. Update trajectory to get next target
128        self.trajectory.update()
129        target_pose = self.trajectory.get_pose()
130        target_velocity = self.trajectory.get_velocity()
131        target_acceleration = self.trajectory.get_acceleration()
132        
133        # 3. Get current end-effector pose
134        current_ee_matrix = self.qontrol_model.get_frame_pose(self.tip_frame)
135        current_se3 = pin.SE3(current_ee_matrix[:3, :3], current_ee_matrix[:3, 3])
136        
137        # 4. Compute SE(3) log6 error
138        cartesian_error = compute_log6_error(target_pose, current_se3)
139        error_norm = np.linalg.norm(cartesian_error)
140        
141        # 5. Get current end-effector velocity
142        ee_velocity = self.qontrol_model.get_frame_velocity(self.tip_frame)
143        
144        # 6. Compute desired Cartesian acceleration (PD control + feedforward)
145        # xdd_star = Kp * err + Kd * (v_des - v_curr) + a_des
146        desired_acceleration = (self.p_gains * cartesian_error + 
147                              self.d_gains * (target_velocity - ee_velocity) + 
148                              target_acceleration)
149        
150        # Set main task target
151        self.main_task.set_target_acceleration(desired_acceleration)
152        
153        # 7. Solve QP for joint velocities
154        t_qp_start = time.perf_counter()
155        self.acceleration_problem.update(self.sim.dt)
156        t_qp = time.perf_counter() - t_qp_start
157        
158        # 8. Apply velocities directly to simulation state (bypass actuators)
159        solution_found = self.acceleration_problem.solution_found()
160        if solution_found:
161            joint_velocities = self.acceleration_problem.get_joint_velocity_command()
162            velocity_norm = np.linalg.norm(joint_velocities)
163            self.sim.set_joint_velocities(joint_velocities)
164        else:
165            print("Warning: No QP solution found!")
166            velocity_norm = 0.0
167            self.sim.set_joint_velocities(np.zeros(self.ndof))
168        
169        # Statistics
170        update_time = time.perf_counter() - start_time
171        self.update_times.append(update_time)
172        self.iteration_count += 1
173        
174        # Verbose diagnostics
175        if self.verbose and self.iteration_count % 100 == 0:
176            print(f"\nIteration {self.iteration_count} (t={self.trajectory.time:.2f}s):")
177            print(f"  QP solve:       {t_qp*1000:.3f} ms")
178            print(f"  Total:          {update_time*1000:.3f} ms")
179            print(f"  Error norm:     {error_norm:.4f}")
180            print(f"  Velocity norm:  {velocity_norm:.2f} rad/s")
181            print(f"  Solution:       {'FOUND' if solution_found else 'NOT FOUND'}")
182            target_pos = target_pose.translation
183            current_pos = current_se3.translation
184            pos_error = np.linalg.norm(target_pos - current_pos)
185            print(f"  Position error: {pos_error*1000:.2f} mm")
186        
187        # Print statistics every second
188        if not self.verbose:
189            current_time = time.perf_counter()
190            if current_time - self.last_stats_time >= 1.0:
191                avg_update_time = np.mean(self.update_times[-1000:])
192                progress = (self.trajectory.time / self.trajectory.duration) * 100
193                print(f"[Acceleration] t={self.trajectory.time:.1f}s/{self.trajectory.duration:.1f}s ({progress:.0f}%) | Compute: {avg_update_time*1000:.2f} ms | QP: {t_qp*1000:.2f} ms | Error: {error_norm:.4f}")
194                self.last_stats_time = current_time
195        
196        # Check if trajectory is finished
197        if self.trajectory.is_finished():
198            return False  # Signal to stop simulation
199        
200        return True  # Continue simulation
201    
202    def print_statistics(self):
203        """Print final statistics"""
204        if self.update_times:
205            avg_time = np.mean(self.update_times)
206            update_rate = 1.0 / avg_time if avg_time > 0 else 0
207            print(f"\nFinal Statistics:")
208            print(f"Average update rate: {update_rate:.1f} Hz")
209            print(f"Average update time: {avg_time*1000:.2f} ms")
210            print(f"Max update time: {np.max(self.update_times)*1000:.2f} ms")
211            print(f"Min update time: {np.min(self.update_times)*1000:.2f} ms")
212            print(f"Total iterations: {self.iteration_count}")
213            print(f"Trajectory duration: {self.trajectory.duration:.2f}s")
214
215
216def main():
217    """Main entry point"""
218    parser = argparse.ArgumentParser(
219        description="Non-interactive acceleration control with trajectory following",
220        formatter_class=argparse.RawDescriptionHelpFormatter,
221        epilog="""
222Examples:
223  %(prog)s panda
224  %(prog)s panda --verbose
225        """
226    )
227    parser.add_argument("robot", type=str, help="Robot name (e.g., panda)")
228    parser.add_argument("--verbose", action="store_true", help="Print detailed diagnostics")
229    
230    args = parser.parse_args()
231    
232    # Get resource path
233    script_dir = os.path.dirname(os.path.abspath(__file__))
234    resource_path = os.path.join(script_dir, "..", "..", "..", "examples", "resources", args.robot)
235    
236    if not os.path.exists(resource_path):
237        print(f"Error: Resource path not found: {resource_path}")
238        return 1
239    
240    # Initialize MuJoCo simulator
241    sim = MujocoSimulator(args.robot, resource_path, "scene.xml")
242    
243    # Initialize controller
244    controller = QontrolAccelerationController(args.robot, resource_path, sim, verbose=args.verbose)
245    
246    print("\n" + "="*60)
247    print("Starting acceleration control with trajectory following")
248    print("Press ESC to exit")
249    print("="*60 + "\n")
250    
251    try:
252        # Run simulation with controller (will stop when trajectory finishes)
253        sim.run(controller.update)
254    except KeyboardInterrupt:
255        print("\nInterrupted by user")
256    except Exception as e:
257        print(f"\nError during simulation: {e}")
258        import traceback
259        traceback.print_exc()
260        return 1
261    finally:
262        # Print final statistics
263        controller.print_statistics()
264    
265    return 0
266
267
268if __name__ == "__main__":
269    exit(main())

Explanation of the code

Declaration

First we declare all the objects that will be used to define our problem.

        self.qontrol_model = qontrol.RobotModel.load_from_file(urdf_path)

We use pinocchio for our model library.

        self.acceleration_problem = qontrol.JointAccelerationProblem(self.qontrol_model, self.solver)

The output of our qp controller is at the acceleration level.

        self.main_task = self.acceleration_problem.task_set.add_cartesian_acceleration("MainTask", 1.0)

The main task is expressed as a Cartesian acceleration task.

Initialization

    def __init__(self, robot_name: str, resource_path: str, mujoco_sim: MujocoSimulator, verbose: bool = False):
        """
        Initialize the Qontrol acceleration controller
        
        Args:
            robot_name: Name of the robot (e.g., 'panda')
            resource_path: Path to robot resources folder
            mujoco_sim: MuJoCo simulator instance
            verbose: Print detailed diagnostics each iteration (default: False)
        """
        self.robot_name = robot_name
        self.resource_path = resource_path
        self.sim = mujoco_sim
        
        # =================================================================
        # Qontrol Setup - Robot Model and QP Problem
        # =================================================================
        urdf_path = os.path.join(resource_path, "robot.urdf")
        if not os.path.exists(urdf_path):
            raise FileNotFoundError(f"URDF file not found: {urdf_path}")
        
        print(f"Loading Qontrol model from: {urdf_path}")
        self.qontrol_model = qontrol.RobotModel.load_from_file(urdf_path)

During initialization we instantiate the model with the robot urdf.

        self.solver = qontrol.create_qpmad_solver()
        
        # Create acceleration-level problem
        self.acceleration_problem = qontrol.JointAccelerationProblem(self.qontrol_model, self.solver)

We initialize the problem by giving it the model. By default, the qpmad library is used.

        self.regularization_task = self.acceleration_problem.task_set.add_joint_acceleration("RegularizationTask", 1e-5)

We then fill the task set of acceleration_problem with the main task and the regularisation task. Each tasks is given a name and a relative weight \(\omega\). This weight can be modified at any time. In this example, the regularisation tasks is defined as a joint acceleration task. Its means that this task will minimize the overall robot joint acceleration.

        self.joint_config_constraint = self.acceleration_problem.constraint_set.add_joint_configuration("JointConfigurationConstraint")
        self.joint_vel_constraint = self.acceleration_problem.constraint_set.add_joint_velocity("JointVelocityConstraint")

We then fill the constraint set of acceleration_problem with the two pre-implemented constraints. Each constraint is given a name. These constraints will automatically be updated during the update of Qontrol.

        self.robot_state = qontrol.RobotState()
        self.robot_state.resize(self.ndof)

We also create the robot state and resize it.

        self.trajectory = TrajectoryGeneration(trajectory_path, self.sim.dt)

We create a simple trajectory that has been precalculated and store in a csv file. This trajectory start at the robot current Cartesian pose and does a translation of (-0.1, -0.1, -0.1) m.

Update

    def update(self):
        """Main control update: compute and apply joint velocity commands"""
        start_time = time.perf_counter()
        
        # 1. Update Qontrol model with current state from simulation
        qpos, qvel = self.sim.get_joint_state()
        self.robot_state.joint_position = qpos
        self.robot_state.joint_velocity = qvel
        self.qontrol_model.set_robot_state(self.robot_state)
        
        # 2. Update trajectory to get next target
        self.trajectory.update()

The update function is called every milliseconds. At the beginning of each update we fill the new robot state according to the simulated robot. We also update the trajectory so that it gives the next Cartesian pose to reach in 1 ms.

        desired_acceleration = (self.p_gains * cartesian_error + 
                              self.d_gains * (target_velocity - ee_velocity) + 
                              target_acceleration)

We then compute the desired Cartesian acceleration using a simple derivate proportionnal controller. Pinocchio is used to compute the error between the desired Cartesian pose and the current Cartesian pose. This is done by the log6 function. The p_gains are the proportionnal gains of the controller and the d_gains are the derivate gains of the controller.

        self.main_task.set_target_acceleration(desired_acceleration)

The desired Cartesian acceleration is then fed to the main task.

        self.acceleration_problem.update(self.sim.dt)
        t_qp = time.perf_counter() - t_qp_start
        
        # 8. Apply velocities directly to simulation state (bypass actuators)
        solution_found = self.acceleration_problem.solution_found()
        if solution_found:
            joint_velocities = self.acceleration_problem.get_joint_velocity_command()
            velocity_norm = np.linalg.norm(joint_velocities)
            self.sim.set_joint_velocities(joint_velocities)

Once we updated the necassary tasks and constraints we can update the whole problem. If a solution to the problem exist we can then get it and send it to the simulated robot.

Main function

def main():
    """Main entry point"""
    parser = argparse.ArgumentParser(
        description="Non-interactive acceleration control with trajectory following",
        formatter_class=argparse.RawDescriptionHelpFormatter,
        epilog="""
Examples:
  %(prog)s panda
  %(prog)s panda --verbose
        """
    )
    parser.add_argument("robot", type=str, help="Robot name (e.g., panda)")
    parser.add_argument("--verbose", action="store_true", help="Print detailed diagnostics")
    
    args = parser.parse_args()
    
    # Get resource path
    script_dir = os.path.dirname(os.path.abspath(__file__))
    resource_path = os.path.join(script_dir, "..", "..", "..", "examples", "resources", args.robot)
    
    if not os.path.exists(resource_path):
        print(f"Error: Resource path not found: {resource_path}")
        return 1
    
    # Initialize MuJoCo simulator
    sim = MujocoSimulator(args.robot, resource_path, "scene.xml")
    
    # Initialize controller
    controller = QontrolAccelerationController(args.robot, resource_path, sim, verbose=args.verbose)
    
    print("\n" + "="*60)
    print("Starting acceleration control with trajectory following")
    print("Press ESC to exit")
    print("="*60 + "\n")
    
    try:
        # Run simulation with controller (will stop when trajectory finishes)
        sim.run(controller.update)

The main function fetches the robot name via argparse and starts the Mujoco simulation.