I am trying to factorize the following Laplacian matrix in terms of $ D^TD$, D is the first derivative matrix. The tridiagonal form of the secon derivative matrix using Neumann boundary condition is given by,

$\begin{bmatrix} -1 & 1 & 0 & 0 & 0 \\ 1 & -2 & 1 & 0 &0\\ 0 & 1 & -2 &1&0\\0&0&1&-2&1\\0&0&0&1&-1\end{bmatrix}$.

In order to write the above in terms of the gradient( first derivative) operator, I tried factorizing the Laplacian using Cholesky factorization. But, I couldn't succeed in factorizing.

Any suggestions on how to find a decomposition of the form $D^TD$?

  • 1
    $\begingroup$ naïve Cholesky won't work, because your Laplacian matrix is not positive definite. there are many square roots (=choices of $D$) to select from. to recover the "upwind" discretization of the first order derivative operator, you'll want to use the LDL decomposition. see here for a more thorough discussion. $\endgroup$ – GoHokies Nov 14 '18 at 20:34
  • 3
    $\begingroup$ Isn't it just D=[-1 +1 0 0 0; 0 -1 +1 0 0; 0 0 -1 +1 0; 0 0 0 -1 +1]? I think you can also extend to higher dimensions, using the (signed) edge-to-vertex adjacency matrix. (Essentially a discrete version of the gradient operator). $\endgroup$ – rchilton1980 Nov 14 '18 at 20:49
  • $\begingroup$ @rchilton1980 yes, and that's exactly what the LDL decomposition produces. $\endgroup$ – GoHokies Nov 14 '18 at 21:17
  • $\begingroup$ But the Laplace matrix is not just $D^T D$. That's because the Laplacian operator is the divergence of the gradient, so the two operators are transposes of each other. But the Laplace matrix is not of this form. $\endgroup$ – Wolfgang Bangerth Nov 15 '18 at 4:26
  • 1
    $\begingroup$ Even if the inverse does not exist, take a pseudo-inverse. It will also not be sparse. $\endgroup$ – Wolfgang Bangerth Nov 15 '18 at 14:17

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Browse other questions tagged or ask your own question.