Re: PDE with RecurrenceTable
- To: mathgroup at smc.vnet.net
- Subject: [mg124832] Re: PDE with RecurrenceTable
- From: Alexei Boulbitch <Alexei.Boulbitch at iee.lu>
- Date: Thu, 9 Feb 2012 05:34:10 -0500 (EST)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
Am 03.02.2012 08:13, schrieb Alexei Boulbitch:
> Dear Community,
>
>
>
> I am trying to make a simple numeric FEM solution of parabolic PDE. It seems that RecurrenceTable function is designed exactly for such a job. Indeed, in the Help/RecurrenceTable/Scope/Partial Difference Equations one finds an example. My problem is that this example is only one, and I would say, it is not basic enough.
>
>
>
> My question: do you know some other examples of the use of the RecurrenceTable for this sort of equations??
>
>
>
> I would like to explain: I already went through the MathGroup archive and have seen numerous posts recommending various sophisticated FEM packages. My question is not about them. I want to learn to make simple programs of this sort myself to fast test an equation at hand.
>
>
>
> For example, here is a classical equation of temperature conductivity:
>
>
>
> (\[PartialD]u(x,t))/\[PartialD]t=a^2*(\[PartialD]^2u(x,t))/\[PartialD]x^2
>
>
>
> with a=0.1, the boundary conditions: u[0,t]==1 and u[1,t]==0 and the initial condition u[x,0]== Cos[3 \[Pi]*x/2];
>
> using the explicit method on the rectangular lattice, taken from a textbook:
>
>
>
> u[j,k+1]=\[Sigma]*u[j+1,k]+(1-2*\[Sigma])*u[j,k]+\[Sigma]*u[j-1,k];
>
> \[Sigma]=(a^2*\[Tau])/h^2;
>
>
>
> Tau and h are temporal and spatial step sizes.
>
>
>
> This is the code:
>
>
>
> a = 0.1;
>
> h = 0.1;
>
> \[Tau] = 0.0001;
>
> \[Sigma] = a^2*\[Tau]/h^2;
>
> lst2 = RecurrenceTable[{u[j,
>
> k + 1] == \[Sigma]*u[j + 1, k] + (1 - 2*\[Sigma])*
>
> u[j, k] + \[Sigma]*u[j - 1, k], u[j, 0] == Cos[3 Pi*j/20.],
>
> u[0, k] == 1, u[10, k] == 0}, u, {j, 0, 10}, {k, 0, 100}];
>
>
>
>
>
> Show[{
>
> ListPlot3D[lst2,
>
> AxesLabel -> {Style["t", 16, Italic], Style["x", 16, Italic],
>
> Style["u", 16, Italic]}, PlotStyle -> Blue],
>
> Plot3D[0, {t, 0, 100}, {x, 0, 10},
>
> PlotStyle -> {Yellow, Opacity[0.4]}]
>
> }]
>
>
>
> The solution obtained this way, however, does not show evolution. It is clear that with increasing t the temperature, u, should forget its initial form and approach to a straight line. What is wrong?
>
>
>
> Thank you, Alexei
>
>
>
>
>
> Alexei BOULBITCH, Dr., habil.
>
> IEE S.A.
>
> ZAE Weiergewan,
>
> 11, rue Edmond Reuter,
>
> L-5326 Contern, LUXEMBOURG
>
>
>
> Office phone : +352-2454-2566
>
> Office fax: +352-2454-3566
>
> mobile phone: +49 151 52 40 66 44
>
>
>
> e-mail: alexei.boulbitch at iee.lu<mailto:alexei.boulbitch at iee.lu>
Hi Alexei,
it seems to me that your tau (delta-t) and/or maximum k is far too small
to observe the flattening.
I'll use ListConvolve, because it is much faster; the algorithm is still
the same:
a=1/10; h=1/10; \[Tau]=1/10000;
\[Sigma]=a^2*\[Tau]/h^2;
kmax=121000; (* just to get the zero of u[3/4,t] *)
To avoid memory problems when plotting the result, I'll take ~ 100 data
lines:
lst2 = NestList[
Flatten[{1, ListConvolve[{\[Sigma], 1 - 2*\[Sigma], \[Sigma]}, #1],
0}] & ,
Table[Cos[(3.*Pi*j)/20], {j, 0, 10}],
kmax][[1 ;; All ;; Floor[kmax/100]]];
Now the flattening is visible "to the naked eye":
Show[ListPlot3D[lst2, AxesLabel -> {Style["x", 16, Italic],
Style["t", 16, Italic], Style["u", 16, Italic]},
PlotStyle -> {Lighter[Blue], Specularity[White, 22]},
DataRange -> {{0, 1}, {0, kmax*\[Tau]}}, PlotRange -> All],
Plot3D[0, {x, 0, 1}, {t, 0, kmax*\[Tau]},
PlotStyle -> {Yellow, Opacity[0.4]}]]
Peter
Thank you, Peter,
Works really faster than my code and the boundary conditions are correctly forced.
Thank you once more, Alexei
Alexei BOULBITCH, Dr., habil.
IEE S.A.
ZAE Weiergewan,
11, rue Edmond Reuter,
L-5326 Contern, LUXEMBOURG
Office phone : +352-2454-2566
Office fax: +352-2454-3566
mobile phone: +49 151 52 40 66 44
e-mail: alexei.boulbitch at iee.lu<mailto:alexei.boulbitch at iee.lu>