1
$\begingroup$

I recently asked a question pertaining to the appliciation of Jacobi's method to a semilinear elliptic PDE (Poisson's equation)

$$ \nabla^2u = -\rho~e^{-u} $$

A more efficient method like the Bi conjugate gradient stabilised method was recommended. I have tested this method out and it is indeed much faster. But I am unsure of what the matrix representation of a semilinear system would look like. For an ordinary linear PDE like

$$ \nabla^2u=-\rho $$

it looks like $$ \frac{1}{h^2}\left( \begin{array}{ccc} 2 & -1 & & & \\ -1 & 2 & -1 & & \\ & -1 & 2 & -1 & \\ & & -1 & 2 & -1 \\ & & & -1 & 2 \\ \end{array} \right)\left( \begin{array}{c} u_1 \\ u_2 \\ u_3 \\ u_4 \\ u_5\end{array} \right) = \left( \begin{array}{c} \rho_1+g \\ \rho_2 \\ \rho_3 \\ \rho_4 \\ \rho_5+g\end{array} \right) $$ where $g$ is the Dirichlet boundary condition.

My question: What would the corresponding matrix representation of the set of simultaneous equations for the semilinear case look like? I'm guessing something like

$$ \frac{1}{h^2}\left( \begin{array}{ccc} 2 & -1 & & & \\ -1 & 2 & -1 & & \\ & -1 & 2 & -1 & \\ & & -1 & 2 & -1 \\ & & & -1 & 2 \\ \end{array} \right)\left( \begin{array}{c} u_1 \\ u_2 \\ u_3 \\ u_4 \\ u_5\end{array} \right) = \left( \begin{array}{c} \rho_1e^{-u_1}+g \\ \rho_2e^{-u_2} \\ \rho_3e^{-u_3} \\ \rho_4e^{-u_4} \\ \rho_5e^{-u_5}+g\end{array} \right) $$

But this doesn't leave me with all $u$ values in a single vector.

Would it make sense to do something like:

1) Solve the linear case $$ \nabla^2u = -\rho $$ 2) Use the resultant $u$ to construct a new linear case $$ \nabla^2u_i = -\rho~C $$ where $$ C = e^{-u_{old}} $$ 3) Repeat step 2 until self-consistency is reached.

$\endgroup$
2
$\begingroup$

Since you don't end up with a linear system, you have to use a nonlinear root finder. The most common choice is Newton's method, because it converges quite quickly (if your initial guess is good enough). You've probably seen Newton's method for a single variable in a single equation. There is an analogous method for systems of equations.

Given a vector function $f(\mathbf{u}^*) = 0$, and an initial guess $\mathbf{u}^0$, we can approximate the root $\mathbf{u}^*$ by solving the sequence of linear systems $$ J(\mathbf{u}^k)\Delta\mathbf{u}^k = -f(\mathbf{u}^k),$$ where $\Delta \mathbf{u}^k = \mathbf{u}^{k+1} - \mathbf{u}^{k}$ and $J(\mathbf{u}^k)$, called the Jacobian, is defined as

$$ J(\mathbf{u}^k) = \begin{pmatrix} \frac{\partial f_1}{\partial u_1}(\mathbf{u}^k) & \cdots & \frac{\partial f_1}{\partial u_n}(\mathbf{u}^k)\\ \vdots & \ddots & \vdots\\ \frac{\partial f_n}{\partial u_1}(\mathbf{u}^k) & \cdots &\frac{\partial f_n}{\partial u_n}(\mathbf{u}^k)\end{pmatrix}.$$

In your case $f(\mathbf{u})$ is given by

$$ f(\mathbf{u}) = \frac{1}{h^2}A\mathbf{u} - \mathbf{b}(\mathbf{u})$$ where $A$ is your finite difference matrix and $\mathbf{b}(\mathbf{u})$ is your right hand side. The Jacobian is given by

$$ J(\mathbf{u}) = \frac{1}{h^2}A + \begin{pmatrix} \rho_1 e^{-u_1}\\ \mathbf{0} & \ddots & \mathbf{0} \\ & & \rho_n e^{-u_n}\end{pmatrix}.$$

Your initial guess could come from solving the homogeneous system for example. Then given $\mathbf{u}^0$, we can evaluate $J(\mathbf{u}^0)$ and $f(\mathbf{u}^0)$ and solve for $\Delta \mathbf{u}^0$ using $$ J(\mathbf{u}^0)\Delta\mathbf{u}^0 = -f(\mathbf{u}^0).$$

We then find $\mathbf{u}^1$ from the definition of $\Delta \mathbf{u}^0$. We can repeat the process to find $\mathbf{u}^2, \mathbf{u}^3, \cdots$, stopping once $|\Delta \mathbf{u}^k|$ is below some tolerance.

$\endgroup$
1
$\begingroup$

@H.H. already gives a good answer on how to implement things for the Newton iteration. If you want to get an overview of different approaches to address nonlinear problems, take a look at lectures 31.5 and following here: http://www.math.colostate.edu/~bangerth/videos.html .

$\endgroup$

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.