# 2D-line mesh coupling with PETSc

**URL:** <https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254>\
**Category:** General Discussion\
**Created:** [April 15, 2026, 9:53pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254 "2026-04-15T21:53:51Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![aszaboa](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/aszaboa/32/2918_2.png) [@aszaboa](https://community.freefem.org/u/aszaboa)\
**Post date:** [April 15, 2026, 9:53pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/1 "2026-04-15T21:53:51Z")

</div>

Hi everyone,

I am trying to do the following. I am solving a fluid dynamics problem with FreeFEM + PETSc, and then, I would like to calculate the derivatives of the velocity field with a simple finite element projection method (mass matrix unknown derivative = varf(dx(u)v)). In general this seems to work, but near the boundaries (that correspond to walls) the derivative seems ‘noisy’, the accuracy of the derivative calculation is reduced. To enforce boundary conditions in the derivative calculation, I would like to prescribe that the tangential component of the velocity derivative is zero, and enforcing this with a Lagrange-multiplier seems to be the most straightforward method. I tried to implement the method using examples found online, but for some reason the results seem to be worse than before (see the image below). I am a bit stuck finding the issue, any help is appreciated. I attach a MWEish script, which is a modification of the [FreeFem-sources/examples/hpddm/navier-stokes-2d-PETSc.edp at develop · FreeFem/FreeFem-sources · GitHub](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/navier-stokes-2d-PETSc.edp) example. Any help is appreciated!

[navier-stokes-2d\_BF\_couplingTest\_v5.edp](https://community.freefem.org/uploads/short-url/4yW5M9vzwiFzMW7TxtGiW6LHFTH.edp) (7.0 KB)

 ![Screenshot from 2026-04-15 17-53-23](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/2X/3/372768486b4c2f5924f4ea44ed25ba205d2e0d9f.png)

---

<div class="post-metadata">

**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\
**Post date:** [April 16, 2026, 1:49pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/2 "2026-04-16T13:49:10Z")

</div>

Your code does not run in debug mode, which makes it difficult to debug.

```auto
Rank 0/2 has th2L: 240x47822
Assertion failed: (hmin <= longmini_box), function SameVertex, file GenericMesh.hpp, line 1348.

```

I’m guessing you are calling `trunc()` or a similar function with an empty mesh (or an always-false condition). Could you fix that?

---

<div class="post-metadata">

**Author:** ![aszaboa](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/aszaboa/32/2918_2.png) [@aszaboa](https://community.freefem.org/u/aszaboa)\
**Post date:** [April 16, 2026, 6:20pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/3 "2026-04-16T18:20:18Z")

</div>

I am not sure if I was calling the `trunc()` or `extract()`commands on an empty mesh. Maybe the problem is related to the fact that since I use a distributed mesh, some processes do not intersect with the boundary line that I want to extract so maybe.

Nevertheless, I attach a newer version of the code in which the `extract()`and `trunc()` commands should be definitely good. The error still persist. Do you have any further guesses why?

[navier-stokes-2d\_BF\_couplingTest\_v6.edp](https://community.freefem.org/uploads/short-url/4YQ3tKY8tb56Aef8yY5QV3GbUf5.edp) (6.9 KB)

```auto
int[int] lab2Extract = [1]; // label of the surface
meshL thLG = extract(ThGlobal, label=lab2Extract);
meshL ThL;
{
    fespace PhL(thLG, P0);
    PhL u = chi(Th);
    ThL = trunc(thLG,abs(u-1.0)<1.0e-2);
}

```

---

<div class="post-metadata">

**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\
**Post date:** [April 16, 2026, 6:34pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/4 "2026-04-16T18:34:25Z")

</div>

No, the error is still there.

```auto
Rank 0/2 has th2L: 240x47822
Assertion failed: (hmin <= longmini_box), function SameVertex, file GenericMesh.hpp, line 1348.

```

You are likely missing a `if (MeshName.nt)` or `if (functionUsedInTrunc[].max > 1.0E-10)`.

---

<div class="post-metadata">

**Author:** ![aszaboa](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/aszaboa/32/2918_2.png) [@aszaboa](https://community.freefem.org/u/aszaboa)\
**Post date:** [April 16, 2026, 6:52pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/5 "2026-04-16T18:52:17Z")

</div>

Ok, I think I got it.

```auto
{
    fespace PhL(thLG, P0);
    PhL u = chi(Th);
    if(u[].max > 1.0e-2){
        cout << "Rank " << mpirank << "/" << mpisize << " has thLG: " << thLG.nt << " triangles, max chi: " << u[].max << endl;
        ThL = trunc(thLG,abs(u-1.0)<1.0e-2);
    }
}

```

[navier-stokes-2d\_BF\_couplingTest\_v6.edp](https://community.freefem.org/uploads/short-url/ef4wF8GfI9cJs2MlxM5wolSvxb.edp) (7.1 KB)

---

<div class="post-metadata">

**Author:** ![aszaboa](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/aszaboa/32/2918_2.png) [@aszaboa](https://community.freefem.org/u/aszaboa)\
**Post date:** [April 16, 2026, 9:33pm UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/6 "2026-04-16T21:33:14Z")

</div>

I think I figured it out. What works for me is calculating the constraints on the line mesh first, then multiplying it with the 2D-1D interpolation matrix. I attach the corredponsing code segment and the full code. Thanks @prj for the quick replies!

```auto
/* -------------------------------------------------------------- */
/* scalar FE space for the derivative calculation */
varf vMv(u1,v1) = int2d(Th)(u1 * v1);
Mv = vMv(Vh,Vh,tgv=-1);

/* line mesh extraction */
int[int] lab2Extract = [1]; // label of the surface
meshL thLG = extract(ThGlobal, label=lab2Extract);
meshL ThL;
{
    fespace PhL(thLG, P0);
    PhL u = chi(Th);
    if(u[].max > 1.0e-2){
        cout << "Rank " << mpirank << "/" << mpisize << " has thLG: " << thLG.nt << " triangles, max chi: " << u[].max << endl;
        ThL = trunc(thLG,abs(u-1.0)<1.0e-2);
    }
}

/* restriction matrix for line FE space */
fespace VhL(ThL,PkV);
matrix th2L = interpolate(VhL,Vh);
cout << "Rank " << mpirank << "/" << mpisize << " has th2L: " << th2L.n << "x" << th2L.m << endl;
Mat B(Mv, restriction = th2L);
Mat restrictionMatrix(B, Mv, th2L);

/* constraints for the tangential derivative */
varf vMvDer1(u1, lambda) = int1d(ThL)(lambda * u1*Tl.x);
varf vMvDer2(u1, lambda) = int1d(ThL)(lambda * u1*Tl.y);
matrix vMvDer1temp = vMvDer1(VhL,VhL);
matrix vMvDer2temp = vMvDer2(VhL,VhL);
Mat MvDer1(B, vMvDer1temp);
Mat MvDer2(B, vMvDer2temp);

/* Calculating the actual constrains with matrix - matrix multiplication */
Mat MvDer1F2L, MvDer2F2L;
MatMatMult(MvDer1, restrictionMatrix, MvDer1F2L);
MatMatMult(MvDer2, restrictionMatrix, MvDer2F2L);

Mat Ablock = [[Mv, 0, MvDer1F2L'], 
              [0, Mv, MvDer2F2L'], 
             [MvDer1F2L, MvDer2F2L, 0]];

```

[navier-stokes-2d\_BF\_couplingTest\_v7.edp](https://community.freefem.org/uploads/short-url/pRYOXV2EqTXFQhoPRGTjuIixOuB.edp) (7.3 KB)

---

<div class="post-metadata">

**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\
**Post date:** [April 17, 2026, 4:28am UTC](https://community.freefem.org/t/2d-line-mesh-coupling-with-petsc/4254/7 "2026-04-17T04:28:25Z")

</div>

I wasn’t of much help, but glad you got it working.