> # sfb 31 Oct 2002
> # numerical approximations to u(x,t)
> #what is a matrix?
> A:=Matrix(3,2,[[1, 2],[3, 4],[5, 6]]);
> B:=Matrix(2,3,[[1, 2, 3],[4, 5, 6]]);
> C:=Matrix(3,2,[[1, 1],[2, 2],[3, 3]]);
> A+B;
> A+C,A-C;
> # We need a matrix to hold the numbers, first some sizing
> numX:=20;deltaX:=1/numX;numT:=101;deltaT:=deltaX;
> M:=Matrix(numX+1,numT+1);# x = (i-1)*deltaX, t = (j-1)*deltaT
> for i from 1 to numX+1 do
> for j from 1 to numT+1 do
>  M[i,j]:=evalf(sin((i-1)*Pi/20)*sin((j-1)*Pi/50));
> od;
> od:
> # Filling the matrix with values
> M[3,34];
> with(plots);
> matrixplot(M,axes=boxed,labels=["x","time","u"]);
> L:=convert(M,listlist):
> listplot3d(L,axes=boxed,shading=zhue,title="Graphing a Matrix");
> #
> # The Wave Equation example u_xx = u_tt
> #
> solve((u(x+h,t)-2*u(x,t)+u(x-h,t))/h^2=(u(x,t+k)-2*u(x,t)+u(x,t-k))/k^2,u(x,t+k));
> # initial values
> for i from 1 to numX+1 do
> M[i,1]:=0;
> M[i,2]:=0;
> od:
> #Start the wave.
> M[2,1]:=1;M[3,1]:=1;M[3,2]:=1;M[4,2]:=1;
> # boundry values -- ``walls to bounce and invert''
> for j from 1 to numT+1 do
> M[1,j]:=0; M[numX+1,j]:=0;
> od:
> # internal values
> for j from 3 to numT + 1 do
> for i from 2 to numX do
> M[i,j]:=deltaT^2*(M[i+1,j-1]-2*M[i,j-1]+M[i-1,j-1])/deltaX^2+2*M[i,j-1]-M[i,j-2];
> od;
> od:
> L:=convert(M,listlist):
> listplot3d(L,axes=boxed,shading=zhue,title="Bouncing Wave");
> f:=(i,j)->M[i,j];
> M2:=Matrix(numX+1,numX/2,f);
> L2:=convert(M2,listlist):
> listplot3d(L2,axes=boxed,shading=zhue,title="Close up",labels=["x/delta x+1","time/delta t+1","u(x,t)"]);
> M[2,1]:=1;M[3,1]:=1;M[3,2]:=1;M[4,2]:=1;M[20,1]:=2;M[19,1]:=2;M[19,2]:=2;M[18,2]:=2;
> for j from 3 to numT + 1 do
> for i from 2 to numX do
> M[i,j]:=deltaT^2*(M[i+1,j-1]-2*M[i,j-1]+M[i-1,j-1])/deltaX^2+2*M[i,j-1]-M[i,j-2];
> od;
> od:
> L:=convert(M,listlist):
> listplot3d(L,axes=boxed,shading=zhue,title="Waves meet");
>