Message Boards Message Boards

1
|
6922 Views
|
5 Replies
|
10 Total Likes
View groups...
Share
Share this post:

Boundary condition error with a partial differential equation

Slavi T
Posted 11 years ago
Hallo everyone emoticon ,

I'm new to Mathematica so I need some help solving this steady state mass transfer problem in a pipe.
Due to simplisity I set all factors to one and the R as well. Here is my partial differential equation
NDSolve[{2*(1 - (r/R)^2)*
    D[\[Rho][r, z], {z}] == (r*D[\[Rho][r, z], {r, 2}] +
      D[\[Rho][r, z], {r}])/r,
  D[\[Rho][0, z], {r}] == 0, \[Rho][R, z] == 0, \[Rho][r < R, 0] ==
   1}, \[Rho], {r, 0, R}, {z, 0, 10}, Reals]
I get the error that my furst boundary condition is aways true and I can't fix this error.
Do you have an Idea how to fix it. I want to plot the solution as well.
Thanks
Slavi
POSTED BY: Slavi T
5 Replies
Sort By:
To solve the PDE using NDSolve, you will have to formulate the problem as follows:
R = 1 + 10^-15;
eqns = (2 (1 - (r/R)^2)) D[rho[r, z], z] == (r D[rho[r, z], r, r] + D[rho[r, z], r])/r ;
bcs = {rho[R, z] == 0, Derivative[1, 0][rho][10^-19, z] == 0, rho[r, 0] == 1};

sol = NDSolveValue[{eqns, bcs}, rho[r, z], {r, 10^-19, R}, {z, 0, 10}, Method -> {"PDEDiscretization" -> {"MethodOfLines",
"DifferentiateBoundaryConditions" -> {True, "ScaleFactor" -> 10^4}}}]
A couple of points need to be made here. There are singularities at r=0 and (r/R) = 1. So, in order to avoid the left singulariy, the boundary conditions are offset by a tiny bit. The second point is associated with the initial condition rho[z,0] = 1. This is inconsistent with the BC rho[R,z]=0. This is why the additional options have been added in NDSolve to give additional weight to the BCs.

The solution can be plotted as:
Plot3D[sol, {r, 0, R}, {z, 0, 0.5}, PlotRange -> All]
POSTED BY: Paritosh Mokhasi
The usual cause of that warning is that something (probably ? in your case) is already defined.
You might try
??
?z
?r
to make sure none has a value. Or to just remove any definitions,
Clear[?, z, r]
POSTED BY: Bruce Miller
If you are using V8, then you would have to change the syntax a little bit as follows:
R = 1 + 10^-15;
eqns = (2 (1 - (r/R)^2)) D[rho[r, z],z] == (r D[rho[r, z], r, r] + D[rho[r, z], r])/r;
bcs = {rho[R, z] == 0, Derivative[1, 0][rho][10^-19, z] == 0,rho[r, 0] == 1};

sol = Quiet@NDSolve[{eqns, bcs}, rho, {r, 10^-19, R}, {z, 0, 10},Method -> {"MethodOfLines",
"DifferentiateBoundaryConditions" -> {True,"ScaleFactor" -> 10^4}}] ;
If you want to plot the solution from {-R, R}, then you can just mirror the solution.
Show[Plot3D[Evaluate[rho[r, z] /. sol], {r, 0, R}, {z, 0, 0.5},PlotRange -> All] ,
Plot3D[Evaluate[rho[-r, z] /. sol], {r, -R, 0}, {z, 0, 0.5}]]
POSTED BY: Paritosh Mokhasi
Paritosh, this is very interesting. Is it possible to plot this solution?
POSTED BY: Vitaliy Kaurov
Slavi T
Slavi T
Posted 11 years ago
 Paritosh, thank you very much. You probably use Mathematica 9 that's why it didn't work with the plot. But I plot it this way.
R = 1 + 10^-15;
eqns = (2 (1 - (r/R)^2)) D[rho[r, z],
     z] == (r D[rho[r, z], r, r] + D[rho[r, z], r])/r;
bcs = {rho[R, z] == 0, Derivative[1, 0][rho][10^-19, z] == 0,
   rho[r, 0] == 1};

sol = NDSolve[{eqns, bcs}, rho[r, z], {r, 10^-19, R}, {z, 0, 10}]  (*NDSolve instead of NDSolveValue *)
Plot3D[rho[r, z] /. sol, {r, 0, R}, {z, 0, 0.5}, PlotRange -> All] (* I added rho[r, z] /.sol here *)
I wonder how can I solve the problem, if I have to plot this over the hole pipe. The radius is not going from 0 to R, but from -R to R. It sould be like this plot but with additional mirrored plot of this across the z axis. In this case I'll have somewere problem whit the zero.
POSTED BY: Slavi T
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract