(* ::Package:: *)

(************************************************************************)
(* This file was generated automatically by the Mathematica front end.  *)
(* It contains Initialization cells from a Notebook file, which         *)
(* typically will have the same name as this file except ending in      *)
(* ".nb" instead of ".m".                                               *)
(*                                                                      *)
(* This file is intended to be loaded into the Mathematica kernel using *)
(* the package loading commands Get or Needs.  Doing so is equivalent   *)
(* to using the Evaluate Initialization Cells menu command in the front *)
(* end.                                                                 *)
(*                                                                      *)
(* DO NOT EDIT THIS FILE.  This entire file is regenerated              *)
(* automatically each time the parent Notebook file is saved in the     *)
(* Mathematica front end.  Any changes you make to this file will be    *)
(* overwritten.                                                         *)
(************************************************************************)



InternalEq:= (
k=i+1;
A[[k,k-1]]=  1.0;
A[[k,k      ]]=-2.0;
A[[k,k+1]]=   1.0;
)
InternalRHS:=(
k=i+1;
F[[k]]=f[a+i*\[CapitalDelta]x,n*\[CapitalDelta]t];
)
DirichletBoundaryEq:=(
k=i+1;
Switch[ i,
0,   mult=  1; x=a,
nx, mult=-1; x=b
,_, Print["Error"]
];
A[[k,k]]=-1/\[Nu]x;
)
DirichletRHS:=(
k=i+1;
Switch[ i,
0,   mult=  1; x=a,
nx, mult=-1; x=b
,_, Print["Error"]
];
F[[k]]=g[x,t];
)
RobinBoundaryEq:=(
k=i+1;
Switch[ i,
0,   mult=  1; x=a,
nx, mult=-1; x=b
,_, Print["Error"]
];
A[[k,k]]=-2.0+(1 + mult \[Alpha][x] \[CapitalDelta]x);
A[[k,k+mult]]= 1.0;
)
RobinRHS:=(
k=i+1;
Switch[ i,
0,   mult=   1; x=a,
nx, mult=-1; x =b,
_, Print["Error"]
];
F[[k]]=f[x,n*\[CapitalDelta]t]- mult  \[CapitalDelta]x g[x,n*\[CapitalDelta]t];
)
PeriodicBoundaryEq:=(
k=i+1;
Switch[ i,
0, x=a,
nx,x=b,
_, Print["Error"]
];
A[[k,k]]=-2.0;
A[[k,2]]=  1.0;
A[[k,nx]]=1.0;
)
PeriodicRHS:=(
k=i+1;
Switch[ i,
0,   x=a,
nx, x =b,
_, Print["Error"]
];
F[[k]]=f[x,n*\[CapitalDelta]t];
)
IC:=(
k=i+1;
U0[[k]]=u0[a+i*\[CapitalDelta]x];
)


InternalA[nx_]:= Module[
{i,A=SparseArray[{},{nx+1,nx+1}]},
Do[
k=i+1;
A[[k,k-1]]=  1.0;A[[k,k      ]]=-2.0;A[[k,k+1]]=   1.0,
{i,1,nx-1}];
A
]
PeriodicA[nx_]:=Module[
{A=InternalA[nx],m},
A[[nx+1,nx+1]]=-2.0; A[[nx+1,2]]->1.0;A[[nx+1,nx]]->1.0;
A[[1,1]]=-2.0; A[[1,2]]->1.0;A[[1,nx]]->1.0;
A
]
RNDA[nx_, {\[CapitalDelta]x_,{\[Alpha]0_,\[Beta]0_,g0_},{\[Alpha]n_,\[Beta]n_,gn_}}]:=Module[
{A=InternalA[nx],k},
k=nx+1;
A[[k,k-1]]=1.0;A[[k,k]]=-2.0 + (\[Alpha]n-\[Beta]n*\[CapitalDelta]x)/\[Alpha]n ;
k=1;
A[[k,k+1]]=1.0;A[[k,k]]=-2.0 +(\[Alpha]0+\[Beta]0*\[CapitalDelta]x)/\[Alpha]0 ;
A
]


rhs[f_,t_,{a_,b_,\[CapitalDelta]x_}]:=Module[
{x},
Table[f[x,t],{x,a,b,\[CapitalDelta]x}]
]
Periodicrhs[f_,t_,{a_,b_,\[CapitalDelta]x_}]:= rhs[f,t,{a,b,\[CapitalDelta]x}]
RNDrhs[f_,t_,{a_,b_,\[CapitalDelta]x_},{\[Kappa]_,{\[Alpha]0_,\[Beta]0_,g0_},{\[Alpha]n_,\[Beta]n_,gn_}}]:=Module[
{F=rhs[f,t,{a,b,\[CapitalDelta]x}]},
F[[   1]]-=(\[Kappa] g0 )/(\[Alpha]0 \[CapitalDelta]x);
F[[-1]]+= (\[Kappa] gn )/(\[Alpha]n \[CapitalDelta]x);
F
]



