[bddcml-users] Complex number support

Wed Jan 27 17:35:20 CET 2021

Hi Jakub,

Thanks very much for that clarification, that all makes sense!

I've just pushed the new complex_numbers branch to the Github repository.
If you're interested in what's changed, I'd suggest using
the --ignore-space-at-eol option for git diff since my editor was setup to
trim trailing whitespace so without it you'll get a lot of redundant lines
in the diff. In general, the changes are really just in two categories:

   1. Type changes from real to complex
   2. Library subroutine call changes for Lapack, Mumps, and Blopex

I've tested it on the Poisson equation example in examples/ ensuring that
the real and complex branches yield identical results for both direct
solves and the various Krylov methods, in both serial and parallel
execution. There may well still be lingering bugs so by all means let me
know if you find something that is suspect.

This was really fun and I'm excited to now integrate the library with my
code to test further and see what happens. Cheers!

Sebastian

On Tue, Jan 26, 2021 at 4:28 PM Jakub Sistek <sistek at math.cas.cz> wrote:

> Hi Sebastian,
>
> wow, your progress with the complex branch is impressive!
>
> As for 1), the main focus of BDDCML were symmetric problems for a long
> time. Also now, most of the applications are focused on symmetric problems.
> When I wanted to extend the method for nonsymmetric problems, I just
> thought it would be simpler to use BiCGstab than GMRES. Also, there were
> probably no functions for keeping the and orthogonalizing the basis in
> GMRES. In the meantime, such components were added to the PCG method to
> allow Krylov subspace recycling following the early works of Fahrat et al.
> Thinking of that, with these things prepared, adding GMRES now should not
> be a big deal, and in fact a nice addition to BDDCML. We shall look on this!
>
> As for 2), yes, your understanding is correct and this should work.
> However, it would lead to very large yet almost empty matrices. If this was
> desirable, I would probably tweak the BDDCML code and make it change all
> the edge nodes to corners. If I am not mistaken, this should be in
> module_dd.f90, changing the line 9172 to
>
> globtypes(jnodi) = 3
>
> and
> lines 9179-9195 to
>                if (nsubn.eq.1) then
>                   ! it is a face
>                   ifaces = ifaces + 1
>                   globtypes(inodi) = 1
>                else
>                   !if (nglobn.gt.1) then
>                   !   ! it is an edge
>                   !   iedges = iedges + 1
>                   !   globtypes(inodi) = 2
>                   !else
>                   !   ! it is a vertex, add it as a corner
>                   !   icorners = icorners + 1
>                   !   globtypes(inodi) = 3
>                   !end if
>                   icorners = icorners + 1
>                   globtypes(inodi) = 3
>                end if
>
> I would need to test it a bit more but it may force the code to identify
> all edge nodes as corners. In any case, one should think of changing the
> classification of interface nodes rather than the huge user constraints.
>
> I am excited to see your progress,
>
> Jakub
>
>
> On 26/01/2021 00:13, Sebastian Grimberg wrote:
>
> Hi Jakub,
>
> Sure thing, I won't make any changes to master of course. I'll look
> forward to hearing what Juan has to say on this topic, I am also by no
> means an expert here and happy to learn best approaches.
>
> For the weights, I was speaking purely hypothetically as I don't know if I
> will be using the diagonal stiffness. That makes sense about the weights
> for the inner product being a totally different thing than the weights from
> the BDDC point of view, so as long as they are constructed to add up to one
> then you can use them to compute dot products of parallel vectors, since
> the redundant dofs will be weighted to add to one. This is all good and
> makes sense for the Krylov solvers.
>
> Sorry again, but two last quick questions out of curiosity:
> 1. Why BiCGSTAB over GMRES? Would it be substantial to implement GMRES
> given what already exists, and/or would it provide any benefit?
> 2. From my understanding, the "use_arithmetic_constraints" add for each
> subdomain edge a single coarse constraint for the average of the fine dofs
> on that edge (is this correct?). If you wanted to add all edge dofs to the
> coarse space, so that you are actually imposing continuity at each edge dof
> separately, how would you do this? I would imagine the input to
> bddcml_setup_preconditioner would be a matrix of size n_edge_nodes x
> n_nodes where n_edge_nodes is the number of subdomain nodes lying on an
> edge. This matrix would just have a one in each row corresponding to the
> edge node. Is my understanding here correct? Obviously this would lead to a
> very large coarse space but just curious if I understand what is going on.
>
> Thanks, and I've already got most of the complex branch changes compiling
> so just need to test them out to see if I broke anything.
>
> Sebastian
>
> On Mon, Jan 25, 2021 at 1:55 PM Jakub Sistek <sistek at math.cas.cz> wrote:
>
>> Hi Sebastian,
>>
>> I have added you as a developer of BDDCML. Just do not push to the
>> master, please. Another branch, say "complex", is fine!
>>
>> I will add Juan Calvo into the discussion on Nedelec elements vs.
>> BDDC(ML). He has studied these things. I have never worked with these
>> elements too closely really, but would be happy to learn.
>>
>> I think that I now understand your question about the inner products and
>> weights! It is a slight abuse of the concept of partition of unity. In DD
>> methods, the weighting operator plays a crucial role and there is a number
>> of strategies available like the diagonal stiffness etc. Do you actually
>> plan to use the diagonal stiffness at this side?
>>
>> In the Krylov method, on the other hand, the operator is just reused to
>> get the right outcome from the duplicated entries at the interface. There
>> are multiple ways how to do this, but I decided to go for this weighted
>> inner products to avoid the need to assign unknowns to processes in a
>> unique way as this is not needed anywhere else in the code. However, ANY
>> partition of unity would do this job without affecting the performance of
>> the preconditioner. So, the simple weights by cardinality would help here
>> just fine. If you plan to use the diagonal stiffness in the preconditioner,
>> then we will perhaps need to use separate weights here, it is not a problem
>> although we would need to think how to do that. Otherwise, you could just
>> use the
>> wights_type = 0
>> straightaway.
>> Even if you needed the diagonal stiffness weights and thus the new
>> partition of unity for inner products, that should be a rather simple thing
>> to be done. I think I would not put the complex numbers as weights into
>> these inner products.
>>
>> Let us keep in touch on this! I will contact the other developers to get
>> their points of view and put you in the loop.
>>
>> Do not hesitate to ask! The code can be quite hard to read at many places
>> so a quick question can be worth hours of decyphering the code. I am more
>> than happy to help!
>>
>> Best wishes,
>>
>> Jakub
>>
>> --
>> Jakub Sistek, Ph.D.
>>
>> Researcher
>> Institute of Mathematics
>> Czech Academy of Sciences
>> sistek at math.cas.czhttp://www.math.cas.cz/~sistek
>> +420 222 090 710
>>
>>
>
> --
> Jakub Sistek, Ph.D.
>
> Researcher
> Institute of Mathematics
> Czech Academy of Sciences
> sistek at math.cas.czhttp://www.math.cas.cz/~sistek
> +420 222 090 710
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://list.math.cas.cz/pipermail/bddcml-users/attachments/20210127/8ab0ea1e/attachment.html>