processes are involved the impact time(s) provide additional time constants, which
may be much shorter than those of the unimpacting system. The same applies when
other discontinuous or abruptly changing processes are involved, e.g. as with dry
friction or intermittent contact problems.
Stiff solvers operate with large integration steps, at the expense of introducing
possibly large errors on the shortest time scales. If this is unacceptable, do not
choose a stiff solver, despite the span in time scales.
With stiff solvers the cost of each integration step is usually large. So only
choose a stiff solver if the Jacobian of the system does not need to be recomputed
by numerical finite differencing at every time step, that is, if the system is linear or
weakly linear, or if you can provide the Jacobian in analytical form to the solver.
Go through all options of the solver and choose or experiment with these
carefully.
Integration order For smooth systems higher-order algorithms (like MATLAB’s
ode45 Runge-Kutta fourth/fifth order), are usually superior to the corresponding
lower-order ones (like MATLAB’s ode23 Runge-Kutta second/third order). With
non-smooth/discontinuous systems or when accuracy requirements are modest,
lower-order algorithms can be more efficient.
A.2. Computing Chaos-Related Quantities
For studying chaotic vibrations one need to compute phase plane orbits, Poincaré
maps, Lyapunov exponents, and perhaps Lyapunov dimensions and frequency
spectra. Transients should be discarded by allowing sufficient simulation time to
pass before displaying data.
Phase plane graphs are obtained by displaying two state variables of the system
versus each other. For the Duffing equation above, e.g., one displays v(t) versus u
(t) where the two quantities are output from the ODE-solver.
Poincaré maps of periodically driven systems are obtained by solving the ODEs
with a time-step equal to the period of the driving force. Most ODE-solvers then
internally iterate in small time-steps, to output only the solution at the end t k of the
user-specified time interval. To display a Poincaré map for the above Duffing
equation, e.g., one displays v(t k ) versus u(t k ) using small dots.
Lyapunov exponents and dimensions are calculated as described in, e.g., Wolf
et al. (1985), Parker and Chua (1989), Moon (1987), or Kantz (1994). Freely
downloadable software exists for this (see below).
Frequency spectra are obtained by employing an FFT algorithm (Fast Fourier
Transform) to the time data output from the ODE-solver. Usually one takes the
logarithm and displays the result in dB as a bar graph. Most commercial software
libraries contain FFT routines. Consult, e.g., Press et al. (2002), should you consider coding yourself. For a readable treatise on digital signal processing of simulated or experimental vibration data see Brandt (2011).
Appendix A: Performing Numerical Simulation
451
may be much shorter than those of the unimpacting system. The same applies when
other discontinuous or abruptly changing processes are involved, e.g. as with dry
friction or intermittent contact problems.
Stiff solvers operate with large integration steps, at the expense of introducing
possibly large errors on the shortest time scales. If this is unacceptable, do not
choose a stiff solver, despite the span in time scales.
With stiff solvers the cost of each integration step is usually large. So only
choose a stiff solver if the Jacobian of the system does not need to be recomputed
by numerical finite differencing at every time step, that is, if the system is linear or
weakly linear, or if you can provide the Jacobian in analytical form to the solver.
Go through all options of the solver and choose or experiment with these
carefully.
Integration order For smooth systems higher-order algorithms (like MATLAB’s
ode45 Runge-Kutta fourth/fifth order), are usually superior to the corresponding
lower-order ones (like MATLAB’s ode23 Runge-Kutta second/third order). With
non-smooth/discontinuous systems or when accuracy requirements are modest,
lower-order algorithms can be more efficient.
A.2. Computing Chaos-Related Quantities
For studying chaotic vibrations one need to compute phase plane orbits, Poincaré
maps, Lyapunov exponents, and perhaps Lyapunov dimensions and frequency
spectra. Transients should be discarded by allowing sufficient simulation time to
pass before displaying data.
Phase plane graphs are obtained by displaying two state variables of the system
versus each other. For the Duffing equation above, e.g., one displays v(t) versus u
(t) where the two quantities are output from the ODE-solver.
Poincaré maps of periodically driven systems are obtained by solving the ODEs
with a time-step equal to the period of the driving force. Most ODE-solvers then
internally iterate in small time-steps, to output only the solution at the end t k of the
user-specified time interval. To display a Poincaré map for the above Duffing
equation, e.g., one displays v(t k ) versus u(t k ) using small dots.
Lyapunov exponents and dimensions are calculated as described in, e.g., Wolf
et al. (1985), Parker and Chua (1989), Moon (1987), or Kantz (1994). Freely
downloadable software exists for this (see below).
Frequency spectra are obtained by employing an FFT algorithm (Fast Fourier
Transform) to the time data output from the ODE-solver. Usually one takes the
logarithm and displays the result in dB as a bar graph. Most commercial software
libraries contain FFT routines. Consult, e.g., Press et al. (2002), should you consider coding yourself. For a readable treatise on digital signal processing of simulated or experimental vibration data see Brandt (2011).
Appendix A: Performing Numerical Simulation
451