Fit Params to a run

In this tutorial we find the mass, the stiffness and the damping of a robot from a log of its motion, as one does when a model has to match a real robot.

The run

The robot is a mass on a spring and a damper that slides along \(x\). A controller moves it with a virtual spring to a goal that wanders, which gives motions at several speeds, and we log the run. This robot is the truth, and the fit does not see its numbers:

import casadi as ca
import numpy as np
import matplotlib.pyplot as plt
import virtualmodelcontrol as vmc


def slider(mass, stiffness, damping):
    line = vmc.models.JointSpace(1, unit="m")
    robot = vmc.Mechanism("slider", model=line)
    x = robot.joint(0)
    robot.add("m", vmc.Inertance(x, mass))  # [kg]
    robot.add("spring", vmc.LinearSpring(x, stiffness))  # [N/m]
    robot.add("damper", vmc.LinearDamper(x, damping))  # [N·s/m]
    return robot

truth = slider(2.0, 40.0, 1.5)
goal = vmc.Custom(lambda t: 0.1 * ca.sin(2 * t) + 0.05 * ca.sin(9 * t),
                  [vmc.Time()], dim=1, unit="m")
ctrl = vmc.Mechanism("ctrl")
ctrl.add("pull", vmc.LinearSpring(truth.joint(0) - goal, 60.0))
ctrl.add("damp", vmc.LinearDamper(truth.joint(0), 3.0))
controller = vmc.VMCController(
    vmc.compile(vmc.VirtualMechanismSystem(truth, ctrl)))
plant = vmc.sim.ModelPlant(truth, max_step=1e-3)
log = vmc.sim.run(plant, controller, vmc.sim.SimClock(1 / 200), T=10.0)

The fit

We start from a model that has the right structure and the wrong numbers, and name the Params to find. A Param is named after its component and its own name: the mass of the Inertance added as m is m.inertance, the stiffness of spring is spring.stiffness (robot.params lists them all; see Parameters). fit_params takes the robot, the names (or globs) of its Params and the runs: each a log from vmc.sim.run, or a dict of arrays with one row per sample, t, q, v and motor_torque, and a too if you have the accelerations. It asks which values make the robot’s own dynamics, at the logged motion, agree with the torques that were sent. It returns the fitted values and their standard errors std, each by name, and the root-mean-square rms of what is left over:

from virtualmodelcontrol.identification import fit_params

guess = slider(1.0, 10.0, 0.2)
names = ["m.inertance", "spring.stiffness",
         "damper.damping"]  # component.Param
fit = fit_params(guess, names, [log], smoothing=11)
for name in names:
    value, std = float(fit.values[name]), float(fit.std[name])
    print(f"{name:18s} {value:8.4f} ± {std:.1e}")
m.inertance          2.0009 ± 3.3e-05
spring.stiffness    39.9983 ± 1.4e-03
damper.damping       1.5191 ± 2.2e-04

The true values are 2 kg, 40 N/m and 1.5 N·s/m. The log holds no accelerations, so the fit gets them by smoothing the velocity and differentiating it, which is where its error comes from (smoothing is the number of samples of the smoothing). A run that has an a of its own, as a simulation can give, is used as it is. After the fit, 0.0027 N is left over per sample. The numbers after ± are the standard errors: how far the values may be off if what is left over is noise. Here it is not noise, it is the smoothing, so the errors are much larger than the ± says (the damping is 1.3 % off, 85 times its standard error). Read the ± to compare how well the run determines one Param against another, not as a bound on the error.

What it asks of the data

A Param shows only in motion that depends on it. The mass shows in the accelerations, the stiffness in the position, the damping in the velocity, so a run that only sits still fits none of them: every Param keeps the value it started from, with a standard error of inf. A Param of a nonlinear component, such as the saturating damper TanhDamper, is fitted the same way, by iterating from the values the robot has, and the fit keeps every Param within its bounds.

To use a result, put its values into the robot: guess.params[name].value = fit.values[name] for each name. The robot then simulates, plans and is controlled with the fitted numbers. The hanging arm fits the stiffness and damping of a real kind of arm from logged runs, and Tuning and Optimizing a virtual mechanism use such a model.