For a film loop to make sense to a viewer,
the subject of the **movie** must be periodic,
and organized so that the last frame leads naturally into the first.
In the **movie** created by `wavemovie()`
there is a parameter
`lambda`
that controls the basic repetition rate of wave pulses fired
onto the screen from the top.
When a wavelet travels one-quarter of the way down the frame,
another is sent in.
This is defined by the line

Take any point in (*x*,*z*)-space.
The signal there will be a superposition of sinusoids
of various frequencies, .We can choose what frequencies
we will use in the calculation
and what amplitudes and phases we will attach to the initial conditions
at those frequencies.
Here we will simply take uniformly spaced sinusoids of unit amplitude
and no phase.
The `nw` frequencies are
,,...,
`nw`.The lowest frequency `dw` = must be inversely proportional to the wavelength
`lambda`

Finally, the time duration of the film loop must equal the period of the lowest-frequency sinusoid

This latter equation defines the time interval on the line If you use more frequencies, you might like the result better because the wave pulses will be shorter, and the number of wavelengths between the pulses will increase. Thus the quiet zones between the pulses will get quieter. The frequency components can be weighted differently--but this becomes a digression into simple Fourier analysis.# from par: integer n3:nt=12, n2:nx=48, n1:nz=96, nw=2, nlam=4 # from par: real dx=2, dz=1, v=1 # subroutine wavemovie( nt, nx, nz, nw, nlam, dx,dz,v, p, cd, q) integer it,nt,ix,nx,iz,nz,iw,nw, nlam real dx,dz,v, phase,pi2,z0,x0,dt,dw,lambda,w,wov,x, p(nz,nx,nt) complex aa,a,b,c,cshift, cd(nx),q(nx) lambda=nz*dz/nlam; pi2=2.*3.141592; dw=v*pi2/lambda; dt=pi2/(nt*dw) x0 = nx*dx/3; z0 = nz*dz/3 call null( p, nz*nx*nt) do iw = 1,nw { # superimpose nw frequencies w = iw*dw; wov = w/v # frequency / velocity do ix = 1,nx { # initial conditions for a x = ix*dx-x0; # collapsing spherical wave phase = -wov*sqrt( z0**2+x**2) q(ix) = cexp( cmplx( 0.,phase)) } aa = (0.,1.)*dz/(4.*dx**2*wov) # tridiagonal matrix coefficients a = -aa; b = 1.+2.*aa; c = -aa do iz = 1,nz { # extrapolation in depth do ix = 2,nx-1 # diffraction term cd(ix) = aa*q(ix+1) + (1.-2.*aa)*q(ix) + aa*q(ix-1) cd(1) = 0.; cd(nx) = 0. call ctris( nx,-a,a,b,c,-c,cd,q) # Solves complex tridiagonal equations cshift = cexp( cmplx( 0.,wov*dz)) do ix = 1,nx # shifting term q(ix) = q(ix) * cshift do it=1,nt { # evolution in time cshift = cexp( cmplx( 0.,-w*it*dt)) do ix = 1,nx p(iz,ix,it) = p(iz,ix,it) + q(ix)*cshift } } } return; end

12/26/2000