8 years, 28 days

## @helix Don't worry about it, thanks...

@helix

Don't worry about it, thanks for your help. I will keep looking for a work around and will post it here if I find one.

## @helix Hmm is (* *) not a multiple ...

@helix
Hmm is (* *) not a multiple line comment in your version of maple?

## @helix Hmm is (* *) not a multiple ...

@helix
Hmm is (* *) not a multiple line comment in your version of maple?

## @helix Here is the code. restart;Co...

@helix

Here is the code.

restart;
Code Packages
with(plots):
with(CodeGeneration):
with(LinearAlgebra):
Important Numbers and Functions
1 Square Numbers and Functions
phi:=proc(n,m,l1,l2)

RETURN((2/sqrt(l1*l2))*sin(((n*Pi)/l1)*x)*sin(((m*Pi)/l2)*y));

end proc:

lamda:=proc(n,m,l1,l2)

RETURN(((n*Pi)/l1)^2 + ((m*Pi)/l2)^2);

end proc:

omega:=proc(n,m,l1,l2)

RETURN(sqrt(lamda(n,m,l1,l2)));

end proc:

wave:=proc(n,m,l1,l2)

RETURN(phi(n,m,l1,l2) * cos(omega(n,m,l1,l2) * t));

end proc:

wave2:=proc(n,m,l1,l2)

RETURN(phi(n,m,l1,l2) * sin(omega(n,m,l1,l2) * t));

end proc:

2 Squares Numbers and Functions
HeavyStepLeft:= piecewise(x<0,1);
HeavyStepRight:= piecewise(x>=0,1);
aPrime:= proc(n1,n2,a,m,l2)

local AP;

AP:= sqrt(((n2^2/n1^2)*a^2) + ((m^2/l2^2)*Pi^2)*((n2^2/n1^2) -1 ));

RETURN(AP);

end proc:

interval_array_ten:= Array([[1.25,1.75,2,2.25,2.75,3,3.5,3.75,4,4.5],[1.25,1.75,2,2.5,2.75,3.25,3.5,4,4.25,4.75],[1.25,2,2.5,2.75,3.25,3.75,4,4.25,4.75,5],[1.75,2.25,2.75,3.25,3.5,4,4.25,4.75,5,5.5],[1.25,2,2.5,3,3.5,4,4.25,4.75,5,5.5],[1.25,2.25,2.75,3.25,3.75,4.25,4.5,5,5.5,6],[1.75,2.5,3,3.5,4,4.5,5,5.25,5.75,6.25],[1.25,2,2.5,3.25,3.75,4.25,4.75,5.25,5.75,6.25],[1.25,2.25,2.75,3.5,4,4.5,5,5.5,6,6.5],[1.5,2.5,3,3.75,4.25,4.75,5.25,5.75,6.25,6.75]]);

A_norm:= proc(n1,n2,l1,l2,a,m)

local A;

if(sin(aPrime(n1,n2,a,m,l2)*l1) = 0)
then

A:= aPrime(n1,n2,a,m,l2)*cos(aPrime(n1,n2,a,m,l2)* l1) * (1/(sqrt(((aPrime(n1,n2,a,m,l2))^2 * cos(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + a^2 *(cos(a*l1))^2 *(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

else

A:= sin(aPrime(n1,n2,a,m,l2)* l1) * (1/(sqrt( (sin(aPrime(n1,n2,a,m,l2)*l1))^2* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + (sin(a*l1))^2*(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

fi;

RETURN(evalf(A));

end proc:

APrime_norm:= proc(n1,n2,l1,l2,a,m)

local A;

if(sin(a*l1) = 0)
then

A:= -a*cos(a*l1) * (1/(sqrt(((aPrime(n1,n2,a,m,l2))^2 * cos(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + a^2 *(cos(a*l1)^2)*(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

else

A:= sin(a*l1) * (1/(sqrt((sin(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + (sin(a*l1))^2 *(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

fi;

RETURN(evalf(A));

end proc:

lambda_Two_squares:= proc(n_value,a_value,m,l2)

RETURN((1/n_value^2)*(a_value^2 + (m^2/l2^2)*Pi^2));

end proc:
omega_Two_squares:= proc(n_value,a_value,m,l2)

RETURN(sqrt(lambda_Two_squares(n_value,a_value,m,l2)));

end proc:
Phi_left:= proc(n1,n2,l1,l2,a,m)

local A;

A:= A_norm(n1,n2,l1,l2,a,m);

RETURN(A*sin(a*(x+l1))*sin(((m*Pi)/l2)*y));

end proc:

Phi_right:= proc(n1,n2,l1,l2,a,m)

local A;
local ap;

A:= APrime_norm(n1,n2,l1,l2,a,m);
ap:= aPrime(n1,n2,a,m,l2);

RETURN((A*sin(ap*(l1-x))*sin(((m*Pi)/l2)*y)));

end proc:

PhiFunc:= proc(n1,n2,l1,l2,a,m)

RETURN(Phi_left(n1,n2,l1,l2,a,m) *HeavyStepLeft + Phi_right(n1,n2,l1,l2,a,m) *HeavyStepRight);

end proc:
n_of_x:= piecewise(x<0,n1,x>0,n2);
wave_Two_squares:= proc(n1,n2,l1,l2,a,m)

RETURN(PhiFunc(n1,n2,l1,l2,a,m)* cos(omega_Two_squares(n1,a,m,l2) * t));

end proc:
wave_Two_squares2:= proc(n1,n2,l1,l2,a,m)

RETURN(PhiFunc(n1,n2,l1,l2,a,m)* sin(omega_Two_squares(n1,a,m,l2) * t));

end proc:

Both
FMaker := proc(x,y,Vx,Vy,Sigma)

local f;

f:= exp(-Sigma*((x-(Pi/2))^2 + ((y-(Pi/2))^2)))*cos(Vx*x + Vy*y);

RETURN (f);

end proc:
GMaker := proc(f, Vx, Vy, Sigma)

local g;
local DX;
local DY;

local ExPart := (simplify(f/exp(-Sigma*(y-(Pi/2))^2)))/cos(Vx*x + Vy*y);
local EyPart := (simplify(f/exp(-Sigma*(x-(Pi/2))^2)))/cos(Vx*x + Vy*y);

DX := ((EyPart*cos(Vy*y))*diff(ExPart*cos(Vx*x),x)-((EyPart*sin(Vy*y))*diff(ExPart*sin(Vx*x),x)));
DY := ((ExPart*cos(Vx*x))*diff(EyPart*cos(Vy*y),y)-((ExPart*sin(Vx*x))*diff(EyPart*sin(Vy*y),y)));

g:= Vx*DX + Vy*DY;

RETURN (g);

end proc:

1 Square Functions
CoeffFinder:=proc(f,Sigma,sum,l1,l2,Vx,Vy)

local n;
local m;
local F;
local SinPartx;
local SinParty;

local C := 2/sqrt(l1*l2);
local ExPart := (simplify(f/exp(-Sigma*(y-(Pi/2))^2)))/cos(Vx*x + Vy*y);
local EyPart := (simplify(f/exp(-Sigma*(x-(Pi/2))^2)))/cos(Vx*x + Vy*y);

F:=Array(1..sum,1..sum);

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

SinPartx:= sin(((n*Pi)/l1)*x);
SinParty:= sin(((m*Pi)/l2)*y);

F[n,m]:= evalf(C*((evalf(int(ExPart*SinPartx*cos(Vx*x),x=0..l1,numeric=true))*evalf((int(EyPart*SinParty*cos(Vy*y),y=0..l2,numeric=true)))-(evalf((int(ExPart*SinPartx*sin(Vx*x),x=0..l1,numeric=true)))*evalf((int(EyPart*SinParty*sin(Vy*y),y=0..l2,numeric=true)))))));

end do;

end do;

RETURN(F);

end proc:

CoeffFinder2:=proc(g,sum,l1,l2,Vx,Vy)

local n;
local m;
local G;

G:=Array(1..sum,1..sum);

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

G[n,m]:= (evalf(int(int(g*phi(n,m,l1,l2),x=0..l1,numeric=true),y=0..l2,numeric=true))/omega(n,m,l1,l2));

end do;

end do;

RETURN(G);

end proc:

WaveMaker:=proc(F,G,sum,l1,l2)

local n;
local m;
local totalWave;
totalWave := 0;

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

totalWave := totalWave + (F[n,m] * wave(n,m,l1,l2)) + (G[n,m] * wave2(n,m,l1,l2));

end do;

end do;
RETURN(totalWave);

end proc:
AnimationMaker:= proc(Sigma,Vx,Vy,Nodes,l1,l2,Dimensional)

local f;
local g;
local F;
local G;

f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);

F:=CoeffFinder(f,Sigma,Nodes,l1,l2,Vx,Vy);
G:=CoeffFinder2(g,Nodes,l1,l2,Vx,Vy);

if(Dimensional)
then
else
(RETURN(animate(densityplot,[WaveMaker(F,G,Nodes,l1,l2),x=0..l1,y=0..l2,contrast=0.5],t=0..10,frames=150,scaling=constrained)));
end if;

end proc:

1 Square Output
#AnimationMaker(25,1,2,10,Pi,Pi,false);
#Check initial function against f.
(*
Sigma:=10;
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
Nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=CoeffFinder(f,Sigma,Nodes,l1,l2,Vx,Vy);
G:=CoeffFinder2(g,Nodes,l1,l2,Vx,Vy);
Waveee:=WaveMaker(F,G,Nodes,l1,l2):
plot3d(f,x=0..l1,y=0..l2);
t:=0;
plot3d(Waveee,x=0..l1,y=0..l2);
plot3d((Waveee-f),x=0..l1,y=0..l2,scaling=constrained);
*)
2 Squares with Boundary Functions
a_root_finder:= proc(n1,n2,m,l1,l2,interval1,interval2)

local M;
local a;

M:= Matrix([[sin(a*l1),-sin(aPrime(n1,n2,a,m,l2)*l1)],[a*cos(a*l1),aPrime(n1,n2,a,m,l2)*cos(aPrime(n1,n2,a,m,l2)*l1)]]);

RETURN(fsolve(Determinant(M) = 0, a=interval1..interval2));

end proc:

a_array_assembler:= proc(n1,n2,l1,l2,no_of_roots,no_of_ms,interval_array)

local i;
local m_local;
local start_point;
local a_Array;
local interval1;
local interval2;

start_point:= 0.1;
a_Array:= Array(1..no_of_roots,1..no_of_roots);

for m_local from 1 by 1 to no_of_ms
do

start_point := 0.1;

for i from 1 by 1 to no_of_roots
do

interval1:= start_point;
interval2:= interval_array[m_local,i];

a_Array[m_local,i] := a_root_finder(n1,n2,m_local,l1,l2,interval1,interval2);
start_point := interval2;

od;

od;

RETURN(a_Array);

end proc:
#This generates the values of 'a' which will be used in the left and right Eigenfunctions. Each row corresponds to m=1..10.
a_array:= a_array_assembler(1,2,Pi,Pi,10,10,interval_array_ten);

Two_squares_CoeffFinder:= proc(f,nodes,n1,n2,l1,l2, a_Array)

local F;
local count_a;
local m;
local a;

F:= Array(1..nodes,1..nodes);

for m from 1 by 1 to nodes
do

for count_a from 1 by 1 to nodes
do

a:= a_array[m,count_a];

#F[m,count_a]:= ((n1^2)*(int(int((Phi_left(n1,n2,l1,l2,a,m))*f,x=-l1..0,numeric=true),y=0..l2,numeric=true))) + ((n2^2)*(int(int((Phi_right(n1,n2,l1,l2,a,m))*f,x=0..l1,numeric=true),y=0..l2,numeric=true)));

F[m,count_a]:= int(int((n_of_x)^2 * PhiFunc(n1,n2,l1,l2,a,m) * f,x=-l1..l1,numeric=true),y=0..l1,numeric=true);

od;

od;

RETURN(F);

end proc:
Two_squares_CoeffFinder2:= proc(g,nodes,n1,n2,l1,l2,a_Array)

local G;
local count_a;
local m;
local a;

G:= Array(1..nodes,1..nodes);

for m from 1 by 1 to nodes
do

for count_a from 1 by 1 to nodes
do

a:= a_array[m,count_a];

#G[m,count_a]:= (1/omega_Two_squares(n1,a,m,l2))*((n1^2)*(int(int((Phi_left(n1,n2,l1,l2,a,m))*g,x=-l1..0,numeric=true),y=0..l2,numeric=true))) + ((n2^2)*(int(int((Phi_right(n1,n2,l1,l2,a,m))*g,x=0..l1,numeric=true),y=0..l2,numeric=true)));

G[m,count_a]:= (1/omega_Two_squares(n1,a,m,l2)) * int(int((n_of_x)^2 * PhiFunc(n1,n2,l1,l2,a,m) * g,x=-l1..l1,numeric=true),y=0..l1,numeric=true);
od;

od;

RETURN(G);

end proc:

WaveMaker_Two_squares:=proc(F,G,nodes,l1,l2,a_array,n1,n2)

local a;
local a_count;
local m;
local totalWave;
totalWave := 0;

for m from 1 by 1 to nodes
do

for a_count from 1 by 1 to nodes
do

a:= a_array[m,a_count];
totalWave := totalWave + (F[m,a_count] * wave_Two_squares(n1,n2,l1,l2,a,m)) + (G[m,a_count] * wave_Two_squares2(n1,n2,l1,l2,a,m));

end do;

end do;
RETURN(totalWave);

end proc:
AnimationMaker_Two_squares:= proc(Sigma,Vx,Vy,nodes,l1,l2,Dimensional,n1,n2,a_array,time)

local f;
local g;
local F;
local G;
local wavee;
f:= FMaker(x,y,Vx,Vy,Sigma); g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array); wavee:= WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2);
if(Dimensional)
then
else
(RETURN(animate(densityplot,[wavee,x=-l1..l1,y=0..l2],t=0..time,frames=150,scaling=constrained))); end if;

end proc:
2 Squares with Boundary Output

#AnimationMaker_Two_squares(10,1,0,10,Pi,Pi,false,1,2,a_array,10);

#Check initial function against f

(*
Sigma:=10;
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array);
wavee:=WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2):
plot3d(f,x=-Pi..Pi,y=0..Pi,scaling=constrained);
t:=0;
plot3d(wavee,x=-l1..l1,y=0..l1,scaling=constrained);
*)

#Check phi normalisation
(*
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
m:=2;
a:=1;
n_of_x:= piecewise(x<0,n1,x>0,n2);
PhiFunc:= Phi_left(n1,n2,l1,l2,a,m) *HeavyStepLeft + Phi_right(n1,n2,l1,l2,a,m) *HeavyStepRight;
int((n_of_x)^2 *(PhiFunc)^2,x=-l1..l1,y=0..l2);
*)

#Check WaveMaker

Vx:=1;
Vy:=0;
Sigma:=10;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array);
wavee:= WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2);
animate(densityplot,[wavee,x=-l1..l1,y=0..l2],t=0..10,frames=150,scaling=constrained);

## @helix Here is the code. restart;Co...

@helix

Here is the code.

restart;
Code Packages
with(plots):
with(CodeGeneration):
with(LinearAlgebra):
Important Numbers and Functions
1 Square Numbers and Functions
phi:=proc(n,m,l1,l2)

RETURN((2/sqrt(l1*l2))*sin(((n*Pi)/l1)*x)*sin(((m*Pi)/l2)*y));

end proc:

lamda:=proc(n,m,l1,l2)

RETURN(((n*Pi)/l1)^2 + ((m*Pi)/l2)^2);

end proc:

omega:=proc(n,m,l1,l2)

RETURN(sqrt(lamda(n,m,l1,l2)));

end proc:

wave:=proc(n,m,l1,l2)

RETURN(phi(n,m,l1,l2) * cos(omega(n,m,l1,l2) * t));

end proc:

wave2:=proc(n,m,l1,l2)

RETURN(phi(n,m,l1,l2) * sin(omega(n,m,l1,l2) * t));

end proc:

2 Squares Numbers and Functions
HeavyStepLeft:= piecewise(x<0,1);
HeavyStepRight:= piecewise(x>=0,1);
aPrime:= proc(n1,n2,a,m,l2)

local AP;

AP:= sqrt(((n2^2/n1^2)*a^2) + ((m^2/l2^2)*Pi^2)*((n2^2/n1^2) -1 ));

RETURN(AP);

end proc:

interval_array_ten:= Array([[1.25,1.75,2,2.25,2.75,3,3.5,3.75,4,4.5],[1.25,1.75,2,2.5,2.75,3.25,3.5,4,4.25,4.75],[1.25,2,2.5,2.75,3.25,3.75,4,4.25,4.75,5],[1.75,2.25,2.75,3.25,3.5,4,4.25,4.75,5,5.5],[1.25,2,2.5,3,3.5,4,4.25,4.75,5,5.5],[1.25,2.25,2.75,3.25,3.75,4.25,4.5,5,5.5,6],[1.75,2.5,3,3.5,4,4.5,5,5.25,5.75,6.25],[1.25,2,2.5,3.25,3.75,4.25,4.75,5.25,5.75,6.25],[1.25,2.25,2.75,3.5,4,4.5,5,5.5,6,6.5],[1.5,2.5,3,3.75,4.25,4.75,5.25,5.75,6.25,6.75]]);

A_norm:= proc(n1,n2,l1,l2,a,m)

local A;

if(sin(aPrime(n1,n2,a,m,l2)*l1) = 0)
then

A:= aPrime(n1,n2,a,m,l2)*cos(aPrime(n1,n2,a,m,l2)* l1) * (1/(sqrt(((aPrime(n1,n2,a,m,l2))^2 * cos(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + a^2 *(cos(a*l1))^2 *(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

else

A:= sin(aPrime(n1,n2,a,m,l2)* l1) * (1/(sqrt( (sin(aPrime(n1,n2,a,m,l2)*l1))^2* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + (sin(a*l1))^2*(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

fi;

RETURN(evalf(A));

end proc:

APrime_norm:= proc(n1,n2,l1,l2,a,m)

local A;

if(sin(a*l1) = 0)
then

A:= -a*cos(a*l1) * (1/(sqrt(((aPrime(n1,n2,a,m,l2))^2 * cos(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + a^2 *(cos(a*l1)^2)*(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

else

A:= sin(a*l1) * (1/(sqrt((sin(aPrime(n1,n2,a,m,l2) * l1)^2)* (l2/2) * (n1^2) * ((-sin(2*a*l1) + 2*a*l1)/(4*a)) + (sin(a*l1))^2 *(l2/2) * (n2^2) * ((-sin(2*aPrime(n1,n2,a,m,l2)*l1) + 2*aPrime(n1,n2,a,m,l2)*l1)/(4*aPrime(n1,n2,a,m,l2))))));

fi;

RETURN(evalf(A));

end proc:

lambda_Two_squares:= proc(n_value,a_value,m,l2)

RETURN((1/n_value^2)*(a_value^2 + (m^2/l2^2)*Pi^2));

end proc:
omega_Two_squares:= proc(n_value,a_value,m,l2)

RETURN(sqrt(lambda_Two_squares(n_value,a_value,m,l2)));

end proc:
Phi_left:= proc(n1,n2,l1,l2,a,m)

local A;

A:= A_norm(n1,n2,l1,l2,a,m);

RETURN(A*sin(a*(x+l1))*sin(((m*Pi)/l2)*y));

end proc:

Phi_right:= proc(n1,n2,l1,l2,a,m)

local A;
local ap;

A:= APrime_norm(n1,n2,l1,l2,a,m);
ap:= aPrime(n1,n2,a,m,l2);

RETURN((A*sin(ap*(l1-x))*sin(((m*Pi)/l2)*y)));

end proc:

PhiFunc:= proc(n1,n2,l1,l2,a,m)

RETURN(Phi_left(n1,n2,l1,l2,a,m) *HeavyStepLeft + Phi_right(n1,n2,l1,l2,a,m) *HeavyStepRight);

end proc:
n_of_x:= piecewise(x<0,n1,x>0,n2);
wave_Two_squares:= proc(n1,n2,l1,l2,a,m)

RETURN(PhiFunc(n1,n2,l1,l2,a,m)* cos(omega_Two_squares(n1,a,m,l2) * t));

end proc:
wave_Two_squares2:= proc(n1,n2,l1,l2,a,m)

RETURN(PhiFunc(n1,n2,l1,l2,a,m)* sin(omega_Two_squares(n1,a,m,l2) * t));

end proc:

Both
FMaker := proc(x,y,Vx,Vy,Sigma)

local f;

f:= exp(-Sigma*((x-(Pi/2))^2 + ((y-(Pi/2))^2)))*cos(Vx*x + Vy*y);

RETURN (f);

end proc:
GMaker := proc(f, Vx, Vy, Sigma)

local g;
local DX;
local DY;

local ExPart := (simplify(f/exp(-Sigma*(y-(Pi/2))^2)))/cos(Vx*x + Vy*y);
local EyPart := (simplify(f/exp(-Sigma*(x-(Pi/2))^2)))/cos(Vx*x + Vy*y);

DX := ((EyPart*cos(Vy*y))*diff(ExPart*cos(Vx*x),x)-((EyPart*sin(Vy*y))*diff(ExPart*sin(Vx*x),x)));
DY := ((ExPart*cos(Vx*x))*diff(EyPart*cos(Vy*y),y)-((ExPart*sin(Vx*x))*diff(EyPart*sin(Vy*y),y)));

g:= Vx*DX + Vy*DY;

RETURN (g);

end proc:

1 Square Functions
CoeffFinder:=proc(f,Sigma,sum,l1,l2,Vx,Vy)

local n;
local m;
local F;
local SinPartx;
local SinParty;

local C := 2/sqrt(l1*l2);
local ExPart := (simplify(f/exp(-Sigma*(y-(Pi/2))^2)))/cos(Vx*x + Vy*y);
local EyPart := (simplify(f/exp(-Sigma*(x-(Pi/2))^2)))/cos(Vx*x + Vy*y);

F:=Array(1..sum,1..sum);

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

SinPartx:= sin(((n*Pi)/l1)*x);
SinParty:= sin(((m*Pi)/l2)*y);

F[n,m]:= evalf(C*((evalf(int(ExPart*SinPartx*cos(Vx*x),x=0..l1,numeric=true))*evalf((int(EyPart*SinParty*cos(Vy*y),y=0..l2,numeric=true)))-(evalf((int(ExPart*SinPartx*sin(Vx*x),x=0..l1,numeric=true)))*evalf((int(EyPart*SinParty*sin(Vy*y),y=0..l2,numeric=true)))))));

end do;

end do;

RETURN(F);

end proc:

CoeffFinder2:=proc(g,sum,l1,l2,Vx,Vy)

local n;
local m;
local G;

G:=Array(1..sum,1..sum);

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

G[n,m]:= (evalf(int(int(g*phi(n,m,l1,l2),x=0..l1,numeric=true),y=0..l2,numeric=true))/omega(n,m,l1,l2));

end do;

end do;

RETURN(G);

end proc:

WaveMaker:=proc(F,G,sum,l1,l2)

local n;
local m;
local totalWave;
totalWave := 0;

for n from 1 by 1 to sum
do

for m from 1 by 1 to sum
do

totalWave := totalWave + (F[n,m] * wave(n,m,l1,l2)) + (G[n,m] * wave2(n,m,l1,l2));

end do;

end do;
RETURN(totalWave);

end proc:
AnimationMaker:= proc(Sigma,Vx,Vy,Nodes,l1,l2,Dimensional)

local f;
local g;
local F;
local G;

f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);

F:=CoeffFinder(f,Sigma,Nodes,l1,l2,Vx,Vy);
G:=CoeffFinder2(g,Nodes,l1,l2,Vx,Vy);

if(Dimensional)
then
else
(RETURN(animate(densityplot,[WaveMaker(F,G,Nodes,l1,l2),x=0..l1,y=0..l2,contrast=0.5],t=0..10,frames=150,scaling=constrained)));
end if;

end proc:

1 Square Output
#AnimationMaker(25,1,2,10,Pi,Pi,false);
#Check initial function against f.
(*
Sigma:=10;
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
Nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=CoeffFinder(f,Sigma,Nodes,l1,l2,Vx,Vy);
G:=CoeffFinder2(g,Nodes,l1,l2,Vx,Vy);
Waveee:=WaveMaker(F,G,Nodes,l1,l2):
plot3d(f,x=0..l1,y=0..l2);
t:=0;
plot3d(Waveee,x=0..l1,y=0..l2);
plot3d((Waveee-f),x=0..l1,y=0..l2,scaling=constrained);
*)
2 Squares with Boundary Functions
a_root_finder:= proc(n1,n2,m,l1,l2,interval1,interval2)

local M;
local a;

M:= Matrix([[sin(a*l1),-sin(aPrime(n1,n2,a,m,l2)*l1)],[a*cos(a*l1),aPrime(n1,n2,a,m,l2)*cos(aPrime(n1,n2,a,m,l2)*l1)]]);

RETURN(fsolve(Determinant(M) = 0, a=interval1..interval2));

end proc:

a_array_assembler:= proc(n1,n2,l1,l2,no_of_roots,no_of_ms,interval_array)

local i;
local m_local;
local start_point;
local a_Array;
local interval1;
local interval2;

start_point:= 0.1;
a_Array:= Array(1..no_of_roots,1..no_of_roots);

for m_local from 1 by 1 to no_of_ms
do

start_point := 0.1;

for i from 1 by 1 to no_of_roots
do

interval1:= start_point;
interval2:= interval_array[m_local,i];

a_Array[m_local,i] := a_root_finder(n1,n2,m_local,l1,l2,interval1,interval2);
start_point := interval2;

od;

od;

RETURN(a_Array);

end proc:
#This generates the values of 'a' which will be used in the left and right Eigenfunctions. Each row corresponds to m=1..10.
a_array:= a_array_assembler(1,2,Pi,Pi,10,10,interval_array_ten);

Two_squares_CoeffFinder:= proc(f,nodes,n1,n2,l1,l2, a_Array)

local F;
local count_a;
local m;
local a;

F:= Array(1..nodes,1..nodes);

for m from 1 by 1 to nodes
do

for count_a from 1 by 1 to nodes
do

a:= a_array[m,count_a];

#F[m,count_a]:= ((n1^2)*(int(int((Phi_left(n1,n2,l1,l2,a,m))*f,x=-l1..0,numeric=true),y=0..l2,numeric=true))) + ((n2^2)*(int(int((Phi_right(n1,n2,l1,l2,a,m))*f,x=0..l1,numeric=true),y=0..l2,numeric=true)));

F[m,count_a]:= int(int((n_of_x)^2 * PhiFunc(n1,n2,l1,l2,a,m) * f,x=-l1..l1,numeric=true),y=0..l1,numeric=true);

od;

od;

RETURN(F);

end proc:
Two_squares_CoeffFinder2:= proc(g,nodes,n1,n2,l1,l2,a_Array)

local G;
local count_a;
local m;
local a;

G:= Array(1..nodes,1..nodes);

for m from 1 by 1 to nodes
do

for count_a from 1 by 1 to nodes
do

a:= a_array[m,count_a];

#G[m,count_a]:= (1/omega_Two_squares(n1,a,m,l2))*((n1^2)*(int(int((Phi_left(n1,n2,l1,l2,a,m))*g,x=-l1..0,numeric=true),y=0..l2,numeric=true))) + ((n2^2)*(int(int((Phi_right(n1,n2,l1,l2,a,m))*g,x=0..l1,numeric=true),y=0..l2,numeric=true)));

G[m,count_a]:= (1/omega_Two_squares(n1,a,m,l2)) * int(int((n_of_x)^2 * PhiFunc(n1,n2,l1,l2,a,m) * g,x=-l1..l1,numeric=true),y=0..l1,numeric=true);
od;

od;

RETURN(G);

end proc:

WaveMaker_Two_squares:=proc(F,G,nodes,l1,l2,a_array,n1,n2)

local a;
local a_count;
local m;
local totalWave;
totalWave := 0;

for m from 1 by 1 to nodes
do

for a_count from 1 by 1 to nodes
do

a:= a_array[m,a_count];
totalWave := totalWave + (F[m,a_count] * wave_Two_squares(n1,n2,l1,l2,a,m)) + (G[m,a_count] * wave_Two_squares2(n1,n2,l1,l2,a,m));

end do;

end do;
RETURN(totalWave);

end proc:
AnimationMaker_Two_squares:= proc(Sigma,Vx,Vy,nodes,l1,l2,Dimensional,n1,n2,a_array,time)

local f;
local g;
local F;
local G;
local wavee;
f:= FMaker(x,y,Vx,Vy,Sigma); g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array); wavee:= WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2);
if(Dimensional)
then
else
(RETURN(animate(densityplot,[wavee,x=-l1..l1,y=0..l2],t=0..time,frames=150,scaling=constrained))); end if;

end proc:
2 Squares with Boundary Output

#AnimationMaker_Two_squares(10,1,0,10,Pi,Pi,false,1,2,a_array,10);

#Check initial function against f

(*
Sigma:=10;
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array);
wavee:=WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2):
plot3d(f,x=-Pi..Pi,y=0..Pi,scaling=constrained);
t:=0;
plot3d(wavee,x=-l1..l1,y=0..l1,scaling=constrained);
*)

#Check phi normalisation
(*
Vx:=1;
Vy:=2;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
m:=2;
a:=1;
n_of_x:= piecewise(x<0,n1,x>0,n2);
PhiFunc:= Phi_left(n1,n2,l1,l2,a,m) *HeavyStepLeft + Phi_right(n1,n2,l1,l2,a,m) *HeavyStepRight;
int((n_of_x)^2 *(PhiFunc)^2,x=-l1..l1,y=0..l2);
*)

#Check WaveMaker

Vx:=1;
Vy:=0;
Sigma:=10;
l1:=Pi;
l2:=Pi;
n1:=1;
n2:=2;
nodes:=10;
f:= FMaker(x,y,Vx,Vy,Sigma);
g:= GMaker(f,Vx,Vy,Sigma);
F:=Two_squares_CoeffFinder(f,nodes,n1,n2,l1,l2,a_array);
G:=Two_squares_CoeffFinder2(g,nodes,n1,n2,l1,l2,a_array);
wavee:= WaveMaker_Two_squares(F,G,nodes,l1,l2,a_array,n1,n2);
animate(densityplot,[wavee,x=-l1..l1,y=0..l2],t=0..10,frames=150,scaling=constrained);

## @helix The code I'm trying to execu...

@helix
The code I'm trying to execute and save the animation from is quite long, is there a way to upload a maple file here?

## @helix The code I'm trying to execu...

@helix
The code I'm trying to execute and save the animation from is quite long, is there a way to upload a maple file here?

## K5ky...

I have tried waiting an awful long time and still no luck, this version of maple seems to have a lot of problems and crashes frequently.

## K5ky...

I have tried waiting an awful long time and still no luck, this version of maple seems to have a lot of problems and crashes frequently.

 1 2 Page 2 of 2
﻿