0
|
2788 Views
|
2 Replies
|
0 Total Likes
View groups...
Share
GROUPS:

# Num. Solution of Partial Differential Eq.--The Numerical Method of Lines

Posted 11 years ago
 Hi All,I am trying to solve partial diff. eq. by using method of Numerical Method of LinesThis is the linkhttps://reference.wolfram.com/mathematica/tutorial/NDSolvePDE.htmlI've got some ideaBut it does not give me the right graph.Can anyone check my Mathematica code?This is my PDESubscript[h, t] = 0= Subscript[\[Alpha], 1] Subscript[h, xxxx] +   Subscript[\[Alpha], 2] Subscript[h, xx] +   Subscript[\[Alpha], 5] Subscript[h, xx] Subscript[h, x]^2I really appriciate it n = 100;(*the number of steps*) L = 2*\[Pi];(*the length of the period*) h = L/n;(*the size of the steps*) k = 1;(*frequency*) \[Delta] = 0.5 U[t_] = Table[Subscript[u, i][t], {i, 0, n}]; initc = Thread[U[0] == Table[0, {n + 1}]]; Subscript[\[Alpha], 1] = -8; Subscript[\[Alpha], 2] = -216; \ Subscript[\[Alpha], 5] = -3 Subscript[\[Alpha], 2]/2;wavenum =   Sqrt[(Subscript[\[Alpha], 2]/Subscript[\[Alpha], 1])]/Sqrt[2];count1 = 1;While[count1 <= n + 1,   initc[[count1]] =    Thread[U[0][[count1]] == 1 + \[Delta] Cos[k (count1 - 1) h]];   count1++];count1 = 1;N[initc]count2 = 1;eqnsnes = Table[0, {n + 1}];While[count2 <= n + 1,   eqnsnes[[count2]] = Thread[D[U[t][[count2]], t] ==           If[count2 == 1,       Subscript[\[Alpha],         1] (1/h^4 U[t][[count2 + 2]] - 4 U[t][[count2 + 1]] +           6 U[t][[count2]] - 4 U[t][[n]] + U[t][[n - 1]]) +        Subscript[\[Alpha],         2] ((U[t][[count2 + 1]] - 2 U[t][[count2]] + U[t][[n]])/h^2)/         wavenum^2 +        Subscript[\[Alpha],         5] ((U[t][[count2 + 1]] - 2 U[t][[count2]] + U[t][[n]])/h^2)        *((U[t][[count2 + 1]] - U[t][[n]])/(2 h))^2,              If[count2 == 2,        Subscript[\[Alpha],          1] (1/h^4 (U[t][[count2 + 2]] - 4 U[t][[count2 + 1]] +              6 U[t][[count2]] - 4 U[t][[count2 - 1]] + U[t][[n]]) ) +         Subscript[\[Alpha],          2] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +              U[t][[count2 - 1]])/h^2)/wavenum^2 +                Subscript[\[Alpha],          5] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +             U[t][[count2 - 1]])/          h^2)*((U[t][[count2 + 1]] - U[t][[count2 - 1]])/(2 h))^2,              If[count2 == n,        Subscript[\[Alpha],           1] (1/h^4 ((U[t][[2]] - 4 U[t][[count2 + 1]] +                6 U[t][[count2]] - 4 U[t][[count2 - 1]] +                U[t][[count2 - 2]]) )) +          Subscript[\[Alpha],           2] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +               U[t][[count2 - 1]])/h^2)/wavenum^2 +          Subscript[\[Alpha],           5] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +              U[t][[count2 - 1]])/           h^2)*((U[t][[count2 + 1]] - U[t][[count2 - 1]])/(2 h))^2,                If[count2 == n + 1,          Subscript[\[Alpha],            1] (1/h^4 (U[t][[3]] - 4 U[t][[2]] + 6 U[t][[count2]] -                4 U[t][[count2 - 1]] + U[t][[count2 - 2]]) ) +           Subscript[\[Alpha],            2] ((U[t][[2]] - 2 U[t][[count2]] + U[t][[count2 - 1]])/             h^2)/wavenum^2 +           Subscript[\[Alpha],            5] ((U[t][[2]] - 2 U[t][[count2]] + U[t][[count2 - 1]])/            h^2)*((U[t][[2]] - U[t][[count2 - 1]])/(2 h))^2,          Subscript[\[Alpha],            1] (1/h^4 (U[t][[count2 + 2]] - 4 U[t][[count2 + 1]] +                6 U[t][[count2]] - 4 U[t][[count2 - 1]] +                U[t][[count2 - 2]]) ) +           Subscript[\[Alpha],            2] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +                U[t][[count2 - 1]])/h^2)/wavenum^2 +           Subscript[\[Alpha],            5] ((U[t][[count2 + 1]] - 2 U[t][[count2]] +               U[t][[count2 - 1]])/            h^2)*((U[t][[count2 + 1]] - U[t][[count2 - 1]])/(            2 h))^2]]]]]; count2++];eqnsneslines = NDSolve[{eqnsnes, initc}, U[t], {t, 0, 1}, MaxSteps -> 20000]ParametricPlot3D[ Evaluate[Table[{i*h, t, First[Subscript[u, i][t] /. lines]}, {i, 0,     n}]], {t, 0, 1}, PlotRange -> All, AxesLabel -> {"Spatial(x)", "Time(t)", "Height(h)"}]
2 Replies
Sort By:
Posted 11 years ago
 yes, it is second derivative. I dont why it shows that way. Also, in my code variable is t instead of x..how can I copy and paste my code as it is orijinal?
Posted 11 years ago
 Are you thinking that Subscript[h, xx] performs a differentiation?  It is just a static typeset thing.