Hello Roland,

thank you for your quick reply..

I have just talked to co-workers in the office (that are programming using MPI) and they told me i could just reduce the values of rho obtained "locally" on every processor that fall in the shell and that should give me the right result..

this seems to be kind of similar to what you proposed..

ie (written in pseudocode)

for(int iradii=1; iradii < nradii; iradii ++) {

LC_LOOP3 //local on each processor { CCTK_INT i3D = CCTK_GFINDEX3D(cctkGH, i, j, k); r[i3D] = sqrt(x^2+y^2+z^2); if (shell_radius[iradii-1] < r[i3D] && r[i3D] < shell_radius[iradii]) { rho_in_shell[iradii] += rho[i3D]; } } }

and then reduce over all radial shells

for(int iradii=0; iradii < nradii; iradii ++) { CCTK_Reduce(cctkGH,-1, &rho_total_shell[iradii], CCTK_VAR = rho_in_shell[iradii]) }

would that be an option? or even work?

best wishes,

Vassili

On Fri, Nov 15, 2013 at 5:07 PM, Roland Haas <roland.haas@physics.gatech.edu

wrote:

-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1

Hello Vassilios,

...

Is there a way to use the global Cactus reduction functions on fractions of the entire grid?

I am interested in summing grid scalars (say rho) in spherical shells in the grid. Is there a way to perform the "sum" reduction on those shells specified by a global radius variable? Or do I need to write my own MPI function that sums the contributions from all processors that have grid point that fall within the shell I am interested in?

There is no such facility build into Cactus unfortunately. The easiest way to achieve a reduction using only points satisfying some criterion (having coordinates inside the shells in your case) is to set up a helper grid function which is set to:

      / rho : if condition is true

helper = { \ 0 : otherwise

For multiple shells you will have to set up multiple helper grid functions (or re-use the same one several times). One way of doing this is:

schedule.ccl: STORAGE: helper[1] schedule setup_integrand in CCTK_ANALYSIS { OPTIONS: global loop-local LANG: C } "compute 'helper' grid function for integrand" schedule integrate in CCTK_ANALYSIS AFTER setup_integrand { OPTIONS: global LANG: C } "reduce 'helper' grid function on all processes"

You can have a look at the thorn ADMMass (in AEIThorns) and its ADMMass_Volume routines to see how this is done in practise.

Doing the reduction by hand can be tricky since you will have to properly handle mesh refinement boundaries. Roughly speaking you have to multiply your integrand by the CarpetReduce::weight grid function and only include the non-ghost zones points (ie start cctk_nghostzones points from the edge of each component). Outer boundaries require special treatment (detectable via cctk_bbox) and you have to include all points there. I am right now not sure how symmetry boundaries need to be handled (ie whether for vertex centering weight is set to 1/2 already or of you'll have to query the symmetry thorn whether an outer boundary is a symmetry boundary).

Yours, Roland

---

My email is as private as my paper mail. I therefore support encrypting and signing email messages. Get my PGP key from http://keys.gnupg.net. -----BEGIN PGP SIGNATURE----- Version: GnuPG v1.4.15 (GNU/Linux) Comment: Using GnuPG with Icedove - http://www.enigmail.net/

iEYEARECAAYFAlKGRr4ACgkQTiFSTN7SboWpsQCfUMEzLBth6C+s+lrqNz3Anb59 6eYAn3SdDPotHueI34ho6dm0rfVyKi/J =VQEr -----END PGP SIGNATURE----- _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users

Thanks a lot, I will try this!

best wishes,

Vassili

On Fri, Nov 15, 2013 at 5:39 PM, Roland Haas rhaas@tapir.caltech.eduwrote:

-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1

Hello Vassilios,

...

and then reduce over all radial shells

for(int iradii=0; iradii < nradii; iradii ++) { CCTK_Reduce(cctkGH,-1, &rho_total_shell[iradii], CCTK_VAR = rho_in_shell[iradii]) }

would that be an option? or even work?

Yes, that is essentially option 2 that I had outlined. You need to modify your first loop to be CCTK_LOOP3_INTBND (described in the user guide, be careful when trying to use OpenMP since rho_in_shell[iradii] might see contributions from different threads) so that you exclude ghost zones. You also have to multiply rhop[i3D] by CarpetReduce's weight function (easiest is to inherit from CarpetReduce which will give you access to the grid function). You can do the Reduction all in one MPI call if you'd like:

const int op_sum = CCTK_ReductionArrayHandle("sum"); const int ierr = CCTK_ReduceArray(cctkGH,-1, op_sum, nradii, CCTK_VARIABLE_REAL, rho_total_shell, 1/*num_dims*/, 1/*num_in_arrays*/, CCTK_VARIALBE_REAL, nradii, rho_in_shell);

or one of of the 1D ArrayToArray helper routines (CCTK_ReduceLocArrayToArray1D) described in the user guide.

Yours, Roland

---

My email is as private as my paper mail. I therefore support encrypting and signing email messages. Get my PGP key from http://keys.gnupg.net. -----BEGIN PGP SIGNATURE----- Version: GnuPG v1.4.15 (GNU/Linux) Comment: Using GnuPG with Icedove - http://www.enigmail.net/

iEYEARECAAYFAlKGTk8ACgkQTiFSTN7SboWrIgCgzzgU3hS/FzvXTzyzup9b0Tt0 3WYAoNFx6ZROsng8/gHVysQ68zp7os1r =nH9a -----END PGP SIGNATURE----- _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users

-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1

Hello Vassili,

how can i now perform the reduction on this "processor-local" quantity?

You can use the same as for the first example that I gave:

schedule.ccl: schedule integrate in CCTK_ANALYSIS AFTER setup_integrand { OPTIONS: global LANG: C } "reduce 'helper' grid scalar on all processes"

interface.ccl: CCTK_REAL rho_total_shell[nradii] TYPE=SCALAR "blah" CCTK_REAL rho_in_shell[nradii] TYPE=SCALAR "blah"

C-code: void integrate(CCTK_ARGUMENTS) { ... CCTK_ReduceArray(cctkGH,-1, op_sum, nradii, CCTK_VARIABLE_REAL, rho_total_shell, 1/*num_dims*/, 1/*num_in_arrays*/, CCTK_VARIALBE_REAL, nradii, rho_in_shell); }

Once intergrate ran then on all processors rho_total_shell will contain the sums of rho_in_shell (one sum per radius).

if i write a new function, will it have access to all the instances every processor creates?

I don't understand. Any scheduled function always only has access to processor local data (since you cannot access remote memory using MPI).

if i calculate the reduction inside the "processor-local" function, will it work? because the function is global loop-local and i will need a loop to perform the integration on every shell i specify..

Do you mean loop-local function? You cannot call the reduction inside of the loop-local routine. It will most likely cause an MPI hang. Instead you have to schedule a separate GLOBAL (no loop-local) routine to call do the reduction.

i tried to manually change modes using the carpet c++ macros, but that has led to segmentation faults only so far...

any idea how i could do that?

The thorn LSUThorns Refluxing does that (in its correct.cc file as far as I remember). I'd check there. Otherwise Erik has at times explained how this is done. Most recently here:

http://lists.einsteintoolkit.org/pipermail/users/2013-May/003004.html

Yours, Roland

- -- My email is as private as my paper mail. I therefore support encrypting and signing email messages. Get my PGP key from http://keys.gnupg.net.

On Thu, Nov 21, 2013 at 1:59 AM, Roland Haas <roland.haas@physics.gatech.edu

wrote:

-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1

Hello Vassili,

...

...

how can i now perform the reduction on this "processor-local" quantity?

You can use the same as for the first example that I gave:

Also, before I forget again: the reduction result will only be correct at "coarse" timesteps ie when all refinement levels are aligned in time. For the other steps you will mix current and past values from coarse and fine levels (you are doing a 0th order interpolation in time). Carpet will to 2nd order interpolation in time for reductions if not all data is aligned in time.

yes, i already noticed that when trying things out..

i also got the switching between levels sorted out (i hope, al least there are no more assertions faults...)

basically my question is whether the local integration of the scalar, ie

LOOP{ rho_in_shell[iradii] += BLAH

}

that has to be run as LOOP-LOCAL (or using the macros to manually switch the modes) has to be run in the SAME function as the later reduction, which has to be run as GLOBAL...

is that true? do they have to be in the same function? because then i have to switch modes again, because i want to do the reduction in a FOR loop across all iradii, and that would be local again, because the entire function integrate was scheduled as GLOBAL LOOP_LOCAL...

maybe i am missing things as well..

in the meantime, i have implemented the way used ADMBase using a loop_counter....

while i trust this method more, it will also be slower than the scalar way, because i will have to calculate the value of the now gridfunction rho_in_shell nradii times for every grid point..

best wishes,

Vassili

Yours, Roland

---

My email is as private as my paper mail. I therefore support encrypting and signing email messages. Get my PGP key from http://keys.gnupg.net. -----BEGIN PGP SIGNATURE----- Version: GnuPG v1.4.15 (GNU/Linux) Comment: Using GnuPG with Icedove - http://www.enigmail.net/

iEYEARECAAYFAlKNWwEACgkQTiFSTN7SboWXlgCgvkqmENaO7XBDup6cxl2ORXmN LBYAoJPKHV23EKuQhqX5iPXc8pMuxPYF =ynMM -----END PGP SIGNATURE----- _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users