Velocity Qontrol#

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

\[\begin{split}\begin{array}{ccc}\boldsymbol{\dot{q}}^{opt} = & \underset{\boldsymbol{\dot{q}}}{\mathrm{argmin}} & ||J(\boldsymbol{q})\boldsymbol{\dot{q}} - \boldsymbol{v}^{target} || + \omega || \boldsymbol{\dot{q}} ||^2\\& \textrm{s.t.} & \boldsymbol{\dot{q}^{min}} \leq \boldsymbol{\dot{q}} \leq \boldsymbol{\dot{q}^{max}}. \\ & & \boldsymbol{q}^{min} \leq \boldsymbol{q}(\boldsymbol{\dot{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
 19using namespace Qontrol;
 20
 21class MujocoQontrol : public MujocoSim
 22{
 23public:
 24  //------------------------------------------- simulation -------------------------------------------
 25  std::shared_ptr<Qontrol::Model::RobotModel> model;
 26  std::shared_ptr<Qontrol::JointVelocityProblem> velocity_problem;
 27  std::shared_ptr<Qontrol::Task::CartesianVelocity<Qontrol::ControlOutput::JointVelocity>> main_task;
 28
 29  Qontrol::RobotState robot_state;
 30
 31  TrajectoryGeneration *traj;
 32  std::string resource_path;
 33
 34  void initController() override
 35  {
 36    model =
 37        Model::RobotModel::loadModelFromFile(resource_path + "robot.urdf");
 38
 39    const int ndof = model->getNrOfDegreesOfFreedom();
 40
 41    velocity_problem = std::make_shared<Qontrol::JointVelocityProblem>(model);
 42    main_task = velocity_problem->task_set->add<Task::CartesianVelocity>("MainTask");
 43    auto regularisation_task = velocity_problem->task_set->add<Task::JointVelocity>("RegularisationTask", 1e-5);
 44
 45    auto joint_configuration_constraint = velocity_problem->constraint_set->add<Constraint::JointConfiguration>("JointConfigurationConstraint");
 46    auto joint_velocity_constraint = velocity_problem->constraint_set->add<Constraint::JointVelocity>("JointVelocityConstraint");
 47    joint_configuration_constraint->setHorizon(15);
 48    mju_copy(d->qpos, m->key_qpos, m->nu);
 49    robot_state.joint_position.resize(ndof);
 50    robot_state.joint_velocity.resize(ndof);
 51
 52    traj = new TrajectoryGeneration(resource_path + "trajectory.csv",
 53                                    m->opt.timestep);
 54  }
 55
 56  void updateController() override
 57  {
 58    const int ndof = model->getNrOfDegreesOfFreedom();
 59
 60    for (int i = 0; i < ndof; ++i)
 61    {
 62      robot_state.joint_position[i] = d->qpos[i];
 63      robot_state.joint_velocity[i] = d->qvel[i];
 64    }
 65    model->setRobotState(robot_state);
 66
 67    traj->update();
 68    pinocchio::SE3 traj_pose(traj->pose.matrix());
 69
 70    pinocchio::SE3 current_pose(model->getFramePose(model->getTipFrameName()).matrix());
 71    const pinocchio::SE3 tipMdes = current_pose.actInv(traj_pose);
 72    auto err = pinocchio::log6(tipMdes).toVector();
 73    Eigen::Matrix<double, 6, 1> p_gains;
 74    p_gains << 10, 10, 10, 10, 10, 10;
 75    Eigen::Matrix<double, 6, 1> xd_star = p_gains.cwiseProduct(err);
 76
 77    main_task->setTargetVelocity(xd_star);
 78
 79    velocity_problem->update(m->opt.timestep);
 80
 81    if (velocity_problem->solutionFound())
 82    {
 83      sendJointVelocity(velocity_problem->getJointVelocityCommand());
 84    }
 85  }
 86};
 87
 88int main(int argc, const char **argv)
 89{
 90  MujocoQontrol mujoco_qontrol;
 91  Qontrol::Log::Logger::parseArgv(argc, argv);
 92
 93  mjvCamera cam;
 94  mjv_defaultCamera(&cam);
 95
 96  mjvOption opt;
 97  mjv_defaultOption(&opt);
 98
 99  mjvPerturb pert;
100  mjv_defaultPerturb(&pert);
101
102  // simulate object encapsulates the UI
103  auto sim = std::make_unique<mj::Simulate>(
104      std::make_unique<mj::GlfwAdapter>(),
105      &cam, &opt, &pert, /* is_passive = */ false);
106
107  std::string robot = argv[1];
108  std::string mujoco_scene = "./resources/" + robot + "/scene.xml";
109  mujoco_qontrol.resource_path = "./resources/" + robot + "/";
110
111  // start physics thread
112  std::thread physicsthreadhandle(&MujocoQontrol::PhysicsThread, mujoco_qontrol, sim.get(), mujoco_scene.c_str());
113
114  // start simulation UI loop (blocking call)
115  sim->RenderLoop();
116  physicsthreadhandle.join();
117
118  return 0;
119}

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::JointVelocityProblem> velocity_problem;

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

  std::shared_ptr<Qontrol::Task::CartesianVelocity<Qontrol::ControlOutput::JointVelocity>> main_task;

The main task is expressed as a Cartesian velocity 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::JointVelocityProblem> velocity_problem;

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

  std::shared_ptr<Qontrol::Task::CartesianVelocity<Qontrol::ControlOutput::JointVelocity>> main_task;

  Qontrol::RobotState robot_state;

  TrajectoryGeneration *traj;
  std::string resource_path;

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

    const int ndof = model->getNrOfDegreesOfFreedom();

    velocity_problem = std::make_shared<Qontrol::JointVelocityProblem>(model);
    main_task = velocity_problem->task_set->add<Task::CartesianVelocity>("MainTask");
    auto regularisation_task = velocity_problem->task_set->add<Task::JointVelocity>("RegularisationTask", 1e-5);

We then fill the task set of velocity_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 veloicty task. Its means that this task will minimize the overall robot joint veloicty.

    auto joint_configuration_constraint = velocity_problem->constraint_set->add<Constraint::JointConfiguration>("JointConfigurationConstraint");
    auto joint_velocity_constraint = velocity_problem->constraint_set->add<Constraint::JointVelocity>("JointVelocityConstraint");

We then fill the constraint set of velocity_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_position.resize(ndof);
    robot_state.joint_velocity.resize(ndof);

We also create the robot state and resize it.

    traj = new TrajectoryGeneration(resource_path + "trajectory.csv",

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
  {
    const int ndof = model->getNrOfDegreesOfFreedom();

    for (int i = 0; i < ndof; ++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.

    pinocchio::SE3 traj_pose(traj->pose.matrix());

    pinocchio::SE3 current_pose(model->getFramePose(model->getTipFrameName()).matrix());
    const pinocchio::SE3 tipMdes = current_pose.actInv(traj_pose);
    auto err = pinocchio::log6(tipMdes).toVector();
    Eigen::Matrix<double, 6, 1> p_gains;
    p_gains << 10, 10, 10, 10, 10, 10;
    Eigen::Matrix<double, 6, 1> xd_star = p_gains.cwiseProduct(err);

We then compute the desired Cartesian velocity using a simple 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.

    main_task->setTargetVelocity(xd_star);

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

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

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

Explanation of the code

Note: The following explanation refers to the logic implemented in C++, but the Python API directly mirrors this structure.

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.velocity_problem = qontrol.JointVelocityProblem(self.qontrol_model, self.solver)

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

        self.main_task = self.velocity_problem.task_set.add_cartesian_velocity("MainTask", 1.0)

The main task is expressed as a Cartesian velocity task.

Initialization

During initialization we instantiate the model with the robot urdf.

        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)

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

        self.solver = qontrol.create_qpmad_solver()
        
        # Create velocity-level problem
        self.velocity_problem = qontrol.JointVelocityProblem(self.qontrol_model, self.solver)

We then fill the task set of velocity_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 veloicty task. Its means that this task will minimize the overall robot joint veloicty.

        self.main_task = self.velocity_problem.task_set.add_cartesian_velocity("MainTask", 1.0)
        
        # Add joint velocity regularization (minimize joint velocities)
        self.regularization_task = self.velocity_problem.task_set.add_joint_velocity("RegularizationTask", 1e-5)

We then fill the constraint set of velocity_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.joint_config_constraint = self.velocity_problem.constraint_set.add_joint_configuration("JointConfigurationConstraint")
        self.joint_vel_constraint = self.velocity_problem.constraint_set.add_joint_velocity("JointVelocityConstraint")

We also create the robot state and resize it.

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

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.

        trajectory_path = os.path.join(resource_path, "trajectory.csv")
        if not os.path.exists(trajectory_path):
            raise FileNotFoundError(f"Trajectory file not found: {trajectory_path}")
        
        print(f"Loading trajectory from: {trajectory_path}")
        self.trajectory = TrajectoryGeneration(trajectory_path, self.sim.dt)

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.

    def update(self):
        """Main control update: compute and apply joint velocities"""
        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)

We then compute the desired Cartesian velocity using a simple 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.

        self.trajectory.update()
        target_pose = self.trajectory.get_pose()
        
        # 3. Get current end-effector pose
        current_ee_matrix = self.qontrol_model.get_frame_pose(self.tip_frame)
        current_se3 = pin.SE3(current_ee_matrix[:3, :3], current_ee_matrix[:3, 3])
        
        # 4. Compute SE(3) log6 error
        cartesian_error = compute_log6_error(target_pose, current_se3)
        error_norm = np.linalg.norm(cartesian_error)
        
        # 5. Compute desired Cartesian velocity (P control)
        desired_velocity = self.p_gains * cartesian_error

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

        self.main_task.set_target_velocity(desired_velocity)

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.

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

Main function

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

def main():
    """Main entry point"""
    parser = argparse.ArgumentParser(
        description="Non-interactive velocity 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 = QontrolVelocityController(args.robot, resource_path, sim, verbose=args.verbose)
    
    print("\n" + "="*60)
    print("Starting velocity 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)
    except KeyboardInterrupt:
        print("\nInterrupted by user")
    except Exception as e:
        print(f"\nError during simulation: {e}")
        import traceback
        traceback.print_exc()
        return 1
    finally:
        # Print final statistics
        controller.print_statistics()
    
    return 0