# Finite differences and Neumann boundary conditions

I am dealing with a highly nonlinear system of two PDEs. I already have a code to solve the system in case of Dirichlet boundary conditions. The explicit system is:

$$\begin{eqnarray*} \partial_{t}u & = & D_{11}\partial_{x}^{2}u+\partial_{x}D_{11}\cdot\partial_{x}u+D_{12}\partial_{x}^{2}v+\partial_{x}D_{12}\cdot\partial_{x}v\\ \partial_{t}v & = & D_{21}\partial_{x}^{2}u+\partial_{x}D_{21}\cdot\partial_{x}u+D_{22}\partial_{x}^{2}v+\partial_{x}D_{22}\cdot\partial_{x}v \end{eqnarray*}$$

with boundary conditions:

$$u(0,t)=0 \ (\partial_x v)(0,t) = 0$$

And the discretized form of the scheme are:

$$\begin{eqnarray*} -D_{21}(x_{j})u_{n+1,j-1}-D_{22}(x_{j})v_{n+1,j-1}\\ {}[\mu+2D_{11}(x_{j})]u_{n+1,j}+2D_{12}(x_{j})v_{n+1,j}\\ -D_{11}(x_{j})u_{n+1,j+1}-D_{12}(x_{j})v_{n+1,j+1} & = & \frac{1}{4}[4\mu-D_{11}(x_{j+1})+D_{11}(x_{j-1})]u_{n,j}\\ & & -[D_{12}(x_{j+1})-D_{12}(x_{j})]v_{n,j}\\ & & +\frac{1}{4}[D_{11}(x_{j+1})-D_{11}(x_{j-1})]u_{n,j+1}\\ & & +\frac{1}{4}[D_{12}(x_{j+1})-D_{12}(x_{j-1})]v_{n,j+1} \end{eqnarray*}$$

$$\begin{eqnarray*} -D_{21}(x_{j})u_{n+1,j-1}-D_{22}(x_{j})v_{n+1,j-1}\\ 2D_{21}(x_{j})u_{n+1,j}+[\mu+2D_{22}(x_{j})]v_{n+1,j}\\ -D_{21}(x_{j})u_{n+1,j+1}-D_{22}(x_{j})v_{n+1,j+1} & = & \frac{1}{4}[4\mu-D_{21}(x_{j+1})+D_{21}(x_{j-1})]u_{n,j}\\ & & -[D_{22}(x_{j+1})-D_{22}(x_{j})]v_{n,j}\\ & & +\frac{1}{4}[D_{21}(x_{j+1})-D_{21}(x_{j-1})]u_{n,j+1}\\ & & +\frac{1}{4}[D_{22}(x_{j+1})-D_{22}(x_{j-1})]v_{n,j+1} \end{eqnarray*}$$

That is, I am solving a linear system of the form:

$$A\boldsymbol{x}=\boldsymbol{b}$$

where $A$ has $2J\times 2J$ entries, according to the structure:

$$\left(\begin{array}{cccccccccc} \mbox{1st equation}\rightarrow & | & c_{u}^{j} & c_{v}^{j} & c_{u}^{j+1} & c_{v}^{j+1} & 0 & & ... & 0\\ \mbox{2nd equation}\rightarrow & | & d_{u}^{j} & d_{v}^{j} & d_{u}^{j+1} & d_{v}^{j+1} & 0 & & ... & 0\\ \vdots & | & . & . & . & . & . & . & . & .\\ \mbox{1st}\rightarrow & | & 0 & 0 & c_{u}^{j-1} & c_{v}^{j-1} & c_{u}^{j} & c_{v}^{j} & c_{u}^{j+1} & c_{v}^{j+1}\\ \mbox{2nd }\rightarrow & | & 0 & 0 & d_{u}^{j-1} & d_{v}^{j-1} & d_{u}^{j} & d_{v}^{j} & d_{u}^{j+1} & d_{v}^{j+1}\vdots \end{array}\right)$$

where the $c$'s and $d$'s are the the coefficients corresponding to the scheme above.

Now, I am trying to implement the Neumann boundary condition (the one on $v$), using the ghost point method as discussed, for instance, here for a simpler equation.

Thus, I would like to change my code in order to add some extra terms to the array, but I can't work out what to put in them. Also, in the respective components of the matrix $A$. Because when I write the condition:

$$\frac{u_1-u_{-1}}{(\Delta x)^2}=0$$ (and same for $v$)

I need to eliminate $u_{-1}$ and $v_{-1}$, using also the third and fourth equations I wrote in this post. If I inverted for $u_{-1}$ and $v_{-1}$ (assuming I know how to), I do not see how I could write down a code which works for any kind of scheme I have, without putting the explicit inversion in the matrix.

• Your question does not state which part of the implementation is unclear, what other references you have considered, what you have already tried and didn't work, etc. A question "how can I do it" is too unspecific to elicit any answers. Please rewrite your question to be more specific. Aug 2 '14 at 11:05
• @WolfgangBangerth You're right, the fact is that I am quite stuck before starting, I tried another way to include the BC before considering the ghost point, but that would be off-topic. Anyway, I edited the question a little, I hope this gives a clearer idea of what I am struggling with. Aug 2 '14 at 18:20
• Have you implemented Neumann Boundary Conditions for a simpler PDE, namely Laplace Equation? If the question is regarding the BC you should try them in a easier scenario first. Aug 2 '14 at 18:59
• How is your neumann condition stated? Please write out the entire condition for the equation.
– Paul
Aug 3 '14 at 14:04
• Have you seen this post? scicomp.stackexchange.com/questions/5355/…
– Paul
Aug 3 '14 at 15:10