# Saving and loading global/distibuted arrays **URL:** \ **Category:** General Discussion\ **Created:** [July 6, 2020, 4:27pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493 "2020-07-06T16:27:47Z")\ **Posts on this page:** 20\ **Page:** 1
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 6, 2020, 4:27pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/1 "2020-07-06T16:27:47Z")
Hello FreeFEM developers! Is there an easy way to gather/scatter distributed arrays on local partitions to/from a global array for file I/O under the HPDDM framework developed by @prj? In the examples I’ve seen (e.g. [navier-stokes-2d-PETSc.edp](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/navier-stokes-2d-PETSc.edp)), the file I/O is handled in an efficient distributed fashion where the local solution on each mpirank is saved/loaded in a separate file. How could I perform file I/O via a single file on a single process? In essence, I’d like to gather a distributed solution into a single array on mpirank 0 and save it as a single file for post-processing purposes. Then, I could similarly load that single solution file on mpirank 0 and scatter it to the various other processors for other computations. I imagine that macro\_ddm.idp already has a convenient way to do this, but I don’t know of any examples or documentation which could give some clues. Any suggestions would be helpful! Thank you, Chris ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 6, 2020, 4:32pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/2 "2020-07-06T16:32:12Z")
Sure thing, see [this](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/maxwell-3d-PETSc.edp#L27-L29) followed by [that](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/maxwell-3d-PETSc.edp#L73-L84), or [this](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/laplace-adapt-3d-PETSc.edp#L23-L25) followed by [that](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/laplace-adapt-3d-PETSc.edp#L42-L53). Both examples do exactly what you are asking for: from an initial global mesh, distribute it, keep the information to go from global to local, and once the parallel payload is over, go back from local to global and gather the solution (either on process 0 or everywhere). ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 9, 2020, 10:02pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/3 "2020-07-09T22:02:38Z")
What is the purpose of [this](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/maxwell-3d-PETSc.edp#L73-L75)? (pasted below) ``` real[int] sol; changeNumbering(A, solx[], sol); changeNumbering(A, solx[], sol, inverse = true); ``` As I see it, this 1. Copies the fe array `solx[]` in FreeFEM numbering to `sol` in PETSc numbering. 2. Copies `sol` back into FreeFEM numbering. I feel as if I’m missing a subtle trick here… ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 10, 2020, 6:05am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/4 "2020-07-10T06:05:34Z")
You have to always keep in mind two things when using PETSc in FreeFEM: 1. FreeFEM numbering has ghost unknowns (overlapping decomposition) 2. PETSc numbering has no ghost unknown (each unknown in PETSc numbering is assigned to a single process) So what this little piece of code does is basically: 1. go from FreeFEM to PETSc, which is very cheap, basically create a vector `sol` of dimension lower than `solx[]` and only keep the PETSc unknowns assigned to this local process 2. go from PETSc to FreeFEM. The trick here is that the additional parameter `exchange` is not used, false by default, so the value of the ghost unknowns are not updated, and set to 0 Why all this is needed? It’s just a very fancy way of making sure duplicated unknowns (because of ghost elements) are not summed multiple times in the [follow-up MPI reduction](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/maxwell-3d-PETSc.edp#L82). You can also simply remove those two `changeNumbering`, and do something like `solx[] .*= A.D;`, which will scale the solution by the partition of unity. Thus, duplicated unknowns (because of ghost elements) will be summed _but_ weighted, so you’ll still get a consistent solution w.r.t. a run on 1 process. See [this](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/#pf8c) tutorial and the attached example for more. [exchange\_changeNumbering.edp](https://community.freefem.org/uploads/short-url/5X3002C44AqlljsQeOpXthqajt9.edp) (1.9 KB) ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 10, 2020, 11:09am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/5 "2020-07-10T11:09:13Z")
Very enlightening. Thank you for the detailed explanation! 😌 ---
**Author:** ![cyrillebreard](https://avatars.discourse-cdn.com/v4/letter/c/90ced4/32.png) [@cyrillebreard](https://community.freefem.org/u/cyrillebreard)\ **Post date:** [July 12, 2020, 10:53pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/6 "2020-07-12T22:53:16Z")
I am having the same issue. From the reply, I m understanding that the information from global to local is given from “see this followed by that”. In That it requires the matrix A in “changeNumbering(A, solx[], sol);” so it requires “CreateMat” , does it ? Am I understanding that this should include { macro def(i)[i, i#y, i#z]// EOM // vector field definition macro init(i)[i, i, i]// EOM // vector field initialization createMat(Th, A, Pk) } Thanks for your help. Cyrille ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 13, 2020, 4:53am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/7 "2020-07-13T04:53:46Z")
That plus the definition of the macro `YourMeshName#N2O` so that you can later on use `int[int] rest = restrict(VhLocal, VhGlobal, YourMeshName#N2O);` to get the interpolation from local to global. ---
**Author:** ![cyrillebreard](https://avatars.discourse-cdn.com/v4/letter/c/90ced4/32.png) [@cyrillebreard](https://community.freefem.org/u/cyrillebreard)\ **Post date:** [July 13, 2020, 6:55am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/8 "2020-07-13T06:55:15Z")
Dear Chris, Are you willing to share your experience code specially from global to local? As you, I am trying to load solution (mean flow from NavierStokes Equation) on mpirank 0 and scatter it to all processors for eigenvalue analysis. There are some concepts that I can’t understand. Thanks Cyrille ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 13, 2020, 12:17pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/9 "2020-07-13T12:17:21Z")
@cyrillebreard here is how I do it: ``` load "PETSc" macro dimension()2// include "macro_ddm.idp" string workdir = getARGV("-dir", ""); string meshfile = getARGV("-msh", ""); string filein = getARGV("-fi", ""); string fileout = getARGV("-fo", ""); func Pk = []; // define FE space mesh Th = readmesh(workdir + meshfile + ".msh"); fespace Vh(Th, Pk); Vh def(u); // define u globally u[] = loadfile(workdir, filein, meshfile); // function to load file Mat A; int[int] n2o; macro ThN2O()n2o// EOM buildDmesh(Th); createMat(Th, A, Pk); def(u) = def(u); // define u on local partitions /*Main part of code ... */ real[int] q; changeNumbering(A, u[], q); if(fileout !=""){ changeNumbering(A, u[], q, inverse = true); mesh Thg = readmesh(workdir + meshfile + ".msh"); fespace Vhg(Thg, Pk); Vhg def(ug), def(uo); int[int] rest = restrict(Vh, Vhg, n2o); for[i, v : rest] ug[][v] = u[][i]; mpiReduce(ug[], uo[], processor(0, mpiCommWorld), mpiSUM); real[int] qo = uo[]; if (mpirank == 0) printfile(workdir, fileout, qo); // function to write file } ``` The load and print functions are basically just `ifstream` and `ofstream`. Note that I made some simplifications from my actual code to (hopefully) make this easier to understand, but this may have also led to a few typos or inconsistencies. I think the important thing here is knowing when you are working with local or global variables. Once you call `buildDmesh` and `createMat`, lots of things change! I spent a good chunk of time figuring out what the `macros_ddm.idp` file is actually doing before I got everything to work. By the way @prj, are there any efforts underway to expand the documentation for the PETSc interface in FreeFEM, or are things still evolving too rapidly to make that effort worthwhile? Obviously I’m still learning, but I’d be happy to help however I could. ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 13, 2020, 12:30pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/10 "2020-07-13T12:30:20Z")
I’ve added some documentation in my tutorial w.r.t.: 1. `restrict` and `new2old` functions, see [this](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/#pf27) slide and [that](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/section_3/example8.edp) example 2. `N2O` for `Dmesh`, see [this](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/#pfb2) slide and [that](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/section_8/example13.edp) example @cyrillebreard, what you want to do is exactly, I think, what is done in [this](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/navier-stokes-2d-PETSc.edp) script which computes a baseflow, followed by [that](https://github.com/FreeFem/FreeFem-sources/blob/develop/examples/hpddm/navier-stokes-2d-SLEPc-complex.edp) script which does linear stability analysis. 3D scripts are available [there](https://github.com/prj-/moulin2019al) as well. @cmd, the interface is almost complete now, so it’s not evolving very fast, just new functionalities here and there, plus bug fixes. What would you want more to be put in the [tutorial](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/)? Would you like to contribute? It’s missing parallel mesh adaptation and parallel interpolation. Probably will add these sometime soon. ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 13, 2020, 12:51pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/11 "2020-07-13T12:51:26Z")
@prj I will message you separately to avoid making this a very long thread 🙂 ---
**Author:** ![cyrillebreard](https://avatars.discourse-cdn.com/v4/letter/c/90ced4/32.png) [@cyrillebreard](https://community.freefem.org/u/cyrillebreard)\ **Post date:** [July 14, 2020, 8:49am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/12 "2020-07-14T08:49:18Z")
Dear Chris, Thanks for sharing the coding, and specially this line: def(u) = def(u); // define u on local partitions I would never find out this is the way to “transfer” the data from global to local. Thanks again. @prj like Chris I want to read one “sol” file and not several files like in " \*\*navier-stokes-2d-SLEPc-complex.edp" Thanks for your prompt reply. Cyrille ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 14, 2020, 11:10am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/13 "2020-07-14T11:10:56Z")
> def(u) = def(u); // define u on local partitions This is actually not the smartest thing to do since @cmd is using `macro ThN2O()n2o// EOM`. Indeed, this line relies on the interpolator, which may be very costly, and even inaccurate for very complex meshes. Much faster and safer to use the `n2o` array + the `restrict` keyword. See the penultimate line of this [example](http://jolivet.perso.enseeiht.fr/FreeFem-tutorial/section_8/example13.edp). ```auto u[] = uGlobal[](R); // notice the [] here, which means there is indeed no interpolation ``` No interpolation, just integer manipulation. ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 14, 2020, 2:44pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/14 "2020-07-14T14:44:01Z")
@prj Thanks for the tip! This is indeed a better solution. ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 15, 2020, 1:10pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/15 "2020-07-15T13:10:49Z")
@prj would it be worth adding a something similar to the `ThName#n2o` macros to automatically get this restriction information using the DDM macros? I’m thinking of something like: ``` int[int] n2o, restriction; macro ThN2O()n2o// EOM macro ThRestriction()restriction// EOM buildDmesh(Th); createMat(Th, A, Pk); ``` where `createMat` will automatically populate the `restriction` array if the `ThRestriction` macro exists. This might be a little bit misleading since the restriction is linked to the FE space on the mesh and not the mesh itself, but I think this would be a nice feature if it doesn’t exist already. 🙂 ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 15, 2020, 1:21pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/16 "2020-07-15T13:21:10Z")
It doesn’t exist, and it’s not that trivial to add, because you also need to pass the original global mesh. Maybe something like ```auto int[int] n2o, restriction; mesh backup = Th; macro ThRestriction()restriction// EOM macro ThOriginal()backup// EOM buildDmesh(Th) createMat(Th, A, Pk) ``` But this is a bit of a mouthful, IMHO… ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 15, 2020, 1:58pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/17 "2020-07-15T13:58:18Z")
Anyway, here you go [f36059d](https://github.com/FreeFem/FreeFem-sources/commit/f36059d03e492500f7f8becc71c932b0e3403881). With the following `diff` for example13.edp. ```auto diff --git a/section_8/example13.edp b/section_8/example13.edp index cda1adc..63fc255 100644 --- a/section_8/example13.edp +++ b/section_8/example13.edp @@ -19,3 +19,3 @@ broadcast(processor(0), ThGlobal); -int[int] n2o; +int[int] n2o, R; mesh Th = ThGlobal; @@ -27,5 +27,8 @@ fespace Vh(Th, Pk); // full global mesh fespace VhGlobal(ThGlobal, Pk); // distributed local mesh +macro ThOriginal()ThGlobal// EOM +macro ThRestriction()R// EOM +Mat A; +createMat(Th, A, Pk) Vh u = x^2 + y^2; -int[int] R = restrict(Vh, VhGlobal, n2o); VhGlobal uGlobal, uReduce; @@ -42,4 +45,2 @@ plot(uReduce, wait = 1, cmm = "u on the global mesh (weighted before summed)"); -Mat A; -createMat(Th, A, Pk) u = x^2 + y^2; ``` ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [July 15, 2020, 4:54pm UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/18 "2020-07-15T16:54:02Z")
Yes, I agree. This is not as elegant as I thought it might be. 😐 Oh well! Hopefully it will be useful to some. ---
**Author:** ![cyrillebreard](https://avatars.discourse-cdn.com/v4/letter/c/90ced4/32.png) [@cyrillebreard](https://community.freefem.org/u/cyrillebreard)\ **Post date:** [July 16, 2020, 7:34am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/19 "2020-07-16T07:34:05Z")
I made huge progress thanks to your help. But I don’t get the same result with 2 npus compared to one cpu. I can see that the problem is at the ghost cell. I am getting the following message: " Warning: – Your set of boundary condition is incompatible with the mesh label." This message does not appear with only one cpu. Any hints/tips? " ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 16, 2020, 9:19am UTC](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493/20 "2020-07-16T09:19:04Z")
The warning is to be expected, because there are “artificial” boundaries created at the interfaces between subdomains due to the domain decomposition. You can suppress it using `-v 0` in the command line. If you have problems at the ghost cells, you are probably missing a synchronization of the values, like, `exchange = true` during the `changeNumbering`. But this is only one of the many possible explanations. If you can, it would be best to extract the part of the code which does not give satisfying results into a [MWE](https://en.wikipedia.org/wiki/Minimal_working_example), and I’ll help you debug. [Next page](https://community.freefem.org/t/saving-and-loading-global-distibuted-arrays/493.md?page=2)