Consider the 1D Poisson equation$$\nabla^2 u = f.$$Using finite difference method on cell corner data and a uniform grid with ghost points, I think we can write the system of equations with Neumann BCs as:$$Au = f-Au_{BC},$$where$$\begin{aligned}A &= \frac{1}{\Delta x^2} \left[\begin{array}{ccccccccc}-2    & 2         &           &           &           &           &     \\1     & -2        & 1         &           &           &           &     \\      &           & \ddots    & \ddots    & \ddots    &           &     \\      &           &           &           & 1         & -2        & 1   \\      &           &           &           &           &  2        & -2  \\\end{array} \right]\\Ax_{BC} &= \frac{1}{\Delta x^2} \left[\begin{array}{ccccccccc}    2 \hat{n} \theta \Delta x   \\    0  \\       \\\vdots \\       \\    0  \\    2 \hat{n} \theta \Delta x   \\\end{array} \right]\end{aligned}$$And $\hat{n},\theta$ is the outward facing normal and the prescribed slope of the derivative at the boundaries.Notably, the system is singular, which can be addressed by removing the mean of the RHS.QuestionI know that $A$ must be symmetric, positive definite. I've done some tests without multiplying the first and last rows by 0.5 and they seem to work fine, so my question is:Must the first and last rows of the left and right hand side be multiplied by 0.5? In other words, are there cases where it will not work without multiplication of 0.5?NotesI imagine that both sides of the equation would be balanced the same with and without this multiplication.As a final note, I know that the code / algorithm may be more clear with the multiplication, since the equation explicitly satisfies symmetry, but I'm more interested in whether the multiplication is necessary or not.Any help is greatly appreciated.
CG question: is symmetry always necessary?
finite difference;numerical analysis;conjugate gradient;symmetry
The CG method works by producing in every iteration a vector $x^k\in\mathbb{R}^n$ that solves the minimization problem (assuming $x^0=0$ for simplicity)$$\min_{x\in \mathcal{K}_k} \frac12(x-x^*)^TA(x-x^*),\tag{1}\label{cg1}$$where $x^*=A^{-1}b$ and $\mathcal{K}_k$ is a Krylov space of (apart from exceptional circumstances) dimension $k$. (The CG algorithm is just a clever way of computing this minimizer efficiently by reusing as much information from the previous iteration as possible.)So, for $k=n$ (at the latest), $x^n$ minimizes the function $f(x) = \frac12(x-x^*)^TA(x-x^*)$ over $\mathbb{R}^n$, which means that the gradient has to vanish. After expanding the product and using the product rule, this yields$$0 = \nabla f(x^n) = \frac12(A+A^T)x^n - \frac12(A+A^T)x^*.\tag{2}\label{cg2}$$If $A$ is symmetric ($A^T=A$) and $Ax^*=b$, this of course implies that $Ax^n=b$; otherwise you only get a solution to $(A+A^T)x = (A+A^T)A^{-1}b$. In your case, $A$ is nearly symmetric and fairly well-conditioned, so the difference between the solutions to these two equations is probably so small that you don't see it.Another question is whether the numerical algorithm breaks down if $A$ is not symmetric. Looking at the algorithm, this can only happen if $x^TA x =0$ for some $x\neq 0$ (which of course is impossible for positive definite matrices by definition). Now note that you can write any matrix $A$ as$$A = \frac12 (A+A^T) + \frac12(A-A^T) =: A_{\text{sym}} + A_{\text{skew}},$$where $A_{\text{sym}}$ is the symmetric part and $A_{\text{skew}}$ the skew-symmetric part of $A$. Furthermore, $x^TA_{\text{skew}}x=0$ and hence $x^TAx = x^TA_{\text{sym}}x$ for all $x\in\mathbb{R}^n$. So the algorithm works as long as the symmetric part $A_{\text{sym}}$ of $A$ is positive definite. Conveniently, this implies that the minimization problem \eqref{cg1} is strictly convex and hence has a unique minimizer, which is the only solution to \eqref{cg2}.TL;DR: CG works for nonsymmetric matrices as long as the symmetric part is positive definite, and will give a reasonable approximation to the solution if the skew-symmetric part is negligible.