finite difference scheme for nonlinear partial differential equations

Question

I have the following second order partial differential equation (PDE) on $[0,T] \times \mathbb{R},~ T >0 $: \begin{equation} \left(1 + \frac{1}{(1 + b f)^2}\right) \frac{\partial f}{\partial t} -\left(\frac{1}{(1 + b f)^2}\right) \frac{\partial f}{\partial x} = \frac{\partial^2 f}{\partial x^2} . \end{equation} supplied with following initial and boundary conditions: \begin{equation} f(x,t=0)=\begin{cases} 1, & x\leq 0\\ 0, & x > 0, ~ \text{the Heaviside unit step function}. \end{cases} \end{equation}

\begin{equation} \frac{\partial f}{\partial x} =0, ~\mbox{as}~|x| \rightarrow \infty, \end{equation}

Here $b$ is a non-negative constant.

My questions:

1. Does there exists any analytic solution to this PDE (e.g., some approximation solution).
2. Does there exists any finite difference scheme or any numerical scheme to solve this PDE.

P.S. I have some idea how to solve non-linear PDEs with constant coefficients for time derivative. Buy I have no guess how to start for stated PDE.

It will be great help, if someone suggest some idea. Thanks!

Answer 1

A little more details as requested (in private email)

Sorry I made a mistake above, of course you want $b \ll 1$. If that is the case you can use whats called "regular perturbation theory" (as opposed to singular perturbation theory which is what you would need if for example the $f_{xx}$ term was multiplied by a small parameter).

Make the following assumption $f(x,t) = f_0(x,t) + b f_1(x,t) + O(b^2)$

Now plug that in and collect terms corresponding to different powers of $b$. You get the following

$$2 \partial_t f_0 - \partial_x f_0 = \partial_x^2 f_0$$

And $f_0$ has the same BCs as you had for $f$ above. This problem is much easier to solve. Then and should give you an approximation when $b \ll 1$. You can even derive a correction by considering the equation you would get for the $f_1$ term. This equation will be forced (i.e. be non-homogenous) by the solution of $f_0$.

Have you looked at any finite difference schemes for parabolic PDEs? Unless you actually try to implement some scheme and run into specific problems, I don't know what you're expceting people on here to tell you about numericas for this thing.

Answer 2

This problem is quite tricky especially considering that it doesn't have a strict space boundary. I assume that you're not very familiar with numerical methods, thus I can give you an advice on where one would start the numerical treatment if it had a more strict space bound, for example $[0, T]\times[a,b], T>0$.

1) We discretize time domain using first order forward difference and space domain using 2nd order central difference:

$$ \left(1 + \frac{1}{(1 + bf_i^n)^2}\right) \frac{f_i^{n+1} - f_i^n}{\Delta t} = \frac{f_{i+1}^n - 2f_i^n + f_{i-1}^n}{\Delta x^2} +\left(\frac{1}{(1 + b f_i^n)^2}\right) \frac{f_{i+1}^n-f_{i-1}^n}{\Delta x} $$

where $f_i^n = f(x_0+i*\Delta x, t_0+n*\Delta t), 0 \leq i \leq H, H=\frac{a-b}{\Delta x}, 0 \leq n \leq N, N = \frac{T}{\Delta t}$. Which means that we discretize space into H+1 discrete intervals and time into N+1 discrete intervals correspondingly.

2) From there we write down an "update" formula for the value of the function at the next time step $n+1$ and obtain an explicit Euler scheme:

$$ f_i^{n+1}=f_i^n +{\Delta t}\left(c\frac{f_{i+1}^n - 2f_i^n + f_{i-1}^n}{\Delta x^2} +(1-c)\frac{f_{i+1}^n-f_{i-1}^n}{\Delta x}\right) $$ where $c=\frac{(1+bf_i^n)^2}{1+(1+bf_i^n)^2}$.

3) So here comes the first tricky part, since normally you get a stencil matrix with constant values of the form

$$ A= \left( \begin{array} -4c & c(2-\Delta x) + \Delta x & 0 & ... & 0 \\ c(2+\Delta x) - \Delta x & -4c & c(2-\Delta x) + \Delta x & ... & 0 \\ 0 & c(2+\Delta x) - \Delta x & -4c & c(2-\Delta x) + \Delta x & ... \\ ... & ... & ... & ... & ... \\ 0 & ... & 0 & c(2+\Delta x) - \Delta x & -4c \\ \end{array} \right) $$

with which its easy to solve a system of ODEs $F^{n+1}=F^n+\Delta t AF^n$, where $F^n=[f_0^n, f_1^n, ..., f_H^n]^T$ Here our $c$ depends on the function value $f(x_i,t_n)$ which makes the whole system of ODEs harder to compute.

4) At this step you would normally use the values of the function defined by the BC and the fact that the domain of our PDE is bounded on $[a,b]$. Thus you would explicitly compute the values of $F^n$.

Now, if to come back to the fact that $\frac{\partial f}{\partial x} =0, ~\mbox{as}~|x| \rightarrow \infty,$ it makes the whole numerical scheme kinda useless, because you don't even know where to start the computation from.

P.S.: Hope that helps. :)