Annulus.mws

Vibrating Annulus

This worksheet allow one to plot modes of vibration for a circular annular membrane with fixed homogeneous boundary conditions and to visualize the vibrations in time. The radii are a<b and the wave speed is given by c.
 

>    restart: a:=2: b:=4: c:=1:   with(linalg): with(plots):

Warning, the protected names norm and trace have been redefined and unprotected

Warning, the name changecoords has been redefined

Set up boundary conditions and plot to see zeros.

>    Eq:=det([[BesselJ(m,a*lambda),BesselY(m,a*lambda)],[BesselJ(m,b*lambda),BesselY(m,b*lambda)]]);

Eq := BesselJ(m,2*lambda)*BesselY(m,4*lambda)-BesselY(m,2*lambda)*BesselJ(m,4*lambda)

>    plot({seq(Eq,m=0..3)},lambda=0.1..5);

[Maple Plot]

The frequency of oscillation of the mn-th mode  is given by

>    for n from 1 to 3 do nu[n]:=[seq(c*fsolve(Eq,lambda=2*(n-1)..2*n),m=0..2)]; od;

nu[1] := [1.561515460, 1.598289190, 1.703460713]

nu[2] := [3.136717857, 3.156174755, 3.213882961]

nu[3] := [4.709103771, 4.722232463, 4.761426135]

>    for m from 0 to 2 do
  for n from 1 to 3 do
    C[m+1,n]:=BesselY(m,nu[n][m+1]/c*b)/BesselJ(m,nu[n][m+1]/c*b):
  od:
od:

The product solutions take the following form for a zero initial velocity, though one can easily add a phase shift in time to simulate other motions.

>    u:=(m,n,r,theta,t)->cos(c*nu[n][m+1]*t)*cos(m*theta)*(C[m+1,n]*BesselJ(m,nu[n][m+1]/c*r)-BesselY(m,nu[n][m+1]/c*r));

u := proc (m, n, r, theta, t) options operator, arrow; cos(c*nu[n][m+1]*t)*cos(m*theta)*(C[m+1,n]*BesselJ(m,nu[n][m+1]/c*r)-BesselY(m,nu[n][m+1]/c*r)) end proc

The initial profile of the membrane:

>    m0:=0: n0:=1: plot3d([r*cos(theta),r*sin(theta),u(m0,n0,r,theta,0)],r=a..b,theta=0..2*Pi,grid=[40,40]);

[Maple Plot]

The evolution of this profile

>    animate3d([r*cos(theta),r*sin(theta),u(m0,n0,r,theta,t)],r=a..b,theta=0..2*Pi,t=0..6/nu[n0][m0+1],frames=20);

[Maple Plot]

The following code shows how to save animations as animated gifs. The first line sets up the plot to be gif and the third line reverts the plotting back oto the screen.

>    plotsetup(gif,plotoutput=`amesh21.gif`,plotoptions=`height=200,width=200`):

>    m0:=2: n0:=1: animate3d([r*cos(theta),r*sin(theta),u(m0,n0,r,theta,t)],r=a..b,theta=0..2*Pi,t=0..6/nu[n0][m0+1],frames=20);

>    plotsetup(inline);

>