Re: bad PDE solution...

*To*: mathgroup at smc.vnet.net*Subject*: [mg129442] Re: bad PDE solution...*From*: "Nasser M. Abbasi" <nma at 12000.org>*Date*: Mon, 14 Jan 2013 23:27:44 -0500 (EST)*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com*Delivered-to*: l-mathgroup@wolfram.com*Delivered-to*: mathgroup-newout@smc.vnet.net*Delivered-to*: mathgroup-newsend@smc.vnet.net*References*: <kd03f3$6e0$1@smc.vnet.net>*Reply-to*: nma at 12000.org

On 01/13/2013 11:01 PM, Fred Bartoli wrote: > Hello, > > As a first simple test case of a more complicated pb I want to solve the > plane diffusion equation. > It seems NDSolve gives a bad answer to the simple step case but I can't > find what's happening... > > (* Set up the equation: *) > > eq1=D[h[t,x],x,x]- k1 D[h[t,x],t]==0 > > const=k1-> 10000 > > (*"solve" it...*) > > sol = NDSolve[ > {D[h[t, x], {x, 2}] - k1 D[h[t, x], t] == 0, h[0, x] == 0, > h[t, 1] == 0, h[t, 0] == 1 - Exp[-10^6*t]} /. const, > h, {t, 0, 10^-3}, {x, 0, 0.01}][[1]] > > Plot3D[(h[t,x])/.mag/.solCart,{x,0,10^-4},{t,0,10^-4},AxesLabel->{"x","t","h"},PlotRange->All,PlotLabel->"Bogus > solution"] > Plot3D[Erfc[x/(2*Sqrt[t/k1])]/.const,{x,0,10^-4},{t,0,10^-4},AxesLabel->{"X","t","h"},PlotRange->All,PlotLabel->"The > right one"] > > (* Test the Erfc[x/(2Sqrt[t/k1])] solution *) > Release[Hold[D[h[t,x],x,x]- k1*D[h[t,x],t]]/. > h[t,x]->Erfc[x/(2*Sqrt[t/k1])]]//Simplify > > > I hope I'm doing something wrong, but I can find what... > First, your code, little cleaned up to make it easy to see and to remove the `mag` and `solCart` which I have no idea what it is. There is no need to write `eq=` but later on, copy the same equation again. Why? Also, it is best to write the IC in one variable, and the equations in another, so things are easier to see. ---------------------------------------------- eq1 = D[h[t, x], x, x] - k1 D[h[t, x], t] == 0 ic = {h[0, x] == 0, h[t, 1] == 0, h[t, 0] == 1 - Exp[-10^6*t]} const = k1 -> 10000; sol = First@ NDSolve[Flatten@{eq1, ic} /. const, h, {t, 0, 1}, {x, 0, 0.01}]; plot[f_, x_, t_, title_] := Plot3D[f, {x, 0, 0.01}, {t, 0, 1}, AxesLabel -> {x, t, h}, PlotRange -> All, PlotLabel -> title] plot[h[t, x] /. sol, x, t, "Bogus solution"] plot[Erfc[x/(2*Sqrt[t/k1])] /. const, x, t, "The right one"] ------------------------------------------ The question though, why do you think that the initial conditions will not make any difference to the solution of the pde? You are comparing a numerical solution for some arbitrary IC you specficed above with one particular solution that you have found. --Nasser