Dealing with "almost symmetric" matrices

Linear Algebra
1

I’m working on a linear time-evolution problem with a symmetric differential operator. By chance, I have discovered that the matrix assembled by assemble_system, denote it A, is not symmetric. Inquiring further however, I have discovered that A is “almost” symmetric, meaning that the entries of the difference A-A.T are very small compared to the entries of A itself (the max(abs(A-A.T)) ca. 7e-15 and max(abs(A)) ca. 40). So it appears that A is a symmetric matrix plus some numerical distortion.

Next, playing around with simple variational forms I have discovered that “almost symmetric” matrices are in fact quite common output of assemble_system. So I assume there must be some good practices to deal with such matrices.

I would like to find out a bit about such good practices. To start, I would like to ask the following questions:

  1. Does “almost symmetry” have a negative impact on the behavior of symmetric solvers? For example, Conjugate Gradient is known to suffer from the loss of conjugacy of subsequent search directions which may result from the accumulation of numerical error. Isn’t even a small distortion to the matrix symmetry an issue likely to amplify such undesired effects and so obstruct the convergence of Conjugate Gradient?

  2. Are there any tricks in FEniCS that help make the assembled matrix symmetric? I have been trying to obtain symmetry by changing the epsilon parameter passed to assemble_system via the argument form_compiler_parameter, but with no success. I have also tried to clip very small entries of the already assembled matrix to zero (since such elements are likely to be just distorted zeros), but the resulting matrix also was non-symmetric.

  3. Would it be a good idea to replace the “almost symmetric” matrix A with the exactly symmetric matrix B=0.5*(A+A.T) and call solve on B rather than A? Would there be any ill side effects of such trick which I should beware of?

If anyone feels wizard enough to shed some light on these issues, please do not hesitate to answer.

FEniCS version: 2019.1.0
Installation method: Docker image

I don’t attach an MWE since I find it not necessary in the context of this topic. If an MWE will turn out helpful in the run of discussion, I will post one.

2

What you mean to say is that A is symmetric to machine precision.

Your worries are unnecessary. Error incurred by iterative methods is much larger than machine precision (unless you set extremely strict tolerances of course).

Use a symbolic algebra package if you really need to represent exact symmetry, and not a discretisation scheme.

1 Like
3

So the good practice is just to treat matrices as they were exactly symmetric (and choose reasonable solver tolerances). This is the answer to my inquiries - thanks!

For completeness, I have tested the convergence of Conjugate Gradient for an almost symmetric matrix A and the exactly symmetric matrix B=0.5(A+A.T). The solver has converged within almost the same number of iterations and the difference between the resulting solutions has occurred very small. This confirms what you say.

2 Likes
4

Hi,

I am adding a question to this because this is very much related to what I want to ask.

I am also solving a linear time dependant problem (Acoustic Wave Propagation) using an own solver with Conjugate Gradient. I need to assemble mass, stiffness and damping matrices seperately and thus couldn’t use assemble_system. I used explicit calls to assemble(M), assemble©, assemble(K) and noticed that they are not completely symmetric. Also, the difference in values are not to machine precision but of significant ranges to other values. How can I resolve this ?

5

Please add a minimal working example illustrating this behavior.