Update Programming/simulation chapter.

PiperOrigin-RevId: 847424298
Change-Id: Icc6efb05dee5e67a5813ed04caa15a5703e59aee

doc/APIreference/functions.rst

@@ -89,7 +89,7 @@ Copy real-valued arrays from model to spec, returns 1 on success.
 Recompile spec to model, preserving the state. Like :ref:`mj_compile`, this function compiles an :ref:`mjSpec` to an
 :ref:`mjModel`, with two differences. First, rather than returning an entirely new model, it will
 reallocate existing :ref:`mjModel` and :ref:`mjData` instances in-place. Second, it will preserve the
-:ref:`integration state<geIntegrationState>`, as given in the provided :ref:`mjData` instance, while accounting for
+:ref:`integration state<siIntegrationState>`, as given in the provided :ref:`mjData` instance, while accounting for
 newly added or removed degrees of freedom. This allows the user to continue simulation with the same model and data
 struct pointers while editing the model programmatically.
 
@@ -2753,7 +2753,7 @@ outputs of derivative functions are the trailing rather than leading arguments.
 
 Compute finite-differenced discrete-time transition matrices.
 
-Letting :math:`x, u` denote the current :ref:`state<gePhysicsState>` and :ref:`control<geInput>`
+Letting :math:`x, u` denote the current :ref:`state<siPhysicsState>` and :ref:`control<siInput>`
 vector in an mjData instance, and letting :math:`y, s` denote the next state and sensor
 values, the top-level :ref:`mj_step` function computes :math:`(x,u) \rightarrow (y,s)`
 :ref:`mjd_transitionFD` computes the four associated Jacobians using finite-differencing.
@@ -2804,7 +2804,7 @@ These matrices and their dimensions are:
 
 Finite differenced continuous-time inverse-dynamics Jacobians.
 
-Letting :math:`x, a` denote the current :ref:`state<gePhysicsState>` and acceleration vectors in an mjData instance, and
+Letting :math:`x, a` denote the current :ref:`state<siPhysicsState>` and acceleration vectors in an mjData instance, and
 letting :math:`f, s` denote the forces computed by the inverse dynamics (``qfrc_inverse``), the function
 :ref:`mj_inverse` computes :math:`(x,a) \rightarrow (f,s)`. :ref:`mjd_inverseFD` computes seven associated Jacobians
 using finite-differencing. These matrices and their dimensions are:

doc/APIreference/functions_override.rst

@@ -48,7 +48,7 @@ If compilation fails, :ref:`mj_compile` returns ``NULL``; the error can be read
 Recompile spec to model, preserving the state. Like :ref:`mj_compile`, this function compiles an :ref:`mjSpec` to an
 :ref:`mjModel`, with two differences. First, rather than returning an entirely new model, it will
 reallocate existing :ref:`mjModel` and :ref:`mjData` instances in-place. Second, it will preserve the
-:ref:`integration state<geIntegrationState>`, as given in the provided :ref:`mjData` instance, while accounting for
+:ref:`integration state<siIntegrationState>`, as given in the provided :ref:`mjData` instance, while accounting for
 newly added or removed degrees of freedom. This allows the user to continue simulation with the same model and data
 struct pointers while editing the model programmatically.
 
@@ -643,7 +643,7 @@ outputs of derivative functions are the trailing rather than leading arguments.
 
 Compute finite-differenced discrete-time transition matrices.
 
-Letting :math:`x, u` denote the current :ref:`state<gePhysicsState>` and :ref:`control<geInput>`
+Letting :math:`x, u` denote the current :ref:`state<siPhysicsState>` and :ref:`control<siInput>`
 vector in an mjData instance, and letting :math:`y, s` denote the next state and sensor
 values, the top-level :ref:`mj_step` function computes :math:`(x,u) \rightarrow (y,s)`
 :ref:`mjd_transitionFD` computes the four associated Jacobians using finite-differencing.
@@ -689,7 +689,7 @@ These matrices and their dimensions are:
 
 Finite differenced continuous-time inverse-dynamics Jacobians.
 
-Letting :math:`x, a` denote the current :ref:`state<gePhysicsState>` and acceleration vectors in an mjData instance, and
+Letting :math:`x, a` denote the current :ref:`state<siPhysicsState>` and acceleration vectors in an mjData instance, and
 letting :math:`f, s` denote the forces computed by the inverse dynamics (``qfrc_inverse``), the function
 :ref:`mj_inverse` computes :math:`(x,a) \rightarrow (f,s)`. :ref:`mjd_inverseFD` computes seven associated Jacobians
 using finite-differencing. These matrices and their dimensions are:

doc/changelog.rst

@@ -29,6 +29,14 @@ MJX
 ^^^
 - Added ``actuator_length``, ``cdof`` and ``cdof_dof`` fields to ``mjx.Data``.
 
+
+Documentation
+^^^^^^^^^^^^^
+- General improvements to the :ref:`Programming/Simulation<Simulation>` chapter. Notably, the main discussion of
+  :ref:`state<siStateControl>` has been moved there, and the section on :ref:`mjModel changes<siChange>` has been
+  expanded.
+
+
 Bug fixes
 ^^^^^^^^^
 
@@ -1081,8 +1089,8 @@ Python bindings
 12. Improved the implementation of the :ref:`rollout<PyRollout>` module. Note the changes below are breaking, dependent
     code will require modification.
 
-    - Uses :ref:`mjSTATE_FULLPHYSICS<geFullPhysics>` as state spec, enabling divergence detection by inspecting time.
+    - Uses :ref:`mjSTATE_FULLPHYSICS<siFullPhysics>` as state spec, enabling divergence detection by inspecting time.
-    - Allows user-defined control spec for any combination of :ref:`user input<geInput>` fields as controls.
+    - Allows user-defined control spec for any combination of :ref:`user input<siInput>` fields as controls.
     - Outputs are no longer squeezed and always have dim=3.
 13. The ``sync`` function for the :ref:`passive viewer<PyViewerPassive>` can now pick up changes to rendering flags in
     ``user_scn``, as requested in :issue:`1190`.

doc/computation/index.rst

@@ -599,126 +599,29 @@ Fast implicit-in-velocity (``implicitfast``)
      ``implicitfast`` yet conserves energy well under ``RK4``. Note that under ``implicit``, this model doesn't diverge
      but rather loses energy.
 
+
+.. Leave links below to sections that were previously here and have moved to the simulation chapter.
+.. _gePhysicsState:
+.. _geFullPhysics:
+.. _geInput:
+.. _geWarmstart:
+.. _geIntegrationState:
+.. _geSimulationState:
+
 .. _geState:
 
 The State
 ~~~~~~~~~
 
-To complete our description of the general framework we will now discuss the notion of *state*. MuJoCo has a compact,
+To complete our description of the general framework we will quickly discuss the notion of *state*. MuJoCo has a
-well-defined internal state which, together with the :ref:`deterministic computational pipeline<piReproducibility>`,
+compact, well-defined internal state which, together with the :ref:`deterministic pipeline<piReproducibility>`, means
-means that operations like resetting the state and computing dynamics derivatives are also well-defined.
+that operations like (re)setting the state and computing dynamics derivatives are also well-defined.
 
 The state is entirely encapsulated in the :ref:`mjData` struct and consists of several components. The components are
-enumerated in :ref:`mjtState` as bit flags, along with several common combinations, corresponding to the groupings
+enumerated in :ref:`mjtState` as bit flags. Concatenated state vectors can be conveniently read from and written into
-below. Concatenated state vectors can be conveniently read from and written into :ref:`mjData` using :ref:`mj_getState`
+:ref:`mjData` using :ref:`mj_getState` and :ref:`mj_setState`, respectively.
-and :ref:`mj_setState`, respectively.
 
-.. _gePhysicsState:
+More detail can be found in the :ref:`State and Control<siStateControl>` section in the Simulation chapter.
-
-Physics state
-^^^^^^^^^^^^^
-The *physics state* (:ref:`mjSTATE_PHYSICS<mjtState>`) contains the main quantities which are time-integrated during
-stepping. These are ``mjData.{qpos, qvel, act}``:
-
-Position: ``qpos``
-  The configuration in generalized coodinates, denoted above as :math:`q`.
-
-Velocity: ``qvel``
-  The generalized velocities, denoted above as :math:`v`.
-
-Actuator activation: ``act``
-  ``mjData.act`` contains the internal states of stateful actuators, denoted above as :math:`w`.
-
-.. _geFullPhysics:
-
-Full physics state
-^^^^^^^^^^^^^^^^^^
-
-The *full physics state* (:ref:`mjSTATE_FULLPHYSICS<mjtState>`) contains the physics state and two additional
-components:
-
-Time: ``time``
-  The simulation time is given by the scalar ``mjData.time``. Since physics is time-invariant, it is
-  excluded from the *physics state*; exceptions include time-dependent user callbacks and plugins (e.g., an open-loop
-  controller), in which case time should be included.
-
-Plugin state: ``plugin_state``
-  ``mjData.plugin_state`` are states declared by :ref:`engine plugins<exPlugin>`. Please see the :ref:`exPluginState`
-  section for more details.
-
-.. _geInput:
-
-User inputs
-^^^^^^^^^^^
-
-These input fields (:ref:`mjSTATE_USER<mjtState>`) are set by the user and affect the physics simulation, but are
-untouched by the simulator. All input fields except for MoCap poses default to 0.
-
-Control: ``ctrl``
-  Controls are defined by the :ref:`actuator<actuator>` section of the XML. ``mjData.ctrl`` values either produce
-  generalized forces directly (stateless actuators), or affect the actuator activations in ``mjData.act``, which then
-  produce forces.
-
-Auxiliary Controls: ``qfrc_applied`` and ``xfrc_applied``
-  | ``mjData.qfrc_applied`` are directly applied generalized forces.
-  | ``mjData.xfrc_applied`` are Cartesian wrenches applied to the CoM of individual bodies. This field is used for
-    example, by the :ref:`native viewer<saSimulate>` to apply mouse perturbations.
-  | Note that the effects of ``qfrc_applied`` and ``xfrc_applied`` can usually be recreated by appropriate actuator
-    definitions.
-
-MoCap poses: ``mocap_pos`` and ``mocap_quat``
-  ``mjData.mocap_pos`` and ``mjData.mocap_quat`` are special optional kinematic states :ref:`described here<CMocap>`,
-  which allow the user to set the positions and orientations of static bodies in real-time, for example when streaming
-  6D poses from a motion-capture device. The default values set by :ref:`mj_resetData` are the poses of the bodies at
-  the default configuration.
-
-Equality constraint toggle: ``eq_active``
-  ``mjData.eq_active`` is a byte-valued array that allows the user to toggle the state of equality constraints at
-  runtime. The initial value of this array is ``mjModel.eq_active0`` which can be set in XML using the
-  :ref:`active<equality-connect-active>` attribute of :ref:`equality constraints<coEquality>`.
-
-User data: ``userdata``
-  ``mjData.userdata`` acts as a user-defined memory space untouched by the engine. For example it can be used by
-  callbacks. This is described in more detail in the :ref:`Programming chapter<siSimulation>`.
-
-.. _geWarmstart:
-
-Warmstart acceleration
-^^^^^^^^^^^^^^^^^^^^^^
-
-``qacc_warmstart``
-  ``mjData.qacc_warmstart`` are accelerations used to warmstart the constraint solver, saved from the previous step.
-  When using a slowly-converging :ref:`constraint solver<Solver>` like PGS, these can speed up simulation by reducing
-  the number of iterations required for convergence. Note however that the default Newton solver converges so quickly
-  (usually 2-3 iterations), that warmstarts often have no effect on speed and can be disabled.
-
-  Different warmstarts have no perceptible effect on the dynamics but should be saved if perfect numerical
-  reproducibility is required when loading a non-initial state. Note that even though their effect on physics is
-  negligible, many physical systems will accumulate small differences  `exponentially
-  <https://en.wikipedia.org/wiki/Lyapunov_exponent>`__ when time-stepping, quickly leading to divergent trajectories
-  for different warmstarts.
-
-.. _geIntegrationState:
-
-Integration state
-^^^^^^^^^^^^^^^^^
-
-The *integration state* (:ref:`mjSTATE_INTEGRATION<mjtState>`) is the union of all the above :ref:`mjData` fields and
-constitutes the entire set of inputs to the *forward dynamics*. In the case of *inverse dynamics*, ``mjData.qacc`` is
-also treated as an input variable. All other :ref:`mjData` fields are functions of the integration state.
-
-Note that the full integration state as given by :ref:`mjSTATE_INTEGRATION<mjtState>` is maximalist and includes fields
-which are often unused. If a small state size is desired, it might be sensible to avoid saving unused fields.
-In particular ``xfrc_applied`` can be quite large (``6 x nbody``) yet is often unused.
-
-.. _geSimulationState:
-
-Simulation state: ``mjData``
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-The *simulation state* is the entirety of the :ref:`mjData` struct and associated memory buffer. This state includes
-all derived quantities computed during dynamics computation. Because the :ref:`mjData` buffers are preallocated for the
-worst case, it is often significantly faster to recompute derived quantities from the *integration state* rather than
-using :ref:`mj_copyData`.
 
 .. _Constraint:
 
@@ -1824,8 +1727,8 @@ Control callback
 
 Force/acceleration
 ''''''''''''''''''
-The stages below compute quantities that depend on :ref:`user inputs<geInput>`. Due to the sequential nature
+The stages below compute quantities that depend on :ref:`user inputs<siInput>`. Due to the sequential nature
-of the pipeline, the actual dependence is on the entire :ref:`integration state<geIntegrationState>`.
+of the pipeline, the actual dependence is on the entire :ref:`integration state<siIntegrationState>`.
 
 20. Compute the actuator forces and activation dynamics if defined: :ref:`mj_fwdActuation`
 21. Compute the joint acceleration resulting from all forces except for the (still unknown) constraint forces:
@@ -1871,8 +1774,8 @@ MuJoCo's simulation pipeline is entirely deterministic and reproducible -- if a
 saved and reloaded and :ref:`mj_step` called again, the resulting next state will be identical. However, there are some
 important caveats:
 
-- Save all the required :ref:`integration state<geIntegrationState>` components. In particular :ref:`warmstart
+- Save all the required :ref:`integration state<siIntegrationState>` components. In particular :ref:`warmstart
-  accelerations<geWarmstart>` have only a very small effect on the next state, but should be saved if bit-wise equality
+  accelerations<siWarmstart>` have only a very small effect on the next state, but should be saved if bit-wise equality
   is required.
 - Any numerical difference between states, no matter how small, will become significant upon integration, especially for
   systems with contact. Contact events have high `Lyapunov exponents

doc/programming/modeledit.rst

@@ -21,15 +21,15 @@ The new API augments the traditional workflow of creating and editing models usi
 *compile* steps. As summarized in the :ref:`Overview chapter<Instance>`, the traditional workflow is:
 
  1. Create an XML model description file (MJCF or URDF) and associated assets. |br|
- 2. Call :ref:`mj_loadXML`, obtain an mjModel instance.
+ 2. Call :ref:`mj_loadXML`, obtain an :ref:`mjModel` instance.
 
 The new workflow using :ref:`mjSpec` is:
 
- 1. :ref:`Create<mj_makeSpec>` an empty mjSpec or :ref:`parse<mj_parseXML>` an existing XML file.
+ 1. Create an empty :ref:`mjSpec` using :ref:`mj_makeSpec` or parse an existing XML file using :ref:`mj_parseXML`.
- 2. Programmatically edit the mjSpec datastructure by adding, modifying and removing elements.
+ 2. Programmatically edit the :ref:`mjSpec` datastructure by adding, modifying and removing elements.
- 3. :ref:`Compile<mj_compile>` the mjSpec to an mjModel instance.
+ 3. Compile the :ref:`mjSpec` to an :ref:`mjModel` instance using :ref:`mj_compile`.
 
- After compilation, the mjSpec remains editable, so steps 2 and 3 are interchangeable.
+ After compilation, the :ref:`mjSpec` remains editable, so steps 2 and 3 are interchangeable.
 
 
 .. _meUsage:

doc/programming/simulation.rst

@@ -8,20 +8,16 @@ Simulation
 Initialization
 ~~~~~~~~~~~~~~
 
-After the :ref:`version <inVersion>` check, the next step is to allocate and initialize the main data structures needed
+:ref:`mjModel` and :ref:`mjData` should never be allocated directly by the user. Instead they are allocated and
-for simulation, namely mjModel and mjData. Additional initialization steps related to visualization and callbacks will
+initialized by the corresponding API functions. These are very elaborate data structures, containing (arrays of) other
-be discussed later.
+structures, preallocated data arrays for all intermediate results, as well as an :ref:`internal stack <siStack>`. Our
-
+strategy is to allocate all necessary heap memory at the beginning of the simulation, and free it after the simulation
-mjModel and mjData should never be allocated directly by the user. Instead they are allocated and initialized by the
+is done, so that we never have to call the C memory allocation and deallocation functions during the simulation. This is
-corresponding API functions. These are very elaborate data structures, containing (arrays of) other structures,
+done for speed, avoidance of memory fragmentation, future GPU portability, and ease of managing the state of the entire
-preallocated data arrays for all intermediate results, as well as an :ref:`internal stack <siStack>`. Our strategy is
+simulator during a reset. It also means however that the maximal variable-memory allocation given by the :at:`memory`
-to allocate all necessary heap memory at the beginning of the simulation, and free it after the simulation is done, so
+attribute in the :ref:`size <size>` MJCF element, which affects the allocation of :ref:`mjData`, must be set to a
-that we never have to call the C memory allocation and deallocation functions during the simulation. This is done for
+sufficiently large value. If this maximal size is exceeded during simulation, it is not increased dynamically, but
-speed, avoidance of memory fragmentation, future GPU portability, and ease of managing the state of the entire
+instead an error is generated. See also :ref:`diagnostics <siDiagnostics>` below.
-simulator during a reset. It also means however that the maximal variable-memory allocation given by the
-:at:`memory` attribute in the :ref:`size <size>` MJCF element, which affects the allocation of ``mjData``, must be
-set to a sufficiently large value. If this maximal size is exceeded during simulation, it is not increased
-dynamically, but instead an error is generated. See also :ref:`diagnostics <siDiagnostics>` below.
 
 First we must call one of the functions that allocates and initializes mjModel and returns a pointer to it. The
 available options are
@@ -43,6 +39,9 @@ available options are
    // option 5: deep copy from existing mjModel
    mjModel* m = mj_copyModel(NULL, mexisting);
 
+   // option 6: compile model from mjSpec
+   mjModel* m = mj_compile(spec, vfs);
+
 All these functions return a NULL pointer if there is an error or warning. In the case of XML parsing and model
 compilation, a description of the error is returned in the string provided as argument. For the remaining functions, the
 low-level :ref:`mju_error` or :ref:`mju_warning` is called with the error/warning message; see :ref:`error handling
@@ -57,7 +56,7 @@ API Reference chapter.
 In addition to mjModel which holds the model description, we also need mjData which is the workspace where all
 computations are performed. Note that mjData is specific to a given mjModel. The API functions generally assume that
 users know what they are doing, and perform minimal argument checking. If the mjModel and mjData passed to any API
-function are incompatible (or NULL) the resulting behavior is unpredictable. mjData is created with
+function are incompatible (or NULL) the resulting behavior is undefined. mjData is created with
 
 .. code-block:: C
 
@@ -80,12 +79,7 @@ done, we can delete them with
 
 The code samples illustrate the complete initialization and termination sequence.
 
-MuJoCo simulations are deterministic with one exception: sensor noise can be generated when this feature is enabled.
+MuJoCo simulations are :ref:`deterministic <piReproducibility>`.
-This is done by calling the C function rand() internally. To generate the same random number sequence, call srand()
-with a desired seed after the model is loaded and before the simulation starts. The model compiler calls srand(123)
-internally, so as to generate random dots for procedural textures. Therefore the noise sequence in the sensor data
-will change if the specification of procedural textures changes, and the user does not call srand() after model
-compilation.
 
 .. _siSimulation:
 
@@ -102,7 +96,7 @@ function :ref:`mj_step` in a loop such as
      mj_step(m, d);
 
 This by itself will simulate the passive dynamics, because we have not provided any control signals or applied forces.
-The default (and recommended) way to control the system is to implement a control callback, for example
+The default way to control the system is to implement a control callback, for example
 
 .. code-block:: C
 
@@ -126,9 +120,7 @@ global control callback pointer :ref:`mjcb_control`:
    mjcb_control = mycontroller;
 
 Now if we call :ref:`mj_step`, our control callback will be executed whenever the control
-signal is needed by the simulation pipeline, and as a result we will end up simulating the controlled dynamics (except
+signal is needed by the simulation pipeline, and as a result we will end up simulating the controlled dynamics.
-damping does not really do justice to the notion of control, and is better implemented as a passive joint property,
-but these are finer points).
 
 Instead of relying on a control callback, we could set the control vector ``mjData.ctrl`` directly. Alternatively we
 could set applied forces as explained in :ref:`state and control <siStateControl>`. If we could compute these control-
@@ -144,7 +136,7 @@ control callback) would become
 
 Why would we not be able to compute the controls before ``mj_step`` is called? After all, isn't this what causality
 means? The answer is subtle but important, and has to do with the fact that we are simulating in discrete time. The top-
-level simulation function ``mj_step`` basically does two things: compute the :ref:`forward dynamics <siForward>` in
+level simulation function ``mj_step`` does two things: compute the :ref:`forward dynamics <siForward>` in
 continuous time, and then integrate over a time period specified by ``mjModel.opt.timestep``. Forward dynamics computes
 the acceleration ``mjData.qacc`` at time ``mjData.time``, given the :ref:`state and control <siStateControl>` at time
 ``mjData.time``. The numerical integrator then advances the state and time to ``mjData.time + mjModel.opt.timestep``.
@@ -153,7 +145,7 @@ controller can be a very complex function, depending on various features of the
 computed by MuJoCo as intermediate results of the simulation. These may include contacts, Jacobians, passive forces.
 None of these quantities are available before ``mj_step`` is called (or rather, they are available but *outdated by one
 time step*). In contrast, when ``mj_step`` calls our control callback, it does so as late in the computation as possible
-- namely after all the intermediate results dependent on the state but not on the control have been computed.
+-- namely after all the intermediate results dependent on the state but not on the control have been computed.
 
 The same effect can be achieved without using a control callback. This is done by breaking ``mj_step`` in two parts:
 before the control is needed, and after the control is needed. The simulation loop now becomes
@@ -166,14 +158,14 @@ before the control is needed, and after the control is needed. The simulation lo
      mj_step2(m, d);
    }
 
-There is one complication however: this only works with Euler integration. The Runge-Kutta integrator (as well as other
+There is one complication however: this only works with the single-step :ref:`integrators<geIntegration>` (Euler,
-advanced integrators we plan to implement) need to evaluate the entire dynamics including the feedback control law
+implicit, implicitfast). The Runge-Kutta integrator needs to evaluate the entire dynamics including the feedback control
-multiple times per step, which can only be done using a control callback. But with Euler integration, the above
+law multiple times per step, which can only be done using a control callback. With the single-step integrators, the
-separation of ``mj_step`` into :ref:`mj_step1` and :ref:`mj_step2` is sufficient to provide the control law with the
+above separation of ``mj_step`` into :ref:`mj_step1` and :ref:`mj_step2` is sufficient to provide the control law with
-intermediate results of the computation.
+the intermediate results of the computation.
 
-To make the above discussion more clear, we provide the internal implementation of mj_step, mj_step1 and mj_step2,
+To make the above discussion clearer, we provide the internal implementation of :ref:`mj_step`, :ref:`mj_step1` and
-omitting some code that computes timing diagnostics. The main simulation function is
+:ref:`mj_step2`, omitting some code that computes timing diagnostics. The main simulation function is
 
 .. code-block:: C
 
@@ -184,15 +176,24 @@ omitting some code that computes timing diagnostics. The main simulation functio
      mj_forward(m, d);
      mj_checkAcc(m, d);
 
-     // compare forward and inverse solutions if enabled
-     if (mjENABLED(mjENBL_FWDINV))
-       mj_compareFwdInv(m, d);
-
      // use selected integrator
-     if (m->opt.integrator == mjINT_RK4)
+     switch ((mjtIntegrator) m->opt.integrator) {
-       mj_RungeKutta(m, d, 4);
+     case mjINT_EULER:
-     else
        mj_Euler(m, d);
+       break;
+
+     case mjINT_RK4:
+       mj_RungeKutta(m, d, 4);
+       break;
+
+     case mjINT_IMPLICIT:
+     case mjINT_IMPLICITFAST:
+       mj_implicit(m, d);
+       break;
+
+     default:
+       mjERROR("invalid integrator");
+     }
    }
 
 The checking functions reset the simulation automatically if any numerical values have become invalid or too large.
@@ -200,8 +201,8 @@ The control callback (if any) is called from within the forward dynamics functio
 
 Next we show the implementation of the two-part stepping approach, although the specifics will make sense only after
 we explain the :ref:`forward dynamics <siForward>` later. Note that the control callback is now called directly, since
-we have essentially unpacked the forward dynamics function. Note also that we always call the Euler integrator in
+we have essentially unpacked the forward dynamics function. Note also that we always call a single-step integrator in
-mj_step2 regardless of the setting of ``mjModel.opt.integrator``.
+:ref:`mj_step2`; if RK4 is selected, the integrator will default to Euler.
 
 .. code-block:: C
 
@@ -227,19 +228,17 @@ mj_step2 regardless of the setting of ``mjModel.opt.integrator``.
      mj_sensorAcc(m, d);
      mj_checkAcc(m, d);
 
-     // compare forward and inverse solutions if enabled
+     // integrate with Euler or implicit; RK4 defaults to Euler
-     if (mjENABLED(mjENBL_FWDINV))
+     if (m->opt.integrator == mjINT_IMPLICIT || m->opt.integrator == mjINT_IMPLICITFAST)
-       mj_compareFwdInv(m, d);
+       mj_implicit(m, d);
-
+     else
-     // integrate with Euler; ignore integrator option
+       mj_Euler(m, d);
-     mj_Euler(m, d);
    }
 
 .. _siStateControl:
 
 State and control
 ~~~~~~~~~~~~~~~~~
-
 MuJoCo has a well-defined state that is easy to set, reset and advance through time. This is closely related to the
 notion of state of a dynamical system. Dynamical systems are usually described in the general form
 
@@ -247,122 +246,201 @@ notion of state of a dynamical system. Dynamical systems are usually described i
 
      dx/dt = f(t, x, u)
 
-where ``t`` is the time, ``x`` is the state vector, ``u`` is the control vector, and ``f`` is the function that
+where ``t`` is the time, ``x`` is the state vector, ``u`` is the control vector, and ``f`` is the function that computes
-computes the time-derivative of the state. This is a continuous-time formulation, and indeed the physics model
+the time-derivative of the state. This is a continuous-time formulation, and indeed the physics model simulated by
-simulated by MuJoCo is defined in continuous time. Even though the numerical integrator operates in discrete time, the
+MuJoCo is defined in continuous time. Even though the numerical integrator operates in discrete time, the main part of
-main part of the computation---namely the function :ref:`mj_forward`---corresponds to the
+the computation---namely the function :ref:`mj_forward`---corresponds to the continuous-time dynamics function
-continuous-time dynamics function ``f(t,x,u)`` above. Here we explain this correspondence.
+``f(t,x,u)`` above. Below we explain this correspondence.
 
-The state vector in MuJoCo is:
+.. _siStateComponents:
 
-.. code-block:: Text
+State components
+^^^^^^^^^^^^^^^^
+The state is composed of distinct components, described in the :ref:`mjtState` bitfield enum, which enumerates both
+individual components and combinations of components. These are:
 
-     x = (mjData.time, mjData.qpos, mjData.qvel, mjData.act)
+.. _siPhysicsState:
 
-For a second-order dynamical system the state contains only position and velocity, however MuJoCo can also model
+Physics state
-actuators (such as cylinders and biological muscles) that have their own activation states assembled in the vector
+"""""""""""""
-``mjData.act``. While the physics model is time-invariant, user-defined control laws may be time-varying; in particular
+The *physics state* (:ref:`mjSTATE_PHYSICS<mjtState>`) contains the main quantities which are time-integrated during
-control laws obtained from trajectory optimizers would normally be indexed by ``mjData.time``.
+stepping. These are ``mjData.{qpos, qvel, act}``:
 
-The reason for the "official" caveat above is because user callbacks may store additional state variables that change
+Position: ``qpos``
-over time and affect the callback outputs; indeed the field ``mjData.userdata`` exists mostly for that purpose. Other
+  The configuration in generalized coodinates, denoted in the :ref:`Numerical Integration<geIntegration>` section as
-state-like quantities that are part of mjData and are treated as inputs by forward dynamics are ``mjData.mocap_pos`` and
+  :math:`q`.
-``mjData.mocap_quat``. These quantities are unusual in that they are meant to change at each time step (normally driven by a
-motion capture device), however this change is implemented by the user, while the simulator treats them as constants. In
-that sense they are no different from all the constants in mjModel, or the function callback pointers set by the user:
-such constants affect the computation, but are not part of the state vector of a dynamical system.
 
-The warm-start mechanism in the constraint solver effectively introduces another state variable. This mechanism uses
+Velocity: ``qvel``
-the output of forward dynamics from the previous time step, namely the acceleration vector ``mjData.qacc``, to estimate
+  The generalized velocities, denoted in the :ref:`Numerical Integration<geIntegration>` section as :math:`v`. In the
-the current constraint forces via inverse dynamics. This estimate then initializes the optimization algorithm in the
+  presence of quaternions (i.e., when free or ball joints are used), the position vector ``mjData.qpos`` has higher
-solver. If this algorithm runs until convergence the warm-start will affect the speed of convergence but not the final
+  dimensionality than the velocity vector ``mjData.qvel`` and so this is not a simple time-derivative in the sense of
-solution (since the underlying optimization problem is convex and does not have local minima), but in practice the
+  scalars, but instead takes quaternion algebra into account.
-algorithm is often terminated early, and so the warm-start has some (usually very small) effect on the solution.
 
-Next we turn to the controls and applied forces. The control vector in MuJoCo is
+Actuator activation: ``act``
+  For a second-order mechanical system, the state contains only position and velocity, but MuJoCo also models stateful
+  actuators (such as biological muscles) that have their own activation states assembled in ``mjData.act``, denoted
+  as :math:`w` in the :ref:`Numerical Integration<geIntegration>` section.
 
-.. code-block:: Text
+.. _siFullPhysics:
 
-     u = (mjData.ctrl, mjData.qfrc_applied, mjData.xfrc_applied)
+Full physics state
+""""""""""""""""""
+`t, x` above correspond to the *full physics state* (:ref:`mjSTATE_FULLPHYSICS<mjtState>`) -- everything which advances
+in time. It is the :ref:`Physics state<siPhysicsState>` and two additional components:
 
-These quantities specify control signals (``mjData.ctrl``) for the actuators defined in the model, or directly apply
+Time: ``time``
-forces and torques specified in joint space (``mjData.qfrc_applied``) or in Cartesian space (``mjData.xfrc_applied``).
+  While mechanics is time-invariant, user-defined control laws may be time-dependent; in particular control laws
+  obtained from trajectories are often time-indexed. The time ``t`` (``mjData.time``) is therefore a state component
+  with ``dt/dt == 1``.
 
-Finally, calling mj_forward which corresponds to the abstract dynamics function ``f(t,x,u)`` computes the
+Plugin state: ``plugin_state``
-time-derivative of the state vector. The corresponding fields of mjData are
+  ``mjData.plugin_state`` are states declared by :ref:`engine plugins<exPlugin>`. Please see the :ref:`exPluginState`
+  section for more details.
 
-::
+.. _siInput:
 
-     dx/dt = f(t,x,u) = (1, mjData.qvel, mjData.qacc, mjData.act_dot)
+User inputs
+"""""""""""
+These input fields (:ref:`mjSTATE_USER<mjtState>`) are set by the user and affect the physics simulation, but are
+untouched by the simulator. All input fields except for MoCap poses default to 0. A general property of all :ref:`User
+input<siInput>` arrays is that they are untouched by the library. Therefore, in the sense that a value written to this
+memory is persistent, they can also be considered stateful.
 
-In the presence of quaternions (i.e., when free or ball joints are used), the position vector ``mjData.qpos`` has higher
+The control vector ``u`` mostly corresponds to the array ``mjData.ctrl``, containing the actuation signal set by the
-dimensionality than the velocity vector ``mjData.qvel`` and so this is not a simple time-derivative in the sense of
+user. "Mostly" because torques and wrenches can also be applied directly using ``mjData.qfrc_applied`` and
-scalars, but instead takes quaternion algebra into account.
+``mjData.xfrc_applied``, respectively. The poses of mocap bodies, which are :ref:`user-controlled static
+bodies<CMocap>`, are also a user-input. The field ``mjData.userdata`` is a fixed-size memory block (allocated by setting
+:ref:`nuserdata<size-nuserdata>`) meant to serve the user for any purpose, and can be used to store various state-like
+and control-like quantities.
 
-To illustrate how the simulation state can be manipulated, suppose we have two mjData pointers src and dst
+Control: ``ctrl``
-corresponding to the same mjModel, and we want to copy the entire simulation state from one to the other (leaving out
+  Controls are defined by the :ref:`actuator<actuator>` section of the XML. ``mjData.ctrl`` values either produce
-internal diagnostics which do not affect the simulation). This can be done as
+  generalized forces directly (stateless actuators), or affect the actuator activations in ``mjData.act``, which then
+  produce forces. Note that while all actuators produce forces, the semantics of ``ctrl`` and ``act`` depend on the
+  specifc parameters of the :ref:`actuation model<geActuation>`.
+
+Auxiliary Controls: ``qfrc_applied`` and ``xfrc_applied``
+  | ``mjData.qfrc_applied`` are directly applied generalized forces.
+  | ``mjData.xfrc_applied`` are Cartesian wrenches applied to the CoM of individual bodies. This field is used for
+    example, by the :ref:`native viewer<saSimulate>` to apply mouse perturbations.
+  | Note that the effects of ``qfrc_applied`` and ``xfrc_applied`` can be recreated by appropriate actuator
+    definitions.
+
+MoCap poses: ``mocap_pos`` and ``mocap_quat``
+  ``mjData.mocap_pos`` and ``mjData.mocap_quat`` are special optional kinematic states :ref:`described here<CMocap>`,
+  which allow the user to set the positions and orientations of static bodies in real-time, for example when streaming
+  6D poses from a motion-capture device. The default values set by :ref:`mj_resetData` are the poses of the bodies at
+  the default configuration.
+
+Equality constraint toggle: ``eq_active``
+  ``mjData.eq_active`` is a byte-valued array that allows the user to toggle the state of equality constraints at
+  runtime. The initial value of this array is ``mjModel.eq_active0`` which can be set in XML using the
+  :ref:`active<equality-connect-active>` attribute of :ref:`equality constraints<coEquality>`.
+
+User data: ``userdata``
+  ``mjData.userdata`` acts as a user-defined memory space untouched by the engine. For example it can be used by
+  callbacks. This is described in more detail in the :ref:`Programming chapter<siSimulation>`.
+
+.. _siWarmstart:
+
+Warmstarts
+""""""""""
+
+Warmstart accelerations: ``qacc_warmstart``
+  ``mjData.qacc_warmstart`` are the previous step's accelerations used to warmstart the constraint solver.
+  Assuming that the current solution is not very different from the previous one, this can speed up simulation by
+  reducing the number of iterations required for convergence. When using a slowly-converging :ref:`constraint
+  solver<Solver>` like PGS, these can speed up simulation by reducing the number of iterations required for convergence.
+  Note however that the default Newton solver converges so quickly (usually 2-3 iterations), that warmstarts often have
+  a negligible effect on speed and can be :ref:`disabled<option-flag-warmstart>`.
+
+  Because our optimization problem is :ref:`strictly convex<Solver>` with a single global minimum, different solver
+  initialization have no perceptible effect on the solution, assuming that convergence was achieved. The effect becomes
+  significant if numerical convergence is not achieved, either due to slow convergence or if
+  :ref:`iterations<option-iterations>` or :ref:`tolerance<option-tolerance>` are capped, as is sometimes done in
+  :ref:`MJX<MjxPerformance>`.
+
+  The other case where warmstarts are critical is if perfect numerical reproducibility is required, when loading a
+  non-initial state (since the initial state is always cold-started). Note that even though their effect on physics is
+  negligible, many physical systems will accumulate small differences `exponentially
+  <https://en.wikipedia.org/wiki/Lyapunov_exponent>`__ when time-stepping, quickly leading to divergent trajectories for
+  different warmstarts. See :ref:`Reproducibility<piReproducibility>` for more details.
+
+.. _siIntegrationState:
+
+Integration state
+"""""""""""""""""
+The *integration state* (:ref:`mjSTATE_INTEGRATION<mjtState>`) is the union of all the above :ref:`mjData` fields and
+constitutes the entire set of inputs to the *forward dynamics*. The pipeline output of two :ref:`mjData` instances with
+the same integration state will be identical. In the case of *inverse dynamics*, ``mjData.qacc`` is also treated as an
+input variable. All other :ref:`mjData` fields are functions of the integration state.
+
+Note that the full integration state as given by :ref:`mjSTATE_INTEGRATION<mjtState>` is maximalist and includes fields
+which are often unused. If a small state size is desired, it might be sensible to avoid saving unused fields.
+In particular ``xfrc_applied`` can be quite large (``nbody x 6``) yet is often unused.
+
+.. _siSimulationState:
+
+Simulation state
+""""""""""""""""
+The *simulation state* is the entirety of the :ref:`mjData` struct and associated memory buffer. This state includes all
+derived quantities computed during dynamics computation. Because the :ref:`mjData` buffers are preallocated for the
+worst case, it is often significantly faster to recompute derived quantities from the :ref:`integration
+state<siIntegrationState>` rather than using :ref:`mj_copyData`. See :ref:`Notes on sleeping <siSleepNotes>` for caveats
+regarding the simulation state when sleeping is enabled.
+
+.. _siStateManipulation:
+
+State manipulation
+^^^^^^^^^^^^^^^^^^
+Manipulation of the state is facilitated by the :ref:`mjtState` bitfield enum, which enumerates the state components
+documented above. Combinations of the components, some of which are available in the enum itself, can be OR'ed together
+to form bitfield values, for example
 
 .. code-block:: C
 
-   // copy simulation state
+   int sig = mjSTATE_TIME | mjSTATE_QPOS | mjSTATE_CTRL;  // custom choice of state components
-   dst->time = src->time;
-   mju_copy(dst->qpos, src->qpos, m->nq);
-   mju_copy(dst->qvel, src->qvel, m->nv);
-   mju_copy(dst->act,  src->act,  m->na);
 
-   // copy mocap body pose and userdata
+The functions using these bitfields are :ref:`mj_getState`, :ref:`mj_setState`, :ref:`mj_copyState` and
-   mju_copy(dst->mocap_pos,  src->mocap_pos,  3*m->nmocap);
+:ref:`mj_extractState`. For example after copying the :ref:`integration state<siIntegrationState>` from an :ref:`mjData`
-   mju_copy(dst->mocap_quat, src->mocap_quat, 4*m->nmocap);
+instance ``src`` into another instance ``dst``:
-   mju_copy(dst->userdata,   src->userdata,   m->nuserdata);
-
-   // copy warm-start acceleration
-   mju_copy(dst->qacc_warmstart, src->qacc_warmstart, m->nv);
-
-Now, assuming the controls are also the same (see below) and that any installed callbacks are not relying on
-user-defined state variables that are different between src and dst, calling mj_forward(m, src) or mj_step(m, src)
-yields the same result as calling mj_forward(m, dst) or mj_step(m, dst) respectively. Similarly, calling mj_inverse(m,
-src) yields the same result as calling mj_inverse(m, dst). More on :ref:`inverse dynamics <siInverse>` later.
-
-The entire mjData can also be copied with the function :ref:`mj_copyData`. This involves
-less code but is much slower. Indeed using the above code to copy the state and then calling mj_forward to recompute
-everything can sometimes be faster than copying mjData. This is because the preallocated buffers in mjData are large
-enough to hold the intermediate results in the worst case where all possible constraints are active, but in practice
-only a small fraction of constraints tend to be active simultaneously.
-
-To illustrate how the control vector can be manipulated, suppose we want to clear all controls and applied forces
-before calling mj_step, so as to make sure we are simulating the passive dynamics (assuming no control callback of
-course). This can be done as
 
 .. code-block:: C
 
-   // clear controls and applied forces
+   mj_copyState(model, src, dst, mjSTATE_INTEGRATION);
-   mju_zero(dst->ctrl, m->nu);
-   mju_zero(dst->qfrc_applied, m->nv);
-   mju_zero(dst->xfrc_applied, 6*m->nbody);
 
-If the user has installed a control callback :ref:`mjcb_control` different from the default callback (which is a NULL
+stepping ``src`` or ``dst`` will produce identical results. States can be retrieved and set from a single
-pointer), the user callback would be expected to set some of the above fields to non-zero. Note that MuJoCo will not
+:ref:`mjtNum` array:
-clear these controls/forces at the end of the time step. This is the responsibility of the user.
+
+.. code-block:: C
+
+   int sig = mjSTATE_TIME | mjSTATE_QPOS | mjSTATE_CTRL;
+   int size = mj_stateSize(model, sig);
+   mjtNum* state = mju_malloc(size * sizeof(mjtNum));
+   mj_getState(model, src, state, sig);  // copy time, qpos and ctrl from src into state
+   mj_setState(model, dst, state, sig);  // copy time, qpos and ctrl from state into dst
+
+The entire mjData can also be copied with the function :ref:`mj_copyData` but is of course much slower than
+:ref:`mj_copyState`.
 
 Also relevant in this context is the function :ref:`mj_resetData`. It sets ``mjData.qpos`` equal to the model reference
 configuration ``mjModel.qpos0``, ``mjData.mocap_pos`` and ``mjData.mocap_quat`` equal to the corresponding fixed body
-poses from mjModel; and all other state and control variables to 0.
+poses from mjModel; and all other state and control variables to 0. When some trees are *initialized asleep*, this
+function does more work, see :ref:`sleeping <siSleepSleeping>` below.
 
 .. _siForward:
 
 Forward dynamics
 ~~~~~~~~~~~~~~~~
-
 The goal of forward dynamics is to compute the time-derivative of the state, namely the acceleration vector
-mjData.qacc and the activation time-derivative ``mjData.act_dot``. Along the way it computes everything else needed to
+``mjData.qacc`` and the activation time-derivative ``mjData.act_dot``. Along the way it computes everything else needed
-simulate the dynamics, including active contacts and other constraints, joint-space inertia and its LTDL
+to simulate the dynamics, including active contacts and other constraints, joint-space inertia and its :math:`L^TDL`
-decomposition, constraint forces, sensor data and so on. All these intermediate results are available in mjData and
+decomposition, constraint forces, sensor data and so on. All these intermediate results are available in :ref:`mjData`
-can be used in custom computations. As illustrated in the :ref:`simulation loop <siSimulation>` section above, the
+and can be used in custom computations. As illustrated in the :ref:`simulation loop <siSimulation>` section above, the
-main stepper function mj_step calls mj_forward to do most of the work, and then calls the numerical integrator to
+main stepper function :ref:`mj_step` calls :ref:`mj_forward` to do most of the work, and then calls the numerical
-advance the simulation state to the next discrete point in time.
+integrator to advance the simulation state to the next discrete point in time.
 
-The forward dynamics function mj_forward internally calls :ref:`mj_forwardSkip` with
+The forward dynamics function :ref:`mj_forward` internally calls :ref:`mj_forwardSkip` with
 skip arguments (mjSTAGE_NONE, 0), where the latter function is implemented as
 
 .. code-block:: C
@@ -396,9 +474,9 @@ skip arguments (mjSTAGE_NONE, 0), where the latter function is implemented as
        mj_sensorAcc(m, d);
    }
 
-Note that this is the same sequence of calls as in mj_step1 and mj_step2 above, except that checking of real values
+Note that this is the same sequence of calls as in :ref:`mj_step1` and :ref:`mj_step2` above, except that checking of
-and computing features such as sensor and energy are omitted. The functions being called are components of the
+real values and computing features such as sensor and energy are omitted. The functions being called are components of
-simulation pipeline. In turn they call sub-components.
+the simulation pipeline. In turn they call sub-components.
 
 The integer argument skipstage determines which parts of the computation will be skipped. The possible skip levels are
 
@@ -414,9 +492,9 @@ mjSTAGE_VEL
    stabilization.
 
 The intermediate result fields of mjData are organized into sections according to which part of the state is needed in
-order to compute them. Calling mj_forwardSkip with mjSTAGE_POS assumes that the fields in the first section (position
+order to compute them. Calling :ref:`mj_forwardSkip` with mjSTAGE_POS assumes that the fields in the first section
-dependent) have already been computed and does not recompute them. Similarly, mjSTAGE_VEL assumes that the fields in
+(position dependent) have already been computed and does not recompute them. Similarly, mjSTAGE_VEL assumes that the
-the first and second sections (position and velocity dependent) have already been computed.
+fields in the first and second sections (position and velocity dependent) have already been computed.
 
 When can we use the above machinery and skip some of the computations? In a regular simulation this is not possible.
 However, MuJoCo is designed not only for simulation but also for more advanced applications such as model-based
@@ -455,12 +533,12 @@ to full convergence, we would have
      mjData.qfrc_inverse = mjData.qfrc_applied + Jacobian'*mjData.xfrc_applied + mjData.qfrc_actuator
 
 where ``mjData.qfrc_actuator`` is the joint-space force produced by the actuators and the Jacobian is the mapping from
-joint to Cartesian space. When the "fwdinv" flag in ``mjModel.opt.enableflags`` is set, the above identity is used to
+joint to Cartesian space. When the :ref:`fwdinv<option-flag-fwdinv>` flag in ``mjModel.opt.enableflags`` is set, the
-monitor the quality of the forward dynamics solution. In particular, the two components of ``mjData.solver_fwdinv`` are
+above identity is used to monitor the quality of the forward dynamics solution. In particular, the two components of
-set to the L2 norm of the difference between the forward and inverse solutions, in terms of joint forces and
+``mjData.solver_fwdinv`` are set to the L2 norm of the difference between the forward and inverse solutions, in terms of
-constraint forces respectively.
+joint forces and constraint forces respectively.
 
-Similar to forward dynamics, ``mj_inverse`` internally calls :ref:`mj_inverseSkip` with skip arguments
+Similar to forward dynamics, :ref:`mj_inverse` internally calls :ref:`mj_inverseSkip` with skip arguments
 ``(mjSTAGE_NONE, 0)``. The skip mechanism is the same as in forward dynamics, and can be used to speed up structured
 sampling. The result ``mjData.qfrc_inverse`` is obtained by using the Recursive Newton-Euler algorithm to compute the
 net force acting on the system, and then subtracting from it all internal forces.
@@ -698,8 +776,6 @@ essential to keep it in mind at all times. All MuJoCo utility functions that ope
 :ref:`mju_mulMatMat`, :ref:`mju_mulMatVec` etc. assume this matrix layout. For vectors there is of course no
 difference between row-major and column-major formats.
 
-.. TODO(tassa): update this section when qM is migrated to CSR.
-
 When possible, MuJoCo exploits sparsity. This can make all the difference between O(N) and O(N^3) scaling. The inertia
 matrix ``mjData.qM`` and its LTDL factorization ``mjData.qLD`` are always represented as sparse. ``qM`` uses a custom
 indexing format designed for matrices that correspond to tree topology, while ``qLD`` uses the standard CSR format.
@@ -1022,6 +1098,8 @@ negative value means a tree is awake, non-negative means asleep. Maximally awake
 island are ready to sleep, they are put to sleep during state advancement and their associated values in ``tree_asleep``
 are set to a (non-negative) index cycle: the "sleeping island". If any tree in the island is woken, all are woken.
 
+.. _siSleepPolicy:
+
 Sleep policy
 ^^^^^^^^^^^^
 
@@ -1030,6 +1108,8 @@ automatically determines the :ref:`policy<mjtSleepPolicy>` to be either "allowed
 overridden using the :ref:`body/sleep <body-sleep>` attribute (see documentation therein). There is also a special
 "init" sleep policy, see next section.
 
+.. _siSleepSleeping:
+
 Sleeping
 ^^^^^^^^
 
@@ -1059,6 +1139,8 @@ Sleeping can happen in one of two ways:
   example XML that will produce a compilation error because this condition is not met. Finally, note that the
   initialized-asleep feature is only available for the default configuration (and not keyframes, see discussion below).
 
+.. _siSleepWaking:
+
 Waking
 ^^^^^^
 
@@ -1084,6 +1166,8 @@ is not always the case. For example, if free bodies on the floor are put to slee
 will remain sleeping in place until woken for another reason. The most extreme example of non-physicality are islands
 which are initialized asleep. These can be placed in mid-air or in deep collisions, but will not move until woken.
 
+.. _siSleepNotes:
+
 Notes
 ^^^^^
 

doc/python.rst

@@ -825,9 +825,9 @@ The basic usage form is
 - ``data`` is either a single instance of MjData or a sequence of compatible MjDatas of length ``nthread``.
 - ``initial_state`` is an ``nbatch x nstate`` array, with ``nbatch`` initial states of size ``nstate``, where
   ``nstate = mj_stateSize(model, mjtState.mjSTATE_FULLPHYSICS)`` is the size of the
-  :ref:`full physics state<geFullPhysics>`.
+  :ref:`full physics state<siFullPhysics>`.
 - ``control`` is a ``nbatch x nstep x ncontrol`` array of controls. Controls are by default the ``mjModel.nu`` standard
-  actuators, but any combination of :ref:`user input<geInput>` arrays can be specified by passing an optional
+  actuators, but any combination of :ref:`user input<siInput>` arrays can be specified by passing an optional
   ``control_spec`` bitflag.
 
 If a rollout diverges, the current state and sensor values are used to fill the remainder of the trajectory.
@@ -1084,7 +1084,7 @@ non-exhaustive list of specific mujoco-py features:
    The MuJoCo library’s computation is deterministic given a specific input, as explained in the :ref:`Programming
    section <Simulation>`. mujoco-py implements methods for getting and setting some of the relevant fields (and
    similarly ``dm_control.Physics`` offers methods that correspond to the flattened case). This functionality is
-   described in the :ref:`state <geState>` section.
+   described in the :ref:`State and Control<siStateControl>` section.
 
 ``sim.model.get_joint_qvel_addr(joint_name)``
    This is a convenience method in mujoco-py that returns a list of contiguous indices corresponding to this joint. The