[petsc-users] Converting complex PDE to real for KNL performance ?

[Sajid Ali]

Hi Jed/PETSc-developers,

My goal is to invert a set of these PDE's to obtain a series of parameters
F_t (with TSSolve and TSAdjoint for function/gradient computation). I was
planning to use TAO for setting up the inverse problem but given that TAO
doesn't support complex scalars, I'm re-thinking about converting this to a
real formulation.

In https://doi.org/10.1007/s00466-006-0047-8, Mark Adams explains that the
K1 formulation of Day/Heroux is well suited to block matrices in PETSc but
the PDE described there has an large SPD operator arising out of the
underlying elliptic PDE. ( Day/Heroux in their paper claim that K2/K3
formulations are problematic for Krylov solvers due to non ideal eigenvalue
spectra created by conversion to real formulation).

Now, the system I have is : u_t = A*(u_xx + u_yy) + F_t*u;  (A is purely
imaginary and F_t is complex, with abs(A/F) ~ 1e-16. The parabolic PDE is
converted to a series of TS solves each being elliptic). I implemented the
K1/K4 approaches by using DMDA to manage the 2-dof grid instead of setting
it up as one large vector of [real,imag] and those didn't converge well
either (at least with simple preconditioners).

Any pointers as to what I could do to make the real formulation well
conditioned ? Or should I not bother with this for now and implement a
first order gradient descent method in PETSc  (while approximating the
regularizer as a cost integrand) ?


Thank You,
Sajid Ali | PhD Candidate
Applied Physics
Northwestern University
s-sajid-ali.github.io

On Wed, Mar 27, 2019 at 9:36 PM Jed Brown <jed at jedbrown.org> wrote:

> When you roll your own equivalent real formulation, PETSc has no way of
> knowing what conjugate transpose might mean, thus symmetry is lost.  I
> would suggest just using the AVX2 implementation for now and putting in
> a request (or contributing a patch) for AVX-512 complex optimizations.
>
> Sajid Ali via petsc-users <petsc-users at mcs.anl.gov> writes:
>
> >  Hi,
> >
> > I'm able to solve the following equation using complex numbers (with
> > ts_type cn and pc_type gamg) :
> >                               u_t = A*u'' + F_t*u;
> > (where A = -1j/(2k) amd u'' refers to u_xx+u_yy implemented with the
> > familiar 5-point stencil)
> >
> > Now, I want to solve the same problem using real numbers. The equivalent
> > equations are:
> > u_t_real   =  1/(2k) * u''_imag + F_real*u_real   - F_imag*u_imag
> > u_t_imag = -1/(2k) * u''_real   + F_imag*u_real - F_real*u_imag
> >
> > Thus, if we now take our new u vector to have twice the length of the
> > problem we're solving, keeping the first half as real and the second half
> > as imaginary, we'd get a matrix that had matrices computing the laplacian
> > via the 5-point stencil in the top-right and bottom-left corners and a
> > diagonal [F_real+F_imag, F_real-F_imag] term.
> >
> > I tried doing this and the gamg preconditioner complains about an
> > unsymmetric matrix. If i use the default preconditioner, I get
> > DIVERGED_NONLINEAR_SOLVE.
> >
> > Is there a way to better organize the matrix ?
> >
> > PS: I'm trying to do this using only real numbers because I realized that
> > the optimized avx-512 kernels for KNL are not implemented for complex
> > numbers. Would that be implemented soon ?
> >
> > Thank You,
> > Sajid Ali
> > Applied Physics
> > Northwestern University
>


-- 
Sajid Ali | PhD Candidate
Applied Physics
Northwestern University
s-sajid-ali.github.io
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200414/67d25f48/attachment.html>