Copy variables between different meshes with PETSc

[tak]

Dear Sir or Madam,

I would like to use different 2 meshes with PETSc.
Additionally, variables want to be copied from each other.
However, as shown following codes, the direct copy, such as “v=u”, does not work.

Would you please let me know how to do it?

Thank you very much in advance.
Best regards,
Tak

// $ ff-mpirun -np 4 script.edp -wg -ne
load “PETSc”
macro dimension()2// EOM
include “macro_ddm.idp”

meshN Sh1=square(30,30,[x,y]);
meshN Sh2=square(50,50,[x,y]);

buildDmesh(Sh1)
buildDmesh(Sh2)

fespace Ph1(Sh1, P2);
Ph1 u=sin(4x)+cos(5y);

fespace Ph2(Sh2, P2);
Ph2 v;
v=u;

macro defPlot(u)u//
plotMPI(Sh2, v, P2, defPlot, real, cmm = “v”)

[prj]

Dear Tak,
You can use the newly implemented routine transfer to help you go between one distributed solution to another. If you want the transfer PETSc Mat, you can use the transferMat function, see this new example PtAP-2d-PETSc.edp.
Here is the function used in your code, and the resulting screenshots.

load "PETSc"
macro dimension()2// EOM
include "macro_ddm.idp"

meshN Sh1=square(30,30,[x,y]);
meshN Sh2=square(50,50,[x,y]);

buildDmesh(Sh1)
buildDmesh(Sh2)

fespace Ph1(Sh1, P2);
Ph1 u=sin(4*x)+cos(5*y);

fespace Ph2(Sh2, P2);
Ph2 v;
// v=u;

macro defPlot(u)u//
transfer(Sh1, P2, u, Sh2, P2, v);
plotMPI(Sh1, u, P2, defPlot, real, cmm = "u")
plotMPI(Sh2, v, P2, defPlot, real, cmm = "v")

[tak]

Dear Professor prj,
Thank you very much for your quick response and your great assistance. Your answer have been of great help to me in my work.

Best,
Tak

[jmorvan (Jorge Morvan)]

Dear Professor Prj,

I’m new to freefem. I’m trying to use the function transfer but I always have some errors. I ran this code that you’ve posted in your reply to Tak and it didn’t work. Is it necessary to have any archive in the same folder of the original code?

Thanks in advance.
Best regards,
Jorge Morvan

[prj]

What is the error, and which FreeFEM version are you using?

[jmorvan (Jorge Morvan)]

Dear professor Prj,

I have a compile error that seems to be associated with “macro_ddm.idp”. This is the message I get and also the line in which the error occurs:

820 @ real timerPartition = mpiWtime The Identifier mpiWtime does not exist

Error line number 820, in file macro: buildOverlapEdgePeriodicRecursive in C:\Program Files (x86)\FreeFem++\idp\macro_ddm.idp, before token mpiWtime

I’m using the version 4.7-1.

Thank you!
Jorge Morvan

[prj]

You are using the FreeFem++ binary, you need to use FreeFem++-mpi.

[jmorvan (Jorge Morvan)]

Ok, I see!

Thank you!
Jorge Morvan

[etj]

Hello,
Does transfer(Sh1, P2, u, Sh2, P2, v); also work with multi-component FE spaces such as [P1,P1] or [P1,P2]?
Thank you.

[prj]

That’s precisely what is done in the FreeFEM example…

[etj]

Indeed, thank you! I missed that example.

Is there any example where the solution to be interpolated is not defined analytically but read from a series of files (previously saved with something like ofstream sol(filename + mpirank + "_" + mpisize + ".sol"); sol << u[];)?

[prj]

That does not change a thing, as long as your solution follows the distribution of the domain decomposition.

[etj]

Ok, thank you. So I guess I’m doing something wrong: the code below

verbosity = 0;
load "msh3"
load "PETSc"
load "gmsh"

int[int] n2oSaved;
int[int] n2oLoaded;

macro dimension()3// EOM            
include "macro_ddm.idp"            
macro def(i)[i, i#B, i#C, i#D]//    
macro init(i)[i, i, i, i]// EOM     

func Pk = [P1b,P1b,P1b, P1];

mesh3 Th0 = readmesh3("mesh0.mesh");
mesh3 Th = readmesh3("mesh.mesh");

Mat A;
macro ThN2O()n2oSaved//
createMat(Th, A, Pk);

Mat A0;
createMat(Th0, A0, Pk);

fespace Wh( Th,  Pk); 
fespace Wh0(Th0, Pk); 

Wh<real>  def(u);
Wh0<real> def(u0);

ifstream sol("sol0_" + mpirank + "_" + mpisize + ".sol");
sol >> u0[];
transfer(Th0, [P1b,P1b,P1b,P1], u0, Th, [P1b,P1b,P1b,P1], u);

returns the following message:
“Out of bound 0 <=0 < 0 array type = P2KNIS_IlEE
current line = 424 mpirank 0 / 8
Exec error : Out of bound in operator []
– number :1
Exec error : Out of bound in operator []
– number :1
err code 8 , mpirank 0”

[prj]

Right, I’ve replaced both readmesh by two cubes, and commented out the ifstream + sol >> u0[] lines (because I don’t have any of the files), and it’s working flawlessly, so the problem is in how you save or load the .sol, I guess.

[etj]

Thank you. I tried to replace the two readmesh by

mesh3 Th0 = cube(20,20,20,[x,y,z]);
mesh3 Th = cube(30,30,30,[x,y,z]);

and commented out the two lines ifstream + sol >> u0[], but I still get the same message.