accurate enough (see below) initial guess for $U^{n,i}$, the solution of (\[E32\]). Then the iterates of Newton’s method (called the [*outer*]{} iterates for reasons that will become clear presently) $U_j^{n,i}$, $j= 1,2,\ldots$ ($U_j^{n,i}$ approximates $U^{n,i}$) satisfy the $2\times 2$ block linear system in $S_h\times S_h$, $$\begin{aligned}
&&\left[\begin{array}{cc}
I+k a_{11} J(U_j^{n,1}) & k a_{12}J(U_j^{n,2})\\
k a_{21}J(U_j^{n,1}) & I+k a_{22} J(U_j^{n,2})
\end{array}\right]\, \left[\begin{array}{c} U_{j+1}^{n,1}\\  U_{j+1}^{n,2} \end{array} \right]=
\left[\begin{array}{c} U^n\\  U^n \end{array} \right]\label{E34}\\
&&-k
\left[\begin{array}{cc}
a_{11} &  a_{12}\\
a_{21} &  a_{22}
\end{array}\right]\, \left[\begin{array}{c} B(U_{j}^{n,1})\\ B( U_{j}^{n,2}) \end{array} \right],\nonumber\end{aligned}$$

where, for $\psi, \phi$ in $S_{h}$ $$J(\phi)\psi=J_1(\phi)\psi+J_2(\phi)\psi,$$ $$J_1(\phi)\psi=-2B(\phi,\psi)-\Theta_1\psi,$$ and $$J_2(\phi)\psi=-\Theta_2\psi.$$

The equations (\[E34\]) represent a $2N\times 2N$ linear system for the coefficients of the new Newton iterates $U_{j+1}^{n,i}$, $i=1,2$, for each $j$, with respect to a basis of $S_{h}$. The two operator equations in (\[E34\]) are uncoupled as follows: We evaluate the entries of the matrix in the left-hand side of (\[E34\]) at a point $U^{\ast}\in S_h$, defined by $$\label{E35}
U^{\ast}=\frac{1}{2}(U_0^{n,1}+U_0^{n,2}),$$ (which makes the operators in the entries of this matrix independent of $j$ and allows them to commute with each other). We may then write (\[E34\]) equivalently as $$\label{E36}
\begin{split}
&\left[\begin{array}{cc}
I+k a_{11} J_1(U^{\ast}) & k a_{12}J_1(U^{\ast})\\
k a_{21}J_1(U^{\ast}) & I+k a_{22} J_1(U^{\ast})
\end{array}\right]\, 