# Plot Generalized Heat Equation

Given the following problem $$u_{t}-\frac{v}{2}\cdot u_{xx}+\frac{v}{2}\cdot x^2 \cdot u(x)=0$$ $$u(x,0)=f(x)$$

Where $$f(x)=y_1(x)+0.2y_4(x)+0.01y_6 (x)$$

$$y_n$$ are the Eigenfunctions that defined as $$y_n(x)=\exp\left(\displaystyle{\frac{-x^2}{2}}\right)H_n(x)$$

By separation we assume $$u(x,t)=X(x)T(t)$$

$$\frac{T'(t)}{T(t)}=v (\frac{1}{2} \frac{X''(x)}{X(x)}-\frac{x^2}{2})=-\lambda$$

So we have

$$T'(x)=-\lambda v T(x)$$

I have found that the general solution of $$T$$ is $$T(x)=c_1e^{-\lambda v t}$$ and for $$\frac{1}{2} \frac{X''(x)}{X(x)}-\frac{x^2}{2}=-\lambda$$ which gives $$-\frac{1}{2}X''(x)+\frac{x^2}{2}X(x)=\lambda X(x)$$

the general solution of $$X$$ is $$X_n(x)=\exp\left(\displaystyle{\frac{-x^2}{2}}\right)H_n(x)$$

So the general solution to our initial problem is $$u_n(x,t)=\sum_{n=1}^{\infty}A_n \exp\left(\displaystyle{\frac{-x^2}{2}}\right)H_n(x) c_1e^{-\lambda v t}$$

Using the I.C we got $$u(x,0)=\sum_{n=1}^{\infty}A_n \exp\left(\displaystyle{\frac{-x^2}{2}}\right)H_n(x)=f(x)$$

I found an analytical solution with the following code.

f[x_]:=Exp[-x^2/2] (1+0.2*HermiteH[4,x]+0.01*HermiteH[6,x])
n=10;
A=Table[1/(Sqrt[Pi] 2^m m!) Integrate[f[x] HermiteH[m,x] Exp[x^2/2],{x,0,1}],{m,1,n}];
u[x_,t_]:=Sum[A[[m]] Exp[-m/2 t] Exp[-x^2/2] HermiteH[m-1,x],{m,1,n}]
Plot[{u[x,0],u[x,0.3],u[x,0.7]},{x,-20,20},PlotRange->All,PlotLegends->{"t=0","t=0.3","t=0.7"}]
fig1=Plot3D[u[x,t], {t,0,10},{x,-21,21},PlotRange->All,ColorFunction->"BlueGreenYellow"]


I tried also to write the following code, thus I could compare which one gives me better results but it does not seem to work

(* Define the Eigenfunctions *)
y[n_, x_] := Exp[-x^2/2] HermiteH[n, x]
(* Define the constant *)
v = 1;
(* Find the general solution for T *)
T[x_, t_] := c1 Exp[-v \[Lambda] t]
(* Find the general solution for X *)
X[n_, x_] := y[n, x]
(* Combine to get the general solution for u *)
u[n_, x_, t_] := X[n, x] T[x, t]
(* Define the initial condition *)
f[x_] := y[1, x] + 0.2 y[4, x] + 0.01 y[6, x]
(* Solve for the coefficients *)
An[n_] := (Integrate[f[x] X[n, x], {x, -\[Infinity], \[Infinity]}]) / (Integrate[X[n, x]^2, {x, -\[Infinity], \[Infinity]}])
(* Define the analytical solution *)
u[x_, t_] := Sum[An[n] X[n, x] T[x, t], {n, 1, \[Infinity]}]
(* Plot the solution for t=0, 0.1, 0.2, 0.3 *)
Plot[{u[x, 0], u[x, 0.1], u[x, 0.2], u[x, 0.3]},{x, -5, 5},PlotRange -> All]

• Your "general solution Xn" doesn't solve the ode Xn''[x]==(x^2-n) Xn[x]I think! Mar 22 at 9:40
• @UlrichNeumann the general solution is $X_n(x)=\exp\left(\displaystyle{\frac{-x^2}{2}}\right)H_n(x)$ Mar 22 at 9:47
• That's understood. But this Xn should solve Xn''[x]==(x^2-n) Xn[x] Mar 22 at 9:57
• Now about your $X_n(x)$ solution. This is an eigenvalue BVP ode. Not a normal ode. i.e. you can't just use DSolve to solve it. It will not work. You need to use DEigensystem to find the eigenvalues and eigenfunctions. As mentioned in the answer to your question [mathematica.stackexchange.com/questions/282594/… it is not possible to analytically solve this eigenvalue BVP to find the eigenvalues $\lambda$. But can be done numerically. Mar 22 at 13:24
• fyi, this page shows how to solve eigenvalue BVP for all possible homogeneous BC's for one standard second order ode. You will have to do the same for the $X(x)$ ode above. But because your ode is much harder, this makes it not possible to do analytically. But it is the same logic that needs to be applied. Mar 23 at 4:49

Here is 3D plot, analytical vs. numerical, using solution given in Analytical Solution in Generalized Heat Equation ## Code

ClearAll[x, t,n];
v = 1;
sol[x_, t_] :=
Exp[-x^2/2]*(Exp[-v/2*3*t]*HermiteH[1, x] +
1/5*Exp[-v/2*9*t]*HermiteH[4, x] +
1/100*Exp[-v/2*13*t]*HermiteH[6, x])
Y[n_] := Exp[-x^2/2]*HermiteH[n, x]
ic = u[x, 0] == Y + 0.2 Y + 0.01 Y
pde = D[u[x, t], t] + v/2*x^2*u[x, t] == v/2*D[u[x, t], {x, 2}];
L = 10; (*for NDSolve only*)
bc = {u[-L, t] == 0, u[L, t] == 0}
nsol = NDSolveValue[{pde, ic, bc}, u, {x, -L, L}, {t, 0, 3}];

Manipulate[
Module[{t},
Grid[{{Plot3D[Evaluate@sol[x, t0], {x, -L, L}, {t, 0, t0},
PlotRange -> {Automatic, {0, 2}, {-2, 2}},
PerformanceGoal -> "Quality",
PlotRange -> {Automatic, {-4, 4}},
ImageSize -> 300,
PlotLabel -> "Analytical",
AxesLabel -> {"x", "time", "u(x,t)"},
BaseStyle -> 14]}
,
{Plot3D[Evaluate@nsol[x, t0], {x, -L, L}, {t, 0, t0},
PlotRange -> {Automatic, {0, 2}, {-2, 2}},
PerformanceGoal -> "Quality",
PlotRange -> {Automatic, {-4, 4}},
ImageSize -> 300,
PlotLabel -> "Numerical",
AxesLabel -> {"x", "time", "u(x,t)"},
BaseStyle -> 14]}
}]
],
{{t0, 0.01, "time"}, 0.01, 2, .01, Appearance -> "Labeled",
ContinuousAction -> True}, TrackedSymbols :> {t0}]


And

Plot3D[sol[x, t], {x, -L, L}, {t, 0, 6},
PlotRange -> {Automatic, {0, 6}, {-4, 4}}] 