[genius] Sun Apr 26 11:57:32 2015 Jiri (George) Lebl <jirka 5z com>

[George Lebl <jirka src gnome org>]

commit f85228a0c8caf95c1d3d9a0eb56fec43f410845d
Author: Jiri (George) Lebl <jiri lebl gmail com>
Date:   Sun Apr 26 11:57:37 2015 -0500

    Sun Apr 26 11:57:32 2015  Jiri (George) Lebl <jirka 5z com>
    
        * examples/Makefile.am, examples/explicit-fdm-heat.gel: Add example
          heat equation solved with explicit FDM.

 ChangeLog                      |    5 +++
 examples/Makefile.am           |    1 +
 examples/explicit-fdm-heat.gel |   69 ++++++++++++++++++++++++++++++++++++++++
 3 files changed, 75 insertions(+), 0 deletions(-)
---
diff --git a/ChangeLog b/ChangeLog
index 2819bf6..fb8878b 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,8 @@
+Sun Apr 26 11:57:32 2015  Jiri (George) Lebl <jirka 5z com>
+
+       * examples/Makefile.am, examples/explicit-fdm-heat.gel: Add example
+         heat equation solved with explicit FDM.
+
 Tue Apr 21 00:35:54 2015  Jiri (George) Lebl <jirka 5z com>
 
        * examples/laplace-fdm.gel: update to add more things to play around
diff --git a/examples/Makefile.am b/examples/Makefile.am
index 7fb7f2f..8d6adfa 100644
--- a/examples/Makefile.am
+++ b/examples/Makefile.am
@@ -5,6 +5,7 @@ example_DATA = \
        dalemb-pulse.gel \
        eulers-method-exp-run.gel \
        eulers-method-graphs-exp.gel \
+       explicit-fdm-heat.gel \
        laplace-fdm.gel \
        linapprox.gel \
        lorenz.gel \
diff --git a/examples/explicit-fdm-heat.gel b/examples/explicit-fdm-heat.gel
new file mode 100644
index 0000000..1987738
--- /dev/null
+++ b/examples/explicit-fdm-heat.gel
@@ -0,0 +1,69 @@
+# Category: Differential Equations
+# Name: Heat equation solved using Explicit Finite Difference Method
+#
+# Solve heat equation using explicit FDM.  We try to animate as different
+# times are solved for, though of course we must skip graphing some times when
+# there are too many steps.
+#
+# The equation is u_t = u_{xx} where 0 < x < 1
+# The boundary conditions are insulated on x=0, u_x(0,t) = 0
+# and Dirichlet condition on x=1, u(1,t) = 0
+
+# The number of intervals on the x axis
+n := 20;
+h := float(1/n);
+
+function initialf(x) = sin(2*pi*x);
+
+# The value of k/h^2 must be less than or equal 0.5 for the method to be
+# stable.  try putting 0.55 in the line below instead of 0.5.  Try setting it
+# to less than 0.5 to what it does to accuracy
+k := 0.5*h^2;
+
+# Maximum t to go up to
+maxt := 0.06;
+
+m := round(maxt/k);
+
+# Never plot more than 4*n steps in the t (that is y) direction
+# rather skip rows.  If you do not want to skip rows, comment out
+# the next line.
+plotevery = max(1,floor(m/(2*n)));
+
+#Set up initial value
+u := null;
+data := null;
+for j=1 to n+1 do (
+  x = (j-1)*h;
+  u@(j) := initialf(x);
+  data := [data;[x,0,u@(j)]]
+);
+
+SurfacePlotDrawLegends = true;
+
+# If you change it from 0 to plotevery, then no plots will be generated
+toplot = 0;
+#toplot = plotevery;
+
+for i=1 to m do (
+  v = u;
+  for j=2 to n do (
+    u@(j) := v@(j) + (k/(h^2))*(v@(j+1)+v@(j-1)-2*v@(j))
+  );
+  u@(1) := u@(2);
+  u@(n+1) := 0;
+
+  increment toplot;
+  if toplot == plotevery then (
+    toplot = 0;
+    for j=2 to n do (
+      data = [data;[(j-1)*h,i*k,u@(j)]]
+    );
+    data = [data;[0,i*k,u@(1)]];
+    data = [data;[1,i*k,u@(n+1)]];
+    SurfacePlotData(data,"u",[0,1,0,maxt])
+  )
+);
+
+#print left hand endpoint
+print ("Approximate value of u(0," + maxt + ") ~~ " + u@(1));