Suggestions for preconditioners of implicit time discretized wave equations

[yy.wayne]

Hello 

I hope to precondition a system similar to step-23, with implicit time discretization (backward euler or Crank Nicolson). For example, the matrix has the form of M+dt^2K, where M is the mass matrix and K is the stiffness matrix.

I'm preconditioning it with Jacobi preconditioner, which is very efficient because the system is dominated by M and applies to matrix-free framework. However, when I increase the number of freedomes, iterations increase with just Jacobi preconditioner. Is there another recommended preconditioner or modification to simple jacobi preconditioners? I cannot find one since jacobi preconditioner is quite efficient for this problem, the lack of independency to h is its only drawback. 

[Wolfgang Bangerth]

Have you tried something like SSOR?

Best
WB

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[yy.wayne]

I tried by SSOR gives same behavior as jacobi. 

For mesh with 0, 1, and 2 times of refinement, it needs (jacobi: 20, 31, 67) and (SSOR: 10, 15, 31) iterations.

Besides, is SSOR not applicable in matrix-free framework?

[Wolfgang Bangerth]

On 10/4/24 01:38, 'yy.wayne' via deal.II User Group wrote:
> **

>
> I tried by SSOR gives same behavior as jacobi.
> For mesh with 0, 1, and 2 times of refinement, it needs (jacobi: 20, 31, 67)
> and (SSOR: 10, 15, 31) iterations.
> Besides, is SSOR not applicable in matrix-free framework?

Yes, SSOR needs access to matrix elements. If these iteration counts are too
large for you, then the only alternative you have is to go with multigrid-type
preconditioners, for which you can use matrix-free methods.

Best
W.

[yy.wayne]

Thank you Prof. Bangerth,

I've considered multigrid preconditioners as the iterations are relatively constant, but isn't it suits diffusion problems more,

as mentioned here (A former discussion)?

I assume multigrid works well for M+dt^2K type problem, but it's expensive for such 'easy' problem(when dt is small). The only drawback

for jacobi preconditioner is the iteration number scales with problem size ...

[Wolfgang Bangerth]

On 10/4/24 08:10, 'yy.wayne' via deal.II User Group wrote:
>
> I've considered multigrid preconditioners as the iterations are relatively
> constant, but isn't it suits diffusion problems more,
> as mentioned here (A former discussion

> <https://groups.google.com/g/dealii/c/det9e4HWGrk/m/q0oj-yQlBAAJ#:~:text=-%20lumping%20of%20mass%20matrices%0A-%20AMG/GMG%20for%20laplace-like%20operators%0A-%20ILU/Jacobi%20for%20mass%20matrix-like%20operators%0A-%20spectrally-equivalent%20matrices%20to%20approximate%20Schur%20complements%0A-%20and%20more>)?

Multigrid works well for diffusion problems. It doesn't work well for
advection problems.

But if your matrix really is M + dt^2 K (i.e., with dt^2 instead of dt), then
you must be solving a problem with second time derivatives -- say the wave
equation. That's not a diffusion problem either.

> I assume multigrid works well for M+dt^2K type problem, but it's expensive for
> such 'easy' problem(when dt is small). The only drawback
> for jacobi preconditioner is the iteration number scales with problem size ...

Well, that's the trade-off you have. Your simple method takes too many
iterations. You need to invest in a more complicated method that perhaps per
iteration is more expensive, but at least results in a smaller number of
iterations. The total time to solution is (time per iteration) * (number of
iterations). Your current approach results in (cheap) * (growing number)
whereas multigrid might be (expensive) * (constant). At some point, the latter
will be better than the former.

[yy.wayne]

Yes, I'm solving a transient wave equation. Should the behavior between M + dt*K and  M + dt^2*K

be similar when solving the linear system? I can only understand them as positive helmholtz equation,

the latter more dominated by mass matrix.

I agree with the trade-off between iteration cost vs iteration numbers. I plan to increase the dofs to 100 million

or even 1 billion, so multigrid (geometric or algebraic) might be the final option. I guess the growing iterations

are mainly contributed by the stiffness matrix?

[Wolfgang Bangerth]

> Yes, I'm solving a transient wave equation. Should the behavior
> between M + dt*K and M + dt^2*K be similar when solving the linear
> system? I can only understand them as positive helmholtz equation,
> the latter more dominated by mass matrix.

The condition number of K is h^{-2}, so if you choose dt=h, for the wave
equation you should generally obtain a matrix
M + dt^2 K
whose condition number is bounded as you refine the mesh. For the heat
equation, you get
M + dt K
which will (if you choose dt=h) have a condition number that grows like
h^{-1}.

As a consequence, I'm a *bit* surprised that your iterations grow by so
much, but it's been a long time since I've solved the wave equation and
so I don't recall whether what you see is expected or not.

Best
W.

[yy.wayne]

Sorry that I didn't state the equation clear.

I'm solving a time-dependent wave equation, and dt is the time step. For backward Euler

time integration, I get M+dt^2*K = (v,u) + (\delta T)^2*(\nabla v, \nabla u). The time step

dt does not scale with mesh size h, so the condition number is not bounded.

[Wolfgang Bangerth]

On 10/4/24 20:39, 'yy.wayne' via deal.II User Group wrote:
>
> I'm solving a time-dependent wave equation, and dt is the time step. For
> backward Euler
> time integration, I get M+dt^2*K = (v,u) + (\delta T)^2*(\nabla v, \nabla u).
> The time step
> dt does not scale with mesh size h, so the condition number is not bounded.

I see. It's of course true that for an implicit discretization you don't
*have* to choose dt proportional to h. You can then get into a situation where
indeed the condition number grows.

But from a practical perspective, for accuracy (not stability) reasons, you
will want to choose dt \approx h/c where c is the wave speed. You just won't
get accurate solutions otherwise.

[yy.wayne]

I didn't know that dt should be propotional to h/c for accuracy results. Thank you for 

pointing it out.  

[yy.wayne]

I decrease dt by half for each mesh refinement. Now the iterations remain 

constant (20, 21, and 23), same as you predicted. I guess in practical I should

keep dt and h constant (based on frequency and wave speed) for a acceptable

accuracy, while increase the dimension of computational space to test the

number of iterations?

[Wolfgang Bangerth]

On 10/4/24 21:44, 'yy.wayne' via deal.II User Group wrote:
> **

>
> I decrease dt by half for each mesh refinement. Now the iterations remain
> constant (20, 21, and 23), same as you predicted. I guess in practical I should
> keep dt and h constant (based on frequency and wave speed) for a acceptable
> accuracy, while increase the dimension of computational space to test the
> number of iterations?

I think you should do plenty of computational experiments to gain the
experience what works and what doesn't, and to assess the accuracy you get! :-)

[yy.wayne]

Thank you Prof. Bangerth for your time and very insightful suggestions. 

I'll definitely do more computational experiments and check literatures .

Best,

Wayne