How to calculate analytical solution to fokker planck equation in finite domain?

Question

I am looking for the solution to the following Fokker-Planck equation.

$$\frac{\partial f(x,t)}{\partial t}= (k_1 x - v)\frac{\partial f(x,t)}{\partial x} + (k_2 x + D) \frac{\partial^2 f(x,t)}{\partial x^2},$$ where $k_1, v, k_2$ and $D$ are constants. In a domain $x \in (0, \infty)$. with boundary conditions that flux at $x=0$ is zero and $f(0,t)=0$ and initial condition, $f(1,0)=1$. The particle has an absorbing site at $x=\bar{x}$. In particular I want to find out the distribution of first passage time for this Fokker-Planck equation.

Answer 1

First of all this is not really a Fokker-Planck equation since in this form the probability flux is not conserved but let us put aside this problem now.

In order to solve the partial differential equation(PDE) we make an ansatz $f(x,t) = X(x) T(t)$. Inserting this into the PDE we get: \begin{eqnarray} (k_2 x+D) X^{''}(x) + (k_1 x-v) X^{'}(x) + E \cdot X(x)=0 \quad (i)\\ T^{'}(t) + E \cdot T(t) = 0 \quad (ii) \end{eqnarray} where $E \ge 0$ is a constant. Now using Yet another example of linear second order ODEs being reduced to hypergeometric functions. we get the solution. We have: \begin{eqnarray} &&T(t) = \exp(-E \cdot t)\\ &&X(x) = e^{-\frac{k_1 x}{k_2}} \cdot \\ &&\left( C_1 \, _1F_1\left(-\frac{D_1 k_1^2+k_2 v k_1+E k_2^2}{k_1 k_2^2};-\frac{D_1 k_1+k_2 v}{k_2^2};\frac{\sqrt{k_1^2} (D_1+k_2 x)}{k_2^2}\right)+ C_2 U\left(-\frac{D_1 k_1^2+k_2 v k_1+E k_2^2}{k_1 k_2^2},-\frac{D_1 k_1+k_2 v}{k_2^2},\frac{\sqrt{k_1^2} (D_1+k_2 x)}{k_2^2}\right)\right) \quad (iii) \end{eqnarray} where $U$ is the confluent hypergeometric function.

In[1470]:= 
k1 =.; k2 =.; DD =.; v =.; EE =.; D1 =.; x =.;

y[x_] := E^(-((k1 x)/ 
    k2)) (C[1] Hypergeometric1F1[-((D1 k1^2 + EE k2^2 + k1 k2 v)/(
        k1 k2^2)), -((D1 k1 + k2 v)/k2^2), (Sqrt[k1^2] (D1 + k2 x))/
       k2^2] + C[
       2] HypergeometricU[-((D1 k1^2 + EE k2^2 + k1 k2 v)/(
        k1 k2^2)), -((D1 k1 + k2 v)/k2^2), (Sqrt[k1^2] (D1 + k2 x))/
       k2^2] );
eX = ((k2 x + D1) D[y[x], {x, 2}] + (k1 x - v) D[y[x], x] + (EE) y[x]);
{k2, D1, k1, v, EE, x} = RandomReal[{0, 5}, 6, WorkingPrecision -> 50];
eX // Expand




Out[1474]= 0.*10^-45 C[1] + 0.*10^-43 C[2]

Now in order to proceed further we reduce equation (i) to a Sturm-Liouville equation using the standard way (see https://en.wikipedia.org/wiki/Sturm%E2%80%93Liouville_theory#The_integrating_factor_for_a_general_second-order_differential_equation). We define: \begin{equation} \mu(x):= e^{\frac{k_1 x}{k_2}} (D+k_2 x)^{-\frac{D k_1+k_2 (k_2+v)}{k_2^2}} \end{equation} and we have: \begin{equation} \frac{d}{d x} \left[\mu(x) (k_2 x+D) X^{'}(x)\right] + \mu(x) \cdot E \cdot X(x) = 0 \quad (iv) \end{equation} Now in order to achieve orthogonality (cmp https://proofwiki.org/wiki/Orthogonality_of_Solutions_to_the_Sturm-Liouville_Equation_with_Distinct_Eigenvalues) we need to make sure that the boundary conditions are satisfied , i.e \begin{equation} X(0)=X(\infty)=0 \quad (v) \end{equation} We use $(iii)$ and we conclude that $C_1=0$ and $C_2=1$ and \begin{eqnarray} E_n := \text{arg} \left. \left\{ U\left(-\frac{D k_1^2+k_1 k_2 v+k_2^2 x}{k_1 k_2^2},-\frac{D k_1+k_2 v}{k_2^2},\frac{D \sqrt{k_1^2}}{k_2^2}\right) = 0 \right| x\in [n, n+1) \right\} \end{eqnarray} where $n=0,1,2,\cdots$.

Below I took $(k_1,k_2,D,v)=(0.819434, 0.440462, 0.378185, 0.562809)$ and I ploted the functions $X_n(x) := \left. X(x) \right|_{E=E_n}$ for $n=0,\cdots,3$. We have:

Now equation $(iv)$ along with boundary conditions $(v)$ ensure orthogonality ,i.e \begin{equation} \int\limits_0^\infty \mu(x) \cdot X_n(x) \cdot X_m(x) dx = 0 \quad \mbox{if $n\neq m$} \end{equation}

Therefore the general solution to our original PDE now reads: \begin{equation} f(x,t) = \sum\limits_{n=0}^\infty C_n \cdot \exp(-E_n \cdot t) \cdot X_n(x) \end{equation} where the coefficients read: \begin{eqnarray} C_n := \frac{\int\limits_0^\infty f(x,0) \cdot \mu(x) \cdot X_n(x) dx}{ \int\limits_0^\infty [X_n(x)]^2 \cdot \mu(x) dx} \end{eqnarray} for $n=0,1,\cdots$.