1
$\begingroup$

I'm in the midst of solving the given problem

Given the setup I tried writing up the following:

ClearAll["Global`*"]
β = 0.05;
sol = NDSolveValue[{
        D[Subscript[U, z][η, τ], τ] == -Sin[τ] + β (1/η D[η*Subscript[U, z][η, τ], η, η]),
        Subscript[U, z][η, 0] == 0, Subscript[U, z][1, τ] == 0,},
        Subscript[U, z], {τ, 0, 10}]

How can I model the second boundary condition accordingly, and what should I do next to model the solution on a 3D plot?

$\endgroup$
  • 3
    $\begingroup$ Have a look at Menu/Help/WolframDocumentation/NeumannValue and there go to Applications/Time-Dependent Problems. There you will find examples which are no exactly equal to yours (that is, not radial), but very close and with comparable type of Neumann and initial conditions. $\endgroup$ – Alexei Boulbitch Nov 19 '19 at 8:48
  • $\begingroup$ The book does not reference that function so I will not use that $\endgroup$ – TexMexDex Nov 19 '19 at 21:22
4
$\begingroup$

A couple things. Since $\eta$ is in the denominator of your pde, your solution will have a problem when $\eta$ is 0. Luckily since zero is one of the end points we can pick a number close to zero without changing the answer much. A numerical solution is an approximation after all. Also your NDSolveValue needs the limits of $\eta$ as well as $\tau$. Fixing that and choosing a very small number for zero, we have:

ClearAll["Global`*"]
β = 0.05;
zero = 10^-12;

sol = NDSolveValue[{D[U[η, τ], τ] == -Sin[τ] + β*(1/η)*D[η*D[U[η, τ], η], η], 
    U[η, 0] == 0, U[1, τ] == 0, Derivative[1, 0][U][zero, τ] == 0}, U, {η, zero, 1}, 
   {τ, 0, 10}]

tp = Table[
   Plot[Evaluate[sol[η, τ]], {η, zero, 1}, 
    PlotRange -> {-2, 1}], {τ, 0, 10, .1}];
ListAnimate[tp]

enter image description here

In this particular case we can include $\eta = 0$ in the solution because the derivative wrt to $\eta$ is zero there. Expanding the pde we have

Derivative[0, 1][U][η, τ] == (β*Derivative[1, 0][U][η, τ])/η + 
   β*Derivative[2, 0][U][η, τ] - Sin[τ]

the term with $\eta$ in the denominator is $0/0$ at $\eta = 0$, so we can use L'Hospital's rule. Rewriting the pde with Piecewise functions to separate the zero case from the general case and solving we get:

Clear["Global`*"]

pde = Derivative[0, 1][U][η, τ] == 
   Piecewise[{{β*Derivative[2, 0][U][η, τ], η == 0}, 
      {(β*Derivative[1, 0][U][η, τ])/η, True}}] + β*Derivative[2, 0][U][η, τ] - 
    Sin[τ]

β = .05

sol = NDSolveValue[{pde, U[η, 0] == 0, U[1, τ] == 0, 
   Derivative[1, 0][U][0, τ] == 0}, 
  U, {η, 0, 1}, {τ, 0, 10}];

The resulting plot looks the same as before, but we can now start the solution at $\eta = 0$.

$\endgroup$
  • $\begingroup$ What is the purpose of the code: β*Derivative[2, 0][U][η, τ] $\endgroup$ – TexMexDex Nov 22 '19 at 0:21
  • $\begingroup$ It comes from L'Hospital's rule. When you have f[x]/g[x] = 0/0 at x = 0. You can often find the limit as x->0 by finding f'[x]/g'[x]. So when you have in this case an expression like u'[x]/x, find the limit at 0 by taking the numerator and denominator derivatives and get u''[x]/1. $\endgroup$ – Bill Watts Nov 22 '19 at 1:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.