Hello,
I’m currently testing the behavior of the “Articulation-Body” feature, which uses PhysX SDK 4.0 included in Unity3D.
To investigate the behavior of a revolute joint a simple test case was used. A single beam is connected to a immovable base structure (Cube) by a revolute joint. At a set time (0.1s) the target angle is changed to a new position and the displacement of the endpoint of the beam is tracked.
As a reference a similar test case was build in Matlab/Simulink. As shown in the figure the damping in Unity3D is significantly greater compared to the matlab results.
Before further testing can take place it would require to know the reason where this difference comes from. Is there any simplification in PhysX that could cause this amplified damping effect?
Best regards,
Niklas
Hi Niklas,
There are multiple damping parameters that may be causing this discrepancy. In addition to the joint-drive damping, there is joint friction, and linear and angular rigid-body damping. Can you please verify that all are set to zero? I.e. specifically, please check:
and
and
Best,
Philipp
Hi Philipp,
thank you for the quick reply. The problem was solved yesterday during further testing. I will give a short summary for similar cases:
Only the stiffness and damping of the joint-drive were set to the desired values. Joint friction, linear and angular damping were set to 0.
Assuming that PhysX uses an explicit Euler method with a fixed timestep the set timestep was reduced greatly to improve the results. Using a fixed timestep of 0.0001s produces comparable results to Matlab. The target angle is not set by a smooth trajectory but a sudden jump, which explains the necessity of the small timestep.
Could you please give me some information about the used integration scheme in PhysX?
In the documentation I only found the PGS and TGS solver, respectively.
https://gameworksdocs.nvidia.com/PhysX/4.0/documentation/PhysXGuide/Manual/RigidBodyDynamics.html#temporal-gauss-seidel
Best regards,
Niklas
Hi Niklas,
To clarify my understanding: In the initial plot that you posted, was there still any nonzero friction or linear/angular damping, or is the difference entirely due to the integration-scheme differences between PhysX and Matlab?
The integration scheme for the joint drives is implicit. In particular, the force applied by the drive is constant over a time step, but it is computed from the projected end-of-timestep joint state given the applied force and response of the joint at the beginning of the timestep.
The effective time-step-size will depend on the solver selection (TGS/PGS). You can find more info on TGS and its difference to PGS here: https://matthias-research.github.io/pages/publications/smallsteps.pdf
Best,
Philipp
Hi Philipp,
the following plot shows the different runs in Unity3D compare to matlab. The settings were always set to the values described in the previous reply. The solver type was set to PGS. Only the timestep was altered, thus the assumption of a explicit euler scheme.
Is their any documentation of all the used algorithms and methods to better understand the function of the physics engine? E.g. the Unity3D Robotics Hub states that the Featherstone algorithm is used to compute an articualation-body. Part of my work is to compare the different approaches in PhysX and Matlab.
Best regards,
Niklas
What are you simulating in Matlab? Is it a continuous-time or discrete-time PD controller? If it is a continuous-time controller, are you discretizing the gains before transferring them to the discrete-time Unity simulation?
Hi Philipp,
the revolute joint and the beam are part of a Simscape Multibody Model. The torque of the joint is calculated using a PD-controller and comparing the desired position and velocity to the actual position and velocity. A fixed-step solver with a step size of 1e-4 is used. The desired angle of the beam is then set at 0.1s.
Your last question is unclear to me, which parameter do you refer to as gains?
I did some research on physics based animation and found the following approach by K. Erleben: https://www.researchgate.net/profile/Kenny-Erleben/publication/247181209_Physics-Based_Animation/links/5e1b2ed04585159aa4cb43d8/Physics-Based-Animation.pdf
There the equations of motion are dicretized using a explicit euler scheme. Including the contact conditions and a friction law results in the linear complementarity problem (LCP) given in equation 7.56. This LCP can then be solved using the PGS method. This approach seems similar to the one used in PhysX but I found no source confirming this.
My question stays the same: Is their a documentation about how exactly the mechanical problem are solved in PhysX?
Best regards,
Niklas
Hi Niklas,
In general, we don’t provide detailed documentation on the PhysX solver. While we document how to use our simulators and what physical interpretation parameters such as, for example, a drive stiffness have, we consider the simulation technology that we develop an implementation detail and want to be free to change and improve that transparently to our users as we develop innovations.
As for gains, I am referring to the drive stiffness and damping - in control theory, these are equivalent to a proportional and derivative gain, respectively. How exactly these gains are used to compute a torque that is applied to the joint will determine the closed-loop dynamics that you observe. My suspicion is that the difference you are seeing in the plot is from the implicit formulation of the PhysX joint drives, which will appear more “damped” compared to an explicit approach that you may have in the Matlab simulation.
So hence my question: Is the Matlab drive/controller using an explicit formulation where the torque is computed from the beginning-of-timestep joint position and velocity and then assumed constant over the timestep?
If you do want to figure out the exact updates for PhysX, I think that should be straightforward for the simple example of a beam on a revolute joint without any constraints. Just:
- Compute the drive torque using the joint position and velocity at the end of a timestep
- Assume that the computed drive torque is constant during the time step to compute the velocity update, i.e. velAfter = velBefore + 1/inertia * dt * driveTorque
- Use an implicit euler update on the joint position, i.e. posAfter = posBefore + dt * velAfter
- Solve for velAfter from combining all equations above.
Note that my comment above about the TGS/PGS effective time step is wrong - they both use the full timestep in the equations above.
Best,
Philipp
Hi Philipp,
Thank you for your answer!
The matlab simulation uses Heun’s method and therefore an explicit scheme.
In the documentation of PhysX 4.1 it states that for Articulation “velocity clamps or damping may be required to ensure stability of the simulation at large angular velocities, which would not be required when using the maximal coordinate articulations”. So I either get the damping effect shown in the figure above or a joint error whan switching to a maximal coordinate articulation.
I did some more research regarding the mathematical background of physics engines. The paper that explains the TGS method speaks of a position based approach explained in the following paper: Müller, Matthias, et al. “Position based dynamics.” Journal of Visual Communication and Image Representation 18.2 (2007): 109-118.
But the position based dynamics are only explained for single point dynamics and not rigid bodies.
In addition to the position based aproach other approaches can be found in the literature, e.g. penalty based, constraint based, impulse based. I’m not necessarily interested in which order the constraints are solved or what additional methods are used to improve the simulation. Right now I know that the featherstone algorithm is used to get the variables of the articulation components and that the PGS or TGS is used to solve a (mixed) linear complementarity problem (LCP). The general approach on how that LCP is derived is the missing piece in my research.
Best,
Niklas