Solver for unsteady N-S equation

General Discussion
1

Hi everyone,
This post is out of such reason: I am not familiar with solvers/preconditioners and it may take some time for me to harness them.
As the mesh scale grows, the computation gets slower. Are there easy adjustments to -pc_type lu -ksp_type preonly to improve the efficiency? Or is it a must to implement a specific preconditioner for a certain problem?

In some posts the mAL preconditioner is mentioned, is mAL also okay with unsteady problem? If so, maybe I can mimic the example natural-convection-fieldsplit-2d-PETSc.edp.

The 2D mesh could be up to 3000*3000(square) or even more with [P2,P2,P1] vectorial elements, will 128 cores be sufficient for such problem scale?

Best regards,
H.Weng

2

If it’s a 2D problem, something not too fancy may work out of the box. You could simply try -ksp_type gmres -pc_type asm -sub_pc_type lu.

3

Thanks for your quick reply,
The solver is now set(TS, sparams = "-ksp_type gmres -pc_type asm -sub_pc_type lu -ksp_view_final_residual");.
Unfortunately I got
TSStep has failed due to DIVERGED_STEP_REJECTED
after first step.

4

Does the linear solver converge?

5

图片
图片765×555 59.4 KB

The terminal is like this, the first line is the computation for the initial field. And then the TS routine for N-S and temperature, actually the thermal equation called another TSSolve in the TSSolve routine of N-S equation.
Should I upload the code?

6

Ah, try to add -sub_pc_factor_mat_solver_type mumps.

7

Now it is set(TS, sparams = "-ksp_type gmres -pc_type asm -sub_pc_type lu -ksp_view_final_residual -sub_pc_factor_mat_solver_type mumps");
and it won’t diverge at first step, but after several steps it will diverge.
I tried to use a smaller dt, and dt is finally turn to 3.8147e-08.

0 TS dt 0.01 time 0.
N-S equation solved, solve the thermal equation,
0 TS dt 0. time 0.
thermal equation solved, solve the next step
KSP final norm of residual 8.16106e-09
KSP final norm of residual 3.86178e-14
KSP final norm of residual 5.17033e-09
KSP final norm of residual 1.71385e-14
1 TS dt 0.01 time 0.01
N-S equation solved, solve the thermal equation,
0 TS dt 0.01 time 0.
KSP final norm of residual 7.9199e-32
1 TS dt 0.01 time 0.01
thermal equation solved, solve the next step
KSP final norm of residual 2.52215e-08
KSP final norm of residual 1.5837e-13
2 TS dt 0.01 time 0.02
N-S equation solved, solve the thermal equation,
0 TS dt 0.01 time 0.
KSP final norm of residual 7.40135e-25
1 TS dt 0.01 time 0.01
thermal equation solved, solve the next step
KSP final norm of residual 3.84135e-08
KSP final norm of residual 4.43181e-13
3 TS dt 0.01 time 0.03
N-S equation solved, solve the thermal equation,
0 TS dt 0.01 time 0.
KSP final norm of residual 1.12881e-24
1 TS dt 0.01 time 0.01
thermal equation solved, solve the next step
KSP final norm of residual 7.68258e-09
KSP final norm of residual 8.80407e-14
4 TS dt 0.01 time 0.04
N-S equation solved, solve the thermal equation,
0 TS dt 0.01 time 0.
KSP final norm of residual 1.34192e-24
1 TS dt 0.01 time 0.01
thermal equation solved, solve the next step
KSP final norm of residual 3.08425e-13
KSP final norm of residual 3.09537e-16
KSP final norm of residual 3.93751e-14
KSP final norm of residual 3.92473e-15
KSP final norm of residual 1.00479e-14
KSP final norm of residual 5.35912e-14
KSP final norm of residual 5.24062e-12
KSP final norm of residual 2.05859e-11
KSP final norm of residual 3.92306e-11
KSP final norm of residual 8.81865e-12
KSP final norm of residual 1.01131e-12
KSP final norm of residual 1.99914e-14
KSP final norm of residual 1.60835e-19
5 TS dt 3.8147e-08 time 0.04
N-S equation solved, solve the thermal equation,
0 TS dt 3.8147e-08 time 0.
KSP final norm of residual 1.46592e-24
1 TS dt 3.8147e-08 time 3.8147e-08
thermal equation solved, solve the next step
KSP final norm of residual 3.85305e-14
KSP final norm of residual 7.55988e-19
6 TS dt 3.8147e-08 time 0.0400001
N-S equation solved, solve the thermal equation,
0 TS dt 3.8147e-08 time 0.
KSP final norm of residual 1.46488e-24
1 TS dt 3.8147e-08 time 3.8147e-08
thermal equation solved, solve the next step
KSP final norm of residual 5.22306e-14
KSP final norm of residual 9.46825e-19
7 TS dt 3.8147e-08 time 0.0400001

I also tried to use such solver for my coupled script, N-S equation and thermal equation are solved simultaneously. It end up with the same divergence problem.

8

Does the linear solver converge?

9

I added -ksp_converged_reason and got

Linear solve did not converge due to DIVERGED_ITS iterations 10000
KSP final norm of residual 1.65849e-05
Linear solve did not converge due to DIVERGED_ITS iterations 10000
KSP final norm of residual 1.65849e-05
Linear solve did not converge due to DIVERGED_ITS iterations 10000
KSP final norm of residual 1.65849e-05
10

Well, then yes, you need to investigate fancier preconditioners.

11

Is mAL a candidate? In my problem, the block structure of N-S equation is similar, with a buoyancy term and a shifted mass matrix added.

12

It is worth the try.

13

Hi prj,
I just modified my code to 3d by buildlayers, with two periodic pairs [[5, x, y], [6, x, y], [2, y, z], [4, y, z]].
Since I have not sort out the preconditioners, the solvers are the same default one -pc_type lu -ksp_type preonly.
For a 3d case (16*16*16), I think it would have the same speed approximately as 64*64 in 2d. But it is much slower actually :slightly_frowning_face:.
Is it normal for a 3d case to be much slower than its 2d counterpart?

14

Yes, direct solvers in 3D have a much higher complexity.

15

Hi, @prj ,
I tried to modify my code, start from 2d problem and the direct solver is unchanged.
Noticed that mAL have a different stabilization term, the first thing is to change the varf.
So eps*p*q is replaced by the grad-div stabilization gammaAL * div(u) * div(v) simply.
Nevertheless, the TSSolve/SNES routine became hard to converge after several steps.

0 SNES Function norm 6.931563109443e-11 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 3.74686e-20
    1 SNES Function norm 6.928264440069e-15 
8 TS dt 0.025 time 0.2
N-S equation is solved, try to solve heat equation

0 TS dt 0.025 time 0.
    0 SNES Function norm 2.301090427222e-11 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 5.12756e-27
    1 SNES Function norm 3.463803651916e-18 
1 TS dt 0.025 time 0.025
Heat equation is solved, try to slove next time step

    0 SNES Function norm 2.403739778848e-11 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 2.23604e-21
    1 SNES Function norm 9.080687237421e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 7.89078e-22
    2 SNES Function norm 7.274418214209e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 1.49822e-21
    3 SNES Function norm 7.318620083181e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 2.14804e-21
    4 SNES Function norm 7.267641604000e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 2.05017e-20
    5 SNES Function norm 7.307178962458e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 2.90308e-21
    6 SNES Function norm 7.287125580495e-15 
 

...
...
49 SNES Function norm 7.242230569222e-15 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 4.54047e-21
   50 SNES Function norm 7.255060411073e-15 
    0 SNES Function norm 1.602268542663e-12 
Linear solve converged due to CONVERGED_ITS iterations 1
KSP final norm of residual 1.61509e-21
    1 SNES Function norm 7.257310637315e-15

The first try seems failed, is such stabilization scheme not fit for direct solver?

16

You are still using -ksp_type preonly, which does not make sense with an inexact preconditioner.

17

Sorry that I am a bit confused.
Yes, the solver is still -pc_type lu -ksp_type preonly.
Just take a small step, I want to see the difference if using a different stabilization term.
There are only changes to varf, so I think it is not relevant to preconditioner issues.
Do you mean the grad-div stabilization itself is a preconditioner, in a sense?

18

Why would you use grad-div stabilisation with PCLU?

19

Because I think LU is the most robust one. Okay, so PCLU is not a good choice. I will try flexible GMRES.

20

Hi @prj,
It has been a while for this post. So far, the coupled NSE system is still too complicated to precondition for me. I would like to use a projection method to solve NS equation in parallel.

The projection method solves two equation, the LHSs are (\nabla p,\nabla q) and \alpha(u,v)+\nu (\nabla u,\nabla v). In F. Hecht’s sequential code, the former is solved by CG, and the latter is solved by UMFPACK.

From the Youtube tutorial, Poisson equation maybe easy to solve:
varf(p,q)=int2d(Th)(grad(p)'*grad(q)).

When it comes to a system like
varf(u,v)=int2d(Th)(alpha*u*v+nu*grad(u)'*grad(v)),
what is the proper way to solve it efficiently?

And, if my time step is fixed, it seems that the matrix is unchanged, will FreeFEM or PETSc reuse the inversion of the matrix? If so, I think only matrix multiplication is needed once the inversion is computed.