diff --git a/source/user/manual/analysis/integrator.rst b/source/user/manual/analysis/integrator.rst
index 5af88231..d1beedfc 100644
--- a/source/user/manual/analysis/integrator.rst
+++ b/source/user/manual/analysis/integrator.rst
@@ -44,6 +44,8 @@ Transient Integrators:
    integrator/GeneralizedAlpha
    integrator/TRBDF2
    integrator/ExplicitDifference
+   integrator/ExplicitBathe
+   integrator/ExplicitDifferenceStatic
 
 
 Utility Integrators:
diff --git a/source/user/manual/analysis/integrator/ExplicitBathe.rst b/source/user/manual/analysis/integrator/ExplicitBathe.rst
new file mode 100644
index 00000000..e92bc407
--- /dev/null
+++ b/source/user/manual/analysis/integrator/ExplicitBathe.rst
@@ -0,0 +1,57 @@
+.. _ExplicitBathe:
+
+ExplicitBathe
+-------------
+.. function:: integrator ExplicitBathe $p $compute_critical_timestep?
+
+.. csv-table:: 
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $p, |float|, Damping parameter (0 < p < 1, typically 0.54)
+   $compute_critical_timestep, |int|, Optional flag to compute critical timestep (0 or 1)
+
+.. note::
+   * The method is second-order accurate and explicit.
+   * Only mass matrix is assembled on RHS (no tangent matrix required).
+   * Critical time step is approximately 2× larger than central difference method.
+   * For stability: :math:`\Delta t \leq \frac{2}{\omega_{max}}`
+   * Typical p values: 0.54-0.95 for optimal numerical damping.
+
+Theory
+^^^^^^
+
+The ExplicitBathe method is a two-step explicit integration scheme with built-in numerical damping. The method performs two sub-steps per time step:
+
+1. First step: :math:`t \rightarrow t + p\Delta t`
+2. Second step: :math:`t + p\Delta t \rightarrow t + \Delta t`
+
+The integration coefficients are computed from the damping parameter p:
+
+.. math::
+
+   q_1 = \frac{1 - 2p}{2p(1-p)}
+
+   q_2 = 0.5 - p \cdot q_1
+
+   q_0 = -q_1 - q_2 + 0.5
+
+The method offers enhanced stability compared to standard central difference, with a critical time step approximately twice as large. The parameter p controls numerical damping, with higher values providing more damping but reduced accuracy.
+
+.. admonition:: Example 
+
+   1. **Tcl Code**
+
+   .. code-block:: tcl
+
+      integrator ExplicitBathe 0.54
+
+   2. **Python Code**
+
+   .. code-block:: python
+
+       integrator('ExplicitBathe', 0.54)
+
+.. [Noh2013] Noh, G., & Bathe, K.J. (2013). "An explicit time integration scheme for the analysis of wave propagations." Computers & Structures, 129, 178-193.
+
+Code Developed by: |jaabell|
diff --git a/source/user/manual/analysis/integrator/ExplicitDifferenceStatic.rst b/source/user/manual/analysis/integrator/ExplicitDifferenceStatic.rst
new file mode 100644
index 00000000..890df40a
--- /dev/null
+++ b/source/user/manual/analysis/integrator/ExplicitDifferenceStatic.rst
@@ -0,0 +1,57 @@
+.. _ExplicitDifferenceStatic:
+
+ExplicitDifferenceStatic
+------------------------
+.. function:: integrator ExplicitDifferenceStatic
+
+.. note::
+   * Uses leap-frog integration with velocities at half time steps.
+   * FLAC-style local non-viscous damping (α = 0.59) for pseudo-static analysis.
+   * Only mass matrix is required (no tangent matrix assembly).
+   * Velocity sign memory prevents chatter near zero velocity.
+   * For stability: :math:`\Delta t \leq \frac{2}{\omega_{max}}`
+   * Suitable for quasi-static and dynamic problems with adaptive damping.
+
+Theory
+^^^^^^
+
+The ExplicitDifferenceStatic method is an explicit difference scheme with FLAC-style local non-viscous damping. The method uses a leap-frog approach where:
+
+* Velocities are defined at half time steps: :math:`v_{n+1/2}`
+* Displacements are defined at full time steps: :math:`d_n`
+
+The damping force uses an adaptive FLAC-style approach:
+
+.. math::
+
+   F_{damping} = -\alpha \cdot |F_{unbalanced}| \cdot sign(v)
+
+where α = 0.59 by default. The method includes velocity sign memory with a deadband (1e-4) to prevent numerical chatter when velocities approach zero.
+
+The central difference equations are:
+
+.. math::
+
+   v_{n+1/2} = v_{n-1/2} + \Delta t \cdot a_n
+
+   d_{n+1} = d_n + \Delta t \cdot v_{n+1/2}
+
+This integrator is particularly effective for quasi-static analysis and problems where adaptive local damping is beneficial, such as rock mechanics and soil-structure interaction problems.
+
+.. admonition:: Example 
+
+   1. **Tcl Code**
+
+   .. code-block:: tcl
+
+      integrator ExplicitDifferenceStatic
+
+   2. **Python Code**
+
+   .. code-block:: python
+
+       integrator('ExplicitDifferenceStatic')
+
+.. [Cundall1987] Cundall, P.A. (1987). "Distinct Element Models of Rock and Soil Structure." In Analytical and Computational Methods in Engineering Rock Mechanics, 129-163.
+
+Code Developed by: |jaabell|
diff --git a/source/user/manual/analysis/system.rst b/source/user/manual/analysis/system.rst
index 47acd1c7..89c15c7e 100644
--- a/source/user/manual/analysis/system.rst
+++ b/source/user/manual/analysis/system.rst
@@ -22,3 +22,4 @@ The following contain information about systemType? and the args required for ea
    system/SparseSYM
    system/Mumps
    system/Cusp
+   system/Diagonal
diff --git a/source/user/manual/analysis/system/Diagonal.rst b/source/user/manual/analysis/system/Diagonal.rst
new file mode 100644
index 00000000..a62de346
--- /dev/null
+++ b/source/user/manual/analysis/system/Diagonal.rst
@@ -0,0 +1,129 @@
+Diagonal & MPIDiagonal System
+-----------------------------
+
+This command is used to construct a Diagonal linear system of equation object. This system stores only the diagonal entries of the coefficient matrix, making it extremely memory efficient for problems where off-diagonal coupling can be neglected or lumped into the diagonal. The system is solved by direct inversion of the diagonal entries. The following command is used to construct such a system:
+
+.. function:: system (Diagonal|MPIDiagonal) <-lumped>
+
+   **Optional Parameter:**
+
+   * **-lumped** - If specified, off-diagonal matrix entries are added (lumped) to the diagonal before solving. This is useful for mass matrix lumping or when converting a coupled system to a diagonal approximation.
+
+A diagonal system stores only the diagonal entries of an n×n matrix **A**, where:
+
+:math:`K = \begin{bmatrix} k_{11} & 0 & \cdots & 0 \\ 0 & k_{22} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & k_{nn} \end{bmatrix}`
+
+The solution for :math:`K u = f` is obtained directly as:
+
+:math:`u_i = \frac{f_i}{k_{ii}} \quad \text{for } i = 1, \ldots, n`
+
+.. note::
+
+   1. This solver requires that all diagonal entries are non-zero.
+   2. The system only stores n values instead of n² for a full matrix, providing significant memory savings
+   3. When using the **-lumped** option, the solver accumulates all off-diagonal entries in each row/column into the corresponding diagonal entry
+
+Mass Lumping
+^^^^^^^^^^^^
+
+When the **-lumped** option is specified, the Diagonal or MPIDiagonal system performs **mass lumping** (also called diagonal lumping or row-sum lumping). This technique converts a consistent (coupled) mass matrix into a diagonal form by summing all entries in each row and placing the total on the diagonal.
+
+**Mathematical Formulation:**
+
+For a consistent mass matrix **M** with entries :math:`m_{ij}`, the lumped diagonal entry is:
+
+:math:`\tilde{m}_{ii} = \sum_{j=1}^{n} m_{ij}`
+
+All off-diagonal entries are set to zero: :math:`\tilde{m}_{ij} = 0` for :math:`i \neq j`
+
+   \end{bmatrix}
+
+
+.. warning::
+
+   1. This system is only appropriate when the problem structure allows diagonal treatment
+   2. Mass lumping is generally only recommended for the mass matrix in dynamic analysis (ExplicitDifference integrator or CentralDifference integrator if no damping matrix is used or only mass proportional Rayleigh damping is used)
+   3. Lumping the stiffness matrix is rarely appropriate and can lead to poor results
+   4. For static analysis or implicit dynamics use full sparse solvers.
+   5. Ideally, lumping should be handled at the element level, especially for high-order elements (like the TenNodeTetrahedron or SixNodeTriangle) where this technique leads to negative mass values (!!). 
+
+**Implementation Details:**
+
+When **-lumped** is specified:
+
+.. code-block:: none
+
+   For each element matrix entry m(i,j) at DOF locations id(i), id(j):
+      A[id(i)] += m(i,i)           // Add diagonal entry
+      if (lumped):
+         for all j ≠ i:
+            A[id(i)] += m(j,i)     // Add all column entries to diagonal
+
+This ensures that the total contribution from each element is preserved while creating a diagonal system.
+
+
+.. admonition:: Example 
+
+   The following examples show how to construct a Diagonal system
+
+   **1. Tcl Code - Basic diagonal system (no lumping)**
+
+   .. code-block:: tcl
+
+      system Diagonal
+
+   **2. Tcl Code - Diagonal system with mass lumping for explicit dynamics**
+
+   .. code-block:: tcl
+
+      # Explicit dynamic analysis with lumped mass
+      system Diagonal -lumped
+      
+      constraints Plain
+      numberer Plain
+      test NormDispIncr 1.0e-6 10 0
+      algorithm Linear
+      integrator CentralDifference
+      analysis Transient
+
+   **3. Python Code - Basic diagonal system**
+
+   .. code-block:: python
+
+      ops.system('Diagonal')
+
+   **4. Python Code - Lumped mass for explicit analysis**
+
+   .. code-block:: python
+
+      # Setup for explicit central difference method
+      ops.system('Diagonal', '-lumped')
+      ops.constraints('Plain')
+      ops.numberer('Plain')
+      ops.integrator('CentralDifference')
+      ops.analysis('Transient')
+
+   **5. Python Code - Complete explicit dynamics example**
+
+   .. code-block:: python
+
+      import openseespy.opensees as ops
+      
+      # ... model definition ...
+      
+      # Analysis setup for blast/impact with lumped mass
+      ops.system('Diagonal', '-lumped')
+      ops.constraints('Plain')
+      ops.numberer('Plain')
+      ops.test('NormDispIncr', 1.0e-6, 10, 0)
+      ops.algorithm('Linear')
+      ops.integrator('CentralDifference')
+      ops.analysis('Transient')
+      
+      # Time step (must satisfy CFL condition)
+      dt = 0.0001
+      ops.analyze(1000, dt)
+
+Code Developed by: |fmk|
+
+
diff --git a/source/user/manual/material/matTestCommands.rst b/source/user/manual/material/matTestCommands.rst
index 99f5179c..2df39b43 100644
--- a/source/user/manual/material/matTestCommands.rst
+++ b/source/user/manual/material/matTestCommands.rst
@@ -82,4 +82,59 @@ The command getResponse returns material or section specific information.
    :widths: 10, 10, 40
 
    $args,  |list|, list of the arguments for the material/section response
+    
+NDTest Commands for Direct Material Testing
+===========================================
+
+For more advanced and direct testing of nDMaterial objects, OpenSees provides the :ref:`NDTest Command <NDTest>` which offers additional functionality beyond the basic material testing environment. 
+
+The NDTest command provides comprehensive testing capabilities including:
+
+- Direct strain control with ``NDTest SetStrain``
+- State management with ``NDTest CommitState``
+- Parameter updates during runtime
+- Enhanced response queries
+- Real-time material state inspection
+
+Key differences between the basic material testing environment and NDTest:
+
+- **Basic Testing**: Uses ``testNDMaterial`` with ``setStrain``/``getStress`` commands
+- **NDTest**: Uses ``NDTest`` with more extensive command set and real-time parameter updates
+
+For complete documentation of NDTest commands, see :ref:`NDTest`.
+
+Example comparing both approaches:
+
+TCL:
+
+.. code-block:: tcl
+
+   # Basic testing approach
+   testNDMaterial 1
+   setStrain 0.001 0.0 0.0 0.0 0.0 0.0
+   set stress [getStress]
+   
+   # NDTest approach  
+   NDTest SetStrain 1 0.001 0.0 0.0 0.0 0.0 0.0
+   set stress [NDTest GetStress 1]
+   NDTest CommitState 1
+   # Additional capabilities:
+   set tangent [NDTest GetTangentStiffness 1]
+   NDTest UpdateDoubleParameter 1 1 2.1e-5
+
+Python:
+
+.. code-block:: python
+
+   # Basic testing approach
+   op.testNDMaterial(1)
+   op.setStrain(0.001, 0.0, 0.0, 0.0, 0.0, 0.0)
+   stress = op.getStress()
    
+   # NDTest approach  
+   op.NDTest("SetStrain", 1, 0.001, 0.0, 0.0, 0.0, 0.0, 0.0)
+   stress = op.NDTest("GetStress", 1)
+   op.NDTest("CommitState", 1)
+   # Additional capabilities:
+   tangent = op.NDTest("GetTangentStiffness", 1)
+   op.NDTest("UpdateDoubleParameter", 1, 1, 2.1e-5)
diff --git a/source/user/manual/material/ndMaterial.rst b/source/user/manual/material/ndMaterial.rst
index d2d7d207..d157ea3a 100644
--- a/source/user/manual/material/ndMaterial.rst
+++ b/source/user/manual/material/ndMaterial.rst
@@ -46,6 +46,7 @@ The following contain information about matType? and the args required for each
    ndMaterials/ASDPlasticMaterial
    ndMaterials/OrthotropicRAConcrete
    ndMaterials/SmearedSteelDoubleLayer
+   ndMaterials/NDTest
 
 Concrete Damage Model
 Plane Stress Material
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial.rst
deleted file mode 100644
index a402421d..00000000
--- a/source/user/manual/material/ndMaterials/ASDPlasticMaterial.rst
+++ /dev/null
@@ -1,278 +0,0 @@
-.. _ASDPlasticMaterial:
-
-ASDPlasticMaterial
-^^^^^^^^^^^^^^^^^^
-
-| This command is used to construct an ``ASDPlasticMaterial`` material object. ``ASDPlasticMaterial`` implements a large family of constitutive models based on the classical theory of elastoplasticity. Users build new constitutive models by selecting the yield function, plastic-flow direction, elasticity law, and hardening models for the internal variables from several possible options for each component. 
-
-To create a new model, specify the Yield Function type (``$YieldFunctionType``), Plastic Flow type (``$PlasticFlowType``), and Elasticity type (``$ElasticityType``). All internal variables and model parameters are set to zero unless specified by the user, which might or not make sense depending on context, and this initialization is printed out to the screen. 
-
-After setting the ``$YieldFunctionType``, ``$PlasticFlowType`` and ``$ElasticityType``, you can give initial values to internal variables within the ``Begin_Internal_Variables`` ... ``End_Internal_Variables`` block. Then, the ``Begin_Model_Parameters`` ... ``End_Model_Parameters`` block is used to provide model parameter values (these can be changed during the analysis with the ``setParameter`` command as expected). Finally, model integration options are set within the ``Begin_Integration_Options`` ... ``End_Integration_Options`` code block. Specification blocks can occur in any order or ommitted. 
-
-The complete command looks as follows::
-
-   nDMaterial ASDPlasticMaterial $tag 
-      $YieldFunctionType 
-      $PlasticFlowType $
-      ElasticityType 
-      $IV_TYPE 
-      Begin_Internal_Variables 
-         $InternalVariable1 $$double_value1 $$double_value2... $$double_valueN1  
-         $InternalVariable2 $$double_value1 $$double_value2... $$double_valueN2 
-         #... (depends on how many internal variables the particular selected model has) 
-      End_Internal_Variables 
-      Begin_Model_Parameters 
-         $ModelParameters1 $$double_value1 
-         $ModelParameters2 $$double_value2 
-         #... (depends on how many model parameters the particular selected model has) 
-      End_Model_Parameters 
-      Begin_Integration_Options 
-         f_relative_tol $double_value 
-         stress_relative_tol $double_value 
-         n_max_iterations $int_value 
-         return_to_yield_surface (0 or 1) 
-         method (string) : Forward_Euler | Runge_Kutta_45_Error_Control 
-      End_Integration_Options
-
-
-Explanation
-
-.. csv-table:: 
-   :header: "Argument", "Type", "Description"
-   :widths: 10, 10, 40
-
-   $tag, |integer|, "Unique tag identifying this material."
-   $YieldFunctionType, |string|, "Mandatory. Yield function to be used -> :ref:`YieldFunctionType`"
-   $PlasticFlowType, |string|, "Mandatory. Plastic flow direction to be used -> :ref:`PlasticFlowType`"
-   $ElasticityType, |string|, "Mandatory. Elastic model to be used -> :ref:`ElasticityType`"
-   $IV_TYPE, |string|, "Mandatory. Hardening model for internal variables. Admitted types depend on YF, PF, and EL chosen.  "
-   Begin_Internal_Variables, |string|, "Optional. Marks the beginning of the code block to set the internal variables. If ommitted, all internal variables are initialized to zero. You can specify as many of the following variables as wanted. Ommitted variables are initialized to zero. Number and name of variables that can be set is model dependent (once YF, PF, and EL are specified)"
-   $InternalVariable1, |list of name/value pairs|, "Initial value of internal variable1. Dimension depends on internal variable type.  "
-   $InternalVariable2, |list of name/value pairs|, "Initial value of internal variable2.  "
-   "--", "--", "... (input as many as model supports)"
-   End_Internal_Variables, |string|, "Mandatory if block started. Marks the end of the code block to set the internal variables"
-   Begin_Model_Parameters, |string|, "Optional. Marks the beginning of the code block to set the model parameters"
-   $ModelParameters, |list of name/value pairs|, "Values for parameters of the models to be used. This depends on the particular choices of ``$YieldFunctionType``, ``$PlasticFlowType``, ``$ElasticityType``, and ``$IV_type``. "
-   "--", "--", "... (input as many as model supports)"
-   End_Model_Parameters, |string|, "Mandatory if block started. Marks the beginning of the code block to set the model parameters"
-   Begin_Integration_Options, |string|, "Optional. Marks the beginning of the code block to set the integration options. You can set any ammount "
-   End_Integration_Options, |string|, "Mandatory if block started. Marks the beginning of the code block to set the model parameters"
-
-|  The :ref:`YieldFunctionType`, :ref:`PlasticFlowType`, and :ref:`ElasticityType` have different variables and parameters to be set.
-
-.. admonition:: Specification of internal variables (``$IV_TYPE``)
-
-	All internal variables are specified in a single string with no spaces. Internal variables specifications are separated with a colon ``:``. After the name of the internal variable, the name of the hardening function must be provided in parenthesis. The specification string must end in a colon. The syntax is, thus:
-
-	.. code-block::
-
-		InternalVariable1(HardeningFunction1):InternalVariable2(HardeningFunction2):
-
-	Internal variables are required for the yield function and plastic flow direction (see speficic components for details), and all internal variables must be provided with their hardening function otherwise instantiation fails. 
-
-.. admonition:: Note
-
-    Any internal variable or parameter that is not specified is set to zero by default.  
-
-
-.. toctree::
-   :caption: Arguments detailed description
-   :maxdepth: 2
-
-   ./ASDPlasticMaterial/YieldFunctions
-   ./ASDPlasticMaterial/PlasticFlowType
-   ./ASDPlasticMaterial/ElasticityType
-   ./ASDPlasticMaterial/HardeningFunctions
-
-.. _ASDPlasticTheory:
-
-ASDPlasticMaterial Theory
-""""""""""""""""""""""""""
-
-The ``ASDPlasticMaterial`` internally uses Voigt notation to represent tensors. 
-
-.. math::
-   \newcommand{\vec}[1]{\boldsymbol{#1}}
-   \newcommand{\state}{\sigma, \left\lbrace iv \right\rbrace, \left\lbrace param \right\rbrace }
-   \newcommand{\matorvec}[2]{
-    \left[\begin{array}{#1}
-        #2
-    \end{array}\right]
-    }
-   \vec{\epsilon} = \matorvec{c}{\epsilon_{11} \\ \epsilon_{22} \\ \epsilon_{33} \\ \gamma_{12} \\ \gamma_{23} \\ \gamma_{13} }
-   \vec{\sigma} = \matorvec{c}{\sigma_{11} \\ \sigma_{22} \\ \sigma_{22} \\ \tau_{12} \\ \tau_{23} \\ \tau_{13} }
-
-Where :math:`\gamma_{ij} = 2 \epsilon_{ij}`. 
-
-The material recieves the total (trial) strain :math:`\vec{\epsilon}^{\text{trial}}` from the finite element that contains it. From this, the trial strain increment is computed by subtracting the previously committed total strain:
-
-.. math::
-    \Delta \vec{\epsilon}^{\text{trial}} = \vec{\epsilon}^{\text{trial}} - \vec{\epsilon}^{\text{commit}}
-
-The elastic trial stress increment is computed using the elasticity law (depending on choice of ``$ElasticityType``), which may depend on stress (:math:`\sigma`) internal variables (represented as the list of internal variables of the model :math:`\left\lbrace iv \right\rbrace`) and parameters (represented as the list of parameters of the model :math:`\left\lbrace param \right\rbrace`):
-
-.. math::
-    \Delta \vec{\sigma}^{\text{trial}} = \vec{E}^{\text{current}}(\state) \Delta \vec{\epsilon}^{\text{trial}}
-
-:math:`\vec{E}^{\text{current}}(\state )` is the 6x6 elastic tangent operator. The trial strain is computed as :math:`\mathbf{\sigma}^{\text{trial}} = \mathbf{\sigma}^{\text{commit}} + \Delta \sigma^{\text{trial}}`. With this predictor, the yield-function (:math:`f(\state )`) is evaluated to see whether the stress point falls inside or outside the yield surface (depending on choice of ``$YieldFunctionType``). If it falls inside, then the step is elastic and the integration stops. Otherwise, we need to do some plasticity for which the explicit Forward-Euler version will be used to present all components. 
-
-The exact intersection point with the yield function is computed iteratively with a robust Brent algorithm. The trial stress increment is advanced to the yield-surface intersection point and plastic integration continues. Starting from :math:`\Delta \epsilon^{\text{trial}}` corresponding to a :math:`\Delta \sigma^{\text{trial}}`  on the yield surface, the strain increment is separated into elastic and plastic components (the "trial" qualifier es removed and implied in what follows):
-
-.. math::
-    \Delta \vec{\epsilon}   = \Delta \epsilon^{el} + \Delta \epsilon^{pl}
-
-The plastic increment comes from:
-
-.. math::
-    \Delta \vec{\epsilon}^{pl}   = \Delta \lambda \cdot \mathbf{m} (\state)
-
-Where the vector :math:`\mathbf{m} (\state)` is the plastic flow direction defined with the ``$PlasticFlowType`` option and :math:`\Delta \lambda\` is the plastic multiplier defined as:
-
-.. math::
-    \Delta \lambda\ = \dfrac{ \vec{n}^T E \Delta \vec{\epsilon} }{ \vec{n}^T E \Delta \vec{m} - H }
-
-Where :math:`\mathbf{n} (\state)` is the outward normal to the yield surface, and :math:`H` is the hardening term due to internal variable evolution. For example, if the model has two internal variables :math:`s` a scalar internal variable and  :math:`\vec{T}` a tensor internal variable (e.g. :math:`\left\lbrace iv \right\rbrace = \left\lbrace s,\, \vec{T} \right\rbrace`), the hardening term is computed with:
-
-.. math::
-    H = \left.\dfrac{\partial f}{\partial s}\right\Vert_{(\state)} h_s(\state) + \left.\dfrac{\partial f}{\partial \vec{T}}\right\Vert_{(\state)} h_s(\state) \cdot \vec{h}_T(\state) 
-
-Here the hardening functions of the scalar internal varialbe (a scalar) :math:`h_s(\state)` and the tensor internal variable :math:`\vec{h}_T(\state)` are introduced. These are specified using the ``$IV_type`` option. 
-
-Finally, the trial stress increment considering plasticity is computed using:
-
-.. math::
-    \Delta \vec{\sigma}^{\text{trial}} = \vec{E}(\state) \cdot \Delta \vec{\epsilon}^{pl}
-
-And the new trial state of the material is computed:
-
-.. math::
-    \vec{\sigma}^{\text{trial}} = \vec{\sigma}^{\text{commit}} + \Delta \vec{\sigma}^{\text{trial}}
-
-    s^{\text{trial}} = s^{\text{commit}} + \Delta \lambda  h_s(\state)
-
-    \vec{T}^{\text{trial}} = \vec{T}^{\text{commit}} + \Delta \lambda  \vec{h}_T(\state)
-
-After this, an aditional algorithm can be called to return the stress point to the yield surface (recommended) as well as specific provisions in case integration is ocurring near an apex of the yield surface. 
-
-
-.. _`Integration Options`:
-Integration Options
-"""""""""""""""""""
-
-.. csv-table:: 
-   :header: "Parameter", "Type", "Description"
-   :widths: 10, 10, 40
-
-   $f_relative_tol, |double|, "Relative tolerance to evaluate the yield function crossing."
-   $stress_relative_tol, |double|, "Tolerance for the convergece of the integration algorithm."
-   $n_max_iterations, |int|, "Maximum number of iterations for constitutive integration."
-   $return_to_yield_surface, |0 or 1|, "Whether to apply a return to yield surface algorithm after integration convergence."
-   $method, |string|, "Constitutive integration method. Options: ``Forward_Euler``, ``Runge_Kutta_45_Error_Control``"
-
-The default integration method is **Runge_Kutta_45_Error_Control** that uses the classical RK45 ODE integration algorithm employing a 4-th order prediction of the stress increment together with a 5-order prediction to estimate the integration error. In this scheme the strain increment provided by the element to the Gauss point is sub-divided into sub-increments, a process which is automated such that the provided **$stress_relative_tol** is met. This method is provided as a robust standard method which is applicable across all possible combinations of components, although there are possibly better approaches for specific cases which might become available in the future. 
-
-The different parameters are activated depending on the integration algorithm selected. The *Forward_Euler* algorithm only uses the **$return_to_yield_surface** parameter, while **Runge_Kutta_45_Error_Control** uses the rest. 
-
-The **$f_relative_tol** parameter comes into play when the yield surface is being crossed, that is, when the previous (committed) stress is within the yield surface and the elastic prediction of the stress increment brings the stress state beyond the yield surface. In that case, an elastic increment occurs until the yield surface is touched which requires iterations with the Brent root finding algorithm. This is used by both currently available integration methods. 
-
-
-.. _`Other Features`:
-
-Other Features
-""""""""""""""""""
-
-.. admonition:: General parameters
-
-    These parameters are defined for all models. 
-
-    .. csv-table:: 
-       :header: "Parameter", "Type", "Description"
-       :widths: 10, 10, 40
-
-       ``MassDensity``, scalar, Defines the material mass density :math:`\rho`. 
-       ``InitialP0``, scalar, Defines the initial mean pressure at which material constants will be evaluated at the first step.  
-
-
-.. admonition:: Responses  (`setResponse` and `getResponse` behavior)
-
-    * *Valid queries for recorders*. `stresses` for stress, `strains` for strains, and `pstrains` for plastic strains. 
-    * You can also request any and all internal variables by their specific name as an output. 
-
-.. admonition:: `setParameter` behavior 
-
-   * Not yet provided. 
-
-Implementation details
-""""""""""""""""""""""
-
-``ASDPlasticMaterial`` is implemented using C++ template metaprogramming, with a header-only design and using the "eigen" C++ library for high-performance array operations. This design provides modularity and the ability to mix and match components to create new models, while also providing high-performance because runtime polymorphism is avoided. 
-
-
-Example
-"""""""
-
-The following example defines an instance of ``ASDPlasticMaterial`` with a Drucker-Prager yield function, a constant dilatancy plastic-flow direction, elastic-isotropic elasticity law and linear hardening for both internal variables. 
-
-.. admonition:: TCL code  
-
-    .. code-block:: tcl
-
-        nDMaterial ASDPlasticMaterial 1 \
-            DruckerPrager_YF \
-            ConstantDilatancy_PF \
-            LinearIsotropic3D_EL \
-            BackStress(TensorLinearHardeningFunction):VonMisesRadius(ScalarLinearHardeningFunction): \
-            Begin_Internal_Variables \
-                VonMisesRadius 1. \
-                BackStress 0. 0. 0. 0. 0. 0. \
-            End_Internal_Variables \
-            Begin_Model_Parameters \
-                YoungsModulus 1. \
-                PoissonsRatio 0. \
-                TensorLinearHardeningParameter 0. \
-                ScalarLinearHardeningParameter 0. \
-                Dilatancy 0.02 \
-                MassDensity 2000. \
-            End_Model_Parameters 
-        
-
-.. admonition:: Python code  
-
-    .. code-block:: python
-
-        ops.nDMaterial("ASDPlasticMaterial", 1, 
-        "DruckerPrager_YF",
-        "ConstantDilatancy_PF",
-        "LinearIsotropic3D_EL",
-        "BackStress(TensorLinearHardeningFunction):VonMisesRadius(ScalarLinearHardeningFunction):",
-        "Begin_Internal_Variables",
-            "VonMisesRadius", 1.,
-            "BackStress", 0., 0., 0., 0., 0., 0.,
-        "End_Internal_Variables",
-        "Begin_Model_Parameters",
-            "YoungsModulus", 1.,
-            "PoissonsRatio", 0.,
-            "TensorLinearHardeningParameter", 0.,
-            "ScalarLinearHardeningParameter", 0.,
-            "Dilatancy", 0.02,
-            "MassDensity", 2000.,
-        "End_Model_Parameters",
-        )
-
-If you subject this material to strain-controlled cyclic loading you get the following response.
-
-.. figure:: ASDPlasticMaterial/example1_fig2.png
-  :width: 50%
-  :alt: stress-strain response of the material
-
-  (Left) strain response with input strain history in blue and plastic strains in orange. (Right) stress response. 
-
-Plotted in principal-stress space you can see the material soften as it dilates. 
-
-.. figure:: ASDPlasticMaterial/example1_fig1.png
-  :width: 50%
-
-  Response of the material in principal-stress space with Drucker-Prager surface for reference. 
-
-
-Code Developed by: **José A. Abell** (UANDES, Chile and ASDEA), **Guido Camata** and **Massimo Petracca**  (ASDEA Software, Italy).
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/ElasticityType.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/ElasticityType.rst
index 525d040c..77cc88f8 100644
--- a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/ElasticityType.rst
+++ b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/ElasticityType.rst
@@ -1,16 +1,17 @@
 .. _ElasticityType:
 
 Elasticity Types
-^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^
 
 These components define the stress-strain relationship for the linear part of the model (see :ref:`ASDPlasticTheory`). 
-
 .. math::
    \newcommand{\vec}[1]{\boldsymbol{#1}}
    \vec{\sigma} = \vec{E} \vec{\epsilon}
-
 Where the operator :math:`\vec{E}` may depend on the material state, parameters, etc. 
 
+.. note::
+   When using :ref:`ASDPlasticMaterial3D`, elastic parameters can be updated during analysis using the ``setParameter`` command, allowing for simulation of materials with time-dependent or state-dependent elastic properties. This is particularly useful for modeling damage effects, temperature-dependent stiffness, or staged construction processes.
+
 
 .. _LinearIsotropic3D_EL:
 LinearIsotropic3D_EL
@@ -28,8 +29,24 @@ A classic!
       :header: "Parameter name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``YoungsModulus``, scalar, :math:`E`, Young's modulus
-      ``PoissonsRatio``, scalar, :math:`\nu`, Poisson's ratio
+       ``YoungsModulus``, scalar, :math:`E`, Young's modulus (Pa or psi)
+       ``PoissonsRatio``, scalar, :math:`\nu`, Poisson's ratio (dimensionless, typically 0.0-0.5)
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   LinearIsotropic3D_EL provides the fundamental elastic behavior for most engineering materials. With :ref:`ASDPlasticMaterial3D`:
+   
+   * Constant elastic properties throughout analysis
+   * Efficient computation with closed-form tangent
+   * Suitable for metals, concrete at small strains, and other linear elastic materials
+   
+   Monitor elastic behavior:
+   
+   .. code-block:: tcl
+   
+       # Update elastic parameters during analysis
+       setParameter 1 YoungsModulus 3e10
+       setParameter 1 PoissonsRatio 0.2
 
 .. _DuncanChang_EL:
 DuncanChang_EL
@@ -48,8 +65,29 @@ Where :math:`E_{ref}` (dimensionless) specifies the Young's modulus at reference
       :header: "Parameter name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``ReferenceYoungsModulus``, scalar, :math:`E_{ref}`, Dimensionless reference Young's modulus at reference pressure :math:`p_{ref}`.
-      ``PoissonsRatio``, scalar, :math:`\nu`, Poisson's ratio
-      ``ReferencePressure``, scalar, :math:`p_{ref}`, Reference pressure for the definition of Young's modulus.
-      ``DuncanChang_MaxSigma3``, scalar, :math:`\sigma_{3max}`, Maximum confinement stress.
-      ``DuncanChang_n``, scalar, :math:`\nu`, Exponent for Duncan-Chang law. 
+       ``ReferenceYoungsModulus``, scalar, :math:`E_{ref}`, Dimensionless reference Young's modulus at reference pressure :math:`p_{ref}`.
+       ``PoissonsRatio``, scalar, :math:`\nu`, Poisson's ratio (dimensionless)
+       ``ReferencePressure``, scalar, :math:`p_{ref}`, Reference pressure for the definition of Young's modulus (Pa).
+       ``DuncanChang_MaxSigma3``, scalar, :math:`\sigma_{3max}`, Maximum confinement stress (Pa).
+       ``DuncanChang_n``, scalar, :math:`n`, Exponent for Duncan-Chang law (typically 0.2-1.0).
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   DuncanChang_EL provides pressure-dependent elastic behavior, particularly useful for:
+   
+   * Soil materials with stiffness increasing with confinement
+   * Rock materials exhibiting non-linear elastic response
+   * Interfaces and joints with pressure-dependent behavior
+   
+   The Young's modulus varies with minimum principal stress :math:`\sigma_3` (compression negative):
+   
+   .. math::
+      E(\sigma_3) = E_{ref} \cdot p_{ref} \cdot \left(\dfrac{\vert \sigma_3 \vert}{ p_{ref}} \right)^n
+   
+   Monitor pressure-dependent response:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file pressure.out -time -node 1 -dof 1 PStress
+       recorder Node -file j2_stress.out -time -node 1 -dof 1 J2Stress
+
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/HardeningFunctions.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/HardeningFunctions.rst
index f37ab15a..0b4780aa 100644
--- a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/HardeningFunctions.rst
+++ b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/HardeningFunctions.rst
@@ -1,10 +1,9 @@
 .. _`HardeningFunctions`:
 
 Hardening Functions
-^^^^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^^^^
 
 Internal variables evolve accoding with their own laws to provide plastic hardening to the constitutive model. These laws are specified as hardening functions that specify an ODE for the evolution of the internal variable. Two types of internal variable are considered: scalar-valued and tensor-valued internal variables. These evolve in the following manner:
-
 .. math::
    \newcommand{\vec}[1]{\boldsymbol{#1}}
    \newcommand{\state}{\sigma, \left\lbrace iv \right\rbrace, \left\lbrace param \right\rbrace }
@@ -13,7 +12,10 @@ Internal variables evolve accoding with their own laws to provide plastic harden
 
    \vec{T}^{\text{trial}} = \vec{T}^{\text{commit}} + \Delta \lambda  \vec{h}_T(\state)
 
-Where :math:`s` is a scalar-valued internal variable and :math:`\vec{T}` is tensor-valued. These functions are used to define the hardening term in the plastic multiplier (see :ref:`ASDPlasticTheory`).
+Where :math:`s` is a scalar-valued internal variable and :math:`\vec{T}` is a tensor-valued. These functions are used to define the hardening term in the plastic multiplier (see :ref:`ASDPlasticTheory`). 
+
+.. note::
+   When using :ref:`ASDPlasticMaterial3D`, hardening function parameters can be updated during analysis using ``setParameter``. Internal variables can be accessed by name for monitoring and recording. This enables real-time tracking of hardening evolution and allows for modeling of complex material behavior including state-dependent hardening, temperature effects, or damage-coupled hardening laws.
 
 Available functions:
 
@@ -36,10 +38,26 @@ Where :math:`\vec{m}` is the plastic flow direction computed at the current stat
 .. admonition:: Parameters required
 
    .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
+      :header: "Parameter name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``ScalarLinearHardeningParameter``, scalar, :math:`H`, Linear hardening constant
+       ``ScalarLinearHardeningParameter``, scalar, :math:`H`, Linear hardening constant (Pa per unit plastic strain)
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   When used with :ref:`ASDPlasticMaterial3D`, ScalarLinearHardeningFunction provides:
+   
+   * Linear increase/decrease of scalar variable with plastic deformation
+   * Isotropic hardening/softening behavior
+   * Simple calibration with single parameter
+   
+   Monitor scalar internal variable evolution:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file radius.out -time -node 1 -dof 1 VonMisesRadius
+       # Update hardening parameter during analysis
+       setParameter 1 ScalarLinearHardeningParameter 5e6
 
 
 .. _TensorValuedHF:
@@ -61,10 +79,28 @@ Where :math:`\vec{m}_{dev}` is the deviatoric part of the plastic flow direction
 .. admonition:: Parameters required
 
    .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
+      :header: "Parameter name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``TensorLinearHardeningFunction``, scalar, :math:`H`, Linear hardening constant
+       ``TensorLinearHardeningParameter``, scalar, :math:`H`, Linear hardening constant (dimensionless)
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   When used with :ref:`ASDPlasticMaterial3D`, TensorLinearHardeningFunction provides:
+   
+   * Linear evolution of backstress tensor
+   * Kinematic hardening behavior
+   * Translation of yield surface in stress space
+   
+   The backstress evolves in the direction of deviatoric plastic flow, providing Bauschinger effect modeling.
+   
+   Monitor tensor internal variable evolution:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file backstress.out -time -node 1 -dof 1 2 3 4 5 6 BackStress
+       # Update hardening parameter during analysis
+       setParameter 1 TensorLinearHardeningParameter 0.01
 
 **ArmstrongFrederickHardeningFunction**
 
@@ -84,8 +120,39 @@ Setting :math:`c_r = 0` results in linear hardening.
 .. admonition:: Parameters required
 
    .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
+      :header: "Parameter name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``AF_ha``, scalar, :math:`h_a`, Model constant for the linear part of the hardening model. Controls the rate of saturation.
-      ``AF_cr``, scalar, :math:`c_r`, Model constant for saturation part of the hardening model. Controls the saturation value. 
\ No newline at end of file
+       ``AF_ha``, scalar, :math:`h_a`, Model constant for the linear part of the hardening model. Controls the rate of saturation (Pa).
+       ``AF_cr``, scalar, :math:`c_r`, Model constant for saturation part of the hardening model. Controls saturation value (dimensionless).
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   When used with :ref:`ASDPlasticMaterial3D`, ArmstrongFrederickHardeningFunction provides:
+   
+   * Non-linear kinematic hardening with saturation
+   * Dynamic recovery term for realistic cyclic response
+   * Bauschinger effect modeling with memory
+   
+   The hardening combines:
+   
+   * **Linear term**: (2/3) :math:`h_a \vec{m}_{dev}` - drives hardening
+   * **Recovery term**: :math:`-c_r \sqrt{(2/3) \vec{m}_{dev} \cdot \vec{m}_{dev}} \cdot \vec{T}` - provides saturation
+   
+   Setting :math:`c_r = 0` reduces to linear kinematic hardening.
+   
+   The saturation value for the internal variable will be:
+   
+   .. math::
+      \Vert \vec{T} \Vert = \sqrt{\dfrac{2}{3}} \dfrac{h_a}{c_r}
+   
+   Monitor non-linear hardening:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file backstress.out -time -node 1 -dof 1 2 3 4 5 6 BackStress
+       recorder Node -file eqpstrain.out -time -node 1 -dof 1 eqpstrain
+       # Update hardening parameters during analysis
+       setParameter 1 AF_ha 1e8
+       setParameter 1 AF_cr 100
+
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/PlasticFlowType.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/PlasticFlowType.rst
index eb65488f..22b18b72 100644
--- a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/PlasticFlowType.rst
+++ b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/PlasticFlowType.rst
@@ -1,13 +1,14 @@
 .. _PlasticFlowType:
 
 Plastic Flow Directions
-^^^^^^^^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^^^^^^^^
 
-Specifies the direction of plastic flow :math:`\mathbf{m}` used to define the evolution of the plastic strain (see :ref:`ASDPlasticTheory`).
+Specifies the direction of plastic flow :math:`\mathbf{m}` used to define evolution of the plastic strain (see :ref:`ASDPlasticTheory`). The plastic flow direction determines how plastic strain components develop and influences both volumetric and deviatoric plastic deformation.
 
+Plastic Flow Directions may define internal variables they need for their specification, as well as paramters. When specifying the ``ASDPlasticMaterial`` or ``ASDPlasticMaterial3D`` instance, one must provide internal variables mentioned below, together with their hardening function, when defining internal variables. 
 
-Plastic Flow Directions may define internal variables they need for their specification, as well as paramters. When specifying the ``ASDPlasticMaterial`` instance, once must provide the internal variables mentioned below, together with their hardening function, when defining the internal variables. 
-
+.. note::
+   :ref:`ASDPlasticMaterial3D` provides enhanced monitoring of plastic flow effects through specialized response queries like ``VolStrain`` and ``J2Strain`` to track volumetric and deviatoric plastic deformation separately.
 
 Available functions:
 
@@ -27,13 +28,25 @@ Defines a plastic flow direction derived from the `Von Mises Yield Criterion <ht
       :header: "IV Name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``BackStress``, Rank-6 Tensor , :math:`\vec{\alpha}`,  "Backstress, definining the location in stress space for the axis of the Von-Mises cylinder. "
+       ``BackStress``, Rank-6 Tensor , :math:`\vec{\alpha}`,  "Backstress, definining location in stress space for axis of Von-Mises cylinder. "
 
 .. admonition:: Parameters required
 
-   .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
-      :widths: 10, 10, 10, 70
+   None required for VonMises plastic flow. Parameters are inherited from associated yield function and hardening functions.
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   When using VonMises_PF with :ref:`ASDPlasticMaterial3D`, the plastic flow direction is automatically computed as the stress derivative of the Von Mises yield function. This provides:
+   
+   * Associated plastic flow (flow rule equals yield function gradient)
+   * No volumetric plastic strain (isochoric plastic deformation)
+   * Purely deviatoric plastic response
+   
+   You can monitor plastic strain components:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file plastic_strain.out -time -node 1 -dof 1 2 3 4 5 6 pstrain
 
 
 .. _`DruckerPrager_PF`:
@@ -58,10 +71,21 @@ Defines a plastic flow direction derived from the `Drucker-Prager Yield Criterio
 
 .. admonition:: Parameters required
 
-   .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
-      :widths: 10, 10, 10, 70
+   None required for DruckerPrager plastic flow. Uses parameters from yield function.
+
+.. admonition:: Usage in ASDPlasticMaterial3D
 
+   DruckerPrager_PF provides non-associated flow when used with pressure-sensitive yield functions. Key characteristics:
+   
+   * Pressure-dependent plastic flow (includes volumetric component)
+   * Associated flow with Drucker-Prager yield surface
+   * Suitable for geomaterials exhibiting dilatancy or contractancy
+   
+   Monitor volumetric plastic strain:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file vol_strain.out -time -node 1 -dof 1 VolStrain
 
 
 .. _`ConstantDilatancy_PF`:
@@ -87,7 +111,24 @@ Von-Mises PF provides no volumetric change (dilatancy) during plasticity. On the
 .. admonition:: Parameters required
 
    .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
+      :header: "Parameter Name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-      ``Dilatancy``, Scalar, :math:`D`, "Defines the rate of dilatancy with plastic flow. Positive values specify negative plastic volumetric change."
+       ``Dilatancy``, Scalar, :math:`D`, "Defines the rate of dilatancy with plastic flow. Positive values specify negative plastic volumetric change."
+
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   ConstantDilatancy_PF provides non-associated flow with controllable dilatancy. This is particularly useful for:
+   
+   * Concrete and rock materials (positive dilatancy)
+   * Dense sands (positive dilatancy initially, then contractancy)
+   * Loose sands (continuous contractancy)
+   
+   Monitor dilatancy effects:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file vol_strain.out -time -node 1 -dof 1 VolStrain
+       recorder Node -file j2_strain.out -time -node 1 -dof 1 J2Strain
+
+
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/YieldFunctions.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/YieldFunctions.rst
index 0963090f..b16e063f 100644
--- a/source/user/manual/material/ndMaterials/ASDPlasticMaterial/YieldFunctions.rst
+++ b/source/user/manual/material/ndMaterials/ASDPlasticMaterial/YieldFunctions.rst
@@ -3,11 +3,14 @@
 Yield Functions
 ^^^^^^^^^^^^^^^
 
-Specifies the shape of the yield surface, i.e. the locus of points in stress-space where the response of the material is linear. Points outside the yield surface are not allowed. Yield surfaces can have parameters and internal variables which define its shape and location. 
+Specifies the shape of the yield surface, i.e. locus of points in stress-space where the response of the material is linear. Points outside of yield surface are not allowed. Yield surfaces can have parameters and internal variables which define its shape and location. 
 
-Yield Functions may define internal variables they need for their specification, as well as paramters. When specifying the ``ASDPlasticMaterial`` instance, once must provide the internal variables mentioned below, together with their hardening function. 
+Yield Functions may define internal variables they need for their specification, as well as paramters. When specifying ``ASDPlasticMaterial`` or ``ASDPlasticMaterial3D`` instance, one must provide internal variables mentioned below, together with their hardening function. 
 
-Yield functions are also responsible for providing their stress-derivative (:math:`\mathbf{n} = \partial f / \partial \boldsymbol{\sigma}`) as well as computing the hardening term :math:`H` in the plastic multiplier (see :ref:`ASDPlasticTheory`). 
+Yield functions are also responsible for providing their stress-derivative (:math:`\mathbf{n} = \partial f / \partial \boldsymbol{\sigma}`) as well as computing hardening term :math:`H` in the plastic multiplier (see :ref:`ASDPlasticTheory`). 
+
+.. note::
+   When using :ref:`ASDPlasticMaterial3D`, yield function internal variables can be accessed directly using enhanced response queries. This provides detailed insight into material state evolution during analysis.
 
 Available functions:
 
@@ -44,10 +47,21 @@ Where :math:`p = -(\sigma_{11} + \sigma_{22} + \sigma_{33})/3` is the mean stres
 
 .. admonition:: Parameters required
 
-   .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
-      :widths: 10, 10, 10, 70
+   None required for the yield function itself. Parameters are defined by the internal variable hardening functions and elasticity model.
 
+.. admonition:: Usage in ASDPlasticMaterial3D
+
+   When using VonMises_YF with :ref:`ASDPlasticMaterial3D`, you can access the internal variables using enhanced response queries:
+   
+   * ``VonMisesRadius`` - Current radius of the yield surface
+   * ``BackStress`` - Current backstress tensor (6 components in Voigt notation)
+   
+   These can be recorded using:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file vonmises_radius.out -time -node 1 -dof 1 VonMisesRadius
+       recorder Node -file backstress.out -time -node 1 -dof 1 BackStress
 
 
 .. _DruckerPrager_YF:
@@ -77,13 +91,23 @@ Where all symbols are as before. In this case note that :math:`k`, defines the s
 .. admonition:: Parameters required
 
    .. csv-table:: 
-      :header: "IV Name", "Type", "Symbol", "Description"
+      :header: "Parameter Name", "Type", "Symbol", "Description"
       :widths: 10, 10, 10, 70
 
-.. _RoundedMohrCoulomb_YF:
+      ``Dilatancy``, Scalar, :math:`D`, "Dilatancy parameter from plastic flow model (if used)"
+
+.. admonition:: Usage in ASDPlasticMaterial3D
 
-RoundedMohrCoulomb_YF
-"""""""""""""""""""""
+   DruckerPrager_YF with :ref:`ASDPlasticMaterial3D` provides access to:
+   
+   * ``VonMisesRadius`` - Current strength parameter
+   * ``BackStress`` - Current backstress tensor
+   
+   Enhanced stress measures are particularly useful:
+   
+   .. code-block:: tcl
+   
+       recorder Node -file p_stress.out -time -node 1 -dof 1 PStress
+       recorder Node -file j2_stress.out -time -node 1 -dof 1 J2Stress
 
-Coming soon
 
diff --git a/source/user/manual/material/ndMaterials/ASDPlasticMaterial3D.rst b/source/user/manual/material/ndMaterials/ASDPlasticMaterial3D.rst
new file mode 100644
index 00000000..23f70f17
--- /dev/null
+++ b/source/user/manual/material/ndMaterials/ASDPlasticMaterial3D.rst
@@ -0,0 +1,379 @@
+.. _ASDPlasticMaterial3D:
+
+ASDPlasticMaterial3D
+^^^^^^^^^^^^^^^^^^^
+
+| This command is used to construct an ``ASDPlasticMaterial3D`` material object. ``ASDPlasticMaterial3D`` implements a large family of constitutive models based on the classical theory of elastoplasticity using advanced template metaprogramming. Users build new constitutive models by selecting the yield function, plastic-flow direction, elasticity law, and hardening models for the internal variables from several possible options for each component. 
+
+| This is the 3D-specific implementation that provides enhanced features compared to the base ``ASDPlasticMaterial`` class, including multiple integration methods, tangent operator options, and extensive response capabilities.
+
+To create a new model, specify the Yield Function type (``$YieldFunctionType``), Plastic Flow type (``$PlasticFlowType``), and Elasticity type (``$ElasticityType``). All internal variables and model parameters are set to zero unless specified by the user, which might or not make sense depending on context, and this initialization is printed out to the screen. 
+
+After setting the ``$YieldFunctionType``, ``$PlasticFlowType`` and ``$ElasticityType``, you can give initial values to internal variables within the ``Begin_Internal_Variables`` ... ``End_Internal_Variables`` block. Then, the ``Begin_Model_Parameters`` ... ``End_Model_Parameters`` block is used to provide model parameter values (these can be changed during the analysis with the ``setParameter`` command as expected). Finally, model integration options are set within the ``Begin_Integration_Options`` ... ``End_Integration_Options`` code block. Specification blocks can occur in any order or ommitted. 
+
+The complete command looks as follows::
+
+   nDMaterial ASDPlasticMaterial3D $tag 
+      $YieldFunctionType 
+      $PlasticFlowType $
+      ElasticityType 
+      $IV_TYPE 
+      Begin_Internal_Variables 
+         $InternalVariable1 $$double_value1 $$double_value2... $$double_valueN1  
+         $InternalVariable2 $$double_value1 $$double_value2... $$double_valueN2 
+         #... (depends on how many internal variables the particular selected model has) 
+      End_Internal_Variables 
+      Begin_Model_Parameters 
+         $ModelParameters1 $$double_value1 
+         $ModelParameters2 $$double_value2 
+         #... (depends on how many model parameters the particular selected model has) 
+      End_Model_Parameters 
+      Begin_Integration_Options 
+         f_relative_tol $double_value 
+         stress_relative_tol $double_value 
+         n_max_iterations $int_value 
+         return_to_yield_surface (0 or 1) 
+         method (string) : Forward_Euler | Forward_Euler_Subincrement | Backward_Euler | Backward_Euler_LineSearch | Modified_Euler_Error_Control | Runge_Kutta_45_Error_Control
+         tangent_operator (string) : Elastic | Numerical_Algorithmic_FirstOrder | Numerical_Algorithmic_SecondOrder | Continuum | Secant
+      End_Integration_Options
+
+
+Explanation
+
+.. csv-table:: 
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $tag, |integer|, "Unique tag identifying this material."
+   $YieldFunctionType, |string|, "Mandatory. Yield function to be used -> :ref:`YieldFunctionType`"
+   $PlasticFlowType, |string|, "Mandatory. Plastic flow direction to be used -> :ref:`PlasticFlowType`"
+   $ElasticityType, |string|, "Mandatory. Elastic model to be used -> :ref:`ElasticityType`"
+   $IV_TYPE, |string|, "Mandatory. Hardening model for internal variables. Admitted types depend on YF, PF, and EL chosen.  "
+   Begin_Internal_Variables, |string|, "Optional. Marks the beginning of the code block to set the internal variables. If ommitted, all internal variables are initialized to zero. You can specify as many of the following variables as wanted. Ommitted variables are initialized to zero. Number and name of variables that can be set is model dependent (once YF, PF, and EL are specified)"
+   $InternalVariable1, |list of name/value pairs|, "Initial value of internal variable1. Dimension depends on internal variable type.  "
+   $InternalVariable2, |list of name/value pairs|, "Initial value of internal variable2.  "
+   "--", "--", "... (input as many as model supports)"
+   End_Internal_Variables, |string|, "Mandatory if block started. Marks the end of the code block to set the internal variables"
+   Begin_Model_Parameters, |string|, "Optional. Marks the beginning of the code block to set the model parameters"
+   $ModelParameters, |list of name/value pairs|, "Values for parameters of the models to be used. This depends on the particular choices of ``$YieldFunctionType``, ``$PlasticFlowType``, ``$ElasticityType``, and ``$IV_type``. "
+   "--", "--", "... (input as many as model supports)"
+   End_Model_Parameters, |string|, "Mandatory if block started. Marks the beginning of the code block to set the model parameters"
+   Begin_Integration_Options, |string|, "Optional. Marks the beginning of the code block to set the integration options. You can set any ammount "
+   End_Integration_Options, |string|, "Mandatory if block started. Marks the beginning of the code block to set the model parameters"
+
+|  The :ref:`YieldFunctionType`, :ref:`PlasticFlowType`, and :ref:`ElasticityType` have different variables and parameters to be set.
+
+.. admonition:: Specification of internal variables (``$IV_TYPE``)
+
+	All internal variables are specified in a single string with no spaces. Internal variables specifications are separated with a colon ``:``. After the name of the internal variable, the name of the hardening function must be provided in parenthesis. The specification string must end in a colon. The syntax is, thus:
+
+	.. code-block::
+
+		InternalVariable1(HardeningFunction1):InternalVariable2(HardeningFunction2):
+
+	Internal variables are required for the yield function and plastic flow direction (see speficic components for details), and all internal variables must be provided with their hardening function otherwise instantiation fails. 
+
+.. admonition:: Note
+
+    Any internal variable or parameter that is not specified is set to zero by default.  
+
+
+.. toctree::
+   :caption: Arguments detailed description
+   :maxdepth: 2
+
+   ../ASDPlasticMaterial/YieldFunctions
+   ../ASDPlasticMaterial/PlasticFlowType
+   ../ASDPlasticMaterial/ElasticityType
+   ../ASDPlasticMaterial/HardeningFunctions
+
+.. _ASDPlasticMaterial3DTheory:
+
+ASDPlasticMaterial3D Theory
+"""""""""""""""""""""""""""""
+
+The theoretical foundation is identical to :ref:`ASDPlasticTheory`. However, ``ASDPlasticMaterial3D`` provides additional implementation features:
+
+* **Template Metaprogramming**: Uses C++ template metaprogramming for compile-time optimization
+* **Multiple Integration Methods**: Offers several integration algorithms beyond the basic methods
+* **Advanced Tangent Operators**: Provides different tangent computation strategies
+* **Enhanced Response System**: Supports additional response queries and parameter updates
+
+The material recieves the total (trial) strain :math:`\vec{\epsilon}^{\text{trial}}` from the finite element that contains it. From this, the trial strain increment is computed by subtracting the previously committed total strain:
+
+.. math::
+    \Delta \vec{\epsilon}^{\text{trial}} = \vec{\epsilon}^{\text{trial}} - \vec{\epsilon}^{\text{commit}}
+
+The integration method selected determines how the stress increment is computed and how the plastic corrections are applied.
+
+.. _`ASDPlasticMaterial3DIntegrationMethods`:
+Advanced Integration Methods
+""""""""""""""""""""""""""""
+
+``ASDPlasticMaterial3D`` supports multiple integration algorithms beyond the basic Forward Euler and Runge-Kutta methods:
+
+.. csv-table:: 
+   :header: "Integration Method", "Description"
+   :widths: 30, 50
+
+   ``Forward_Euler``, "Basic explicit forward Euler method with yield surface crossing detection"
+   ``Forward_Euler_Subincrement``, "Forward Euler with automatic subincrementation for large strain increments"
+   ``Backward_Euler``, "Implicit backward Euler method for better stability"
+   ``Backward_Euler_LineSearch``, "Backward Euler with line search for enhanced convergence"
+   ``Modified_Euler_Error_Control``, "Modified Euler method with adaptive error control"
+   ``Runge_Kutta_45_Error_Control``, "4th/5th order Runge-Kutta with adaptive error control (default)"
+
+.. admonition:: Integration Method Selection
+
+   The integration method is selected in the ``Begin_Integration_Options`` block using the ``method`` parameter. If not specified, ``Runge_Kutta_45_Error_Control`` is used by default.
+
+.. _`ASDPlasticMaterial3DTangentOperators`:
+Tangent Operator Types
+"""""""""""""""""""""
+
+``ASDPlasticMaterial3D`` provides several options for computing the consistent tangent operator:
+
+.. csv-table:: 
+   :header: "Tangent Operator", "Description"
+   :widths: 30, 50
+
+   ``Elastic``, "Returns the elastic tangent matrix (ignores plasticity)"
+   ``Numerical_Algorithmic_FirstOrder``, "Computes algorithmic tangent using first-order finite differences"
+   ``Numerical_Algorithmic_SecondOrder``, "Computes algorithmic tangent using second-order finite differences (more accurate but slower)"
+   ``Continuum``, "Returns the continuum tangent operator"
+   ``Secant``, "Returns the average of continuum and elastic tangents"
+
+.. admonition:: Tangent Operator Selection
+
+   The tangent operator type is selected in the ``Begin_Integration_Options`` block using the ``tangent_operator`` parameter. If not specified, an appropriate default is selected based on the integration method.
+
+.. _`ASDPlasticMaterial3DIntegrationOptions`:
+Integration Options
+"""""""""""""""""""
+
+.. csv-table:: 
+   :header: "Parameter", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $f_relative_tol, |double|, "Relative tolerance to evaluate the yield function crossing."
+   $stress_relative_tol, |double|, "Tolerance for the convergece of the integration algorithm."
+   $n_max_iterations, |int|, "Maximum number of iterations for constitutive integration."
+   $return_to_yield_surface, |0 or 1|, "Whether to apply a return to yield surface algorithm after integration convergence."
+   $method, |string|, "Constitutive integration method. See :ref:`ASDPlasticMaterial3DIntegrationMethods` for options"
+   $tangent_operator, |string|, "Tangent operator computation method. See :ref:`ASDPlasticMaterial3DTangentOperators` for options"
+
+The default integration method is **Runge_Kutta_45_Error_Control** that uses the classical RK45 ODE integration algorithm employing a 4-th order prediction of the stress increment together with a 5-order prediction to estimate the integration error. In this scheme the strain increment provided by the element to the Gauss point is sub-divided into sub-increments, a process which is automated such that the provided **$stress_relative_tol** is met. This method is provided as a robust standard method which is applicable across all possible combinations of components, although there are possibly better approaches for specific cases which might become available in the future. 
+
+The different parameters are activated depending on the integration algorithm selected. The *Forward_Euler* family algorithms only uses the **$return_to_yield_surface** parameter, while **Runge_Kutta_45_Error_Control** uses the rest. 
+
+The **$f_relative_tol** parameter comes into play when the yield surface is being crossed, that is, when the previous (committed) stress is within the yield surface and the elastic prediction of the stress increment brings the stress state beyond the yield surface. In that case, an elastic increment occurs until the yield surface is touched which requires iterations with the Brent root finding algorithm. This is used by all currently available integration methods. 
+
+
+.. _`ASDPlasticMaterial3DOtherFeatures`:
+
+Other Features
+""""""""""""""""""
+
+.. admonition:: General parameters
+
+    These parameters are defined for all models. 
+
+    .. csv-table:: 
+       :header: "Parameter", "Type", "Description"
+       :widths: 10, 10, 40
+
+       ``MassDensity``, scalar, Defines the material mass density :math:`\rho`. 
+       ``InitialP0``, scalar, Defines the initial mean pressure at which material constants will be evaluated at the first step.  
+
+
+.. admonition:: Enhanced Responses  (`setResponse` and `getResponse` behavior)
+
+    * *Basic responses*. ``stresses`` for stress, ``strains`` for strains, and ``pstrains`` for plastic strains. 
+    * *Enhanced stress measures*. ``eqpstrain`` for equivalent plastic strain, ``PStress`` for mean stress (:math:`p = (σ_{11}+σ_{22}+σ_{33})/3`), ``J2Stress`` for J2 invariant of stress.
+    * *Enhanced strain measures*. ``VolStrain`` for volumetric strain (:math:`ε_v = ε_{11}+ε_{22}+ε_{33}`), ``J2Strain`` for J2 invariant of strain.
+    * *Internal variables*. You can request any and all internal variables by their specific name as an output. 
+
+.. admonition:: Advanced `setParameter` behavior 
+
+   * *State variable updates*. Can update stress, strain, plasticStrain and their trial variants directly.
+   * *Parameter updates*. Can update any model parameter by name using the generic parameter system.
+   * *Special K0 operations*. ``K02D`` and ``K03D`` for setting K0 stress states in 2D/3D conditions.
+
+.. admonition:: Internal Variable Management
+
+   ``ASDPlasticMaterial3D`` provides methods to inspect and manage internal variables:
+   
+   * ``getInternalVariablesNames()`` - Returns list of internal variable names
+   * ``getInternalVariableSizeByName(name)`` - Returns size of specified internal variable
+   * ``getInternalVariableIndexByName(name)`` - Returns index of specified internal variable
+   * ``setInternalVariableByName(name, size, values)`` - Sets internal variable values
+   * ``getParameterNames()`` - Returns list of model parameter names
+   * ``setParameterByName(name, value)`` - Sets model parameter value
+
+Implementation details
+""""""""""""""""""""""
+
+``ASDPlasticMaterial3D`` is implemented using C++ template metaprogramming, with a header-only design and using the "eigen" C++ library for high-performance array operations. This design provides modularity and the ability to mix and match components to create new models, while also providing high-performance because runtime polymorphism is avoided. 
+
+Key implementation features:
+* **Template-based architecture**: Compile-time generation of optimized code for each material configuration
+* **Eigen library**: High-performance linear algebra operations
+* **Voigt notation**: Efficient tensor representation using 6-component vectors
+* **Component modularity**: Easy extension with new yield functions, flow rules, elasticity models, and hardening functions
+
+
+Example
+"""""""
+
+The following example defines an instance of ``ASDPlasticMaterial3D`` with a Drucker-Prager yield function, a constant dilatancy plastic-flow direction, elastic-isotropic elasticity law and linear hardening for both internal variables, using advanced integration options.
+
+.. admonition:: TCL code  
+
+    .. code-block:: tcl
+
+        nDMaterial ASDPlasticMaterial3D 1 \
+            DruckerPrager_YF \
+            ConstantDilatancy_PF \
+            LinearIsotropic3D_EL \
+            BackStress(TensorLinearHardeningFunction):VonMisesRadius(ScalarLinearHardeningFunction): \
+            Begin_Internal_Variables \
+                VonMisesRadius 1. \
+                BackStress 0. 0. 0. 0. 0. 0. \
+            End_Internal_Variables \
+            Begin_Model_Parameters \
+                YoungsModulus 1e8 \
+                PoissonsRatio 0.3 \
+                TensorLinearHardeningParameter 1e6 \
+                ScalarLinearHardeningParameter 1e4 \
+                Dilatancy 0.02 \
+                MassDensity 2000. \
+            End_Model_Parameters \
+            Begin_Integration_Options \
+                method Runge_Kutta_45_Error_Control \
+                tangent_operator Continuum \
+                f_relative_tol 1e-8 \
+                stress_relative_tol 1e-6 \
+                n_max_iterations 50 \
+                return_to_yield_surface 1 \
+            End_Integration_Options
+
+.. admonition:: Python code  
+
+    .. code-block:: python
+
+        ops.nDMaterial("ASDPlasticMaterial3D", 1, 
+        "DruckerPrager_YF",
+        "ConstantDilatancy_PF",
+        "LinearIsotropic3D_EL",
+        "BackStress(TensorLinearHardeningFunction):VonMisesRadius(ScalarLinearHardeningFunction):",
+        "Begin_Internal_Variables",
+            "VonMisesRadius", 1.,
+            "BackStress", 0., 0., 0., 0., 0., 0.,
+        "End_Internal_Variables",
+        "Begin_Model_Parameters",
+            "YoungsModulus", 1e8,
+            "PoissonsRatio", 0.3,
+            "TensorLinearHardeningParameter", 1e6,
+            "ScalarLinearHardeningParameter", 1e4,
+            "Dilatancy", 0.02,
+            "MassDensity", 2000.,
+        "End_Model_Parameters",
+        "Begin_Integration_Options",
+            "method", "Runge_Kutta_45_Error_Control",
+            "tangent_operator", "Continuum",
+            "f_relative_tol", 1e-8,
+            "stress_relative_tol", 1e-6,
+            "n_max_iterations", 50,
+            "return_to_yield_surface", 1,
+        "End_Integration_Options",
+        )
+
+Example: Using Enhanced Response Queries
+""""""""""""""""""""""""""""""""""""""""
+
+To access the enhanced response capabilities of ``ASDPlasticMaterial3D``, you can use recorders with the following response queries:
+
+.. admonition:: TCL Recorder Example  
+
+    .. code-block:: tcl
+
+        # Basic responses
+        recorder Node -file stress.out -time -node 1 -dof 1 2 3 stress
+        recorder Node -file strain.out -time -node 1 -dof 1 2 3 strain
+        
+        # Enhanced stress measures  
+        recorder Node -file p_stress.out -time -node 1 -dof 1 PStress
+        recorder Node -file j2_stress.out -time -node 1 -dof 1 J2Stress
+        recorder Node -file eq_pstrain.out -time -node 1 -dof 1 eqpstrain
+        
+        # Enhanced strain measures
+        recorder Node -file vol_strain.out -time -node 1 -dof 1 VolStrain  
+        recorder Node -file j2_strain.out -time -node 1 -dof 1 J2Strain
+        
+        # Internal variable by name
+        recorder Node -file backstress.out -time -node 1 -dof 1 BackStress
+        recorder Node -file radius.out -time -node 1 -dof 1 VonMisesRadius
+
+Example: Parameter Updates During Analysis
+""""""""""""""""""""""""""""""""""""""""
+
+``ASDPlasticMaterial3D`` supports parameter updates during analysis:
+
+.. admonition:: TCL Parameter Update  
+
+    .. code-block:: tcl
+
+        # Update model parameters
+        setParameter 1 YoungsModulus 2e8
+        setParameter 1 PoissonsRatio 0.25
+        setParameter 1 ScalarLinearHardeningParameter 2e4
+        
+        # Update internal variables
+        setParameter 1 VonMisesRadius 1.5
+        
+        # Update stress state (for restart or initialization)
+        setParameter 1 stress [list 1e6 1e6 1e6 0 0 0]
+
+Example: Integration Method Selection
+""""""""""""""""""""""""""""""""""""""""
+
+Different integration methods for specific applications:
+
+.. admonition:: Robust Integration (Recommended for most cases)  
+
+    .. code-block:: tcl
+
+        Begin_Integration_Options
+            method Runge_Kutta_45_Error_Control
+            tangent_operator Continuum
+            f_relative_tol 1e-8
+            stress_relative_tol 1e-6
+            n_max_iterations 50
+            return_to_yield_surface 1
+        End_Integration_Options
+
+.. admonition:: Fast Integration (For large models)  
+
+    .. code-block:: tcl
+
+        Begin_Integration_Options
+            method Forward_Euler_Subincrement  
+            tangent_operator Elastic
+            return_to_yield_surface 1
+        End_Integration_Options
+
+.. admonition:: High Accuracy Integration (For research)  
+
+    .. code-block:: tcl
+
+        Begin_Integration_Options
+            method Backward_Euler_LineSearch
+            tangent_operator Numerical_Algorithmic_SecondOrder
+            f_relative_tol 1e-12
+            stress_relative_tol 1e-10
+            n_max_iterations 100
+            return_to_yield_surface 1
+        End_Integration_Options
+
+
+Code Developed by: **José A. Abell** (UANDES, Chile and ASDEA), **Guido Camata** and **Massimo Petracca**  (ASDEA Software, Italy).
\ No newline at end of file
diff --git a/source/user/manual/material/ndMaterials/NDTest.rst b/source/user/manual/material/ndMaterials/NDTest.rst
new file mode 100644
index 00000000..9dac07d4
--- /dev/null
+++ b/source/user/manual/material/ndMaterials/NDTest.rst
@@ -0,0 +1,627 @@
+NDTest Command
+**************
+
+The NDTest command provides a comprehensive suite of testing functions for nDMaterial objects, allowing users to directly manipulate and query material states for testing, validation, and educational purposes. This functionality is particularly useful for material model development, verification, and understanding material behavior under controlled loading conditions.
+
+.. warning::
+   These commands are intended for material testing and debugging purposes. They should not be used in regular structural analysis models where materials should be tested through elements.
+
+Command Syntax
+==============
+
+The general syntax for NDTest commands is:
+
+.. code-block:: tcl
+
+   NDTest $command $arguments
+
+Where ``$command`` specifies the specific operation to perform on the material, and ``$arguments`` are the command-specific parameters.
+
+Available Commands
+==================
+
+The following commands are available within the NDTest framework:
+
+SetStrain
+----------
+
+Sets trial strain components on a specified nDMaterial.
+
+.. code-block:: tcl
+
+   NDTest SetStrain $matTag $eps11 $eps22 $eps33 $gamma12 $gamma23 $gamma31
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+   $eps11, |double|, Normal strain in direction 1 (ε₁₁)
+   $eps22, |double|, Normal strain in direction 2 (ε₂₂)
+   $eps33, |double|, Normal strain in direction 3 (ε₃₃)
+   $gamma12, |double|, Engineering shear strain γ₁₂
+   $gamma23, |double|, Engineering shear strain γ₂₃
+   $gamma31, |double|, Engineering shear strain γ₃₁
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   # Set a uniaxial strain state
+   NDTest SetStrain 1 0.001 0.0 0.0 0.0 0.0 0.0
+
+Python:
+
+.. code-block:: python
+
+   # Set a uniaxial strain state
+   op.NDTest("SetStrain", 1, 0.001, 0.0, 0.0, 0.0, 0.0, 0.0)
+
+CommitState
+-----------
+
+Commits the current trial state of the material, making it the converged state.
+
+.. code-block:: tcl
+
+   NDTest CommitState $matTag
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   NDTest CommitState 1
+
+Python:
+
+.. code-block:: python
+
+   op.NDTest("CommitState", 1)
+
+PrintStrain
+-----------
+
+Prints the current strain state of the material to the OpenSees output stream.
+
+.. code-block:: tcl
+
+   NDTest PrintStrain $matTag
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+PrintStress
+-----------
+
+Prints the current stress state of the material to the OpenSees output stream.
+
+.. code-block:: tcl
+
+   NDTest PrintStress $matTag
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+GetStrain
+----------
+
+Returns the current strain components as a list that can be assigned to a Tcl variable.
+
+.. code-block:: tcl
+
+   set strainList [NDTest GetStrain $matTag]
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+**Returns:** A list of 6 double values [ε₁₁, ε₂₂, ε₃₃, γ₁₂, γ₂₃, γ₃₁]
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   set strains [NDTest GetStrain 1]
+   puts "Strain components: $strains"
+
+Python:
+
+.. code-block:: python
+
+   strains = op.NDTest("GetStrain", 1)
+   print("Strain components:", strains)
+
+GetStress
+----------
+
+Returns the current stress components as a list that can be assigned to a Tcl variable.
+
+.. code-block:: tcl
+
+   set stressList [NDTest GetStress $matTag]
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+**Returns:** A list of 6 double values [σ₁₁, σ₂₂, σ₃₃, τ₁₂, τ₂₃, τ₃₁]
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   set stresses [NDTest GetStress 1]
+   puts "Stress components: $stresses"
+
+Python:
+
+.. code-block:: python
+
+   stresses = op.NDTest("GetStress", 1)
+   print("Stress components:", stresses)
+
+GetTangentStiffness
+-------------------
+
+Returns the current tangent stiffness matrix as a flat list of 36 values.
+
+.. code-block:: tcl
+
+   set stiffnessList [NDTest GetTangentStiffness $matTag]
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+
+**Returns:** A list of 36 double values representing the 6×6 tangent stiffness matrix in row-major order:
+[C₁₁, C₁₂, C₁₃, C₁₄, C₁₅, C₁₆, C₂₁, C₂₂, ..., C₆₆]
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   set tangent [NDTest GetTangentStiffness 1]
+   # Access C11 component
+   set C11 [lindex $tangent 0]
+   # Access C26 component  
+   set C26 [lindex $tangent 17]
+
+Python:
+
+.. code-block:: python
+
+   tangent = op.NDTest("GetTangentStiffness", 1)
+   # Access C11 component
+   C11 = tangent[0]
+   # Access C26 component  
+   C26 = tangent[17]
+
+GetResponse
+-----------
+
+Queries the material for custom responses specific to the material type.
+
+.. code-block:: tcl
+
+   set responseList [NDTest GetResponse $matTag $responseType]
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+   $responseType, |string|, Type of response to query (material-dependent)
+
+**Returns:** A list of double values containing the requested response data.
+
+**Common Response Types:** (varies by material type)
+- "stress" - stress components
+- "strain" - strain components  
+- "tangent" - tangent stiffness
+- "plasticStrain" - plastic strain components
+- "backstress" - backstress tensor components
+- Custom material-specific responses
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   # Get plastic strain from a plasticity model
+   set plasticStrain [NDTest GetResponse 1 "plasticStrain"]
+   
+   # Get damage variable from a damage model
+   set damage [NDTest GetResponse 1 "damage"]
+
+Python:
+
+.. code-block:: python
+
+   # Get plastic strain from a plasticity model
+   plastic_strain = op.NDTest("GetResponse", 1, "plasticStrain")
+   
+   # Get damage variable from a damage model
+   damage = op.NDTest("GetResponse", 1, "damage")
+
+UpdateIntegerParameter
+-----------------------
+
+Updates an integer parameter in the material during runtime.
+
+.. code-block:: tcl
+
+   NDTest UpdateIntegerParameter $matTag $paramID $newValue
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+   $paramID, |integer|, Integer parameter identifier (material-specific)
+   $newValue, |integer|, New integer value for the parameter
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   NDTest UpdateIntegerParameter 1 1 0
+
+Python:
+
+.. code-block:: python
+
+   op.NDTest("UpdateIntegerParameter", 1, 1, 0)
+
+UpdateDoubleParameter
+----------------------
+
+Updates a double-precision parameter in the material during runtime.
+
+.. code-block:: tcl
+
+   NDTest UpdateDoubleParameter $matTag $paramID $newValue
+
+.. csv-table::
+   :header: "Argument", "Type", "Description"
+   :widths: 10, 10, 40
+
+   $matTag, |integer|, Unique tag identifying the nDMaterial object
+   $paramID, |integer|, Integer parameter identifier (material-specific)
+   $newValue, |double|, New double value for the parameter
+
+**Example:**
+
+TCL:
+
+.. code-block:: tcl
+
+   NDTest UpdateDoubleParameter 1 1 2.1e-5
+
+Python:
+
+.. code-block:: python
+
+   op.NDTest("UpdateDoubleParameter", 1, 1, 2.1e-5)
+
+Usage Examples
+===============
+
+Complete Material Test Procedure
+--------------------------------
+
+Here's a complete example of how to use NDTest commands to test a material:
+
+TCL:
+
+.. code-block:: tcl
+
+   # Create a simple elastic material
+   nDMaterial ElasticIsotropic 1 200000 0.3
+
+   # Test procedure 1: Uniaxial tension
+   puts "=== Uniaxial Tension Test ==="
+   
+   # Apply strain increment
+   NDTest SetStrain 1 0.001 0.0 0.0 0.0 0.0 0.0
+   
+   # Get and display results
+   set strains [NDTest GetStrain 1]
+   set stresses [NDTest GetStress 1]
+   puts "Strain: $strains"
+   puts "Stress: $stresses"
+   
+   # Commit the state
+   NDTest CommitState 1
+
+   # Test procedure 2: Shear test  
+   puts "=== Pure Shear Test ==="
+   NDTest SetStrain 1 0.0 0.0 0.0 0.002 0.0 0.0
+   set stresses [NDTest GetStress 1]
+   puts "Shear stress: $stresses"
+   NDTest CommitState 1
+
+   # Test procedure 3: Get tangent stiffness
+   puts "=== Tangent Stiffness ==="
+   set tangent [NDTest GetTangentStiffness 1]
+   puts "C11 = [lindex $tangent 0]"
+   puts "C12 = [lindex $tangent 1]"
+
+Python:
+
+.. code-block:: python
+
+   import openseespy.opensees as op
+   
+   # Create a simple elastic material
+   op.nDMaterial("ElasticIsotropic", 1, 200000, 0.3)
+
+   # Test procedure 1: Uniaxial tension
+   print("=== Uniaxial Tension Test ===")
+   
+   # Apply strain increment
+   op.NDTest("SetStrain", 1, 0.001, 0.0, 0.0, 0.0, 0.0, 0.0)
+   
+   # Get and display results
+   strains = op.NDTest("GetStrain", 1)
+   stresses = op.NDTest("GetStress", 1)
+   print("Strain:", strains)
+   print("Stress:", stresses)
+   
+   # Commit the state
+   op.NDTest("CommitState", 1)
+
+   # Test procedure 2: Shear test  
+   print("=== Pure Shear Test ===")
+   op.NDTest("SetStrain", 1, 0.0, 0.0, 0.0, 0.002, 0.0, 0.0)
+   stresses = op.NDTest("GetStress", 1)
+   print("Shear stress:", stresses)
+   op.NDTest("CommitState", 1)
+
+   # Test procedure 3: Get tangent stiffness
+   print("=== Tangent Stiffness ===")
+   tangent = op.NDTest("GetTangentStiffness", 1)
+   print("C11 =", tangent[0])
+   print("C12 =", tangent[1])
+
+Cyclic Loading Test
+-------------------
+
+TCL:
+
+.. code-block:: tcl
+
+   # Create a material with cyclic behavior
+   nDMaterial J2Plasticity 1 27 0.3 0.001 300 20.0 0.0
+
+   # Cyclic loading parameters
+   set cycles 5
+   set maxStrain 0.01
+   set strainIncrement 0.001
+
+   # Apply cyclic loading
+   for {set cycle 0} {$cycle < $cycles} {incr cycle} {
+       puts "Cycle $cycle"
+       
+       # Loading
+       for {set i 0} {$i <= $maxStrain/$strainIncrement} {incr i} {
+           set strain [expr {$i * $strainIncrement}]
+           NDTest SetStrain 1 $strain 0.0 0.0 0.0 0.0 0.0
+           set stress [lindex [NDTest GetStress 1] 0]
+           puts "$strain $stress"
+           NDTest CommitState 1
+       }
+       
+       # Unloading
+       for {set i [expr {$maxStrain/$strainIncrement}]} {$i >= 0} {incr i -1} {
+           set strain [expr {$i * $strainIncrement}]
+           NDTest SetStrain 1 $strain 0.0 0.0 0.0 0.0 0.0
+           set stress [lindex [NDTest GetStress 1] 0]
+           puts "$strain $stress"
+           NDTest CommitState 1
+       }
+   }
+
+Python:
+
+.. code-block:: python
+
+   import openseespy.opensees as op
+   
+   # Create a material with cyclic behavior
+   op.nDMaterial("J2Plasticity", 1, 27, 0.3, 0.001, 300, 20.0, 0.0)
+
+   # Cyclic loading parameters
+   cycles = 5
+   max_strain = 0.01
+   strain_increment = 0.001
+
+   # Apply cyclic loading
+   for cycle in range(cycles):
+       print(f"Cycle {cycle}")
+       
+       # Loading
+       num_steps = int(max_strain / strain_increment)
+       for i in range(num_steps + 1):
+           strain = i * strain_increment
+           op.NDTest("SetStrain", 1, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           stress = op.NDTest("GetStress", 1)[0]
+           print(f"{strain} {stress}")
+           op.NDTest("CommitState", 1)
+       
+       # Unloading
+       for i in range(num_steps, -1, -1):
+           strain = i * strain_increment
+           op.NDTest("SetStrain", 1, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           stress = op.NDTest("GetStress", 1)[0]
+           print(f"{strain} {stress}")
+           op.NDTest("CommitState", 1)
+
+Parameter Study
+---------------
+
+TCL:
+
+.. code-block:: tcl
+
+   # Create material
+   nDMaterial ElasticIsotropic 1 200000 0.3
+
+   # Study effect of Poisson's ratio
+   set poissons {0.15 0.20 0.25 0.30 0.35 0.40}
+   
+   foreach nu $poissons {
+       # Update Poisson's ratio
+       NDTest UpdateDoubleParameter 1 2 $nu
+       
+       # Apply uniaxial strain
+       NDTest SetStrain 1 0.001 0.0 0.0 0.0 0.0 0.0
+       
+       # Get lateral strain (should be -nu*axial_strain for elastic material)
+       set strains [NDTest GetStrain 1]
+       set lateralStrain [lindex $strains 1]
+       
+       puts "Poisson's ratio: $nu, Lateral strain: $lateralStrain"
+       
+       NDTest CommitState 1
+   }
+
+Python:
+
+.. code-block:: python
+
+   import openseespy.opensees as op
+   
+   # Create material
+   op.nDMaterial("ElasticIsotropic", 1, 200000, 0.3)
+
+   # Study effect of Poisson's ratio
+   poissons = [0.15, 0.20, 0.25, 0.30, 0.35, 0.40]
+   
+   for nu in poissons:
+       # Update Poisson's ratio
+       op.NDTest("UpdateDoubleParameter", 1, 2, nu)
+       
+       # Apply uniaxial strain
+       op.NDTest("SetStrain", 1, 0.001, 0.0, 0.0, 0.0, 0.0, 0.0)
+       
+       # Get lateral strain (should be -nu*axial_strain for elastic material)
+       strains = op.NDTest("GetStrain", 1)
+       lateral_strain = strains[1]
+       
+       print(f"Poisson's ratio: {nu}, Lateral strain: {lateral_strain}")
+       
+       op.NDTest("CommitState", 1)
+
+Implementation Details
+======================
+
+The NDTest command framework is implemented in ``OpenSeesNDTestCommands.cpp`` and consists of:
+
+Core Components
+---------------
+
+1. **Command Registration System**: Uses a function map to associate command strings with their corresponding handler functions
+2. **Material Interface**: Interfaces directly with the NDMaterial class hierarchy
+3. **Error Handling**: Comprehensive error checking and reporting through the OpenSees opserr stream
+4. **Data Conversion**: Handles conversion between OpenSees internal data structures and Tcl variables
+
+Technical Specifications
+------------------------
+
+**Strain/Stress Ordering**: All 6-component vectors follow the engineering convention:
+- Normal components: (ε₁₁, ε₂₂, ε₃₃) or (σ₁₁, σ₂₂, σ₃₃)
+- Shear components: (γ₁₂, γ₂₃, γ₃₁) or (τ₁₂, τ₂₃, τ₃₁)
+
+**Tangent Stiffness Matrix**: Returned as a flat array in row-major order for the 6×6 matrix C where:
+dσ = C : dε (using tensor notation)
+
+**Parameter Updates**: Parameter IDs are material-specific and typically defined in the material's documentation or source code.
+
+Error Handling
+==============
+
+The NDTest commands include comprehensive error checking:
+
+- **Material Existence**: Verifies the material tag exists before operations
+- **Input Validation**: Checks for correct number and type of arguments
+- **Parameter Validation**: Validates parameter ranges where applicable
+- **Memory Management**: Ensures proper handling of temporary data structures
+
+Common error messages and their meanings:
+
+.. code-block:: text
+
+   WARNING too few arguments: NDTest cmd?
+   # Command requires additional arguments
+
+   WARNING NDTest type UnknownCommand is unknown  
+   # Command not recognized
+
+   material with tag X does not exist
+   # Material tag not found in model
+
+   need 6 components of floating-point strains
+   # SetStrain requires exactly 6 strain components
+
+Best Practices
+==============
+
+When using NDTest commands, consider these best practices:
+
+1. **Material Testing**: Use for isolated material testing, not in structural analyses
+2. **State Management**: Always commit states when appropriate to maintain material history
+3. **Error Checking**: Wrap commands in try-catch blocks when possible
+4. **Documentation**: Document parameter IDs and response types for custom materials
+5. **Performance**: Consider performance implications when using in loops
+
+Integration with Analysis
+=========================
+
+While NDTest commands are primarily for material testing, they can be integrated with analysis workflows:
+
+- **Pre-analysis Calibration**: Calibrate material parameters before main analysis
+- **Post-analysis Verification**: Verify material states after analysis completion
+- **Batch Testing**: Automate material testing procedures
+- **Educational Use**: Demonstrate material behavior concepts
+
+For detailed examples and advanced usage patterns, see :ref:`NDTest_examples`.
+
+See Also
+========
+
+- :ref:`nDMaterial` - Main nDMaterial command documentation
+- :ref:`matTestCommands` - Basic material testing commands
+- :ref:`NDTest_examples` - Comprehensive usage examples
+- Material-specific documentation for available response types and parameter IDs
+- OpenSees examples repository for complete testing procedures
\ No newline at end of file
diff --git a/source/user/manual/material/ndMaterials/NDTest_examples.rst b/source/user/manual/material/ndMaterials/NDTest_examples.rst
new file mode 100644
index 00000000..a338b07c
--- /dev/null
+++ b/source/user/manual/material/ndMaterials/NDTest_examples.rst
@@ -0,0 +1,636 @@
+NDTest Examples
+===============
+
+This document provides comprehensive examples of using the NDTest command for testing nDMaterial objects in OpenSees.
+
+Basic Material Testing Example
+===============================
+
+TCL:
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 1: Basic Elastic Material Test
+   # =============================================
+   
+   # Create an elastic isotropic material
+   # E = 200 GPa, ν = 0.3
+   nDMaterial ElasticIsotropic 1 200000 0.3
+   
+   puts "=== Elastic Material Test ==="
+   puts "Material Properties:"
+   puts "  E = 200 GPa"
+   puts "  ν = 0.3"
+   puts ""
+   
+   # Test 1: Uniaxial tension
+   puts "Test 1: Uniaxial Tension"
+   set strain 0.001
+   NDTest SetStrain 1 $strain 0.0 0.0 0.0 0.0 0.0
+   set stresses [NDTest GetStress 1]
+   set stress11 [lindex $stresses 0]
+   puts "  Applied strain: $strain"
+   puts "  Resulting stress σ₁₁: $stress11 MPa"
+   puts "  Expected stress: [expr {200000 * $strain}] MPa"
+   puts ""
+   
+   # Test 2: Get tangent stiffness
+   puts "Test 2: Tangent Stiffness"
+   set tangent [NDTest GetTangentStiffness 1]
+   set C11 [lindex $tangent 0]
+   set C12 [lindex $tangent 1]
+   puts "  C₁₁: $C11 MPa"
+   puts "  C₁₂: $C12 MPa"
+   puts "  Expected C₁₁: [expr {200000}] MPa"
+   puts "  Expected C₁₂: [expr {200000 * 0.3 / (1 - 0.3)}] MPa"
+   puts ""
+   
+   # Commit the state
+   NDTest CommitState 1
+
+Python:
+
+.. code-block:: python
+
+   # =============================================
+   # Example 1: Basic Elastic Material Test
+   # =============================================
+   
+   import openseespy.opensees as op
+   
+   # Create an elastic isotropic material
+   # E = 200 GPa, ν = 0.3
+   op.nDMaterial("ElasticIsotropic", 1, 200000, 0.3)
+   
+   print("=== Elastic Material Test ===")
+   print("Material Properties:")
+   print("  E = 200 GPa")
+   print("  ν = 0.3")
+   print("")
+   
+   # Test 1: Uniaxial tension
+   print("Test 1: Uniaxial Tension")
+   strain = 0.001
+   op.NDTest("SetStrain", 1, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+   stresses = op.NDTest("GetStress", 1)
+   stress11 = stresses[0]
+   print(f"  Applied strain: {strain}")
+   print(f"  Resulting stress σ₁₁: {stress11} MPa")
+   print(f"  Expected stress: {200000 * strain} MPa")
+   print("")
+   
+   # Test 2: Get tangent stiffness
+   print("Test 2: Tangent Stiffness")
+   tangent = op.NDTest("GetTangentStiffness", 1)
+   C11 = tangent[0]
+   C12 = tangent[1]
+   print(f"  C₁₁: {C11} MPa")
+   print(f"  C₁₂: {C12} MPa")
+   print(f"  Expected C₁₁: {200000} MPa")
+   print(f"  Expected C₁₂: {200000 * 0.3 / (1 - 0.3)} MPa")
+   print("")
+   
+   # Commit the state
+   op.NDTest("CommitState", 1)
+
+Plasticity Material Testing
+============================
+
+TCL:
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 2: J2 Plasticity Material Test
+   # =============================================
+   
+   # Create J2 plasticity material
+   # E = 200 GPa, ν = 0.3, yield stress = 300 MPa
+   nDMaterial J2Plasticity 2 200000 0.3 0.001 300 20.0 0.0
+   
+   puts "=== J2 Plasticity Test ==="
+   puts "Material Properties:"
+   puts "  E = 200 GPa, ν = 0.3"
+   puts "  Yield stress = 300 MPa"
+   puts "  Isotropic hardening = 20.0 MPa"
+   puts ""
+   
+   # Monotonic loading beyond yield
+   puts "Monotonic Loading Test:"
+   set maxStrain 0.005
+   set numSteps 50
+   set file [open "j2_test_output.txt" w]
+   puts $file "strain stress11"
+   
+   for {set i 0} {$i <= $numSteps} {incr i} {
+       set strain [expr {$maxStrain * $i / $numSteps}]
+       NDTest SetStrain 2 $strain 0.0 0.0 0.0 0.0 0.0
+       set stresses [NDTest GetStress 2]
+       set stress11 [lindex $stresses 0]
+       puts $file "$strain $stress11"
+       
+       if {$i % 10 == 0} {
+           puts "  Step $i: ε = $strain, σ₁₁ = $stress11 MPa"
+       }
+   }
+   
+   close $file
+   puts "Results saved to j2_test_output.txt"
+   puts ""
+   
+   # Cyclic loading test
+   puts "Cyclic Loading Test:"
+   set cycles 3
+   set amplitude 0.003
+   
+   for {set cycle 1} {$cycle <= $cycles} {incr cycle} {
+       puts "  Cycle $cycle:"
+       
+       # Loading
+       for {set i 0} {$i <= 20} {incr i} {
+           set strain [expr {$amplitude * $i / 20}]
+           NDTest SetStrain 2 $strain 0.0 0.0 0.0 0.0 0.0
+           NDTest CommitState 2
+       }
+       
+       # Unloading
+       for {set i 20} {$i >= 0} {incr i -1} {
+           set strain [expr {$amplitude * $i / 20}]
+           NDTest SetStrain 2 $strain 0.0 0.0 0.0 0.0 0.0
+           NDTest CommitState 2
+       }
+       
+       set stresses [NDTest GetStress 2]
+       puts "    Final stress: [lindex $stresses 0] MPa"
+   }
+
+Python:
+
+.. code-block:: python
+
+   # =============================================
+   # Example 2: J2 Plasticity Material Test
+   # =============================================
+   
+   import openseespy.opensees as op
+   
+   # Create J2 plasticity material
+   # E = 200 GPa, ν = 0.3, yield stress = 300 MPa
+   op.nDMaterial("J2Plasticity", 2, 200000, 0.3, 0.001, 300, 20.0, 0.0)
+   
+   print("=== J2 Plasticity Test ===")
+   print("Material Properties:")
+   print("  E = 200 GPa, ν = 0.3")
+   print("  Yield stress = 300 MPa")
+   print("  Isotropic hardening = 20.0 MPa")
+   print("")
+   
+   # Monotonic loading beyond yield
+   print("Monotonic Loading Test:")
+   max_strain = 0.005
+   num_steps = 50
+   
+   with open("j2_test_output.txt", "w") as f:
+       f.write("strain stress11\n")
+       
+       for i in range(num_steps + 1):
+           strain = max_strain * i / num_steps
+           op.NDTest("SetStrain", 2, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           stresses = op.NDTest("GetStress", 2)
+           stress11 = stresses[0]
+           f.write(f"{strain} {stress11}\n")
+           
+           if i % 10 == 0:
+               print(f"  Step {i}: ε = {strain}, σ₁₁ = {stress11} MPa")
+   
+   print("Results saved to j2_test_output.txt")
+   print("")
+   
+   # Cyclic loading test
+   print("Cyclic Loading Test:")
+   cycles = 3
+   amplitude = 0.003
+   
+   for cycle in range(1, cycles + 1):
+       print(f"  Cycle {cycle}:")
+       
+       # Loading
+       for i in range(21):
+           strain = amplitude * i / 20
+           op.NDTest("SetStrain", 2, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           op.NDTest("CommitState", 2)
+       
+       # Unloading
+       for i in range(20, -1, -1):
+           strain = amplitude * i / 20
+           op.NDTest("SetStrain", 2, strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           op.NDTest("CommitState", 2)
+       
+       stresses = op.NDTest("GetStress", 2)
+       print(f"    Final stress: {stresses[0]} MPa")
+
+Advanced Parameter Study
+=========================
+
+TCL:
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 3: Parameter Study with NDTest
+   # =============================================
+   
+   # Base material properties
+   set E_base 30000
+   set nu_values {0.15 0.20 0.25 0.30 0.35 0.40}
+   set test_strain 0.001
+   
+   puts "=== Poisson's Ratio Effect Study ==="
+   puts "Applied strain: $test_strain"
+   puts ""
+   
+   # Results file
+   set file [open "poisson_study.txt" w]
+   puts $file "poisson_ratio lateral_strain axial_stress"
+   
+   foreach nu $nu_values {
+       # Create material with specific Poisson's ratio
+       nDMaterial ElasticIsotropic 10 $E_base $nu
+       
+       # Apply uniaxial strain
+       NDTest SetStrain 10 $test_strain 0.0 0.0 0.0 0.0 0.0
+       
+       # Get results
+       set strains [NDTest GetStrain 10]
+       set stresses [NDTest GetStress 10]
+       
+       set lateral_strain [lindex $strains 1]
+       set axial_stress [lindex $stresses 0]
+       
+       puts $file "$nu $lateral_strain $axial_stress"
+       
+       puts "ν = $nu:"
+       puts "  Lateral strain ε₂₂: $lateral_strain"
+       puts "  Expected: [expr {-$nu * $test_strain}]"
+       puts "  Axial stress σ₁₁: $axial_stress MPa"
+       puts ""
+       
+       NDTest CommitState 10
+   }
+   
+   close $file
+   puts "Parameter study results saved to poisson_study.txt"
+
+Python:
+
+.. code-block:: python
+
+   # =============================================
+   # Example 3: Parameter Study with NDTest
+   # =============================================
+   
+   import openseespy.opensees as op
+   
+   # Base material properties
+   E_base = 30000
+   nu_values = [0.15, 0.20, 0.25, 0.30, 0.35, 0.40]
+   test_strain = 0.001
+   
+   print("=== Poisson's Ratio Effect Study ===")
+   print(f"Applied strain: {test_strain}")
+   print("")
+   
+   # Results file
+   with open("poisson_study.txt", "w") as f:
+       f.write("poisson_ratio lateral_strain axial_stress\n")
+       
+       for nu in nu_values:
+           # Create material with specific Poisson's ratio
+           op.nDMaterial("ElasticIsotropic", 10, E_base, nu)
+           
+           # Apply uniaxial strain
+           op.NDTest("SetStrain", 10, test_strain, 0.0, 0.0, 0.0, 0.0, 0.0)
+           
+           # Get results
+           strains = op.NDTest("GetStrain", 10)
+           stresses = op.NDTest("GetStress", 10)
+           
+           lateral_strain = strains[1]
+           axial_stress = stresses[0]
+           
+           f.write(f"{nu} {lateral_strain} {axial_stress}\n")
+           
+           print(f"ν = {nu}:")
+           print(f"  Lateral strain ε₂₂: {lateral_strain}")
+           print(f"  Expected: {-nu * test_strain}")
+           print(f"  Axial stress σ₁₁: {axial_stress} MPa")
+           print("")
+           
+           op.NDTest("CommitState", 10)
+   
+   print("Parameter study results saved to poisson_study.txt")
+
+Shear Testing Example
+=====================
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 4: Pure Shear Test
+   # =============================================
+   
+   # Create material
+   nDMaterial ElasticIsotropic 4 70000 0.2
+   
+   puts "=== Pure Shear Test ==="
+   puts "Material: Concrete-like (E = 70 GPa, ν = 0.2)"
+   puts ""
+   
+   # Apply pure shear strain γ₁₂
+   set shear_strain 0.002
+   NDTest SetStrain 4 0.0 0.0 0.0 $shear_strain 0.0 0.0
+   
+   # Get results
+   set strains [NDTest GetStrain 4]
+   set stresses [NDTest GetStress 4]
+   
+   set shear_stress [lindex $stresses 3]
+   set shear_modulus [expr {$shear_stress / $shear_strain}]
+   
+   puts "Applied shear strain γ₁₂: $shear_strain"
+   puts "Resulting shear stress τ₁₂: $shear_stress MPa"
+   puts "Calculated shear modulus G: $shear_modulus MPa"
+   puts "Expected shear modulus: [expr {70000 / (2 * (1 + 0.2))}] MPa"
+   puts ""
+   
+   # Test with multiple shear components
+   puts "Multiple Shear Components Test:"
+   NDTest SetStrain 4 0.0 0.0 0.0 0.001 0.0015 0.0008
+   set stresses [NDTest GetStress 4]
+   puts "  τ₁₂ = [lindex $stresses 3] MPa"
+   puts "  τ₂₃ = [lindex $stresses 4] MPa"  
+   puts "  τ₃₁ = [lindex $stresses 5] MPa"
+   
+   NDTest CommitState 4
+
+Custom Response Testing
+========================
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 5: Custom Material Responses
+   # =============================================
+   
+   # Create a material with custom responses
+   nDMaterial J2Plasticity 5 200000 0.3 0.001 300 20.0 0.0
+   
+   puts "=== Custom Response Testing ==="
+   
+   # Apply plastic deformation
+   NDTest SetStrain 5 0.002 0.0 0.0 0.0 0.0 0.0
+   NDTest CommitState 5
+   
+   # Get custom responses
+   puts "Available custom responses:"
+   
+   # Try to get plastic strain (if supported)
+   if {[catch {set plasticStrain [NDTest GetResponse 5 "plasticStrain"]}]} {
+       puts "  plasticStrain response not available"
+   } else {
+       puts "  Plastic strain: $plasticStrain"
+   }
+   
+   # Try to get backstress (if supported)  
+   if {[catch {set backstress [NDTest GetResponse 5 "backstress"]}]} {
+       puts "  backstress response not available"
+   } else {
+       puts "  Backstress: $backstress"
+   }
+   
+   # Try to get yield function value
+   if {[catch {set yieldFunc [NDTest GetResponse 5 "yieldFunction"]}]} {
+       puts "  yieldFunction response not available"
+   } else {
+       puts "  Yield function value: $yieldFunc"
+   }
+   
+   puts ""
+   puts "Basic responses always available:"
+   set stress [NDTest GetResponse 5 "stress"]
+   set strain [NDTest GetResponse 5 "strain"] 
+   puts "  Stress: $stress"
+   puts "  Strain: $strain"
+
+Runtime Parameter Updates
+==========================
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 6: Runtime Parameter Updates
+   # =============================================
+   
+   # Create material
+   nDMaterial ElasticIsotropic 6 200000 0.3
+   
+   puts "=== Runtime Parameter Update Test ==="
+   
+   # Initial test
+   NDTest SetStrain 6 0.001 0.0 0.0 0.0 0.0 0.0
+   set initial_stress [lindex [NDTest GetStress 6] 0]
+   puts "Initial stress: $initial_stress MPa"
+   
+   # Update Young's modulus (assuming parameter ID 1 for E)
+   # Note: Parameter IDs are material-specific - check material documentation
+   set new_E 150000
+   puts "Updating Young's modulus to: $new_E MPa"
+   
+   if {[catch {NDTest UpdateDoubleParameter 6 1 $new_E}]} {
+       puts "Failed to update parameter (ID may be incorrect for this material)"
+   } else {
+       # Test with updated parameter
+       NDTest SetStrain 6 0.001 0.0 0.0 0.0 0.0 0.0
+       set updated_stress [lindex [NDTest GetStress 6] 0]
+       puts "Updated stress: $updated_stress MPa"
+       puts "Expected stress: [expr {$new_E * 0.001}] MPa"
+   }
+   
+   # Update integer parameter example (if applicable)
+   puts "Attempting to update integer parameter..."
+   if {[catch {NDTest UpdateIntegerParameter 6 1 0}]} {
+       puts "Integer parameter update failed (not supported or incorrect ID)"
+   } else {
+       puts "Integer parameter updated successfully"
+   }
+   
+   NDTest CommitState 6
+
+Automated Testing Procedure
+===========================
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 7: Automated Material Testing Procedure
+   # =============================================
+   
+   proc testMaterial {matTag testName testProcs} {
+       puts "=== $testName ==="
+       
+       foreach procName $testProcs {
+           puts "Running $procName..."
+           if {[catch {$procName $matTag} result]} {
+               puts "  ERROR: $result"
+           } else {
+               puts "  PASSED"
+           }
+       }
+       puts ""
+   }
+   
+   proc testElasticResponse {matTag} {
+       # Test linear elastic response
+       set testStrain 0.001
+       NDTest SetStrain $matTag $testStrain 0.0 0.0 0.0 0.0 0.0
+       set stresses [NDTest GetStress $matTag]
+       set stress11 [lindex $stresses 0]
+       
+       # Get tangent stiffness
+       set tangent [NDTest GetTangentStiffness $matTag]
+       set C11 [lindex $tangent 0]
+       
+       # Basic sanity checks
+       if {$C11 <= 0} {
+           return "Invalid tangent stiffness: $C11"
+       }
+       
+       if {abs($stress11) < 1e-10} {
+           return "No stress response to applied strain"
+       }
+       
+       return "OK"
+   }
+   
+   proc testSymmetry {matTag} {
+       # Test material symmetry for elastic materials
+       # Apply strain in direction 1, check response in direction 2
+       NDTest SetStrain $matTag 0.001 0.0 0.0 0.0 0.0 0.0
+       set stresses1 [NDTest GetStress $matTag]
+       set stress22_1 [lindex $stresses1 1]
+       
+       # Apply strain in direction 2, check response in direction 1  
+       NDTest SetStrain $matTag 0.0 0.001 0.0 0.0 0.0 0.0
+       set stresses2 [NDTest GetStress $matTag]
+       set stress11_2 [lindex $stresses2 0]
+       
+       # For isotropic materials, these should be equal
+       set tolerance 1e-6
+       if {abs($stress22_1 - $stress11_2) > $tolerance} {
+           return "Material symmetry violated: $stress22_1 != $stress11_2"
+       }
+       
+       return "OK"
+   }
+   
+   proc testEnergyConservation {matTag} {
+       # Simple energy conservation test
+       # Apply strain, get work done, remove strain, check work recovery
+       
+       # Loading
+       NDTest SetStrain $matTag 0.001 0.0 0.0 0.0 0.0 0.0
+       set stresses [NDTest GetStress $matTag]
+       set stress11 [lindex $stresses 0]
+       
+       # Unloading  
+       NDTest SetStrain $matTag 0.0 0.0 0.0 0.0 0.0 0.0
+       set stresses_unload [NDTest GetStress $matTag]
+       set stress11_unload [lindex $stresses_unload 0]
+       
+       # For elastic materials, stress should return to zero
+       if {abs($stress11_unload) > 1e-10} {
+           return "Elastic material not returning to zero stress: $stress11_unload"
+       }
+       
+       return "OK"
+   }
+   
+   # Run automated tests
+   nDMaterial ElasticIsotropic 7 200000 0.3
+   set testProcs [list testElasticResponse testSymmetry testEnergyConservation]
+   testMaterial 7 "Elastic Material Validation" $testProcs
+
+Integration with Analysis Workflow
+===================================
+
+.. code-block:: tcl
+
+   # =============================================
+   # Example 8: Integration with Structural Analysis
+   # =============================================
+   
+   # This example shows how to use NDTest for material calibration
+   # before running a full structural analysis
+   
+   proc calibrateMaterial {targetModulus targetPoisson} {
+       puts "Material Calibration:"
+       puts "Target: E = $targetModulus, ν = $targetPoisson"
+       
+       # Start with initial guess
+       set E_guess [expr {$targetModulus * 0.8}]
+       set nu_guess [expr {$targetPoisson * 0.9}]
+       
+       for {set iteration 1} {$iteration <= 10} {incr iteration} {
+           puts "Iteration $iteration: E = $E_guess, ν = $nu_guess"
+           
+           # Create test material
+           nDMaterial ElasticIsotropic 99 $E_guess $nu_guess
+           
+           # Apply test strain
+           NDTest SetStrain 99 0.001 0.0 0.0 0.0 0.0 0.0
+           
+           # Measure response
+           set stresses [NDTest GetStress 99]
+           set measured_stress [lindex $stresses 0]
+           set measured_modulus [expr {$measured_stress / 0.001}]
+           
+           # Measure Poisson effect
+           set strains [NDTest GetStrain 99]
+           set lateral_strain [lindex $strains 1]
+           set measured_poisson [expr {-$lateral_strain / 0.001}]
+           
+           puts "  Measured: E = $measured_modulus, ν = $measured_poisson"
+           
+           # Update guesses (simple proportional correction)
+           set E_guess [expr {$E_guess * $targetModulus / $measured_modulus}]
+           set nu_guess [expr {$nu_guess * $targetPoisson / $measured_poisson}]
+           
+           # Check convergence
+           if {abs($measured_modulus - $targetModulus) < 100 && \
+               abs($measured_poisson - $targetPoisson) < 0.01} {
+               puts "Calibration converged!"
+               return [list $E_guess $nu_guess]
+           }
+       }
+       
+       puts "Calibration did not fully converge"
+       return [list $E_guess $nu_guess]
+   }
+   
+   # Run calibration
+   set calibratedProps [calibrateMaterial 30000 0.2]
+   set finalE [lindex $calibratedProps 0]
+   set finalNu [lindex $calibratedProps 1]
+   
+   puts "Calibrated material properties:"
+   puts "E = $finalE MPa"
+   puts "ν = $finalNu"
+   
+   # Now use calibrated material in structural model
+   model BasicBuilder -ndm 3 -ndf 3
+   nDMaterial ElasticIsotropic 1 $finalE $finalNu
+   
+   # ... continue with structural analysis ...
+
+These examples demonstrate the versatility of NDTest commands for material testing, calibration, validation, and research applications.
\ No newline at end of file