deal.II and complex numbers

[Lucas P]

Hello,

I'm considering starting to use deal.ii with a strong emphasis on complex numbers equations. I read in the FAQ that it is recommended to split the problem into real and imaginary parts. I am not fully convinced by the argument given on the scalar product (which can be easily redefined, I think), for instance in FENICSx I directly use complex numbers.

I saw in the tutorial 58 that it's still possible not to split. On the Github issue tracker, an user bumped into errors when trying to use the monolithic way. But since then I observe that additional work was done to ensure the compatibility with `std::complex<double>` (see for instance this and that).

So my question is: how mature is now deal.ii with `std::complex`? Are there other arguments that the one on scalar product to justify the splitting (or maybe I missed why this argument is so important)?

Best,

Lucas

PS: I am surprised I did not find any other post dealing with this question. If this is a duplicate, please excuse me.

[Wolfgang Bangerth]

On 1/6/25 08:20, Lucas P wrote:
>
> I'm considering starting to use deal.ii with a strong emphasis on complex
> numbers equations. I read in the FAQ

> <https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#can-i-solve-problems-over-complex-numbers> that it is recommended to split the problem into real and imaginary parts. I am not fully convinced by the argument given on the scalar product (which can be easily redefined, I think), for instance in FENICSx I directly use complex numbers.

> I saw in the tutorial 58

> <https://dealii.org/current/doxygen/deal.II/step_58.html> that it's still possible not to split. On the Github issue tracker <https://github.com/dealii/dealii/issues/16364>, an user bumped into errors when trying to use the monolithic way. But since then I observe that additional work was done to ensure the compatibility with `std::complex<double>` (see for instance this <https://github.com/dealii/dealii/pull/16717>and that <https://github.com/dealii/dealii/pull/16754>).

> So my question is: how mature is now deal.ii with `std::complex`? Are there
> other arguments that the one on scalar product to justify the splitting (or
> maybe I missed why this argument is so important)?

Lucas:
I suspect that "almost everything" works, or could be made to work with
relatively few changes, with the exception of iterative solvers. I don't think
that any of the iterative solvers we have take into account complex numbers
and how one needs to/must not take the complex conjugate in any of the
operations involved. That's why step-58 works: It just bypasses the issue by
using a direct solver.

If you wanted to give it a try to work in a monolithic formulation (i.e.,
using complex numbers rather than real/imaginary parts separately), I would
say go ahead and see where you get stuck, then report back. Most things should
be relatively easy to fix, but prepared for the fact that you might have to
spend significant time in the iterative solvers figuring what a given solver
*should* be doing in the complex case, and then *auditing* what it actually
does. Making the necessary fixes is likely going to be quick; it's knowing
where something needs to be fixed.

Best
W.


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

[Lucas P]

Thanks Wolfgang for your answer. But I'm surprised to see that some work should be done regarding the iterative solvers. I would expect that deal.ii relies on PETSc to use iterative solvers. And since PETSc can handle complex numbers, I would not have expected this to be a problem.

Best,

Lucas

[Wolfgang Bangerth]

On 1/20/25 09:00, Lucas P wrote:
> Thanks Wolfgang for your answer. But I'm surprised to see that some work
> should be done regarding the iterative solvers. I would expect that deal.ii
> relies on PETSc to use iterative solvers. And since PETSc can handle complex
> numbers, I would not have expected this to be a problem.

Oh, sure, if you want to use the PETSc solvers I, too, would expect those to
work. The issue is more that if you want to use the PETSc solvers, you can't
use block solvers, you have to recompile PETSc to support complex numbers,
etc. It's a hassle. But if you know that all you want to do is solve a
non-block-matrix, PETSc may very well be your simplest choice. Everything
else, I believe, should be in place or should be easily fixable.

If you do go down that route, please do report back whether you made it work
and what it took. That might make for a nice point of experience for others
wondering about the same issue.

Best
W.