lucasmyers97 / maier-saupe-lc-hydrodynamics

Work at University of Minnesota using finite element methods to simulate hydrodynamics of liquid crystals with a Maier-Saupe field theory free energy
4 stars 0 forks source link

LagrangeMultiplierEfficient giving incorrect values for Q-vectors with zero entries #18

Closed lucasmyers97 closed 2 years ago

lucasmyers97 commented 2 years ago

For the current values found in efficient_lagrange_multiplier, the Jacobian output for LagrangeMultiplierEfficient is different than the Jacobian for LagrangeMultiplier.

The Q-vector given has zeros in the 3rd and 5th entry. The Jacobian values disagree in the (3, 3), (5, 3), (3, 5), and (5, 5) entries, so I suspect it has something to do with having zero entries in those places.

lucasmyers97 commented 2 years ago

This is related to the issue raised in deal.II located here. It turns out that all of the eigenvector algorithms have shortcuts which basically decouple the output from the input (i.e. they just hardcode a value for some of the eigenvalues or eigenvectors, and so the Jacobians get messed up). For completely diagonal matrices there is a method which perturbs them to give reasonably accurate values based on Ogden hyperelastic materials or neohookean solids or something. I am going to make a pull request for a patch which allows one to force the perturbation, but I suspect that one can fix it by just getting rid of a line in the Jacobi method. I'm not sure it's so easily fixable using the QL implicit shifts method.

lucasmyers97 commented 2 years ago

The deal.II developers suggested passing in a flag to force a perturbation on the matrix. I wrote this and tried it, but the Jacobian had a much higher error than the case where I just made all zero entries 1e-10 or something like that.

In any case, this no longer matters because I have written a program which calculates the derivatives of the eigenvectors analytically, and thus any errors with automatic differentiation are null. Additionally, the automatic differentiation scheme gives infinities (or null, I don't remember) for diagonalization with repeated eigenvalues: the reason is because the mapping from matrices to their eigenvectors is necessarily discontinuous at those points. See math Stack Exchange post here for details. When you write an analytic expression, the discontinuity is removable when you map back to the original frame of the Q-tensor. Hence, we deal with degenerate eigenvalues by just using one formula or the other.