Running failure of an additive Schwarz method in two overlapping non-matching grids

Hi there,
I wrote an additive Schwarz method code to solve a semi-linear elliptic problem. It seems fine but when I come to run it does not work properly. I guess the failure from the function uh that the code seems does not recognize it or something like that. I am looking for someone to adjust it as soon as possible.
Here is the algorithm of additive Schwarz method in the context of two overlapping nonmatching subdomains i=1,2.

-Δu_{i}ⁿ+cu_{i}ⁿ=f(u_{i}ⁿ⁻¹) in Ω_{i}
u_{i}ⁿ=uⁿ⁻¹ on γ_{i}
u_{i}ⁿ=0 on Γ_{i}

uⁿ(x)=maxu_{i}ⁿ(x), x∈Ω.

The freefem code is as follows

Additive Schwarz method for semilinear PDEs.edp (2.6 KB)

You code is wrong because uh is a number not a function

Vh1 uh ; 

and set a each iteration uh ,

the Boundary condition with  max(uold1,uold2); 

So my solution is

Additive Schwarz method for semilinear PDEs.edp (2.6 KB)

This finally solves my problem.
Thanks a lot for your help Dr. Hecht. :clap:

Dr. Hecht
Could I ask you one more question please?
What do i have to do to reduce the number of iterations?
This makes it more reliable.

You get a error because the number of node on central piece in not compatible in angle (not 16)
and in vertical direction !!!
And after be careful on vertical surface mesh , the orientation of the split quad face must be se same, ( this is no easy).

Thank you for your fast reply Dr. Hecht.
However, I am sorry I did not understand your meaning with the saying ‘’ not compatible in angle (not 16)’’. To your information, the discretization of both subdomains is in the context of overlapping non-matching grids. So it is natural that the number of nodes for both subdomains on the intersection area is not the same. Could you explain that much more please? And since you are well-experienced in this field (consider me as your student), I need your solution if possible.

Also, what is the fault regarding the orientation? According to my understanding of the mesh, I applied the same as what is in your FreeFEM documentation in this regard, unless you meant something else.

Your help is appreciated.

See on the graphic you see the point do not coincide on both mesh.

Thank you Prof. Hecht for your valuable notes.
I have one more question if you do not mind. It can be noticed that some elements angles of both subdomains are equal to 90 degrees
What I want is to build the above meshing in the sense that the angles of all elements are less than 90 degrees. This is so important as it satisfies what we call it " Discrete Maximum Principle".
Your help is appreciated.

This problem to build a acute mesh is very hard, you can read work of Jean-Marie Mirebeau from Orsay universty in France.

Great Prof.
For the sake of time, could you tell me which work of Jean-Marie Mirebeau exactly?

It do this in this tests,
and in

Adaptive and anisotropic finite element approximation: Theory and algorithms


This paper addresses the problem of optimal mesh adaptation, as minimizing the approximation error over all (possibly anisotropic) triangulations of prescribed cardinality, given a norm X of interest and a function f to be approximated.Expand



Greedy bisection generates optimally adapted triangulations

Math. Comput.

  • 2012


It is proved that as the algorithm progresses, the triangles tend to adopt an optimal aspect ratio which is dictated by the local hessian of f, which means that the adaptive triangulations satisfy a convergence bound which is known to be asymptotically optimal among all possible triangulation.Expand

Thank you Prof. Hecht for the two attached theoretical works, which are great to read after my graduation. But I am looking for some numerical examples that solve my numerical work at this limited time if possible. I found a numerical reference in the personal page of the author in the bibliography, but it is not working.
So, could you guide me with some acute triangulation examples, if possible? That would be a great help for me.

Dear Prof. Hecht

According to the code, ‘‘uold1’’ and ‘‘uold2’’ are two vectors. But as you know, we can not compare and take the maximum of two vectors unless all similar components of them are well compared in appropriate way. For example, we can say that uold1 < uold2 if uold1(i) < uold2(i) for all i components on the boundary. My question is, Does FreeFem++ automatically compare between all components of the two vectors? If yes, then we are in the safe side. If no, then could you suggest me how to do that please?
I appreciate your help.

Dear professor,

I faced a problem. Maybe only you can help me.
I want to simulate the propagation of elastic wave in a 2D model. By the obtained displacement of a specify point with time (i.e., the waveform), I could judge the correctness of the simulation results. I have successfully obtain the correct results using “solve or problem”. However, the results obtained based on “varf” is obvious wrong. The code in “varf. file” is written according to the examples from your homepage I have tried to debug this code, while the results still is wrong. Therefore, I am asking for your help.
The two codes with “solve” and “varf” are attached to this post, respectively.
Sequential_2D_Solve.edp (2.2 KB)
Sequential_2D_varf.edp (2.3 KB)

Please give me some advice.

Best wishes,

If I have make no mistake the solution is:
Sequential_2D_Solve-FH.edp (3.0 KB)

Dear Professor,

Thank you so much for your considerate help.
However, The values and trends of the results from the code modified by you are not same with the results from the simple code with only Solve (I uploaded in the last post), which had been confirmed right.
On the other hand, I want to solve the elastic problem in a 3D model. So, my goal is to use PETSC to solve this elastic problem in parallel mode.Now, I first tried to switch the code with “solve” of the two-dimensional problems into the code of using “varf”. Therefore, “solve” should be not existence, and only use “varf” to solve the problem. Could you give me more advice?

They is no problem build petsc Mat and do the same !

Ok, I will try according to your help.
However, the key problem is that the values and trends of the results from the code modified by you are not same with the results from the simple code with only Solve (I uploaded in the last post), which had been confirmed right.
Do you have any other advice?