[genius] Tue Apr 07 14:28:45 2015 Jiri (George) Lebl <jirka 5z com>

[George Lebl <jirka src gnome org>]

commit 6b269d6c5321cefc47320c52839dfa5522e302cd
Author: Jiri (George) Lebl <jiri lebl gmail com>
Date:   Tue Apr 7 14:28:49 2015 -0500

    Tue Apr 07 14:28:45 2015  Jiri (George) Lebl <jirka 5z com>
    
        * examples/vibrating-drumhead-modes.gel:  Add example showing
          various modes of a vibrating drumhead (wave equation)

 ChangeLog                             |    5 ++
 examples/Makefile.am                  |    3 +-
 examples/vibrating-drumhead-modes.gel |  114 +++++++++++++++++++++++++++++++++
 3 files changed, 121 insertions(+), 1 deletions(-)
---
diff --git a/ChangeLog b/ChangeLog
index 7a1a19c..cdad414 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,8 @@
+Tue Apr 07 14:28:45 2015  Jiri (George) Lebl <jirka 5z com>
+
+       * examples/vibrating-drumhead-modes.gel:  Add example showing
+         various modes of a vibrating drumhead (wave equation)
+
 Mon Mar 09 17:25:09 2015  Jiri (George) Lebl <jirka 5z com>
 
        * Release 1.0.20
diff --git a/examples/Makefile.am b/examples/Makefile.am
index b93e9e0..7fb7f2f 100644
--- a/examples/Makefile.am
+++ b/examples/Makefile.am
@@ -20,7 +20,8 @@ example_DATA = \
        fourier-series-plotting.gel \
        rsa.gel \
        easy-matrix-examples.gel \
-       peano.gel
+       peano.gel \
+       vibrating-drumhead-modes.gel
 
 EXTRA_DIST = \
        $(example_DATA)
diff --git a/examples/vibrating-drumhead-modes.gel b/examples/vibrating-drumhead-modes.gel
new file mode 100644
index 0000000..6cff7aa
--- /dev/null
+++ b/examples/vibrating-drumhead-modes.gel
@@ -0,0 +1,114 @@
+# Category: Differential Equations
+# Name: Vibrating circular drumhead (wave equation)
+#
+# Make an animation of the various modes of a vibrating circular drumhead.
+# The equation is u_{tt} = \nabla^2 u in polar coordinates and the
+# drum is of radius 1.  The varius modes n,m are of the form
+# BesselJn(n,k*r)*cos(n*theta)*sin(k*t)
+# where k is the mth zero of the BesselJn(n,x)
+# The superposition is a superposition of 3 different modes,
+# with some shifts also applied in the theta direction for
+# added complexity
+
+the_answer = AskButtons("Which mode?  n-m where n is the Bessel J_n to use",
+                                         "0-1", # 1
+                                         "0-2", # 2
+                                         "0-3", # 3
+                                         "1-1", # 4
+                                         "1-2", # 5
+                                         "1-3", # 6
+                                         "2-1", # 7
+                                         "2-2", # 8
+                                         "2-3", # 9
+                                         "Superposition" # 9
+                                         );
+
+# the zeros of the Bessel functions
+kmn = [2.4048  , 3.8317, 5.1356
+       5.5201, 7.0156, 8.4172
+       8.6537, 10.1735, 11.6198];
+
+if the_answer == 1 then (
+  n = 0;
+  k = kmn@(1,n+1);
+) else if the_answer == 2 then (
+  n = 0;
+  k = kmn@(2,n+1);
+) else if the_answer == 3 then (
+  n = 0;
+  k = kmn@(3,n+1);
+) else if the_answer == 4 then (
+  n = 1;
+  k = kmn@(1,n+1);
+) else if the_answer == 5 then (
+  n = 1;
+  k = kmn@(2,n+1);
+) else if the_answer == 6 then (
+  n = 1;
+  k = kmn@(3,n+1);
+) else if the_answer == 7 then (
+  n = 2;
+  k = kmn@(1,n+1);
+) else if the_answer == 8 then (
+  n = 2;
+  k = kmn@(2,n+1);
+) else if the_answer == 9 then (
+  n = 2;
+  k = kmn@(3,n+1);
+) else (
+  # Superposition
+  n = -1; # signals superposition
+
+  # Two put together
+  coeff1 = 0.5;
+  shift1 = 0.0;
+  n1 = 1;
+  k1 = kmn@(2,n1+1);
+
+  coeff2 = 0.3;
+  shift2 = 0.5;
+  n2 = 0;
+  k2 = kmn@(3,n2+1);
+
+  coeff3 = 0.3;
+  shift3 = 2;
+  n3 = 2;
+  k3 = kmn@(1,n3+1);
+);
+
+SurfacePlotDrawLegends = false; # don't draw the legend
+PlotWindowPresent(); # Make sure the window is raised
+
+#One mode
+if n >= 0 then (
+  for t=0.0 to 10.0 by 0.01 do (
+    data = null;
+    for r=0 to 1.0 by 1/10.0 do (
+      for theta=0 to 2*pi by pi/15 do (
+        x = r*cos(theta);
+        y = r*sin(theta);
+        data = [data;[x,y,BesselJn(n,k*r)*cos(n*theta)*sin(k*t)]]
+      )
+    );
+
+    # Plot the data
+    SurfacePlotData(data,[-1,1,-1,1,-1,1])
+  )
+) else (
+  for t=0.0 to 10.0 by 0.01 do (
+    data = null;
+    for r=0 to 1.0 by 1/10.0 do (
+      for theta=0 to 2*pi by pi/15 do (
+        x = r*cos(theta);
+        y = r*sin(theta);
+        val = coeff1*BesselJn(n1,k1*r)*cos(n1*theta+shift1)*sin(k1*t) +
+              coeff2*BesselJn(n2,k2*r)*cos(n2*theta+shift2)*sin(k2*t) +
+              coeff3*BesselJn(n3,k3*r)*cos(n3*theta+shift3)*sin(k3*t);
+        data = [data;[x,y,val]]
+      )
+    );
+
+    # Plot the data
+    SurfacePlotData(data,[-1,1,-1,1,-1,1])
+  )
+);