Approximate computation of invh1 for Maxwell eigenvalues with meshes including prisms

General
1

Hi All,

I am try to compute Maxwell eigenvalues using an approximate invh1 obtained from cg:

            math1smooth = math1.CreateSmoother(freedofs=fesh1.FreeDofs())

            cgmath1 = CGSolver(

                mat=math1,pre=math1smooth,

                maxsteps=300,

                precision=1e-8

                )

            class H1Inverse(BaseMatrix):

                def Mult(self, x, y):

                    y.data = cgmath1 \* x

                def CreateColVector(self):

                    return math1.CreateColVector()

                def CreateRowVector(self):

                    return math1.CreateRowVector()

                def Height(self):

                    return math1.height

                def Width(self):

                    return math1.width    

            invh1 = H1Inverse()

rather than construct the invh1 directly as described here https://docu.ngsolve.org/latest/i-tutorials/unit-2.4-Maxwell/Maxwellevp.html

but getting the following warning if I use a mesh with prismatic elements and order =2 H(curl) elements :

GetTetsLocal not implemented for el with 0 nodes

It appears OK for meshes without prismatic elements. Is there something currently not implemented for prisms for the above?

Thanks Paul.

2

Hi Paul,

can you send the whole example?
It looks like you got a mesh with elements having 0 vertices.
For the example from the i-tutorials, replacing the inverse by the CGSolver (with enough accuracy) is working,

Joachim

3

Hi Joachim,

I attach a Jupyter notebook.

It looks like I may had an issue with only setting some of the bilinear forms to be symmetric and not all of them, which I think may have caused the warning about 0 nodes and I think is now resolved (sorry about that).

This is a non-standard Maxwell eigenvalue problem in that it consistent of small object (unit sphere) placed in a large free space (air) region. Gradients are skipped in the air region and regularisation is added here. It has zero Dirichlet bcs on the outer boundary.

Without the prismatic layers in the sphere, both the direct and iterative approaches converge to similar results using the PINVT solver. But with prismatic layers in the sphere the results show additional small spurious eigenvalues for the iterative approach, but not for the direct approach. So I might still be doing something incorrectly?

regularised_eigen_problem.ipynb (9.7 KB)

4

Hi Paul,
the h1 solver with jacobi doesn’t converge well with boundarylayers. and the c++ cg doesn’t tell you that it takes full maxiter (should emit a warning imo, but I’d recommend switching to the python version from ngsolve.krylovspace).
If I raise the maxiter in the h1 cg to 1000 I get the correct ones.
You can also try a better preconditioner, I think our h1amg does not yet work on prismatic elements, but in the github cerbsim repo we have some bindings to amgx and amgcl that work on the matrix level and should work.
Best
Christopher

5

Thanks Christopher.

6

Hi Christopher,

I have been looking at the ngsolve_amgcl.AMGCLPreconditioner. The usage at

suggests that you can use directly with a matrix, but the demo included uses a bilinear form. I have run the demo without issues. Looking at help(ngsolve_amgcl.AMGCLPreconditioner) it is also expecting a bilinear form rather than a matrix. So I was wondering if the interface has changed and if it is still possible to use this with the math1.mat as part of the eigensolver I am currently looking at?

Thanks

Paul.

7

ah sorry, yes internally it was only working on the matrix but the py wrapper did not expose that. Have added that now, you can now also create with sparse matrix and freedofs.
Best
Christopher

8

Thanks Christopher that is great.

Best wishes

Paul.