MATH761

In[1]:=

x \[Element] Reals; $Assumptions = h > 0;

Dp[f_[x_]] := f[x + h] - f[x];

Dm[f_[x_]] := f[x] - f[x - h];

d[f_[x_]] := f[x+h/2]-f[x-h/2];

Dp[a_ f_[x_]]=a Dp[f[x]]; Dp[a_+b_]=Dp[a]+Dp[b]; Dp[a_ (b_+c_)]=a Dp[b]+a Dp[c]; Dp[a_]=0;

Dm[a_ f_[x_]]=a Dm[f[x]]; Dm[a_+b_]=Dm[a]+Dm[b]; Dm[a_ (b_+c_)]=a Dm[b]+a Dm[c]; Dm[a_]=0;

d[a_ f_[x_]]=a d[f[x]]; d[a_+b_]= d[a]+ d[b]; d[a_ (b_+c_)]=a d[b]+a d[c]; d[a_]=0;

fdSort[expr_]:=Assuming[{h > 0, x > 0}, Sort[Level[expr, 1], Level[#1, {Depth[#1] - 1}] < Level[#2, {Depth[#2] - 1}] &]];

Q[x_]=Subscript[Q,x];

Dp[f_[x_], p_] :=

Together[Normal[1/h Series[Log[1 + z], {z, 0, p}]] /.

Table[z^j -> Nest[Dp, f[x], j], {j, 1, p}]];

Dm[f_[x_], p_] :=

Together[Normal[-1/h Series[Log[1 - z], {z, 0, p}]] /.

Table[z^j -> Nest[Dm, f[x], j], {j, 1, p}]];

d[f_[x_], p_] := Together[Normal[2/h Series[ArcSinh[z/2], {z, 0, 2 p - 1}]] /.

Table[z^j -> Nest[d, f[x], j], {j, 1, 2 p - 1}]];

Dp[f_[x_], p_, k_] :=

Together[Expand[(Normal[1/h Series[Log[1 + z], {z, 0, p}]])^k] /.

Table[z^j -> Nest[Dp, f[x], j], {j, 1, p k}]];

Dm[f_[x_], p_, k_] :=

Together[Expand[(Normal[-1/h Series[Log[1 - z], {z, 0, p}]])^k] /.

Table[z^j -> Nest[Dm, f[x], j], {j, 1, p k}]];

d[f_[x_], p_, k_] := Together[Expand[(Normal[2/h Series[ArcSinh[z/2], {z, 0, 2 p - 1}]])^k]

/. Table[z^j -> Nest[d, f[x], j], {j, 1, (2 p - 1)k}]];

{qt1,qxx2,qxx4}={Dp[f[x],1],d[f[x],1,2],d[f[x],2,2]}

${\frac{f (h + x) - f (x)}{h}, \frac{f (x - h) + f (h + x) - 2 f (x)}{h^{2}}, \frac{f (x - 3 h) - 54 f (x - 2 h) + 783 f (x - h) + 783 f (h + x) - 54 f (2 h + x) + f (3 h + x) - 1460 f (x)}{576 h^{2}}}$

In[2]:=

Series[qxx4,{h,0,4}]

$f^{''} (x) - \frac{3}{320} h^{4} f^{(6)} (x) + O (h^{5})$

Note that centered formulas for approximating the second derivative of comparable accuracy may be obtained by using a smaller stencil.

In[3]:=

dc[f_[x_], p_, k_] := Module[{z,n,zn,ser,serc,fd},

ser = Expand[(Normal[2/h Series[ArcSinh[z/2], {z, 0, 2 p - 1}]])^k];

n = Exponent[ser,z]; zn = Coefficient[ser,z,n] z^n;

serc = ser - zn;

fd = serc /. Table[z^j -> Nest[d, f[x], j], {j, 1, (2 p - 1)k}];

Return[Together[fd]]];

qxx4=dc[f[x],2,2]

$\frac{- f (x - 2 h) + 16 f (x - h) + 16 f (h + x) - f (2 h + x) - 30 f (x)}{12 h^{2}}$

In[4]:=

Series[qxx4,{h,0,4}]

$f^{''} (x) - \frac{1}{90} h^{4} f^{(6)} (x) + O (h^{5})$

In[5]:=

2Eigenvalues of matrix from semi-discretization

is obtained. The eigenvalues of

𝐀

are necessary for stability analysis of various time integration schemes to solve (1). Since

\frac{1}{h^{2}} 𝐀 ≅ \partial_{x}^{2}

, the eigenvectors of

𝐀

are discrete approximations of the eigenfunctions of

\partial_{x}^{2}

2.1Homogeneous Dirichlet boundary conditions

Separation of variables applied to

q_{t} = q_{x x}

q (t, x) = T (t) X (x)

gives

In[31]:=

op = -X”[x];

lftBC = DirichletCondition[X[x]==0,x==0];

rgtBC = DirichletCondition[X[x]==0,x==1];

DEigensystem[{op,lftBC,rgtBC},X[x],{x,0,1},7]

$(\begin{array}{ccccccc} π^{2} & 4 π^{2} & 9 π^{2} & 16 π^{2} & 25 π^{2} & 36 π^{2} & 49 π^{2} \\ \sin () & \sin (2) & \sin (3) & \sin (4) & \sin (5 π x) & \sin (6 π x) & \sin (7 π x) \end{array})$

In[48]:=

qxx2

$\frac{f (x - h) + f (h + x) - 2 f (x)}{h^{2}}$

In[49]:=

xk[x_] = Sin[k Pi x];

Axki = qxx2 /. {f->xk, x-> i h}

$\frac{- 2 \sin (π h i k) + \sin (π k (h i - h)) + \sin (π k (h i + h))}{h^{2}}$

In[55]:=

lambda2[h_]=FullSimplify[Axki/xk[i h]]

$\frac{2 (\cos (π h k) - 1)}{h^{2}}$

In[56]:=

Limit[lambda2[h],h->0]

$π^{2} (- k^{2})$

2.2Semi-homogeneous Dirichlet boundary conditions

Separation of variables applied to

q_{t} = q_{x x}

q (t, x) = T (t) X (x)

gives

with eigenfunctions

X_{k} (x) = \sin [(4 k + 1) π x / 2]

λ_{k} =

(4 k + 1) π / 2

$\frac{f (x - h) + f (h + x) - 2 f (x)}{h^{2}}$

In[58]:=

xk[x_] = Sin[(4k+1) Pi x/2];

Axki = qxx2 /. {f->xk, x-> i h}

$\frac{- 2 \sin (\frac{1}{2} π h i (4 k + 1)) + \sin (\frac{1}{2} π (4 k + 1) (h i - h)) + \sin (\frac{1}{2} π (4 k + 1) (h i + h))}{h^{2}}$

In[59]:=

lambda2[h_]=FullSimplify[Axki/xk[i h]]

$- \frac{4 \sin^{2} (π h k + \frac{π h}{4})}{h^{2}}$

In[60]:=

Limit[lambda2[h],h->0]

$- \frac{1}{4} {(4 π k + π)}^{2}$

1Finite difference approximations

2Eigenvalues of matrix from semi-discretization

2.1Homogeneous Dirichlet boundary conditions

2.2Semi-homogeneous Dirichlet boundary conditions