Question #159192 on DOLFIN changed:
https://answers.launchpad.net/dolfin/+question/159192
N.A. Borggren proposed the following answer:
Hello,
I've been working on related Fokker Planck type problems (the 'adjoint' of
what you wrote down) and have encountered this same phenomena. I have not
managed to find a way to get directly the steady state but perhaps you can
learn from what I've tried. The fastest way I've been able to converge to a
steady state is by still solving in the time dependent way (like the
advection/diffusion demo) but using an initial condition 'as close as possible'
to the steady state. For example obtaining an initial guess by histogramming a
few of the corresponding dynamical system or langevin trajectories and solving
poisson's equation gives an initial distribution that converges quite rapidly.
Also artificially increasing the diffusion (high temperature limit) can speed
up the simulation and then you can turn it back down to the value of interest
when it seems closer to converged.
I know those aren't the answers you were hoping for, but that is the best
I've managed to work out so far. Part of the issue, I believe, is the 'steady
state' is, not necessarily the same for all initial conditions. It is often
the case that multiple initial conditions asymptotically convergence to the
same steady state but that technically takes infinitely long. The equation is
still a linear equation, a unique initial condition gives a unique final
condition, and two different initial conditions that seem 'converged' will
actual differ a bit in the details so some tolerance would need to be defined
and set to consider these two solutions identical.
Also if the boundary conditions are at all absorbing then the steady state
actually is the trivial solution unless there is an explicit driving force.
(take t to infinity in, say, Equation 2.40 Zwanzig, Nonequilibrium statistical
mechanics). So depending on your boundary condition dolfin and your form might
actually be working properly.
I was hoping to avoid all that with a reflecting boundary condition case,
(A_i * u)n_i + (d/dx_j (B^2_ij * u)n_i=0, on the mesh boundary with n_i the
normal to the mesh. This is how I keep the probability from going to zero in
the time dependent case as well since it conserves the probability inside the
mesh domain. For the steady state I have also tried to at least avoid the
trivial solution by forcing a normalization condition of 1 with a lagrange
multiplier but haven't yet succeeded with that approach either. The iterative
method of finding a good guess from the dynamics for the initial condition and
propagating that in time has been the fastest, albeit slightly unsatisfying
method for me.
I think for these steady states, as well as the time dependent case,
Fokker-Planck equations would work a lot better if one were able to take into
account that we need only be concerned with positive functions and could
restrict the finite elements used to insist on positivity of a function. Are
there strictly positive finite elements available or good ways to force the
standard elements to stay positive? It could be very useful here. This has
also been the hold up I've had in applying the fokker-planck equations in three
dimensions which seems to quickly exacerbate this issue.
Thanks,
Nathan
--
You received this question notification because you are a member of
DOLFIN Team, which is an answer contact for DOLFIN.
_______________________________________________
Mailing list: https://launchpad.net/~dolfin
Post to : [email protected]
Unsubscribe : https://launchpad.net/~dolfin
More help : https://help.launchpad.net/ListHelp
