Class Introduction
Here is an example that uses a class for calculating the descriptive statistics of an array: class_demo.py
Classes that act like functions
Python classes with function behavior are useful for representing a function with separate coefficients and independent variables. The coefficients themselves may be variable.
Here is an example: class_function_demo.py
Tom Irvine
Introduction
Python can read, write and play sound files. The following scripts can run under Windows and Ubuntu.
* * *
Utility
Each script requires the utility module tompy.py
* * *
wav Files
Here is an example for writing and playing a wav file: generate_wav.py
Here is an example for reading and playing a wav file and for displaying its FFT magnitude: wav_player.py. This script also has an option to convert the wav file to a time history in ASCII text format.
The audio scripts determine whether the operating system is Windows or Linux.
The sound is played using:
* * *
mp3 Files
Here is an mp3 player: mp3_player.py
* * *
pygame
python-pygame can be installed from the Ubuntu Software Center
There are other methods for playing sound files in Python, including SFML, Snack, ossaudiodev, pyglet, PyMedia and wxPython.
See also: Matlab Fun with Sound Files
* * *
Tom Irvine
Acronyms
ODE = ordinary differential equation
SDOF = single-degree-of-freedom
MDOF = multi-degree-of-freedom
* * *
Supporting Functions
The scripts on this page require the utility modules:
* * *
Runge-Kutta
The Runge-Kutta family of numerical methods may be used to solve ordinary differential equations with initial conditions.
The solution is given in the time domain. Thus, the Runge-Kutta method may be used for modal transient analysis.
Note that the Runge-Kutta method may give unstable results for certain “stiff” systems.
The Newmark-beta ODE Solver is thus preferred in most structural dynamics texts.
Another robust method is the Digital Recursive Filtering Relationship ODE Solver.
* * *
odeint
The odeint function in the scipy.integrate package is a wrapper around lsoda, an ODE integrator developed at Lawrence Livermore National Lab, as part of the odepack package. lsoda/odeint automatically switches between stiff and non-stiff integration routines, depending on the characteristics of the solution, and does adaptive time-stepping to achieve a desired level of solution accuracy.
The function may be accessed via
from scipy.integrate import odeint
* * *
Free Vibration Example
A script for the response of a single-degree-of-freedom (SDOF) system subjected to initial condition excitation is given at: sdof_initial.py. An example is given in: sdof_initial_py.pdf.
A similar multi-degree-of-freedom (MDOF) system script is given at: mdof_initial.py
The MDOF script shows how odeint can be applied to simultaneous second-order ODEs.
* * *
Response to Sinusoidal Base Excitation
A script for the response of an SDOF system subjected to sinusoidal base excitation is given at: sdof_base_sine.py
An example is given in: sdof_base_py.pdf.
* * *
Response to Arbitrary Base Input
The odeint solver is intended for initial condition excitation and for external excitation where the forcing function is a continuous analytical function.
The odeint function can be used for a random or arbitrary excitation as follows. First, take a Fourier transform of the forcing function. Second, pass the Fourier coefficients and frequencies to the ode function. Third, represent the excitation time history via a series of analytical sine and cosine terms inside the ode function. This method requires careful bookkeeping. It is also numerically inefficient. Thus, it is only recommended in academic settings.
A faster approach is to simply use the Runge-Kutta Fourth-order method, as shown in the following script.
A script for the response of an SDOF system subjected to arbitrary base excitation is given at: arbit_rk4.py
* * *
SDOF Response to Arbitrary Applied Force
The script is given at: arbit_force_rk4.py
* * *
MDOF Response to Arbitrary Applied Force via Direct Integration
This approach is intended for the case where the damping is specified in terms of a viscous damping coefficient matrix.
The script is given at: mdof_arbit_force_rk4.py.
* * *
MDOF Response to Arbitrary Applied Force via Modal Transient Solution
This approach is intended for the case where the damping is specified in terms of modal damping ratios. The system is decoupled via the eigenvectors as an intermediate step.
The script is given at: mdof_modal_arbit_force_rk4.py.
* * *
MDOF Response to Base Excitation, Modal Transient
The response to base acceleration or enforced acceleration is given at: mdof_modal_enforced_acceleration_rk4.py. It is intended for the case where damping is applied via a modal damping ratio. The system is first partitioned via a transformation matrix. Next it is decoupled using normal modes. Then the response is calculated using the RK4 method as a modal transient solution. A partition tutorial is given at: modal_enforced_motion.pdf
Acronyms
ODE = ordinary differential equation
SDOF = single-degree-of-freedom
MDOF = multi-degree-of-freedom
* * *
Supporting Functions
The scripts on this page require the utility modules:
Some of the scripts use the generalized eigenvalue function to calculate the natural frequencies. This step is not part of the Newmark calculation except for the modal transient implementations. Rather the natural frequencies are given for reference only for the direct implementation.
* * *
Newmark-beta Function
The Newmark-beta method is used for direct integration of a system of structural dynamics equations.
The system equations are second-order ordinary differential equations. The system may be excited by initial conditions or an external forcing function.
* * *
Free Vibration
The following examples use the Newmark function.
The response of a single-degree-of-freedom system to initial excitation is given at: sdof_initial_nm.py
The response of a multi-degree-of-freedom system to initial excitation is given at: mdof_initial_nm.py
* * *
Response to Applied Force, Direct Integration
The response of a multi-degree-of-freedom system to an applied force or forces is given at: mdof_arbit_force_nm.py. It is intended for the case where damping is applied via a damping coefficient matrix.
* * *
Response to Applied Force, Modal Transient
The response of a multi-degree-of-freedom system to an applied force or forces is given at: mdof_modal_arbit_force_nm.py. It is intended for the case where damping is applied via a modal damping ratio. The system is decoupled via the normal modes as an intermediate step.
* * *
Response to Base Excitation, Modal Transient
The response of a multi-degree-of-freedom system to base acceleration or enforced acceleration is given at: mdof_modal_enforced_acceleration_nm.py. It is intended for the case where damping is applied via a modal damping ratio. The system is first partitioned via a transformation matrix. Next it is decoupled using normal modes. Then the response is calculated using the Newmark method as a modal transient solution. A partition tutorial is given at: modal_enforced_motion.pdf
* * *
– Tom Irvine
Acronyms
ODE = ordinary differential equation
SDOF = single-degree-of-freedom
MDOF = multi-degree-of-freedom
* * *
Utilities
The scripts on this page require the utility modules
* * *
Filter Implementation
The digital recursive filtering relationship can be implemented in Python by importing the lfilter function from scipy.
from scipy.signal import lfilter
* * *
Arbitrary Base Input – Time Domain Response for a given SDOF System
Digital recursive filtering relationship method: arbit.py
* * *
Arbitrary Applied Force – Time Domain Response for a given SDOF System
A script for the SDOF response to an applied force is given at: arbit_force.py
* * *
Arbitrary Applied Force – Time Domain Response for a given MDOF System
A script for the MDOF response to applied forces is given at: mdof_modal_arbit_force_ri.py
* * *
External References
Irvine Impulse Response Function 1
Irvine Impulse Response Function 2
Utility
The script on this page requires the utility module tompy.py
* * *
Python Script
The following script solves a system of equations of the form Ax=b
linear_system.py
The key lines are
from scipy import linalg
x = linalg.solve(A, b)
* * *
Fortran
A Fortran program for this problem is: SLE.F
The program is compiled via
gfortran -o SLE SLE.F -llapack
It uses the Lapack subroutine DGESV.
* * *
Tom Irvine
Utility
The scripts on this page require the utility module tompy.py
* * *
Python Scripts
A script for calculating the inverse of a square matrix is given at: inverse_matrix.py
The script uses the function: from numpy import linalg
Here is a script which calculates the pseudo inverse of a singular matrix: pseudo_inverse_matrix.py
* * *
Fortran Program
Here is a Fortran program which performs matrix inversion using the LU decomposition method: INVERSE_MATRIX.F
It is compiled via:
gfortran -o INVERSE_MATRIX INVERSE_MATRIX.F -llapack -lblas
The program uses the subroutines: DGETRF & DGETRI
It compiles & runs under both Ubuntu & Cygwin.
* * *
See also: http://www.nag.com/numeric/fl/nagdoc_fl23/examples/source/f07ajfe.f90
Tom Irvine
Utility
The scripts on this page require the utility module tompy.py
* * *
Rectangular Plates
A script for calculating the natural frequency of a rectangular plate supported at each corner is given at plate_corners.py
The natural frequency is calculated via the Rayleigh method. The resulting fundamental mode shape is displayed using a Matplotlib surface plot.
* * *
Circular Plates
A script for calculating the natural frequencies of the bending modes of a circular plate, with options for both homogenous and honeycomb-sandwich designs is given at circular_plate.py. This script demonstrates a 3D surface plot for a polar base.
* * *
– Tom Irvine
Utility
The scripts on this page require the utility module tompy.py
* * *
f2py
f2py is a Fortran to Python interface generator.
* * *
Sample Project
The following project was performed on Ubuntu.
Here is a project that solves the generalized eigenvalue problem.
The mass & stiffness matrices are read in by a Python script.
The Python script then passes the arrays to a Fortran program for the number-crunching.
The Fortran program then returns the eigenvalues & vectors to the Python script.
The project includes an option for the standard eigenvalue problem.
Note that this project requires Lapack.
Also the python-dev header & static library must be installed.
* * *
Source Code
Here are the two Fortran scripts and their compilation commands: eigen_general.f
f2py -c -m eigen_general eigen_general.f –fcompiler=gnu95 -llapack
eigen_standard.f
f2py -c -m eigen_general eigen_general.f –fcompiler=gnu95 -llapack
Here is the Python script:
python eigenvalues.py -i
Note that the Python script imports the Fortran subrountines via:
import eigen_standard import eigen_general
Related Posts
Generalized Eigenvalue Problem
Vibrationdata Blog Entry – Instructions for the generalized eigenvalue problem including Lapack installation.
* * *
– Tom Irvine
Utilities
The Python scripts on this page require the utility modules
* * *
Primary Method
A script for calculating the eigenvalues and vectors for a mass and stiffness matrix system is given at: gen_eig.py
The script uses the following functions:
from scipy.linalg import eig, eigh
The references for these functions are:
eig
eigh
They provide an interface to the LAPACK routines dgeev and zgeev. These routines implement the QR method with Householder reflectors.
* * *
Sparse Matrices
Here is a script for find a subset of the lowest eigenvalues and vectors for a sparse system: gen_eig_sparse.py
The scripts uses the functions:
from scipy.sparse.linalg import eigsh
from scipy.sparse.linalg import eigs
It also uses these functions in shift-inverted mode to find the lowest eigenvalues with speed and accuracy.
Reference: Sparse Eigenvalue Problems with ARPACK
* * *
Alternate Method
Another method is to use the Locally Optimal Block Preconditioned Conjugate Gradient Method (LOBPCG).
from scipy.sparse.linalg import lobpcg
But there are anecdotal reports that this method is inaccurate in some cases.
* * *
Fortran
Here is a pure Fortran program for the generalized eigenvalue problem: geigen.f.
The program can be compiled via:
$ gfortran -o geigen geigen.f -lblas -llapack
The program uses LAPACK routines DSYGV and DSYEV.
* * *
C/C++
Here is C/C++ program for the generalized eigenvalue problem: gen_eig.cpp.
The program can be compiled via:
gcc -o gen_eig gen_eig.cpp -lblas -llapack -lstdc++
* * *
Scilab
A Scilab program is given at: generalized_eigen_read.sci
* * *
Related Posts
Mixed Fortran & Python – Gives an example using the generalized eigenvalue problem
Vibrationdata Blog Entry – Instructions for the generalized eigenvalue problem including Lapack installation.
* * *
– Tom Irvine
Vibrationdata Python 