Hi all,
This is a question about interpolation in Carpet. I am using InterpGridArrays in a function scheduled in CCTK_Initial for interpolating grid variables and, when I have a single refinement level (unigrid), the interpolation is done correctly with no error messages. When I have two refinement levels (level 0 and level 1), it happens that I need to use values of the variables at grid points on level 0 to interpolate to grid points on level 1. The points on level 0 I want to interpolate from are present on level 0 only. Refinement level 1 does not extend to that region. When I am running it with these 2 refinement levels, Cactus is giving me an error that some points are outside the grid (outside level 1). However, it should be possible for the InterpGridArrays function to interpolate the value at this point from the coarser refinement level as overall speaking the point is within the base grid of the simulation.
Thus, My query is:-- Is there any way to schedule my routine in global mode so that InterpGridArrays can use values at grid points on level 0 while calculating interpolated values for grid points on refinement level 1? If I simply add OPTION: global or OPTIONS: global loop-local in my schedule.ccl file, the run crashes with a segmentation fault.
Thank you for your time. I really hope you could look into it.
Thanks and regards Aryan Sharma ᐧ
Hello Aryan,
Thus, My query is:-- Is there any way to schedule my routine in global mode so that InterpGridArrays can use values at grid points on level 0 while calculating interpolated values for grid points on refinement level 1? If I simply add OPTION: global or OPTIONS: global loop-local in my schedule.ccl file, the run crashes with a segmentation fault.
Scheduling with "OPTION: global" or "OPTION: global-late" are the correct options to use for interpolation (exotic use cases aside).
To be sure: you do not need to write back the interpolated values onto the grid, right?
A segfault in "global" mode would normally indicate that you were trying to access a grid function, since those are not accessible in global mode. However if you also get the same segfault when using loop-local then this does not seem to be the cause of the segfault.
On top of my head I am not aware of any difference in accessible memory between "local" scheduled routines and "global loop-local" ones. If you compiled with enough debug options (and few enough -O levels) so that you get a useful backtrace, do you know the exact line in your code that causes the segfault?
Would it be possible for you to provide a demo thorn (stripped of physics, just the skeleton of the calls) that can be compiled and run demonstrate the issue (for example as an tar / zip attachment to an email)?
Yours, Roland
Dear Dr. Roland, Thanks for your insight and you are correct that I actually wanted to write back the interpolated values back onto the grid and so I was trying to access a grid function in Global mode and which I assume caused the segmentation fault as you told me. Yes, Sir I will try to debug it more and let you know if there is anything else.
Sir, So my one of the query is that:-- Is there any way, so that InterpGridArrays can use values at grid points on refinement level 0 while calculating interpolated values for grid points on refinement level 1?
Thanks for your help and time.
Thanks and regards Aryan Sharma
ᐧ
On Wed, Mar 4, 2020 at 7:04 AM Roland Haas rhaas@illinois.edu wrote:
Hello Aryan,
Thus, My query is:-- Is there any way to schedule my routine in global mode so that InterpGridArrays can use values at grid points on level 0 while
calculating
interpolated values for grid points on refinement level 1? If I simply
add
OPTION: global or OPTIONS: global loop-local in my schedule.ccl file, the run crashes with a segmentation fault.
Scheduling with "OPTION: global" or "OPTION: global-late" are the correct options to use for interpolation (exotic use cases aside).
To be sure: you do not need to write back the interpolated values onto the grid, right?
A segfault in "global" mode would normally indicate that you were trying to access a grid function, since those are not accessible in global mode. However if you also get the same segfault when using loop-local then this does not seem to be the cause of the segfault.
On top of my head I am not aware of any difference in accessible memory between "local" scheduled routines and "global loop-local" ones. If you compiled with enough debug options (and few enough -O levels) so that you get a useful backtrace, do you know the exact line in your code that causes the segfault?
Would it be possible for you to provide a demo thorn (stripped of physics, just the skeleton of the calls) that can be compiled and run demonstrate the issue (for example as an tar / zip attachment to an email)?
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://pgp.mit.edu .
Hello Aryan,
Thanks for your insight and you are correct that I actually wanted to write back the interpolated values back onto the grid and so I was trying to access a grid function in Global mode and which I assume caused the segmentation fault as you told me. Yes, Sir I will try to debug it more and let you know if there is anything else.
ok.
Sir, So my one of the query is that:-- Is there any way, so that InterpGridArrays can use values at grid points on refinement level 0 while calculating interpolated values for grid points on refinement level 1?
Unfortunately I do not know of a simple way to achieve that. You will have to employ Carpet's low-level mode changing macros (defined in Carpet/src/modes.hh and accessible via the carpet.hh public header).
Something like this:
schedule.ccl: schedule InterFunc in CCTK_INITIAL { LANG: C OPTIONS: local } "interpoalte my data"
interface.ccl: USES INCLUDE HEADER carpet.hh
configuration.ccl: REQUIRES Carpet
#include "carpet.hh" extern "C" InterFunc(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS;
// make a copy of the coordinates we want to interpolate to std::vector<CCTK_REAL> coords[3], outdata; for(int i = 0 ; i < 3 ; i++) { coords[i].resize(cctk->ash[0] * cctk->ash[1] * cctk->ash[2]); } outdata.resize(cctk->ash[0] * cctk->ash[1] * cctk->ash[2]); CCTK_LOOP3_ALL(getcoords,i,j,k) { int idx = CCTK_GFINDEX3D(cctkGH, i,j,k); coords[0][idx] = x[idx]; coords[1][idx] = y[idx]; coords[2][idx] = z[idx]; } CCTK_ENDLOOP3_ALL(getcoords);
BEGIN_GLOBAL_MODE(cctkGH) { outvars[0] = &outdata[0]; InterpGrid3D(...); } END_GLOBAL_MODE;
CCTK_LOOP3_ALL(writeresults,i,j,k) { int idx = CCTK_GFINDEX3D(cctkGH, i,j,k); var[idx] = outdata[idx]; } CCTK_ENDLOOP3_ALL(writeresults); }
The important parts is the BEGIN_GLOBAL_MODE()...END_GLOBAL_MODE() section which lets you escape to GLOBAL mode in a LOCAL scheduled function. You need to (or at least you should, it likely works without but only due to abusing the exact way in which DECLARE_CCTK_ARGUMENTS works) make copies of the coordinates and allocate a temporary buffer to hold the interpolation results due to both the coordinates (x,y,z) and the output grid variables (var) being not (legally) acccessible while you are in GLOBAL mode.
Note that in your case you may actually want to not enter GLOBAL mode but LEVEL mode as well asking explicitly for level 0 by adding an ENTER_LEVEL_MODE(cctkGH, /* rl */ 0)..LEAVE_LEVELMODE into the GLOBAL section ("ENTER" this time and not "BEGIN" since we want to descend into and not escape to). Otherwise (in global mode) the interpolator will use the finest refinement level that overlaps a given coordinate which may well be you refined level that you want to interpolate to.
You may also (but this is a pure optimization) skip the routine if the variable Carpet::reflevel is 0 meaning that you are being called on level 0.
Hope this helps.
Please note that my code is untested and largely just schematic. Please also keep in mind that while what you want to do is possible, it is generally not advisable (for example it may interact in odd ways with Carpet's prolongation and restriction operations or how Carpet will use regular prolongation operators if / when you change the grid latter on) and I am not really sure I can see of a reason why you would need to do this.
Basically: I am handing you rope, hoping you will not hang yourself (https://en.wikipedia.org/wiki/Wikipedia:Give_%27em_enough_rope and http://motd.ambians.com/quotes.php/name/linux_computers/toc_id/1-1-4/s/881).
Yours, Roland
Thanks for your help and time.
Thanks and regards Aryan Sharma
ᐧ
On Wed, Mar 4, 2020 at 7:04 AM Roland Haas rhaas@illinois.edu wrote:
Hello Aryan,
Thus, My query is:-- Is there any way to schedule my routine in global mode so that InterpGridArrays can use values at grid points on level 0 while
calculating
interpolated values for grid points on refinement level 1? If I simply
add
OPTION: global or OPTIONS: global loop-local in my schedule.ccl file, the run crashes with a segmentation fault.
Scheduling with "OPTION: global" or "OPTION: global-late" are the correct options to use for interpolation (exotic use cases aside).
To be sure: you do not need to write back the interpolated values onto the grid, right?
A segfault in "global" mode would normally indicate that you were trying to access a grid function, since those are not accessible in global mode. However if you also get the same segfault when using loop-local then this does not seem to be the cause of the segfault.
On top of my head I am not aware of any difference in accessible memory between "local" scheduled routines and "global loop-local" ones. If you compiled with enough debug options (and few enough -O levels) so that you get a useful backtrace, do you know the exact line in your code that causes the segfault?
Would it be possible for you to provide a demo thorn (stripped of physics, just the skeleton of the calls) that can be compiled and run demonstrate the issue (for example as an tar / zip attachment to an email)?
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://pgp.mit.edu .
users@lists.einsteintoolkit.org