how to provide a proper HypreAMS ams for complex linear system when using HypreGMRES

[GKDyah]

I have been using mfem for a year now, and while I find it cool, I still consider myself to be a new user. I have a question about Experiment 25p, where I would like to use HypreGMRES instead of GMRESSolver. I compared these two solvers in a real linear system at Experiment 31p, and I found that HypreGMRES (32s, iteration 19) is faster than GMRESSolver (40s, iteration 22) and requires fewer iterations. However, I am having difficulty dealing with preconditioning when trying to switch solvers in Experiment 25p.

I have read the implementation method in GMRESSolver in Experiment 25p, and it seems that the method forms a block diagonal preconditioner BlockDB for the complex system using two operators obtained from user-defined preconditioning operator ams00. As written in this reference paper (https://www.research-collection.ethz.ch/bitstream/handle/20.500.11850/103394/geo2015-0013.1.pdf?sequence=2), the hidden matrix in BlockDB is similar to an inverse matrix to PA.

My question concerns the role of preconditioning in gmres.Mult(B, X). To my understanding, what really matters is the hidden matrix in BlockDB and the Mult method itself. However, I am having difficulty retrieving/accessing this hidden matrix directly as a HypreParMatrix but ams00.printmatlab.

In Experiment 25p, you defined
HypreAMS* ams00 = new HypreAMS(*PCOpAh.As(),fespace);
ams00 its initial method is associated with fespace, and the size of hidden Matrix and input Matrix *PCOpAh.As<HypreParMatrix>() is limited.

I would like to obtain a small HypreParMartrix from ams00 and combine the small Matrix to form a big Matrix using HypreParMatrixFromBlocks or ComplexHypreParMatrix. Ideally, I’d liketo use HypreParMatrix as precondion directly, intead of using BlockDP in GMRESSolver.

I am not very familiar with the Hypre lib., so I am wondering if there is something special that I should keep in mind when providing a preconditioner with fespace initial. Is it possible to provide a known HypreParMatrix as a preconditioner instead?

[psocratis]

Hi @Kuloss,

First, I find it a bit strange that HypreGMRES and GMRESSolver have different convergence behavior. They should be equivalent. You may want to double check that the stopping criterion is the same.

Second, just to make sure we are on the same page, the HypreAMS is not a matrix, it's a preconditioner designed for a Nedelec (ND) space (H(curl)) based on auxiliary spaces and BoomerAMG. Please see here for more details. In Hypre, AMS is designed to work for 1 ND space. Since example 25 deals with a complex valued problem there are 2 ND spaces (one for the real and one for the imaginary part of the operator) and hence we use a block diagonal preconditioner whose diagonal blocks represent the action of 2 (though the same) AMS preconditioners to precondition the complex system (written in real block form). In principle you could construct a monolithic matrix from the 2 x 2 system, but then it would be harder to construct an AMS preconditioner which would correspond to 2 ND spaces at once. Harder in the sense that one would have to construct all the discrete operators of AMS (a smoother, discrete gradient, ND interpolations) and combine them appropriately for the monolithic 2x2 system. I also believe that in doing that the benefit would not be significant.

I hope this helps.

Best,

Socratis

[GKDyah]

Thank you for your reply, it was helpful. After considering your advice, I have decided to give up switching between HypreGMRES and GMRESSolver solvers. I realized that in my previous question, I neglected to provide complete information. I would like to add some further clarification. I understand that the HypreAMS solver is not a matrix but rather a solver, similar to GMRESSolver defined in hypre.hpp. Recently, I used the command HypreAMS* ams00 = new HypreAMS(*PCOpAh.As(),fespace); and received a printout of ams00.

if (1)
{
     std::ofstream ofs111(".../.../ams00.dat");
     ams00->PrintMatlab(ofs111);
     ofs111.flush();
     ofs111.close();
}

I had a chance to closely examine the matrix from ams00->PrintMatlab through Octave, as shown in Figure 1. However, I encountered difficulty accessing the corresponding HypreParMatrix version, which made me think of it as a "hidden matrix" in ams00.

(exp25p, refine 0, order 1,l-shape.mesh, freq = 5.0)

fig.1

In contrast, *PCOpAh.As<HypreParMatrix>(), a HypreParMatrix can be directly used for the method HypreParMatrixFromBlocks as shown in Fig. 2.

(exp25p, refine 0, order 1, l-shape.mesh, freq = 5.0)

fig.2

I examined the gmres.Mult(B, X) method and identified a specific code segment to investigate the function of prec. Upon closer inspection, I discovered that all instances of prec within the code are being employed similarly with the prec->Mult.function.

 if (prec)
 {
    if (iterative_mode)
    {
       subtract(b, r, w);
       prec->Mult(w, r);    // r = M (b - A x)
    }
    else
    {
       prec->Mult(b, r);
    }
 }

Finally, the 'r' in prec->Mult(w, r); corresponds to the matrix-vector multiplication shown in Fig.3 below.

fig. 3

fig.4

In my opinion, once I have the *ams00 and BlockDP, there is no significant advantage to merely using them as a matrix for computation purposes. Therefore, I am curious as to why we cannot directly use a HypreParMatrix in place of ams00 or BlockDP. Essentially, I can obtain the underlying matrix from ams00, utilize it to develop a monolithic 2x2 system, and HypreAMS can still function as a temporary intermediary.

Based on your comment, it seems that both methods, HypreGMRES and GMRESSolver, have similar speed and that precondition plays a similar role. If there were an option to give a user-defined HypreParMatrix as a precondition, I might be able to bypass the constraints imposed by the AMS precondition when dealing with a complex linear system problem more easily.

I am not well-versed in the field of high performance computing, and my idea may seem a bit fanciful. However, I hope to bring a bit of cheer to your day. I am intrigued by the study of hypre and would like to delve into it. It may take some time for me to fully acquaint myself with the topic, but I am excited to start this new journey.

Thank you again for your email and for taking the time to share your thoughts.

[tzanio]

Hi @Kuloss,

Maybe I missing something, but while the HypreAMS object is a linear operator, you never want to represent it as a matrix. Such a matrix will be dense, and HypreParMatrix will not be an efficient representation for it. Instead, we can efficiently apply the action of the HypreAMS operator on a vector through it's internal structure, think of it as a number of smaller sparse matrices.

This is similar to an outer product matrix $A = u v^$ where $u$ and $v$ are column vectors. Storing $A$ as a matrix is expensive, but you can still compute the action $A x = u (v^t x)$ efficiently by only storting $u$ and $v$.

Hope this is helpful,
Tzanio