periodic_box.prm fails

[geodynamics]

periodic_box.prm fails

[bangerth]

Isn't it correct, though, that in that case there matrix does and should have a null space? Do we need to impose constraints for the testcase to remove the nullspace?

[tjhei]

What we have done successfully in the past is to:

1. ignore the nullspace during assembly/solve with GMRES staying in the same subspace
2. remove the nullspace as a postprocess

[bangerth]

Good point. So what do we do here?

[tjhei]

Before we mark this as resolved, I have one more data point:

Take tests/depth_average_01.prm, increase global refinement by one and switch back to Q4 -> crash. This can be fixed by either switching to Q1 or by setting cheap solver steps = 0. Crazy.

[bangerth]

The way I read the discussion over on the deal.II testcase I created is that this is really a failure in Trilinos to compute things relative, not absolute. I think we need to do something like this in the Trilinos solver wrappers:

  const double scaling_factor = 1./rhs.l2_norm();
  rhs *= scaling_factor;
  solution *= scaling_factor;
  ... do as before...
  rhs /= scaling_factor;
  solution /= scaling_factor;

I'll have to try something along these lines in the coming week.

We could also fix this in ASPECT, but I think it should go into deal.II.

[tjhei]

@bangerth, all this won't solve our issue here though. If you switch to our CG implementation, you will get

The iterative (top left) solver in BlockSchurPreconditioner::vmult did not converge. It reported the following error:
Iterative method reported convergence failure in step 10000. The residual in the last step was 0.107138.

[bangerth]

Well, the change might fix depth_average, not periodic_box. I think these are separate issues.

I don't really know what to do here :-(

[tjhei]

Ideas:

1. figure out how to do diagonal strengthening or something for the coarse KLU.
2. give ML the constant mode for the displacement
3. fix the displacement to make the system uniquely solvable (somewhat complicated)

[bangerth]

I think you're going to have an extremely hard time to figure out the null space just from the domain and the description of what the boundary conditions are. It's sort of a pattern matching problem that humans are good at but that are difficult to put into rules :-(

[ian-r-rose]

For what it's worth, a PETSc build of Aspect runs this test without crashing.

[gassmoeller]

closed by #328