FreeSRS

Info
Download

The shock response spectrum (SRS) allows to estimate transient acceleration signals in terms of a maximum response of a dynamic vibration system. Its origin probably dates from civil engineering and the problem of seismic excitations, i.e. earthquake response of buildings, see Clough Penzien 1975. The SRS is a spectrum, i.e. a curve giving the maximum acceleration response depending on frequency. It is derived from a 1-mass-spring-damper system being excited by a base excitation, i.e. an acceleration. It makes use of the nice fact, that the acceleration response of a one-mass-spring-damper system with base excitation does not depend on the mass, but only on eigenfrequency and damping.

freesrs calculates the shock response spectrum (SRS) and is a small library with Python and Fortran files. The Fortran files consist of

  • SmallwoodFilter.f and

  • GenShockTimeHistories.f

They are linked to Python by using the comfortable f2py tool. The algorithm used is from

David O. Smallwood, see An improved recursive formula for calculating shock response spectra. Shock and Vibration Bulletin, No. 51, May 1981.

For more background information you may also have a look on http://en.wikipedia.org/wiki/Shock_response_spectrum.

The code was developed on a windows computer and compiled using Python(x,y) and f2py.

There is also a pure C++ version which compiles with gcc. I intended to use the Fortran version together with Python, i.e. as the "fast" library. The program is very fast, all CPU-time-critical parts are written in Fortran 77.

There are 2 modules (or libraries):

  • 'SmallwoodFilter' is auto-generated with f2py (version:2).

    Functions:

    ypp = smallwoodfilter(model,fraction,xi,omega_n,ts,t,xpp,n)

    ypp = smallwoodfilterabsaccelmodelexact(xi,omega_n,ts,t,xpp,n)

    ypp = smallwoodfilterreldisplmodelexact(xi,omega_n,ts,t,xpp,n)

    yppm = recursiveequ(b0,b1,b2,a1s,a2s,xppm,xppm1,xppm2,yppm1,yppm2)

    b0,b1,b2,a1s,a2s = absaccelmodel(xi,omega_n,ts,t)

    b0,b1,b2,a1s,a2s = reldisplmodel(xi,omega_n,ts,t)

    ypp = smallwoodfilterapp(b0,b1,b2,a1s,a2s,xpp,n)

    b0,b1,b2,a1s,a2s = absaccelmodelapprox(xi,omega_n,ts)

    b0,b1,b2,a1s,a2s = reldisplmodelapprox(xi,omega_n,ts)

  • 'GenShockTimeHistories' is auto-generated with f2py (version:2).

    Functions:

    t,a = gensquaredsinusshock(f_shock,fs,ampli,n)

    t,a = genrectangularshock(t_shock,fs,ampli,n)

    t = gentimevector(ts,n)

    noofconstlogfreq1 = noofconstlogfreq1(nfo,f1,f2)

    noofconstlogfreq2 = noofconstlogfreq2(nfo,f1,f2)

    f = compfreqconstlog(nfo,nf,f1,f2)

Folder cpp

  • It contains the cpp and exe files (debug and release version) compiled for windows. You may use the acc.inp file (half-sinus transient) for testing.

I recommend using the f77 / Python version. There, all spectra in the time and the frequency domain are accessible to the user via the Python environment, e.g. when running ipython. This allows easy analysis of all kinds of SRS, i.e. positive, negative or primary and residual SRS. For generating the plots I used matplotlib (which is included in Python(x,y)).

The code should also run under Linux or any other OS! Please, let me know when you compile the code for another OS that I can provide it to the community.



Dependencies

The following python packages are needed

  • numpy

  • PyLab

  • Fortran compiler, e.g. gfortran, mingw32, ...

There might be bugs! Please report bugs to dynamicsmec.



Examples

SRS of sinus shock:

SRS of rectangular shock:

Python code:

#     Copyright (C) 2004-2013 by FreeDynamics

#     http://sourceforge.net/projects/FreeSRS/

#     Licensed under the GNU GPL

 

from numpy import *

from pylab import *

import freesrs

 

 

print freesrs.SmallwoodFilter.__doc__

print freesrs.GenShockTimeHistories.__doc__

 

# half sine shock

f_shock = 400. # Hz shock frequency

fs = 50000. # Hz sampling frequency

amp = 1.

 

t_shock = 1./f_shock/2.; ts = 1./fs; nl = int(t_shock/ts)+1

t,a = freesrs.GenShockTimeHistories.gensquaredsinusshock(f_shock,fs,amp,nl)

figure(1)

plot(t,a)

xlabel('sec')

ylabel('m/s**2')

f1 = 10.; f2 = 4000.; nfo = 30

nf = freesrs.GenShockTimeHistories.noofconstlogfreq1(nfo,f1,f2)

#print nf

frev = freesrs.GenShockTimeHistories.compfreqconstlog(nfo,nf,f1,f2)

#print frev

 

xi = 0.01

srs = arange(float(len(frev)))

for i in range(len(frev)):

    omega_n = 2.*pi*frev[i]

    ypp = freesrs.SmallwoodFilter.smallwoodfilter('a',0.1,xi,omega_n,ts,t,a,nl)

    max_1 = max(ypp)

    # determine ypp(tend) and yp(tend) if x0pp = 0 by num integration

    app = ypp[nl-1]; ap = 0.5*(ypp[nl-1]+ypp[nl-2])*ts

    max_2 = freesrs.AnalysePostVib(ap, app, xi, omega_n)

    #print frev[i], max_1, max_2

    if max_1 > max_2:

        srs[i] = max_1

        #print frev[i]

    else:

        srs[i] = max_2

        #print frev[i], max_1, max_2

 

figure(2)

plot (frev, srs)

xlabel('Hz')

ylabel('m/s**2')

show()


				
	









Copyright (C) 2004-2013 by
dynamicsmec