Simulation¶
Define the axes outside the simulation and make one explicit specification:
import dynaresp as dr
spec = dr.SimulationSpec(
system=system,
load=load,
sensors=tuple(sensor_layout),
time=time,
positions=positions,
)
result = dr.run_simulation(spec, seed=42)
time must be finite, strictly increasing, and uniformly sampled. positions
must be finite and strictly increasing. Neither axis is generated or inferred.
The result stores:
| Field | Shape |
|---|---|
time |
(time,) |
positions |
(position,) |
load_values |
(position, time) |
force |
(DOF, time) |
displacement, velocity, acceleration |
(DOF, time) |
signal, clean_signal |
(sensor, time) |
mode_shapes |
(mode, position) |
The model is linear and time invariant:
[ M\ddot{u}(t) + C\dot{u}(t) + Ku(t) = f(t). ]
Inputs are held constant between time samples. Initial displacement and velocity are zero.
Simulation pipeline¶
One call to run_simulation() performs the following operations:
- Validate the system, load, sensors, time axis, and position axis.
- Generate the diagnostic load field and any exact point-force trajectories.
- Project that field into the system's generalized coordinates.
- Solve for generalized displacement, velocity, and acceleration.
- Sample the requested physical response at every sensor position.
- Apply each sensor's measurement effects in order.
- Return both the measured signals and the clean generating physics.
If (h(x_j)) is the system's observation vector at load position (x_j), the columns form the projection matrix
[ H = \begin{bmatrix}h(x_1) & \cdots & h(x_p)\end{bmatrix}. ]
For a distributed load field (F) with shape (p, time), the generalized
force passed to the solver is (f = HF). Each exact point force (q_k(t)) at
(x_k(t)) adds (h(x_k(t))q_k(t)). The same observation model samples a
response at a sensor position (x_s), for example
(u_s = h(x_s)^T u).
This separation matters for beams: changing the reporting positions grid
does not move a point load or alter its generalized force. load_values
remains the nearest-grid diagnostic, and point_forces preserves the exact
trajectories used by the solver.
State-space formulation¶
The solver uses displacement and velocity as its state:
[ x(t) = \begin{bmatrix} u(t) \ \dot{u}(t) \end{bmatrix}. ]
Solving the equation of motion for acceleration gives
[ \ddot{u}(t) = -M^{-1}K u(t) - M^{-1}C\dot{u}(t) + M^{-1}f(t). ]
Combining (\dot{u}) and (\ddot{u}) produces the first-order system
[ \dot{x}(t) = \underbrace{ \begin{bmatrix} 0 & I \ -M^{-1}K & -M^{-1}C \end{bmatrix}}{A} x(t) + \underbrace{ \begin{bmatrix} 0 \ M^{-1} \end{bmatrix}} f(t). ]
For a system with (n) degrees of freedom:
| Quantity | Shape | Meaning |
|---|---|---|
| (M, C, K) | (n, n) |
Mass, damping, and stiffness matrices |
| (u, \dot{u}, \ddot{u}, f) | (n,) |
Response and force at one time |
| (x) | (2n,) |
Displacement and velocity state |
| (A) | (2n, 2n) |
State matrix |
| (B) | (2n, n) |
Force-input matrix |
The implementation computes (M^{-1}K), (M^{-1}C), and (M^{-1}) with linear solves rather than explicitly inverting (M). The solver output is
[ \begin{bmatrix} u \ \dot{u} \ \ddot{u} \end{bmatrix} = \begin{bmatrix} I & 0 \ 0 & I \ -M^{-1}K & -M^{-1}C \end{bmatrix}x + \begin{bmatrix} 0 \ 0 \ M^{-1} \end{bmatrix}f. ]
This lets SciPy return displacement, velocity, and acceleration together.
With interp=False, each force sample is treated as constant until the next
time sample (a zero-order hold).