Macro ThN2O() myN2O

Dear Pierre,

Following up your talk: Pierre Jolivet introduction to parallel FreeFEM part 1 - YouTube Pierre Jolivet introduction to parallel FreeFEM part 2 - YouTube

fespace VhG(Thglobal, [P2, P2, P1]);
fespace Vh(Th, [P2, P2, P1]);

int[int] myN2O;
macro ThN2O() myN2O //

int[int] subIdx = restrict(Vh, VhG, myN2O);

My question is:

Once I remesh (after 10 time steps), both Th and Thglobal would change, so do the FEM space VhG and Vh. I realise there is an error for: int[int] subIdx = restrict(Vh, VhG, myN2O);

I did cout<<"===="<<Vh.ndof<<"===="<<VhG.ndof<<endl; and found that both Vh.ndof and VhG.ndof were different after remeshing. However, It seems I am not allowed to redefine the following, even inside a bracket {…}
int[int] myN2O;
macro ThN2O() myN2O //

Please let me know what is the correct thing to do after remeshing?

Best,
Yongxing

I don’t have any error with this code.

include "macro_ddm.idp"

mesh Th = square(40, 40);
mesh Thglobal = Th;
fespace VhG(Thglobal, [P2, P2, P1]);
fespace Vh(Th, [P2, P2, P1]);

int[int] myN2O;
macro ThN2O() myN2O //

buildDmesh(Th);
int[int] subIdx = restrict(Vh, VhG, myN2O);

Th = square(80, 80);
Thglobal = Th;

buildDmesh(Th);
subIdx.resize(Vh.ndof);
subIdx = restrict(Vh, VhG, myN2O);

Thank you prj, I have now solved my problem using subIdx.resize(Vh.ndof) as you suggested. Can you please help me to look at another problem about savevtk()

I need to call savevtk() twice to save a velocity filed on the whole mesh, at the same I need to save a displacement field on a truncated solid mesh. First, if I only save results based on the local mesh, i.e. the following code works well:

if (n % saveEach == 0)
	savevtk("velocity.vtu",, Th, ux, uy, p, dataname="ux uy p", order=orderOut,append=true);

However, if I try to save the solid displacement as follows, the result are messed up.

if (n % saveEach == 0){
	savevtk("velocity.vtu", Th, ux, uy, p, dataname="ux uy p", order=orderOutVel,append=true);
	
	Wh [sxscaled,syscaled]=[sx,sy];
	sxscaled[].*= M.D; [wGx,wGy]=[0,0];[wRx,wRy]=[0,0];
	wRx[](subIdxW)=sxscaled[];
	mpiAllReduce(wRx[], wGx[], mpiCommWorld, mpiSUM);	

	if(mpirank == 0){
		mesh Ths=trunc(Thglobal,region==solid);
		fespace Sh(Ths, PV2);
		Sh [sgx,sgy];
		[sgx,sgy]=[wGx,wGy];
		savevtk("displacement.vtu", Ths, sgx, sgy, dataname="dx dy", order=orderOutDisp,append=true);	
	}
}

In the above Th is a local mesh and Ths is a global solid mesh which is truncated from global mesh Thglobal.

So what is the correct idea to save results in my case: one (fluid+solid) mesh which is distributed to different processes; at the same time, I also try to visualise displacement on a solid mesh (truncated from the whole/global mesh)?

Thank you very much!
Yongxing.

However, if I try to save the solid displacement as follows, the result are messed up.

Which results? Is wGx correct after the mpiAllReduce? On which FE space is M defined (Wh or Sh)?

Thank you Pierre, all the results are correct. It seems because the implementation of savevtk().

When I remove if(mpirank == 0){ }, and just run savevtk(“displacement.vtu”, Ths,…) at all the processors, it saves the results files correctly.

Do you expect this? the problem has been solved anyway. Thank you very much!

Best,
Yongxing.

Ah, of course, I missed that. Yes, it is expected. If you want to plot “sequential” solution, you need to use the additional parameter communicator = mpiCommSelf in savevtk.