[Users] Carpet: using baseextents?

[Scott Hawley]

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