0
$\begingroup$

So recently I tried modelling diffusion over the rough part of a potential W. To do so I try numerically solving the Fokker-Planck equation: $$\partial_t P(x,t) = -\nabla\cdot J$$ $$J= -D(x)\nabla P-D(x)(\nabla W)P$$ I assume spherical symmetry and can hence reduce the equation to a one dimensional case. I have managed to solve this equation well enough if I assume the Diffusion D and Potential W to be smoothstep functions, however I now wanted to see what would happen if I introduce noise on the potential. The code I used is:

NeuB = NeumannValue[0, x == 0 || x == 30];
molfem[measure_ : Automatic] := {"MethodOfLines", 
   "SpatialDiscretization" -> {"FiniteElement", 
     "MeshOptions" -> MaxCellMeasure -> measure}}; 
SmoothstepPot[x_] = Piecewise[{{h1, x < x1},
     {(-(h2 - h1)/2)*Cos[(Pi/(x2 - x1))*(x - x1)] + (h1 + h2)/2, 
      x1 <= x < x2},
     {h2, x >= x2}}] /. {x1 -> 10, x2 -> 12, h1 -> 0.1, 
    h2 -> 3};
SmoothstepDiffus[x_] =  Piecewise[{{h1, x < x1},
     {(-(h2 - h1)/2)*Cos[(Pi/(x2 - x1))*(x - x1)] + (h1 + h2)/2, 
      x1 <= x < x2},
     {h2, x >= x2}}] /. {x1 -> 10, x2 -> 12, h1 -> 1, h2 -> 30};
J3d[x_] =  -(SmoothstepDiffus[x]*D[P3d[x, t], x]) - (SmoothstepDiffus[x]*D[SmoothstepPot[x], x]*P3d[x, t]);
pde3d = D[P3d[x, t], t] + 1/x^2 D[(x^2*J3d[x]), x] == NeuB;
Evals3d = 
  NDSolveValue[{pde3d, P3d[x, 0] == UnitStep[x - 12]/18}, 
   P3d, {x, 0, 30}, {t, 0, 2000}, Method -> molfem[0.01]];
Manipulate[Plot[Evals3d[x, t], {x, 0, 30}], {t, 0, 500}]

and this works well enough. However, when I try to replace my smooth potential SmoothstepPot with one, that is rough in the transition between the two potential planes:

RoughPot[x_] = Piecewise[{{h1, x < x1},
     {(-((h2 - h1)/2)*Cos[\[Pi]/(x2 - x1)*(x - x1)]) + (h1 + h2)/2 + 
       0.1*(h2 - h1)*
        PDF[NormalDistribution[(x2 + x1)/2, (x2 - x1)/4], 
         x]*(Sin[173 x] + Cos[90 x]), x1 <= x < x2},
     {h2, x >= x2}}] /. {x1 -> 10, x2 -> 12, h1 -> 0.1, h2 -> 3};

several interesting and unoptimal things occur: First the solution becomes highly dependent on the measure of the molfem method. Ergo, changing its parameter from 0.01 to 0.001 can yield very differing results. Second the results are not correct anymore. One would expect the smooth and rough case to behave somewhat similar, especially after a long time, but this is not the case here. I suspect, that this is an issue with how I use the finite element method, but I would appreciate any tipps on how to handle problems such as these!

$\endgroup$

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.