[Users] Carpet: using baseextents?

Erik Schnetter schnetter at cct.lsu.edu
Mon Sep 16 12:02:25 CDT 2013


Scott

I notice that you access the variable baseextents. Instead, consider calling the function "baseextent", as in baseextent(ml,rl).

To get the extent of a particular component, write

extent(ml, rl, c)

where "c" is the component number. The "local" grid has component number "Carpet::component". This would be

extent(mglevel, reflevel, component)

which returns a bbox where lower(), upper(), and stride() define the current component.

-erik

On 2013-09-13, at 0:12 , Scott Hawley <scott.hawley at belmont.edu> wrote:

> Erik,
> 
> I notice there's not "AT(c)" in the lines for the lower() calls, which would tell it which component grid I want the extents of. 
> 
> The result is that both of those lower() lines return [0,0,0].   This is for the finest levels of multiply-nested grids.   Similarly, if I replace "lower" with "upper", I get [768,768,768] -- ie it returns the bounds for the entire domain, not the grid at that level.
> 
> How can we get the integer-coordinate values for the boundaries of the local grid?
> 
> The stride() lines work just fine, by the way.  
> 
> Thanks,
> Scott
> 
> On Sep 12, 2013, at 8:17 PM, "Erik Schnetter" <schnetter at cct.lsu.edu> wrote:
> 
>> On 2013-09-12, at 20:57 , Scott Hawley <scott.hawley at belmont.edu> wrote:
>> 
>>> Goal: I'm trying to determine if a given (fine) grid vertex coincides
>>> with a coarser grid point.
>>> 
>>> Frank pointed me to baseextents, and I've come up with the following
>>> snippet of code...
>>> Just posting this for a sanity check: Am I on the right track?
>>> Is there a more robust way one should do this?
>>> 
>>> My code relies on "baseextents" being the "full" extent of a given
>>> refinement level, apart from any domain decomposition that occurs when
>>> running with multiple MPI processors.   (Well, and it relies on the fact
>>> that fine grid boundaries must lie along coarse grid points.)
>> 
>> The parenthesized remark doesn't sound quite right. See below, coarse_ilo.
>> 
>>> BEGIN_MAP_LOOP (cctkGH, CCTK_GF) {
>>>       BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
>>> 
>>>          DECLARE_CCTK_ARGUMENTS;
>>> 
>>>          assert (maps == 1);
>>>          int const m = 0;            // <--- Really?  m is never > 0?
>> 
>> You are using "m" for two things. One of these should rather be called "ml", a "multigrid-level" that is currently unused, and is thus always 0. The other is the "map", i.e. the patch number; this is 0 unless you have multiple patches.
>> 
>>>          assert (mglevels == 1);
>>>          int const rl = reflevel;
>>> 
>>>          ivect const ilo  = vhh.at(m)->baseextents.at(m).at(rl).lower();
>> 
>> The second .at(m) should be .at(ml) instead.
>> 
>>>          ivect const istr = vhh.at(m)->baseextents.at(m).at(rl).stride();
>> 
>> Same here.
>> 
>>>          int ref_ratio = 2;         // <-- Do we support refinement ratios > 2?
>> 
>> Nothing uses this at the moment; it is untested and thus likely not working.
>> 
>>>          ivect coarse_str = istr * ref_ratio;
>> 
>> Better: vhh.at(ml)->baseextents.at(m).at(rl-1).stride();
>> 
>> Also: coarse_ilo = vhh.at(m)->baseextents.at(ml).at(rl-1).lower();
>> 
>>>          for (CCTK_INT k = 0 ; k < cctk_lsh[2] ; k=k+1) {
>>>              for (CCTK_INT j = 0 ; j < cctk_lsh[1] ; j=j+1) {
>>>                  for (CCTK_INT i = 0 ; i < cctk_lsh[0] ; i=i+1) {
>>> 
>>>                     // does this vertex coincide with a coarser grid vertex?
>>>                     if (  ( 0 ==  (i - ilo[0]) % coarse_str[0]) &&
>>>                           ( 0 ==  (j - ilo[1]) % coarse_str[1]) &&
>>>                           ( 0 ==  (k - ilo[2]) % coarse_str[2])) {
>> 
>> This doesn't seem quite right. ilo[0] is not in the same "units" as i; instead, it is in the same "units" as str.
>> 
>> Alternatively: Check whether the equation   ilo[0] + i * str[0] == coarse_ilo[0] + coarse_i * coarse_str[0]   can be solved for coarse_i, i.e. check whether
>> 
>> (ilo[0] + i * str[0] - coarse_ilo[0]) % coarse_str[0] == 0
>> 
>> and similarly for the other two directions.
>> 
>>>                             CCTK_VInfo(CCTK_THORNSTRING,"   We are over
>>> a CG vertex!");
>>>                     }
>>> 
>>>                  } // end loop over i
>>>               } // end loop over j
>>>           } // end loop over k
>>> 
>>>       } END_LOCAL_COMPONENT_LOOP;
>>> } END_MAP_LOOP;
>>> 
>>> 
>>> Thank you very much,
>> 
>> 
>> To test this, you can output the coordinate of the current grid point, and then check (manually) whether this coordinate exists on the next coarser grid.
>> 
>> -erik
>> 
>> -- 
>> Erik Schnetter <schnetter at cct.lsu.edu>
>> http://www.perimeterinstitute.ca/personal/eschnetter/
>> 
>> 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/.
>> 
> 

-- 
Erik Schnetter <schnetter at cct.lsu.edu>
http://www.perimeterinstitute.ca/personal/eschnetter/

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/.



More information about the Users mailing list