Estimating GRFs while states tracking

[Azarang Asadi]

well hi again,

After months of trying, I'm still having a hard time getting moco to work.
Initially I had some data of stroke gait which included GRFs,IKs, and scaled models (gait 2392). I used these for many studies and then I decided to use moco on the same data, but have been facing some challenges when doing track, inverse is fine.

what am I trying to do : find new motion (new joint trajectories and GRF forces) using experimental data. So basically I generate the muscle activations using moco inverse with inputs of motion,GRF and model. Then I make very small changes to the found activations using inverse so that it represents electrical stimulation. then I use stateTracking and ControlTracking to generate the new motion. Here's the problem: when I add the control Tracking task, it never converges. Statetracking alone works if I use torque driven model, but muscle driven no. The only thing is that I have used foot ground contact models. I even used the muscle driven example in 3Dwalking folder in moco resources to see if the simulation works. That example has the grf forces as input. When I run the original code with grf everything seems fine. But when I use contact model, it doesn't work. So now my question is that if I use the code in 3Dwalking with GRF as input, how can I get the new GRFs? Cause track deviates from the base motion right? so that GRF is not same as the original GRF. How can I get this new GRF? Also Here's how I modeled the contacts:

Code: Select all

		<components>
			<SmoothSphereHalfSpaceForce name="Foot_Ground_R1">
				<!--Path to a Component that satisfies the Socket 'sphere' of type ContactSphere (description: The sphere participating in this contact.).-->
				<socket_sphere>/contactgeometryset/right_one</socket_sphere>
				<!--Path to a Component that satisfies the Socket 'half_space' of type ContactHalfSpace (description: The half-space participating in this contact.).-->
				<socket_half_space>/contactgeometryset/ground</socket_half_space>
				<!--The stiffness constant (i.e., plain strain modulus), default is 1 (N/m^2)-->
				<stiffness>1000000</stiffness>
				<!--The dissipation coefficient, default is 0 (s/m).-->
				<dissipation>1.5</dissipation>
				<!--The coefficient of static friction, default is 0.-->
				<static_friction>0.80000000000000004</static_friction>
				<!--The coefficient of dynamic friction, default is 0.-->
				<dynamic_friction>0.80000000000000004</dynamic_friction>
				<!--The coefficient of viscous friction, default is 0.-->
				<viscous_friction>0.0</viscous_friction>
				<!--The transition velocity, default is 0.01 (m/s).-->
				<transition_velocity>0.20000000000000001</transition_velocity>
				<!--The parameter that determines the smoothness of the transition of the tanh used to smooth the Hertz force. The larger the steeper the transition but the worse for optimization, default is 300.-->
				<hertz_smoothing>300</hertz_smoothing>
				<!--The parameter that determines the smoothness of the transition of the tanh used to smooth the Hunt-Crossley force. The larger the steeper the transition but the worse for optimization, default is 50.-->
				<hunt_crossley_smoothing>50</hunt_crossley_smoothing>
			</SmoothSphereHalfSpaceForce>

Code: Select all

				<ContactSphere name="right_one">
					<!--Path to a Component that satisfies the Socket 'frame' of type PhysicalFrame (description: The frame to which this geometry is attached.).-->
					<socket_frame>/bodyset/calcn_r</socket_frame>
					<!--Location of geometry center in the PhysicalFrame.-->
					<location>0.0019011578840796601 -0.021859 -0.00382630379623308</location>
					<!--Orientation of geometry in the PhysicalFrame (body-fixed XYZ Euler angles).-->
					<orientation>1 0 0</orientation>
					<!--Default appearance for this Geometry-->
					<Appearance>
						<!--The color, (red, green, blue), [0, 1], used to display the geometry. -->
						<color>0 1 1</color>
						<!--Visuals applied to surfaces associated with this Appearance.-->
						<SurfaceProperties>
							<!--The representation (1:Points, 2:Wire, 3:Shaded) used to display the object.-->
							<representation>1</representation>
						</SurfaceProperties>
					</Appearance>
					<!--Radius of the sphere (default: 0).-->
					<radius>0.032000000000000001</radius>
				</ContactSphere>

I have six contacts per foot.

[Aaron Fox]

Hi Azarang,

I don't have a specific solution for your problem right now, but thought I would quickly recommend looking at the example2DWalking code and simulations that come with Moco. These are a much simpler model and problem, but the examples progress in a similar sort of way to what you're proposing. There is an initial tracking problem in this example where the states and GRFs are tracked (using contact spheres), followed by a predictive gait problem where the experimental data is taken out and a combination of goals is used to generate a walking gait in the predictive problem. In newer versions of Moco I think there are even some other predictive simulation approaches (e.g. walking asymmetry). This example should therefore show you how to move from a tracking to a predictive problem, and once you have this you can start experimenting with changing the predictive problem (e.g. adding a controller to simulate electrical stimulation).

Aaron

[Azarang Asadi]

Hi Aaron,

Thanks for your reply. So I've looked at that 2D problem and tried to expand it to 3D. But still it takes years to solve and I just stop the optimization. That's why I started with 3D example in resources folder. My question is this: In the muscle driven states tracking of moco resources, GRFs are provided to the model. Now as moco track can deviate from the motion data as opposed to the inverse that kinematics are prescribed, then GRFs are changes as well after states tracking example, right? IF so, how can I estimate the new GRFs using the states tracking solution?

[Aaron Fox]

Hi Azarang,

My understanding is that if you provide the external loads in the problem, the solution then needs to be consistent with those applied loads (i.e. dynamically consistent). The kinematics can change as there are likely some inconsistencies between experimentally measured kinematics and GRFs (i.e. this is a similar process to what happens in the RRA tool), or the solution can be made dynamically consistent via residual actuators acting on the free body of the model (typically the pelvis). I think residual pelvis actuators are included in the 3D examples but their use heavily penalised.

If you want to estimate the GRFs from a tracking or predictive solution, then you need to take away the experimental external loads and provide a method in the model that estimates the model-ground interaction. This is typically done by foot-ground contact spheres in the Moco examples.

Aaron

[Azarang Asadi]

Thanks Aaron for the response.
That's what I tried to do, foot-contact model but just in 3D example, if you add six contacts per foot ( I used the paper mentioned in moco paper for this) and remove the GRFs, it takes forever to solve and that then using createExternalLoad to get GRFs, they are very weird and different. Is there any resources out there who have used contact 3D models with mocoStateTracking and got reasonable solutions?

[Ross Miller]

Hi Azarang,

Apologies if I'm missing some details here, I wasn't sure if I fully followed what your workflow is:

If you are tracking kinematic data using MocoTrack or MocoStateTrackingGoal, the GRF produced by contact elements on the feet will not typically resemble experimentally measured GRF very closely. You can use MocoContactTrackingGoal to ask the model to also track target GRF data. Here the optimizer will try to minimize the difference between the GRF produced by contact elements vs. the GRF specified as targets in the problem. Is that helpful?

Large 3-D models with muscles and non-rigid tendons take a long time to solve in Moco in my experience. We run those kinds of simulations here on a brand-new Mac Studio computer and they can take 2-3 days sometimes. Are your simulations just taking a long time to convergence but appear to be producing reasonable results (e.g. the cost function score reported by iterations in IPOPT is a reasonable value, constraint violations are small, log(mu) is negative), or are the simulations failing to converge e.g. reporting errors?

We have some 3-D models with contact and related data posted here that might help if you need a starting point that's more detailed than example2Dwalk. I can't remember which of these packages includes codes for actually running the simulations but the one dated November 12, 2020 should have all that I think:

https://simtk.org/projects/umocod

Ross

[Nicholas Bianco]

Hi Azarang,

To add to Ross and Aaron's excellent tips, I will say that when using foot-ground contact models (whether tracking data or not), creating a good initial guess where the contact spheres produce somewhat realistic forces is absolutely essential. You should modify the model and problem in any way you see fit to produce this initial guess before running your actual problem.

My recommendations:

• Even if you aim to predict contact forces, start with MocoContactTrackingGoal if possible.

• Remove muscles from the model and apply very strong actuators to all joints.

• Modify the problem so that the primary objective term is to track GRFs. You may still need light cost weights on effort and kinematic tracking terms.

• Adjust tolerances or mesh interval sizes to speed up the problem. You just need a reasonably good solution with realistic contact forces.

• Finally, you may need to iterate on the previous two bullet points to get an initial guess that is suitable for your final problem. If you are creating muscle-driven simulations, you will need to insert muscle state initial guesses into the initial guess from the previous steps.

When creating tracking simulations using Moco for my last study, I set up my code to make the above workflow possible. You can find that code here: https://github.com/stanfordnmbl/balance ... problem.py.

Best,
Nick

[Azarang Asadi]

Hi all,

thank you so much for your comments. I used the code you provided Nick, and the comments you gave Ross and I ended up getting the simulation working with a muscle driven model, state tracking, contact tracking and control tracking. For control references, I used only one single muscle activation (i.e. soleus) provided from Moco Inverse solution with a little changes to the magnitude of the signal in terminal stance (so I can have better push off force, that is what I'm trying to simulate). Now after two days running still it's running.The tracking weights I have used are :
controlEffortWeight = 0.001;
stateTrackingWeight = 1;
GRFTrackingWeight = 1;
controlWeight = 5;

and also I'm attaching the first iteration and the last ones I have for now. I'd appreciate any comment or insight that might help with the running time.

[Nicholas Bianco]

Hi Azarang,

I would recommend modifying the cost term weights such that the unscaled objective value (the "objective" column in the Ipopt output) is roughly in the range [0.1, 10], as that tends to lead to faster convergence. Right now your objective values are on the order of 1e3.

Nothing else jumps out as being problematic just from the output. You can always set looser tolerances when creating initial guesses to get intermediate solutions that you can investigate and make adjustments for.

Best,
Nick