Macro ThN2O for mixed vectorial and scalar finite element space

Hello everyone!

I’ve recently watched the tutorials from professor Pierre Jolivet on Youtube about parallel computing in FreeFem. I have a code that works fine sequentially and I intend to swicht it to parallel. In the sequential script I have one scalar and one vectorial finite element space, both associated with an initial mesh Th. I need to go from local to global for both vectorial and scalar fields. Can it be obtained by using macro ThN2O for both cases all in the same code?

Thanks!

Yes, why would that matter?

As I’m using macro def(i) [i i#y] and macro init(i) [i i] for the vectorial field, the reduction from local to global works fine for the vectorial field. However, I couldn’t do the same for the scalar field.

Please share a minimal working example. N2O is a mesh information, it has nothing to do with a finite element space whatsoever, so as long as you use it appropriately afterwards, you can do whatever you want with it.

Ok, I see! Thanks for your reply and for the information. I’ll try to make a short code those ideas to share with you. Thanks again!

Here is an example. Of course, you are free to use any kind of vectorial finite elements, I just used Morley here.

load "Morley"
include "macro_ddm.idp"
load "PETSc"

mesh ThGlobal = square(100, 100);
mesh Th = square(100, 100);
fespace WhGlobal(ThGlobal, P2Morley);
fespace VhGlobal(ThGlobal, P1);
WhGlobal [uG, uGx, uGy] = [x^2, y^2, x^2+y^2];
VhGlobal pG = x + y;
int[int] n2o;
macro ThN2O()n2o//
buildDmesh(Th);
fespace Wh(Th, P2Morley);
fespace Vh(Th, P1);
Wh [u, ux, uy] = [x^2, y^2, x^2+y^2];
Vh p = x + y;
int[int] restWh = restrict(Wh, WhGlobal, n2o);
int[int] restVh = restrict(Vh, VhGlobal, n2o);
{
    Mat A;
    createMat(Th, A, P2); // trick => use P2, not P2Morley!
    real[int] tmp;
    ChangeNumbering(A, u[], tmp);
    ChangeNumbering(A, u[], tmp, inverse = true);
    WhGlobal [reducex, reducey, reducez];
    WhGlobal [beforex, beforey, beforez];
    for[i, v : restWh] beforex[][v] = u[][i];
    mpiReduce(beforex[], reducex[], processor(0, mpiCommWorld), mpiSUM);
    uG[] -= reducex[];
    if(mpirank == 0) assert(uG[].l2 < 1.0e-12);
}
{
    Mat A;
    createMat(Th, A, P1);
    real[int] tmp;
    ChangeNumbering(A, p[], tmp);
    ChangeNumbering(A, p[], tmp, inverse = true);
    VhGlobal reduce;
    VhGlobal before;
    before[] = 0;
    for[i, v : restVh] before[][v] = p[][i];
    mpiReduce(before[], reduce[], processor(0, mpiCommWorld), mpiSUM);
    pG[] -= reduce[];
    if(mpirank == 0) assert(pG[].l2 < 1.0e-12);
}

Thank you so much for sharing this example. It’ll certainly help me!

Best regards!