Macro ThN2O for mixed vectorial and scalar finite element space

General Discussion
#1

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!

#2

Yes, why would that matter?

#3

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.

#4

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.

#5

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!

#6

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);
}
#7

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

Best regards!