chenlin/Mujoco_WASM
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Add documentation of ellipsoid-based fluid drag model.

Browse Source 
PiperOrigin-RevId: 566935233
Change-Id: I45980ea2def06d5ffd9091c49543b737b864701d
This commit is contained in:
Yuval Tassa
2023-09-20 05:40:57 -07:00
committed by Copybara-Service
parent dfd48dd821
commit 35e7bb354d
13 changed files with 675 additions and 97 deletions
Download Patch File Download Diff File Expand all files Collapse all files
doc/XMLreference.rst
+33 -16
View File
@@ -2838,31 +2838,48 @@ helps clarify the role of bodies and geoms in MuJoCo.
.. _body-geom-fluidshape:
:at:`fluidshape`: :at-val:`[none, ellipsoid], "none"`
"ellipsoid" Activates geom-level stateless fluid interaction model based on an ellipsoidal approximation of the geom
shape. When active, the model based on :ref:`body inertia sizes <gePassive>` is disabled for the parent body.
"ellipsoid" activates the geom-level fluid interaction model based on an ellipsoidal approximation of the geom
shape. When active, the model based on :ref:`body inertia sizes <flInertia>` is disabled for the parent body.
See section on :ref:`ellipsoid-based<flEllipsoid>` fluid interaction model for details.
.. _body-geom-fluidcoef:
:at:`fluidcoef`: :at-val:`real(5), "0.5 0.25 1.5 1.0 1.0"`
Dimensionless coefficients of fluid interaction model, as follows.
See section on :ref:`ellipsoid-based<flEllipsoid>` fluid interaction model for details.
.. table::
.. list-table::
:width: 60%
:align: left
:widths: 1 5 2 1
:header-rows: 1
* - Index
- Description
- Symbol
- Default
* - 0
- Blunt drag coefficient
- :math:`C_{D, \text{blunt}}`
- 0.5
* - 1
- Slender drag coeficient
- :math:`C_{D, \text{slender}}`
- 0.25
* - 2
- Angular drag coeficient
- :math:`C_{D, \text{angular}}`
- 1.5
* - 3
- Kutta lift coeficient
- :math:`C_K`
- 1.0
* - 4
- Magnus lift coeficient
- :math:`C_M`
- 1.0
+--------+-----------------------------+----------+
| Index | Description | Default |
+========+=============================+==========+
| 0 | Blunt drag coefficient. | 0.5 |
+--------+-----------------------------+----------+
| 1 | Slender drag coeficient. | 0.25 |
+--------+-----------------------------+----------+
| 2 | Angular drag coefficient. | 1.5 |
+--------+-----------------------------+----------+
| 3 | Kutta lift coeficient. | 1.0 |
+--------+-----------------------------+----------+
| 4 | Magnus lift coeficient. | 1.0 |
+--------+-----------------------------+----------+
.. _body-geom-user:
doc/computation/fluid.rst
+436
View File
@@ -0,0 +1,436 @@
Fluid forces
============
Proper simulation of fluid dynamics is beyond the scope of MuJoCo, and would be too slow for the applications we aim to
facilitate. Nevertheless we provide two phenomenological models which are sufficient for simulating behaviors
such as flying and swimming. These models are *stateless*, in the sense that no additional states are assigned to the
surrounding fluid, yet are able to capture the salient features of rigid bodies moving through a fluid medium.
Both models are enabled by setting the :ref:`density<option-density>` and :ref:`viscosity<option-viscosity>` attributes
to positive values. These parameters correspond to the density :math:`\rho` and viscosity :math:`\beta` of the medium.
1. The :ref:`Inertia-based model<flInertia>`, uses only viscosity and density, inferring geometry from body
equivalent-inertia boxes.
2. The :ref:`Ellipsoid-based model <flEllipsoid>` is more elaborate, using an ellipsoid approximation of geoms.
In addition to the global viscosity and density of the medium, this model exposes 5 tunable parameters per
interacting geom.
.. tip::
As detailed in the :ref:`Numerical Integration<geIntegration>` section, implicit integration significantly improves
simulation stability in the presence of velocity-dependent forces. Both of the fluid-force models described below
exhibit this property, so the ``implicit`` or ``implicitfast`` :ref:`intergrators<option-integrator>` are
recommended when using fluid forces. The required analytic derivatives for both models are fully implemented.
.. _flInertia:
Inertia model
-------------
In this model, the shape of each body, for fluid dynamics purposes, is assumed to be the *equivalent inertia box*,
which can also be visualized. Each forward-facing (relative to the linear velocity) face of the box experiences force
along its normal direction. All faces also experience torque due to the angular velocity; this torque is obtained by
integrating the force resulting from the rotation over the surface area. In this sub-section, let :math:`v` and
:math:`\omega` denote the linear and angular body velocity in the body local frame (aligned with the equivalent
inertia box), and :math:`s` the 3D vector of box sizes. When the contributions from all faces are added, the resulting
force and torque applied to the body by a fluid of density :math:`\rho`, in local body coordinates, have the
:math:`i`-th component
.. math::
\begin{aligned}
\text{density force} : \quad &- {1 \over 2} \rho s_j s_k |v_i| v_i \\
\text{density torque} : \quad &- {1 \over 64} \rho s_i \left(s_j^4 + s_k^4 \right) |\omega_i| \omega_i \\
\end{aligned}
This model implicitly assumes high Reynolds numbers, with lift-to-drag ratio equal to the tangent of the angle of
attack. One can also specify a non-zero :ref:`wind<option-wind>`, which is a 3D vector subtracted from the body linear
velocity in the fluid dynamics computation.
Each body also experiences a force and a torque proportional to the viscosity :math:`\beta` and opposite to its linear and
angular velocity. Note that viscosity can be used independent of density, to make the simulation more damped. We use the
formulas for a sphere at low Reynolds numbers, with diameter :math:`d` equal to the average of the equivalent inertia
box sizes. The resulting 3D force and torque in local body coordinates are
.. math::
\begin{aligned}
\text{viscosity force} : \quad &- 3 \beta \pi d v \\
\text{viscosity torque} : \quad &- \beta \pi d^3 \omega \\
\end{aligned}
.. _flEllipsoid:
Ellipsoid model
---------------
.. cssclass:: caption-small
.. figure:: ../images/computation/fruitfly.png
:figwidth: 50%
:align: right
The flight-capable Drosophila Melanogaster model in this figure will be described in a
forthcoming publication.
In this section we describe and derive a stateless model of the forces exerted onto a moving rigid body by the
surrounding fluid, based on an ellipsoidal approximation of geom shape. This model provides finer-grained control of the
different types of fluid forces than the inertia-based model of the previous section. The motivating use-case for this
model is insect flight, see figure on the right.
Summary
~~~~~~~
The model is activated per-geom by setting the :ref:`fluidshape<body-geom-fluidshape>` attribute to ``ellipsoid``, which
also disables the inertia-based model for the parent body. The
5 numbers in the :ref:`fluidcoef<body-geom-fluidcoef>` attribute correspond to the following semantics
.. list-table::
:width: 60%
:align: left
:widths: 1 5 2 1
:header-rows: 1
* - Index
- Description
- Symbol
- Default
* - 0
- Blunt drag coefficient
- :math:`C_{D, \text{blunt}}`
- 0.5
* - 1
- Slender drag coeficient
- :math:`C_{D, \text{slender}}`
- 0.25
* - 2
- Angular drag coeficient
- :math:`C_{D, \text{angular}}`
- 1.5
* - 3
- Kutta lift coeficient
- :math:`C_K`
- 1.0
* - 4
- Magnus lift coeficient
- :math:`C_M`
- 1.0
Elements of the model are a generalization of :cite:t:`andersen2005b` to 3 dimensions.
The force :math:`\mathbf{f}_{\text{fluid}\rightarrow \text{solid}}` and torque
:math:`\mathbf{g}_{\text{fluid} \rightarrow \text{solid}}` exerted by the fluid onto the solid are
the sum of of the terms
.. math::
\begin{align*}
\mathbf{f}_{\text{fluid} \rightarrow \text{solid}} &= \mathbf{f}_A + \mathbf{f}_D + \mathbf{f}_M + \mathbf{f}_K \\
\mathbf{g}_{\text{fluid} \rightarrow \text{solid}} &= \mathbf{g}_A + \mathbf{g}_D
\end{align*}
Where subscripts :math:`A`, :math:`D`, :math:`M` and :math:`K`, denote Added mass, viscous Drag, Magnus lift and
Kutta lift, respectively. The :math:`D`, :math:`M` and :math:`K` terms are scaled by the respective
:math:`C_D`, :math:`C_M` and :math:`C_K` coefficients above, while the added mass term cannot be scaled.
Notation
~~~~~~~~
We describe the motion of the object in an inviscid, incompressible quiescent fluid of density :math:`\rho`. The
arbitrarily-shaped object is described in the model as the equivalent ellipsoid of semi-axes
:math:`\mathbf{d} = \{d_x, d_y, d_z\}`.
The problem is described in a reference frame aligned with the sides of the ellipsoid and moving with it. The
body has velocity :math:`\mathbf{v} = \{v_x, v_y, v_z\}` and angular velocity
:math:`\boldsymbol{\omega} = \{\omega_x, \omega_y, \omega_z\}`. We will also use
.. math::
\begin{align*}
d_\text{max} &= \max(d_x, d_y, d_z) \\
d_\text{min} &= \min(d_x, d_y, d_z) \\
d_\text{mid} &= d_x + d_y + d_z - d_\text{max} - d_\text{min}
\end{align*}
The Reynolds number is the ratio between inertial and viscous forces within a flow and is defined as :math:`Re=u~l/\beta`, where
:math:`\beta` is the kinematic viscosity of the fluid, :math:`u` is the characteristic speed of the flow (or, by change of frame, the
speed of the body), and :math:`l` is a characteristic size of the flow or the body.
We will use :math:`\Gamma` to denote circulation, which is the line integral of the velocity field around a closed curve
:math:`\Gamma = \oint \mathbf{v} \cdot \textrm{d} \mathbf{l}` and, due to Stokes' Theorem,
:math:`\Gamma = \int_S \nabla \times \mathbf{v} \cdot \textrm{d}\mathbf{s}`.
In fluid dynamics notation the symbol :math:`\boldsymbol{\omega}` is often used for the
vorticity, defined as :math:`\nabla \times \mathbf{v}`, rather than the angular velocity. For a rigid-body motion, the
vorticity is twice the angular velocity.
Finally, we use the subscripts :math:`i, j, k` to denote triplets of equations that apply symmetrically to
:math:`x, y, z`. For example :math:`a_i = b_j + b_k` is shorthand for the 3 equations
.. math::
\begin{align*}
a_x &= b_y + b_z \\
a_y &= b_x + b_z \\
a_z &= b_x + b_y
\end{align*}
.. _flProjection:
Ellipsoid projection
~~~~~~~~~~~~~~~~~~~~
We provide the following result without proof. For the derivation, contact the development team.
.. admonition:: Lemma
:class: note
Given an ellipsoid at the origin with semi-axes :math:`(d_x, d_y, d_z)` aligned
with the coordinate axes :math:`(x, y, z)`, and a unit vector :math:`\mathbf{u} = (u_x, u_y, u_z)`,
the area projected by the ellipsoid onto the plane normal to :math:`\mathbf{u}` is
.. math::
A^{\mathrm{proj}}_{\mathbf{u}} = \pi \sqrt{\frac{d_y^4 d_z^4 u_x^2 + d_x^4 d_z^4 u_y^2 + d_x^4 d_y^4 u_z^2}{d_y^2 d_z^2 u_x^2 + d_x^2 d_z^2 u_y^2 + d_x^2 d_y^2 u_z^2}}
Added mass
~~~~~~~~~~
For a body moving in a fluid, added mass or virtual mass measures the inertia of the fluid that is moved due to the
body's motion. It can be derived from potential flow theory (i.e. it is present also for inviscid flows).
Following Chapter 5 of :cite:t:`lamb1932`, the forces :math:`\mathbf{f}_{V}` and torques :math:`\mathbf{g}_{V}` exerted
onto a moving body due to generation of motion in the fluid from rest can be written as:
.. math::
\begin{align*}
\mathbf{f}_{A} &= - \frac{\textrm{d}}{\textrm{d} t} \nabla_{\mathbf{v}} \mathcal{T} + \nabla_{\mathbf{v}} \mathcal{T} \times \boldsymbol{\omega} \\
\mathbf{g}_{A} &= - \frac{\textrm{d}}{\textrm{d} t} \nabla_{\boldsymbol{\omega}} \mathcal{T} + \nabla_{\mathbf{v}} \mathcal{T} \times \mathbf{v} + \boldsymbol{\omega} \times \nabla_{\boldsymbol{\omega}} \mathcal{T}
\end{align*}
where :math:`\mathcal{T}` is the kinetic energy of the fluid alone. These forces are often described as added or
virtual mass because they are due to the inertia of the fluid that is to moved or deflected by the accelerating body. In
fact, for a body with constant linear velocity these forces reduce to zero. We consider the body as having three planes
of symmetry because under this assumption the kinetic energy greatly simplifies and can be written as:
.. math::
2 \mathcal{T} = m_{A, x} v_x^2 + m_{A, y} v_y^2 + m_{A, z} v_z^2 +
I_{A, x} \omega_x^2 + I_ {A, y} \omega_y^2 + I_{A, y} \omega_z^2
For convenience we introduce the added-mass vector :math:`\mathbf{m}_A = \{m_{A, x}, m_{A, y}, m_{A, z}\}` and added-moment of
inertia vector :math:`\mathbf{I}_A = \{I_{A, x}, I_{A, y}, I_{A, z}\}`. Each of these quantities should estimate the inertia
of the moved fluid due the motion of the body in the corresponding direction and can be derived from potential flow
theory for some simple geometries.
For a body with three planes of symmetry, we can write in compact form the forces and torques due to added inertia:
.. math::
\begin{align*}
\mathbf{f}_{A} &= - \mathbf{m}_A \circ \dot{\mathbf{v}} + \left(\mathbf{m}_A \circ \mathbf{v} \right) \times \boldsymbol{\omega} \\
\mathbf{g}_{A} &= - \mathbf{I}_A \circ \dot{\boldsymbol{\omega}} + \left(\mathbf{m}_A \circ \mathbf{v} \right) \times \mathbf{v} + \left(\mathbf{I}_A \circ \boldsymbol{\omega} \right) \times \boldsymbol{\omega}
\end{align*}
Here :math:`\circ` denotes an element-wise product, :math:`\dot{\mathbf{v}}` is the linear acceleration and
:math:`\dot{\boldsymbol{\omega}}` is the angular acceleration. :math:`\mathbf{m}_A \circ \mathbf{v}` and
:math:`\mathbf{I}_A \circ \boldsymbol{\omega}` are the virtual linear and angular momentum respectively.
For an ellipsoid of semi-axis :math:`\mathbf{d} = \{d_x, d_y, d_z\}` and volume :math:`V = 4 \pi d_x d_y d_z / 3`, the
virtual inertia coefficients were derived by :cite:t:`tuckerman1925`. Let:
.. math::
\kappa_i = \int_0^\infty \frac{d_i d_j d_k}{\sqrt{(d_i^2 + \lambda)^3 (d_j^2 + \lambda) (d_k^2 + \lambda)}} \textrm{d} \lambda
It should be noted that these coefficients are non-dimensional (i.e. if all semi-axes are multiplied by the same scalar
the coefficients remain the same). The virtual masses of the ellipsoid are:
.. math::
m_{A, i} = \rho V \frac{\kappa_i}{2 - \kappa_i}
And the virtual moments of inertia are:
.. math::
I_{A, i} = \frac{\rho V}{5} \frac{(d_j^2 - d_k^2)^2 (\kappa_k-\kappa_j)}{2(d_j^2 - d_k^2) + (d_j^2 + d_k^2) (\kappa_j-\kappa_k)}
Viscous drag
~~~~~~~~~~~~
The drag force acts to oppose the motion of the body relative to the surrounding flow. We found that viscous forces
serve also to reduce the stiffness of the equations of motion extended with the fluid dynamic terms. For this reason, we
opted to err on the conservative side and chose approximations of the viscous terms that may overestimate dissipation.
Despite being ultimately caused by viscous dissipation, for high Reynolds numbers the drag is independent of the
viscosity and scales with the second power of the velocity. It can be written as:
.. math::
\begin{align*}
\mathbf{f}_\text{D} = - C_D~\rho~ A_D ~ \|\mathbf{v}\|~ \mathbf{v}\\
\mathbf{g}_\text{D} = - C_D \rho~ I_D ~ \|\boldsymbol{\omega}\| ~ \boldsymbol{\omega}
\end{align*}
Where :math:`C_D` is a drag coefficient, and :math:`A_D` is a reference surface area (e.g. a measure of the projected
area on the plane normal to the flow), and :math:`I_D` a reference moment of inertia.
Even for simple shapes, the terms :math:`C_D`, :math:`A_D` and :math:`I_D` need to be tuned to the problem-specific
physics and dynamical scales :cite:p:`duan2015`. For example, the drag coefficient :math:`C_D` generally decreases with
increasing Reynolds numbers, and a single reference area :math:`A_D` may not be sufficient to account for the skin
drag for highly irregular or slender bodies. For example, experimental fits are derived from problems ranging from
falling playing cards :cite:p:`wang2004,andersen2005a,andersen2005b` to particle transport :cite:p:`loth2008,
bagheri2016`.
We derive a formula for :math:`\mathbf{f}_\text{D}` based on two surfaces :math:`A^\text{proj}_\mathbf{v}` and
:math:`A_\text{max}`. The first, :math:`A^\text{proj}_\mathbf{v}`, is the cylindrical projection of the body onto a
plane normal to the velocity :math:`\mathbf{v}`. The second is the maximum projected surface
:math:`A_\text{max} = 4 \pi d_{max} d_{min}`.
.. math::
\mathbf{f}_\text{D} = - \rho~ \big[ C_{D, \text{blunt}} ~ A^\text{proj}_\mathbf{v} ~ +
C_{D, \text{slender}}\left(A_\text{max} - A^\text{proj}_\mathbf{v} \right) \big] ~ \|\mathbf{v}\|~ \mathbf{v}
The formula and derivation for :math:`A^\text{proj}_\mathbf{v}` is given in the :ref:`lemma<flProjection>` above.
We propose an analogous model for the angular drag. For each Cartesian axis we consider the moment of inertia of the
maximum swept ellipsoid obtained by the rotation of the body around the axis. The resulting diagonal entries of the
moment of inertia are:
.. math::
\mathbf{I}_{D,ii} = \frac{8\pi}{15} ~d_i ~\max(d_j, ~d_k)^4 .
Given this reference moment of inertia, the angular drag torque is computed as:
.. math::
\mathbf{g}_\text{D} = - \rho ~ \boldsymbol{\omega} ~ \Big( \big[ C_{D, \text{angular}} ~ \mathbf{I}_D ~ +
C_{D, \text{slender}} \left(\mathbf{I}_\text{max} - \mathbf{I}_D \right) \big] \cdot \boldsymbol{\omega} \Big)
Here :math:`\mathbf{I}_\text{max}` is a vector with each entry equal to the maximal component of :math:`\mathbf{I}_D`.
The viscosity :math:`\beta`
For Reynolds numbers around or below :math:`O(10)`, the drag is best approximated as linear in the flow velocity
(e.g. Stokes' law). For example, for a sphere the drag force :cite:p:`stokes1850` and torque :cite:p:`lamb1932` are:
.. math::
\begin{align*}
\mathbf{f}_\text{S} &= - 6 \pi r_D \rho ~ \beta \mathbf{v}\\
\mathbf{g}_\text{S} &= - 8 \pi r_D^3 \rho ~ \beta \boldsymbol{\omega}
\end{align*}
Here, :math:`r_D` is the radius of the sphere and :math:`\beta` is the kinematic viscosity of the medium (e.g.
:math:`1.48~\times 10^{-5}~m^2/s` for ambient-temperature air and :math:`0.89 \times 10^{-4}~m^2/s` for water). Here,
for simplicity, we estimate the radius of the equivalent sphere as :math:`r_D = (d_x + d_y + d_z)/3`. To make a
quantitative example, Stokes' law become accurate for room-temperature air if
:math:`u\cdot l \lesssim 2 \times 10^{-4}~m^2/s`, where :math:`u` is the speed and :math:`l` a characteristic length of
the body.
Viscous lift
~~~~~~~~~~~~
The Kutta-Joukowski theorem calculates the lift :math:`L` of a two-dimensional body translating in a uniform flow with
speed `u` as :math:`L = \rho u \Gamma`. Here, :math:`\Gamma` is the circulation around the body. In the next
subsections we define two sources of circulation and the resulting lift forces.
Magnus force
^^^^^^^^^^^^
.. cssclass:: caption-small
.. figure:: ../images/computation/magnus.png
:figwidth: 45%
:align: right
Smoke flow visualization of the flow past a rotating cylinder (WikiMedia Commons, CC BY-SA 4.0). Due to viscosity,
the rotating cylinder deflects the incoming flow upward and receives a downwards force (red arrow).
The Magnus effect describes the motion of a rotating object moving through a fluid. Through viscous effects, a spinning
object induces rotation in the surrounding fluid. This rotation deflects the trajectory of the fluid past the object
(i.e. it causes linear acceleration), and the object receives an equal an opposite reaction. For a cylinder, the Magnus
force per unit length of the cylinder can be computed as :math:`F_\text{M} / L = \rho v \Gamma`, where :math:`\Gamma`
is the circulation of the flow caused by the rotation and :math:`v` the velocity of the object. We estimate this force
for an arbitrary body as:
.. math::
\mathbf{f}_{\text{M}} = C_M ~\rho~ V~ \boldsymbol{\omega}\times\mathbf{v} ,
where :math:`V` is the volume of the body and :math:`C_M` is a coefficient for the force, typically set to 1.
It's worth making an example. To reduce the number of variables, suppose a body rotating in only one direction, e.g.
:math:`\boldsymbol{\omega} = \{0, 0, \omega_z\}`, translating along the other two, e.g. :math:`\mathbf{v} = \{v_x, v_y, 0\}`. The
sum of the force due to added mass and the force due to the Magnus effect along, for example, :math:`x` is:
.. math::
\frac{f}{\pi \rho d_z} = v_y \omega_z \left(2 d_x \min\{d_x, d_z\} - (d_x + d_z)^2\right)
Note that the two terms have opposite signs.
Kutta condition
^^^^^^^^^^^^^^^
A stagnation point is a location in the flow field where the velocity is zero. For a body moving in a flow (in 2D, in
the frame moving with the body) there are two stagnation points: in the front, where the stream-lines separate to either
sides of the body, and in the rear, where they reconnect. A moving body with a sharp trailing (rear) edge will generate
in the surrounding flow a circulation of sufficient strength to hold the rear stagnation point at the trailing edge.
This is the Kutta condition, a fluid dynamic phenomenon that can be observed for solid bodies with sharp corners, such
as slender bodies or the trailing edges of airfoils.
.. cssclass:: caption-small
.. figure:: ../images/computation/kutta_cond_plate.svg
:figwidth: 95%
:align: left
Sketch of the Kutta condition. Blue lines are streamlines and the two magenta points are the stagnation points. The
dividing streamline, which connects the two stagnation points, is marked in green. The dividing streamline and the
body inscribe an area where the flow is said to be "separated" and recirculates within. This circulation produces an
upward force acting on the plate.
For a two-dimensional flow sketched in the figure above, the circulation due to the Kutta condition can be estimated as:
:math:`\Gamma_\text{K} = C_K ~ d_x ~ \| \mathbf{v}\| ~ \sin(2\alpha)`,
where :math:`C_K` is a lift coefficient, and :math:`\alpha` is the angle between the velocity vector and its projection
onto the surface. The lift force per unit length can be computed with the Kutta–Joukowski theorem as
:math:`\mathbf{f}_K / L = \rho \Gamma_\text{K} \times \mathbf{v}`.
In order to extend the lift force equation to three-dimensional motions, we consider the normal
:math:`\mathbf{n}_{s, \mathbf{v}} = \{\frac{d_y d_z}{d_x}v_x, \frac{d_z d_x}{d_y}v_y, \frac{d_x d_x}{d_z}v_z\}`
to the cross-section of the body which generates the body's projection :math:`A^\text{proj}_\mathbf{v}` onto a plane
normal to the velocity given in the :ref:`lemma<flProjection>` above and the corresponding unit vector
:math:`\hat{\mathbf{n}}_{s, \mathbf{v}}`.
We use this direction to decompose :math:`\mathbf{v} = \mathbf{v}_\parallel ~+~ \mathbf{v}_\perp` with
:math:`\mathbf{v}_\perp = \left(\mathbf{v} \cdot \hat{\mathbf{n}}_{s, \mathbf{v}}\right) \hat{\mathbf{n}}_{s, \mathbf{v}}`.
We write the lift force as:
.. math::
\begin{align*}
\mathbf{f}_\text{K} &= \frac{C_K~\rho~ A^\text{proj}_\mathbf{v}}{\|\mathbf{v}\|}
\left( \mathbf{v} \times \mathbf{v}_\parallel\right)\times \mathbf{v} \\
&= C_K~\rho~ A^\text{proj}_\mathbf{v} \left(\hat{\mathbf{v}} \cdot \hat{\mathbf{n}}_{s, \mathbf{v}}\right)
\left( \hat{\mathbf{n}}_{s, \mathbf{v}} \times \mathbf{v} \right)\times \mathbf{v}
\end{align*}
Here, :math:`\hat{\mathbf{v}}` is the unit-normal along :math:`\mathbf{v}`. Note that the direction of :math:`\hat{\mathbf{n}}_{s,
\mathbf{v}}` differs from :math:`\hat{\mathbf{v}}` only on the planes where the semi-axes of the body are unequal. So for
example, for spherical bodies :math:`\hat{\mathbf{n}}_{s, \mathbf{v}} \equiv \hat{\mathbf{v}}` and by construction
:math:`\mathbf{f}_\text{K} = 0`.
Let's unpack the relation with an example. Suppose a body with :math:`d_x = d_y` and :math:`d_z \ll d_x`. Note that the vector
:math:`\hat{\mathbf{n}}_{s, \mathbf{v}} \times \hat{\mathbf{v}}` gives the direction of the circulation induced by the
deflection of the flow by the solid body. Along :math:`z`, the circulation will be proportional to :math:`\frac{d_y d_z}{d_x}v_x v_y
- \frac{d_z d_x}{d_y}v_x v_y = 0` (due to :math:`d_x = d_y`). Therefore, on the plane where the solid is blunt, the motion
produces no circulation.
Now, for simplicity, let :math:`v_x = 0`. In this case also the circulation along :math:`y`, proportional
to :math:`\frac{d_y d_z}{d_x}v_x v_z - \frac{d_y d_x}{d_y}v_x v_z`, is zero. The only non-zero component of the circulation
will be along :math:`x` and be proportional to :math:`\left(\frac{d_x d_z}{d_y} - \frac{d_x d_y}{d_z}\right) v_y v_z \approx
\frac{d_x^2}{d_z} v_y v_z`.
We would have :math:`\mathbf{v}_\parallel = \{v_x, 0, v_z\}` and
:math:`\Gamma \propto \{d_z v_y v_z, ~ 0,~ - d_x v_x v_y \} / \|\mathbf{v}\|`.
The motion produces no circulation on the plane where the solid is blunt, and on the other two planes
the circulation is
:math:`\Gamma \propto r_\Gamma ~ \|\mathbf{v}\|~ \sin(2 \alpha) ~ = ~2 r_\Gamma ~\|\mathbf{v}\| ~\sin(\alpha)~\cos(\alpha)`
with :math:`\alpha` the angle between the velocity and its projection on the body on the plane (e.g. on the plane
orthogonal to :math:`x` we have :math:`\sin(\alpha) = v_y/\|\mathbf{v}\|` and
:math:`\cos(\alpha) = v_z/\|\mathbf{v}\|`), and :math:`r_\Gamma`, the lift surface on the plane (e.g. :math:`d_z` for
the plane orthogonal to :math:`x`). Furthermore, the direction of the circulation is given by the cross product (because
the solid boundary "rotates" the incoming flow velocity towards its projection on the body).
Acknowledgements
~~~~~~~~~~~~~~~~
The design and implementation of the model in this section are the work of Guido Novati.
References
~~~~~~~~~~
.. bibliography::
doc/computation.rst → doc/computation/index.rst
+35 -66
View File
@@ -1,6 +1,12 @@
Computation
===========
.. toctree::
:hidden:
fluid
Introduction
------------
@@ -389,59 +395,16 @@ by MuJoCo are also passive in the sense of physics, i.e., they do not increase e
callback :ref:`mjcb_passive` and add forces to ``mjData.qfrc_passive`` that may increase energy. This will not interfere
with MuJoCo's operation as long as such user forces depend only on position and velocity.
MuJoCo can compute three types of passive forces: spring-dampers in joints and tendons, gravity compensation forces, and
fluid dynamics.
MuJoCo can compute three types of passive forces:
When Euler or the implicit integator are used, joint damping is integrated implicitly which significantly increases
stability. Thus, even though damping can be modeled as an actuator property, it is better to model it as a joint
property. Note also the XML :ref:`joint <body-joint>` attribute springdamper which automates the creation of mass-
spring-dampers with desired time constants and damping ratios; in that case the compiler computes the stiffness and
damping coefficients of the joint by taking the joint inertia into account.
Gravity compensation is a force applied to a body's center of mass opposing gravity, see :ref:`body gravcomp<body>` for
details.
Proper simulation of fluid dynamics is beyond the scope of MuJoCo, and would be too slow for the applications we aim to
facilitate. Nevertheless we provide a phenomenological model which is sufficient for simulating behaviors such as flying
and swimming. It is enabled by setting ``mjModel.opt.viscosity`` and ``mjModel.opt.density`` to positive values (they
are zero by default). These parameters specify the viscosity :math:`\beta` and density :math:`\rho` of the medium and
apply to all bodies. The shape of each body for fluid dynamics purposes is assumed to be the equivalent inertia box,
which can also be visualized. Each forward-facing (relative to the linear velocity) face of the box experiences force
along its normal direction. All faces also experience torque due to the angular velocity; this torque is obtained by
integrating the force resulting from the rotation over the surface area. In this sub-section only, let :math:`v` and
:math:`\omega` denote the linear and angular body velocity in its local frame (which is aligned with the equivalent
inertia box), and :math:`s` the 3D vector of box sizes. When the contributions from all faces are added, the resulting
force and torque in local body coordinates have :math:`i`-th component
.. math::
\begin{aligned}
\text{density force} : \quad &- {1 \over 2} \rho s_j s_k |v_i| v_i \\
\text{density torque} : \quad &- {1 \over 64} \rho s_i \left(s_j^4 + s_k^4 \right) |\omega_i| \omega_i \\
\end{aligned}
This model implicitly assumes high Reynolds numbers, with lift-to-drag ratio equal to the tangent of the angle of
attack. One can also specify a non-zero ``mjModel.opt.wind``, which is a 3D vector subtracted from the body linear
velocity in the fluid dynamics computation.
Each body also experiences a force and a torque proportional to viscosity and opposite to its linear and angular
velocity. Note that viscosity can be used independent of density, to make the simulation more damped. We use the
formulas for a sphere at low Reynolds numbers, with diameter :math:`d` equal to the average of the equivalent inertia
box sizes. The resulting 3D force and torque in local body coordinates are
.. math::
\begin{aligned}
\text{viscosity force} : \quad &- 3 \beta \pi d v \\
\text{viscosity torque} : \quad &- \beta \pi d^3 \omega \\
\end{aligned}
Existing phenomenological models of fluid dynamics tend to be valid for one regime (e.g. flat vs. spherical objects) and
it is difficult to construct a model which can be simulated efficiently and is broadly valid. For example, drag forces
are known to transition from being linear in velocity at low Reynolds numbers, to being quadratic in velocity at high
Reynolds numbers. However it is not clear how this transition should occur for lift forces. This motivated the above
separation, using density to compute all quadratic forces and viscosity to compute all linear forces. Users who need a
more detailed simulation of fluid dynamics should leave the density and viscosity parameters set to zero so as to
disable the built-in mechanism. A custom model of fluid dynamics (or any other force field that depends only on position
and velocity) can be implemented in the callback :ref:`mjcb_passive`.
- Spring-dampers in joints and tendons. See the following attribues for details.
|br| **Joints:**
:ref:`stiffness<body-joint-stiffness>`, :ref:`springref<body-joint-springref>`,
:ref:`damping<body-joint-damping>`, :ref:`springdamper<body-joint-springdamper>`.
|br| **Tendons:** :ref:`stiffness<tendon-spatial-stiffness>`,
:ref:`springlength<tendon-spatial-springlength>`, :ref:`damping<tendon-spatial-damping>`.
- Gravity compensation forces. See the body :ref:`gravcomp<body-gravcomp>` attribute for details.
- Fluid forces exerted by the surrounding medium. See the :doc:`Fluid forces <fluid>` chapter for details.
.. _geIntegration:
@@ -786,7 +749,6 @@ dimensionality of the constraint residual in each case.
computed by the compiler in the initial model configuration ``mjModel.qpos0``.
``distance`` : 1
.. attention::
Distance equality constraints were removed in MuJoCo version 2.2.2. If you are using an earlier version, please
switch to the corresponding version of the documentation.
@@ -968,7 +930,7 @@ each component of :math:`f`. The basis vectors are spatial vectors: 3D force fol
vectors into the columns of a matrix :math:`E`, the force/torque generated by the constraint force in the contact frame
is :math:`E f`. The matrix of basis vectors is constructed as follows.
.. image:: images/computation/contact_frame.svg
.. image:: ../images/computation/contact_frame.svg
:width: 700px
:align: center
@@ -1330,7 +1292,7 @@ representations of the constraint Jacobian and related matrices.
the computational complexity of one matrix-vector multiplication (although the constants are larger). It has
first-order convergence but nevertheless makes rapid progress in a few iterations.
.. image:: images/computation/gPGS.svg
.. image:: ../images/computation/gPGS.svg
:width: 500px
:align: center
@@ -1435,7 +1397,7 @@ coefficients varying from left to right. Mathematically, the penalty in the pyra
the penalty in the elliptic case contains pieces that are quadratics minus square roots of quadratics - allowing
circular contours around the tip of the cone.
.. image:: images/computation/softcontact.png
.. image:: ../images/computation/softcontact.png
:width: 600px
:align: center
@@ -1466,18 +1428,25 @@ others can be pruned quickly without a detailed check. MuJoCo has flexible mecha
checked in detail. The decision process involves two stages: generation and filtering.
Generation
First we generate a list of candidate geom pairs in two ways:
First we generate a list of candidate geom pairs in one of two ways: "pair" or "dynamic". The user can also specify
"all" which merges both sources (and is the default). This is done via the setting ``mjModel.opt.collision``. "Pair"
refers to an explicit list of geom pairs defined with the :ref:`pair <contact-pair>` element in MJCF. It gives the
user full control, however it is a static mechanism (independent of the spatial arrangement of the geoms at runtime)
and can be tedious for large models. It is normally used to supplement the output of the "dynamic" mechanism. Dynamic
generation works with bodies rather than geoms; when a body pair is included this means that all geoms attached to
one body can collide with all geoms attached to the other body.
- Using predefined geom pairs, specified with the :ref:`contact-pair <contact-pair>` element. This mechanism is
independent of the spatial arrangement of the geoms at runtime.
- Generating candidate geom pairs by examining the runtime locations of bodies. First, a sweep-and-prune broadphase
collision stage prunes bodies that cannot collide. Then, a mid-phase collision stage further prunes geom pairs from
remaining body pairs using a static Bounding Volume Hierarchy (a binary tree) of axis-aligned bounding boxes.
Each body is equipped with an AABB tree of its geoms, aligned with the body inertial or geom frames for all inner
or leaf nodes, respectively.
- The user can also explicitly exclude certain body pairs using the :ref:`exclude <contact-exclude>` element.
The body pairs are generated via broad-phase collision detection based on a modified sweep-and-prune algorithm. The
modification is that the axis for sorting is chosen as the principal eigenvector of the covariance matrix of all geom
centers - which maximizes the spread. Then, for each body pair, a mid-phase collision detection using a static
bounding volume hierarchy (a BVH binary tree) of axis-aligned bounding boxes (AABB) is performed. Each body is
equipped with an AABB tree of its geoms, aligned with the body inertial or geom frames for all inner or leaf nodes,
respectively.
These two lists are then merged.
Finally, the user can explicitly exclude certain body pairs using the :ref:`exclude <contact-exclude>` element
in MJCF. Exclusion is applied when "dynamic" or "all" are selected, but not when "pair" is selected. At the end of
this step we have a list of geoms pairs that is typically much smaller than :math:`n (n-1)/2`, but can still be
pruned further before detailed collision checking.
Filtering
Next we apply four filters to the list generated in the previous step. Filters 1 and 2 are applied to all geom pairs.
doc/conf.py
+4
View File
@@ -42,6 +42,7 @@ master_doc = 'index'
# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
# ones.
extensions = [
'sphinxcontrib.bibtex',
'sphinxcontrib.katex',
'sphinxcontrib.youtube',
'sphinx_copybutton',
@@ -51,6 +52,9 @@ extensions = [
'mujoco_include',
]
# Bibtex references for sphinxcontrib.bibtex
bibtex_bibfiles = ['references.bib']
# Add any paths that contain templates here, relative to this directory.
templates_path = ['templates']
doc/css/theme_overrides.css
+4
View File
@@ -22,6 +22,10 @@ body {
font-size: 93%;
}
.caption-small {
font-size: 85%;
}
.sidebar-brand-text {
display: none;
}
doc/images/computation/fruitfly.png
BIN
View File
Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 1.6 MiB

doc/images/computation/kutta_cond_plate.svg
+56
View File
@@ -0,0 +1,56 @@
<svg xmlns="http://www.w3.org/2000/svg" xmlns:xlink="http://www.w3.org/1999/xlink" version="1.1" width="771.87997pt" height="238.96292pt" viewBox="0 0 771.87997 238.96292">
<defs>
<clipPath id="clip_0">
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M0 0H1920V1080H0Z"/>
</clipPath>
<clipPath id="clip_1">
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1106.002 576.4439H1122.823V591.0093H1106.002Z"/>
</clipPath>
</defs>
<g clip-path="url(#clip_0)">
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M0 0H1920V1080H0Z" fill="#ffffff"/>
<path transform="matrix(.9580353,-.2866502,.2866502,.9580353,137.0238,198.89972)" stroke-width="16" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#000000" d="M0 0H528.2787"/>
<path transform="matrix(.9999944,-.003351985,.003351985,.9999944,135.55829,26.180298)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#007b76" d="M501.7182 14.89341C432.7024-50.42863 108.474 117.5608 0 163.1673"/>
<path transform="matrix(1,0,0,1,-16.150513,6.1448976)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M792.5561 47.68689C753.3962 47.0543 713.0286 46.27591 678.7394 27.16879 667.3651 20.83055 657.0775 12.581 645.1292 7.376457 624.4572-1.628081 602.7061-.8732544 580.6685 1.71949 410.6574 21.72136 262.6907 139.3051 96.37235 177.6597 64.72086 184.9589 32.47417 189.2054 0 189.9195"/>
<path transform="matrix(1,0,0,1,-11.400818,73.43042)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M783.3193 .7884597C748.5002-1.851793 713.668 2.487631 679.3967 9.075733 653.5863 14.03738 628.0432 20.27122 602.6815 27.16604 448.6335 69.04565 300.6064 135.3603 141.6754 152.7233 94.62809 157.8632 47.24892 158.5936 0 161.1963"/>
<path transform="matrix(-1,-0,0,-1,766.6359,53.89514)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H6.864602-1.735398"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1367.619 641.5311 1358.019 636.7311 1367.619 631.9311Z" fill="#0076ba"/>
<path transform="matrix(-1,-0,0,-1,762.2249,73.24945)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H6.864602-1.735398"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1363.209 622.1768 1353.609 617.3768 1363.209 612.5768Z" fill="#0076ba"/>
<path transform="matrix(-1,-0,0,-1,12.190674,194.73282)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H8.184024-.4159761"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M611.8548 500.6934 602.2548 495.8934 611.8548 491.0934Z" fill="#0076ba"/>
<path transform="matrix(-1,-0,0,-1,5.006714,233.62622)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H1-7.6"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M611.8548 461.8 602.2548 457 611.8548 452.2Z" fill="#0076ba"/>
<path transform="matrix(-1,-0,0,-1,597.72689,6.1013185)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H8.184024-.4159761"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1197.391 689.3249 1187.791 684.5249 1197.391 679.7249Z" fill="#0076ba"/>
<path transform="matrix(-.9126669,.4087042,-.4087042,-.9126669,335.68159,82.305728)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H1.684309 2.684309"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M930.4307 612.0129 923.6309 603.7085 934.3543 603.2513Z" fill="#0076ba"/>
<path transform="matrix(-.9448149,.3276048,-.3276048,-.9448149,461.72688,31.125794)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H1.684309 2.684309"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1056.811 663.4838 1049.313 655.8036 1059.956 654.4135Z" fill="#0076ba"/>
<path transform="matrix(-.9389535,.3440442,-.3440442,-.9389535,178.89917,152.30274)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H1.684309 2.684309"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M773.9144 542.251 766.5518 534.4412 777.2172 533.2371Z" fill="#0076ba"/>
<path transform="matrix(-.9669135,.2551046,-.2551046,-.9669135,597.5909,99.05774)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H8.184024-.4159761"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1195.984 596.5709 1187.926 589.4807 1198.433 587.2885Z" fill="#0076ba"/>
<path transform="matrix(-.9669135,.2551046,-.2551046,-.9669135,427.16389,151.98743)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H8.184024-.4159761"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1025.556 543.6412 1017.498 536.551 1028.005 534.3588Z" fill="#0076ba"/>
<path transform="matrix(-.9823533,.1870347,-.1870347,-.9823533,195.22657,216.43653)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M0 0H8.184024-.4159761"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M793.968 479.1698 785.4351 472.659 795.7635 469.7392Z" fill="#0076ba"/>
<path transform="matrix(1,0,0,1,389.18428,39.2666)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#0076ba" d="M147.6963 21.83868C158.9935 15.38295 183.9471 9.066525 170.0279 1.104305 152.7453-8.781958-11.57757 50.70998 .6491816 60.18663 9.170442 66.79125 47.41993 54.16165 72.9151 45.28126L73.86098 44.95004"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1061.936 601.5488 1069.41 609.2519 1058.763 610.6093Z" fill="#0076ba"/>
<path transform="matrix(1,0,0,1,389.7519,123.6391)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#000000" d="M0 0H218.4522 219.4522"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1206.452 562.1871 1216.052 566.9871 1206.452 571.7871Z"/>
<text xml:space="preserve" transform="matrix(1 0 -0 1 620.9029 132.77613)" font-size="32" font-family="LatinModernMath"><tspan y="0" x="0">v</tspan></text>
<path transform="matrix(.9578225,-.2873605,.2873605,.9578225,125.16272,158.26215)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#000000" d="M0 0H40.47911 41.47911"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M763.5619 539.3986 771.3777 546.7548 760.8033 548.5937Z"/>
<path transform="matrix(-.2873605,-.9578225,.9578225,-.2873605,125.14148,159.08972)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#000000" d="M0 0H40.47911 41.47911"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M716.3551 571.6877 708.9989 579.5034 707.16 568.929Z"/>
<text xml:space="preserve" transform="matrix(1 0 -0 1 118.88617 108.28851)" font-size="32" font-family="LatinModernMath"><tspan y="0" x="0">y</tspan></text>
<text xml:space="preserve" transform="matrix(1 0 -0 1 165.3255 132.77613)" font-size="32" font-family="LatinModernMath"><tspan y="0" x="0">x</tspan></text>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M1254.748 634.4438C1256.779 632.4717 1256.779 629.2742 1254.748 627.3021 1252.717 625.33 1249.425 625.33 1247.394 627.3021 1245.364 629.2742 1245.364 632.4717 1247.394 634.4438 1249.425 636.416 1252.717 636.416 1254.748 634.4438ZM1254.748 634.4438" fill="#99195e"/>
<path transform="matrix(1,0,0,-1,-598.2481,690.6262)" d="M729.2592 492.4438C731.2899 490.4717 731.2899 487.2742 729.2592 485.3021 727.2285 483.33 723.936 483.33 721.9053 485.3021 719.8746 487.2742 719.8746 490.4717 721.9053 492.4438 723.936 494.416 727.2285 494.416 729.2592 492.4438ZM729.2592 492.4438" fill="#99195e"/>
<path transform="matrix(1,0,0,1,478.9959,91.54443)" stroke-width="2" stroke-linecap="butt" stroke-miterlimit="4" stroke-linejoin="miter" fill="none" stroke="#000000" d="M12.88579 31.99126C13.61055 12.25828 9.315282 1.59453 0 0"/>
</g>
<g clip-path="url(#clip_1)">
<text xml:space="preserve" transform="matrix(1 0 -0 1 506.88386 139.8006)" font-size="32.2243" font-family="STIXGeneral" font-style="italic"><tspan y="-25.97284" x="0">&#x03b1;</tspan></text>
</g>
</svg>
Side by Side

After

Width:  |  Height:  |  Size: 8.6 KiB

doc/images/computation/magnus.png
BIN
View File
Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 565 KiB

doc/index.rst
+1 -1
View File
@@ -8,7 +8,7 @@
:hidden:
overview
computation
computation/index.rst
modeling
XMLreference
programming/index.rst
doc/overview.rst
+13 -12
View File
@@ -48,11 +48,12 @@ Soft, convex and analytically-invertible contact dynamics
In the modern approach to contact dynamics, the forces or impulses caused by frictional contacts are usually defined
as the solution to a linear or non-linear complementarity problem (LCP or NCP), both of which are NP-hard. MuJoCo is
based on a different formulation of the physics of contact which reduces to a convex optimization problem, as
explained in detail in the :doc:`computation` chapter. Our model allows soft contacts and other constraints, and has
a uniquely-defined inverse facilitating data analysis and control applications. There is a choice of optimization
algorithms, including a generalization to the projected Gauss-Seidel method that can handle elliptic friction cones.
The solver provides unified treatment of frictional contacts including torsional and rolling friction, frictionless
contacts, joint and tendon limits, dry friction in joints and tendons, as well as a variety of equality constraints.
explained in detail in the :doc:`computation/index` chapter. Our model allows soft contacts and other constraints,
and has a uniquely-defined inverse facilitating data analysis and control applications. There is a choice of
optimization algorithms, including a generalization to the projected Gauss-Seidel method that can handle elliptic
friction cones. The solver provides unified treatment of frictional contacts including torsional and rolling
friction, frictionless contacts, joint and tendon limits, dry friction in joints and tendons, as well as a variety of
equality constraints.
Tendon geometry
MuJoCo can model the 3D geometry of tendons - which are minimum-path-length strings obeying wrapping and via-point
@@ -746,13 +747,13 @@ multi-threading, but we do so in a way that is different from how object-oriente
Softness and slip
~~~~~~~~~~~~~~~~~
As we will explain at length in the :doc:`computation` chapter, MuJoCo is based on a mathematical model of the physics
of contact and other constraints. This model is inherently soft, in the sense that pushing harder against a constraint
will always result in larger acceleration, and so the inverse dynamics can be uniquely defined. This is desirable
because it yields a convex optimization problem and enables analyses that rely on inverse dynamics, and furthermore,
most contacts that we need to model in practice have some softness. However once we allow soft constraints, we are
effectively creating a new type of dynamics -- namely deformation dynamics -- and now we must specify how these dynamics
behave. This calls for elaborate parameterization of contacts and other constraints, involving the attributes
As we will explain at length in the :doc:`computation/index` chapter, MuJoCo is based on a mathematical model of the
physics of contact and other constraints. This model is inherently soft, in the sense that pushing harder against a
constraint will always result in larger acceleration, and so the inverse dynamics can be uniquely defined. This is
desirable because it yields a convex optimization problem and enables analyses that rely on inverse dynamics, and
furthermore, most contacts that we need to model in practice have some softness. However once we allow soft constraints,
we are effectively creating a new type of dynamics -- namely deformation dynamics -- and now we must specify how these
dynamics behave. This calls for elaborate parameterization of contacts and other constraints, involving the attributes
:at:`solref` and :at:`solimp` that can be set per constraints and will be described later.
An often confusing aspect of this soft model is that gradual contact slip cannot be avoided. Similarly, frictional
doc/programming/simulation.rst
+2 -2
View File
@@ -888,8 +888,8 @@ function is essentially free in terms of CPU cost; so do not hesitate to use thi
Contacts
~~~~~~~~
Collision detection and solving for contact forces were explained in detail in the :doc:`../computation` chapter. Here
we further clarify contact processing from a programming perspective.
Collision detection and solving for contact forces were explained in detail in the :doc:`../computation/index` chapter.
Here we further clarify contact processing from a programming perspective.
The collision detection stage finds contacts between geoms, and records them in the array ``mjData.contact`` of
:ref:`mjContact` data structures. They are sorted such that multiple contacts between the same pair of bodies are
doc/references.bib
+90
View File
@@ -0,0 +1,90 @@
@book{lamb1932,
title={Hydrodynamics. Sixth edition.},
author={Lamb, Horace},
year={1932},
publisher={Cambridge University Press}
}
@book{kennard1967,
title={Irrotational Flow of Frictionless Fluids: Mostly of Invariable Density},
author={Kennard, Earle Hesse},
volume={2299},
year={1967},
publisher={US Government Printing Office}
}
@book{tuckerman1925,
title={Inertia factors of ellipsoids for use in airship design},
author={Tuckerman, LB},
year={1925},
publisher={US Government Printing Office}
}
@article{wang2004,
title={Unsteady forces and flows in low Reynolds number hovering flight: two-dimensional computations vs robotic wing experiments},
author={Wang, Z Jane and Birch, James M and Dickinson, Michael H},
journal={Journal of Experimental Biology},
volume={207},
number={3},
pages={449--460},
year={2004},
publisher={Company of Biologists}
}
@article{andersen2005a,
title={Unsteady aerodynamics of fluttering and tumbling plates},
author={Andersen, A and Pesavento, U and Wang, Z Jane},
journal={Journal of Fluid Mechanics},
volume={541},
pages={65--90},
year={2005},
publisher={Cambridge University Press}
}
@article{andersen2005b,
title={Analysis of transitions between fluttering, tumbling and steady descent of falling cards},
author={Andersen, Anders and Pesavento, Umberto and Wang, Z Jane},
journal={Journal of Fluid Mechanics},
volume={541},
pages={91--104},
year={2005},
publisher={Cambridge University Press}
}
@article{stokes1850,
title={On the effect of internal friction of fluids on the motion of pendulums},
author={Stokes, GG},
journal={Trans. Camb. phi1. S0c},
volume={9},
number={8},
pages={106},
year={1850}
}
@article{loth2008,
title={Drag of non-spherical solid particles of regular and irregular shape},
author={Loth, E},
journal={Powder Technology},
volume={182},
number={3},
pages={342--353},
year={2008},
publisher={Elsevier}
}
@article{bagheri2016,
title={On the drag of freely falling non-spherical particles},
author={Bagheri, Gholamhossein and Bonadonna, Costanza},
journal={Powder Technology},
volume={301},
pages={526--544},
year={2016},
publisher={Elsevier}
}
@article{duan2015,
title={Sphere drag and heat transfer},
author={Duan, Zhipeng and He, Boshu and Duan, Yuanyuan},
journal={Scientific reports},
volume={5},
number={1},
pages={1--7},
year={2015},
publisher={Nature Publishing Group}
}
doc/requirements.txt
+1
View File
@@ -1,5 +1,6 @@
Sphinx==4.5.0
furo==2022.9.29
sphinxcontrib-bibtex==2.6.1
sphinxcontrib-katex==0.9.4
sphinxcontrib-youtube==1.2.0
sphinx-copybutton==0.5.2