# Problems using P3 FE and periodic BC

**URL:** <https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568>\
**Category:** General Discussion\
**Created:** [June 29, 2023, 1:24pm UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568 "2023-06-29T13:24:05Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![Elissa](https://avatars.discourse-cdn.com/v4/letter/e/5e9695/32.png) [@Elissa](https://community.freefem.org/u/Elissa)\
**Post date:** [June 29, 2023, 1:24pm UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/1 "2023-06-29T13:24:05Z")

</div>

I consider the following situation.  
I have a domain comprising of a free flow part Z^+=(0,1)x(0,4) and a porous-medium part Z^-= (0,1)x(0,4)\P, where solid inclusions are present and P denotes the set of solid inclusions.  
Both subdomains are divided by the interface S=(0,1)x{0}.  
The overall flow domain is denoted by Z=Z^+ \cup S \cup Z^-.

Now I consider two problems defined on the domain Z.  
First a Stokes problem with a jump condition across the interface.  
The solution to this problem is {u,p}.  
What I know: u stabilizes exponentially fast to a constant C for x\_2 → 4.

The second problem is then defined as follows:  
problem eta([ue], [uhe], init=reuseMatrix, solver=sparsesolver) =  
int2d(Th)(  
dx(ue)\*dx(uhe)  
+ dy(ue)_dy(uhe)  
)  
+ int2d(Th)((u - C_(y\>0))_uhe)  
+ int1d(Th,interface)( -int2d(Th)((u)_(y\<0)) \* uhe )  
/// we set the following constraint to define eta uniquely  
/// but the exact constant for us constant does not play a role  
+ int1d(Th, bottom)(ue \* uhe)  
;

Since for later calculations I need that the third derivative of ue exists I need to use P3 FEM at least.  
However, when I use P3 FEM I observe that when I compute dx(ue) on the left and right boundary there appear some sort of numerical singularities (see picture attached).  
But indeed, for the second problem itself it is sufficient to use P2 FEM.  
So I did this and using P2 the singularities for dx(ue) do not appear.

Has somebody an idea where these problems come from?

 ![singularities](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/9/94ca4ac289463c23fdcdd6826daefb03b5118d71.jpeg)

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [June 29, 2023, 4:45pm UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/2 "2023-06-29T16:45:49Z")

</div>

Can you post the full code as an attachment? I guess if you just want to  
get around the 3rd derivative issue, you can introduce another  
variable as in mixed formulation, y2=dx(y1) and this may have some  
benefits or not.

---

<div class="post-metadata">

**Author:** ![Elissa](https://avatars.discourse-cdn.com/v4/letter/e/5e9695/32.png) [@Elissa](https://community.freefem.org/u/Elissa)\
**Post date:** [June 29, 2023, 7:34pm UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/4 "2023-06-29T19:34:56Z")

</div>

Hi,

thanks for your fast answer. Unfortunately, I was not allowed to upload the code.  
Attached you find the version I am working with.

Depending on the mesh and the radius of the solid inclusions the singularities appear at different locations.

Attached you find again a plot of dx(ue).

Best,  
Elissa

[test.edp](https://community.freefem.org/uploads/short-url/4uvaohHtl7IjtPZqIMUNsm09Gbj.edp) (4.15 KB)

 ![singularities-r025.png](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/3/33ecaf142a5725cf8686adc0ef3a6802e6ee648c.png)

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [June 29, 2023, 9:13pm UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/5 "2023-06-29T21:13:06Z")

</div>

I would guess that has something to do with the periodic boundary conditions.  
I have not used those enough to comment but I guess there are some  
issues there. I can’t run your code with so many nodes in the adaptmesh  
but may play with it later. In any case the FF experts may have more interest if  
you mention the periodic thing. I take it your porous medium them is  
periodic and not random. Do you expect physics due to this ? Resonances under  
some flow conditions? I guess analytically it would be suitable for FFT analysis  
probably. I guess it is possible under some conditions you get singularities.  
Also 'sum rules" may be important as I found out playing with a coloumb example  
with a net charge per cell 🙂 In any case, if you change the forcing functions curious if  
the peaks move or go away.

---

<div class="post-metadata">

**Author:** ![Elissa](https://avatars.discourse-cdn.com/v4/letter/e/5e9695/32.png) [@Elissa](https://community.freefem.org/u/Elissa)\
**Post date:** [June 30, 2023, 8:06am UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/6 "2023-06-30T08:06:56Z")

</div>

Hi,

thanks again.

I attached a version of the code with less nodes since I am using a smaller domain.  
Nevertheless, the issues stay.

Best,  
Elissa

![blocked.gif](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/2/23a8d201efa0a2f0176ba5b9aed6ae01a1173452.gif)

[test.edp](https://community.freefem.org/uploads/short-url/mkevlyk8WN9TLp6dFIDI7Lq32jh.edp) (4.15 KB)

---

<div class="post-metadata">

**Author:** ![Elissa](https://avatars.discourse-cdn.com/v4/letter/e/5e9695/32.png) [@Elissa](https://community.freefem.org/u/Elissa)\
**Post date:** [June 30, 2023, 8:08am UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/7 "2023-06-30T08:08:47Z")

</div>

That might be the case. But I do not see why the periodic BC should influence the computations only when using P3 FEM and not if I use P2…

The porous medium is periodic, correct, therefore I use the periodic BC.

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [June 30, 2023, 11:55am UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/8 "2023-06-30T11:55:25Z")

</div>

I’m not sure of the issues with adaptmesh and periodic boundary condtions  
but  
I guess given the physics I would try to include 2 periods in one cell,  
I think this is called a supercell, and see if you get a singularity inside the  
compsite cell. If you suddenly reduce the flow rate that fluid has to go  
somewhere and the singularities look like that could be where 🙂 What  
happens with the flow in the x direction which I guess is v?  
You have u and v constrained to zero at the bottom.

I had to reduce the vertex count in adaptmesh and could not  
reproduce the peaks but it sounds like they make some sense.  
I guess you could play with a position dependent height/depth or  
density and viscosity. If the equations don’t have a meaningful  
steady state result the output may vary with details of the model.

I’m not sure if there are any "“long range forces”’ that allow distant  
cells to interact as you see in electromagnetics but maybe you don’t  
need them here 🙂

---

<div class="post-metadata">

**Author:** ![Elissa](https://avatars.discourse-cdn.com/v4/letter/e/5e9695/32.png) [@Elissa](https://community.freefem.org/u/Elissa)\
**Post date:** [July 21, 2023, 7:47am UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/9 "2023-07-21T07:47:31Z")

</div>

Thanks for your answer.  
I run the code now without using adaptmesh but I had the same problems with the singularities on the left and right boundary (they just appeared in different positions).

I did not get the idea with the supercell you meant. Could you please explain again?

For the first problem (called stokes in my code), the problem is well-posed and the constraint setting u=v=0 on the bottom is also physically consistent. For this problem also the solution looks fine.  
Singularities appear only when considering the gradient of ue.  
But as stated before, this is only the case when I use P3 elements (that I need for later purposes) to solve the problem eta. When I use P2 elements the gradient of ue looks fine.

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [July 22, 2023, 1:21am UTC](https://community.freefem.org/t/problems-using-p3-fe-and-periodic-bc/2568/10 "2023-07-22T01:21:10Z")

</div>

If the derivative is positive on one side and negative on the other I guess  
you get a cusp at the boundary. I don’t remember the details but depending on  
the boundary conditions it looked like fluid could be “reflect” from the obstacles  
and the location of singularity looked like a good place for it to pile up.  
IIRC, you had two obstacles one above the other with periodic boundary conditions  
in the X direction. You could create 4 obstacles in a a square pattern and see if you  
get the singularity on the middle of the mesh as well as the periodic boundaries.  
This may also let you perturb their locations slightly if you are worried about  
stuff like that.