/usr/share/genius/examples/laplace-fdm.gel is in genius-common 1.0.23-3.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 | # Category: Differential Equations
# Name: Laplace equation solution using Finite Difference Method
#
# Solve Laplace equation using FDM with Liebmann's iterative method
# to solve the resulting linear equations. The animation is really
# about the iterations of Liebmann's method and how it converges.
#
# The equation is Lu = f(x,y) on the square [0,pi]^2. Normally f(x,y)
# is zero but you can change it.
# The number of intervals on each side. There will be (intervals+1)^2
# total points to solve for and hence that many equations.
intervals = 20;
# zero initial guess
u = zeros(intervals+1,intervals+1);
# how about random initial guess (interval [-1,1])
#u = 2*rand(intervals+1,intervals+1)-ones(intervals+1,intervals+1);
# initial guess of -1s, that's actually somewhat slow to converge if f=0
#u = -ones(intervals+1,intervals+1);
# The side conditions
function bottomside(x) = -sin(x);
function topside(x) = sin(2*x);
function leftside(y) = 0.5*sin(y);
function rightside(y) = 0;
# nonhomogeneity, the f in Lu=f
function f(x,y) = 0;
#function f(x,y) = 2;
#function f(x,y) = -2
#function f(x,y) = 6*exp(-(x-pi/2)^2-(y-pi/2)^2);
#function f(x,y) = 10*sin((x-pi/2)*(y-pi/2));
# set up the side conditions
for n=1 to intervals+1 do
u@(n,1) = bottomside(pi*(n-1)/intervals);
for n=1 to intervals+1 do
u@(n,intervals+1) = topside(pi*(n-1)/intervals);
for n=1 to intervals+1 do
u@(1,n) = leftside(pi*(n-1)/intervals);
for n=1 to intervals+1 do
u@(intervals+1,n) = rightside(pi*(n-1)/intervals);
# don't draw the legend
SurfacePlotDrawLegends = false;
# use standard variable names (in case they got reset)
SurfacePlotVariableNames = ["x","y","z"];
PlotWindowPresent(); # Make sure the window is raised
# plot the data
SurfacePlotDataGrid(u,[0,pi,0,pi]);
# If you want to export the animation to a sequence of .png
# images uncomment here and below
#ExportPlot ("animation" + 0 + ".png");
for n = 1 to 300 do (
# to slow down the animation
wait (0.1);
maxch = 0;
for i=2 to intervals do (
for j=2 to intervals do (
old = u@(i,j);
# This is the equation coming from FDM
u@(i,j) = (1/4)*(u@(i-1,j)+u@(i+1,j)+u@(i,j-1)+u@(i,j+1))
- (1/4)*f(pi*(i-1)/intervals,pi*(j-1)/intervals)*(pi/intervals^2);
# Keep track of maximim change
if |u@(i,j)-old| > maxch then
maxch = |u@(i,j)-old|
)
);
print("Maximum change: " + maxch + " (iteration " + n + ")");
# plot the data
SurfacePlotDataGrid(u,[0,pi,0,pi]);
# Animation saving
#ExportPlot ("animation" + n + ".png");
);
|