Re: i want help in mathematica programing
- To: mathgroup at smc.vnet.net
- Subject: [mg42887] Re: i want help in mathematica programing
- From: jf.alcover at bdpme.fr (jf alcover)
- Date: Fri, 1 Aug 2003 01:25:56 -0400 (EDT)
- References: <bfrhop$qsr$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Farkhanda <farkhanda_yusaf at yahoo.com> wrote in message news:<bfrhop$qsr$1 at smc.vnet.net>...
> solve the heat equation using explicit difference approximation
> with the help of programming in mathematica
> let u= u(x,t)
> with the given boundary conditions u[0, t] =u[1, t] = 0 for 0 <= t <= 1
> and the initial condition u(x,0)=x when 0<=x<=1/2
> u(x,0)=1-x when 1/2 <=x<=1
> The explicit difference scheme is given as
> [u(x,t+delta t)-u(x,t)]/(delta t) =[u(x+delta x,t)-2u(x,t)+u(x-delta x,t)]/(delta x )^2
> with a condition
> (delta t) / (delta x )^2 <=(1/2)
Here is an example (not thoroughly tested,
not the best way, just my way ! )
of what could be done with mma :
The subscripts used are k for x[k] and n for t[n].
ClearAll[u, x, t];
deltat = 0.01;
deltax = 0.2;
(* we check the stability condition : *)
(deltat) / (deltax )^2
0.25
kmax = Floor[1/deltax]
5
nmax = Floor[1/deltat]
100
(* boundary conditions : *)
u[0, n_Integer /; 0 <= n <= nmax] = 0;
u[kmax, n_Integer /; 0 <= n <= nmax] = 0;
u[k_Integer /; 0 <= k <= kmax/2, 0] := k deltax;
u[k_Integer /; kmax/2 <= k <= kmax, 0] := 1 - k deltax;
(* difference scheme : *)
eq = (u[k, n + 1] - u[k, n])/(deltat) == (u[k + 1, n] - 2u[k, n] +
u[k - 1, n])/(deltax )^2;
(* solve for u[k,n+1] : *)
sol = Solve[eq, u[k, n + 1]]
{{u[k, 1 + n] ->
0.01 (100. u[k, n] + 25. (u[-1 + k, n] - 2. u[k, n] + u[1 + k, n]))}}
(* we need another subscript because on left hand side n+1 is not allowed : *)
sol2 = sol[[1]] /. n -> m - 1
{u[k, m] ->
0.01 (100. u[k, -1 + m] +
25. (u[-1 + k, -1 + m] - 2. u[k, -1 + m] + u[1 + k, -1 + m]))}
sol2[[1, 2]] // Simplify
0.25 u[-1 + k, -1 + m] + 0.5 u[k, -1 + m] + 0.25 u[1 + k, -1 + m]
(* we plug this formula into u[k,m] : *)
u[k_Integer, m_Integer /; m < 0] = 0;
u[k_Integer /; k < 0 || k > kmax] = 0;
u[k_Integer /; -kmax <= k <= 2kmax, m_Integer /; m > 0] :=
(* cache is useful for speed *)
u[k, m] = 0.25 u[-1 + k, -1 + m] + 0.50 u[k, -1 + m] + 0.25 u[1 + k, -1 + m];
(* results for the 5 first layers : *)
Table[u[k, n], {k, 0, kmax}, {n, 0, 5}]
{{0., 0., 0., 0., 0., 0.},
{0.2, 0.2, 0.1875, 0.171875, 0.15625, 0.141602},
{0.4, 0.35, 0.3125, 0.28125, 0.253906, 0.229492},
{0.4, 0.35, 0.3125, 0.28125, 0.253906, 0.229492},
{0.2, 0.2, 0.1875, 0.171875, 0.15625, 0.141602},
{0., 0., 0., 0., 0., 0.}}
I hope this will help you to build your own program.