[Carpet] getting dt for the finest grid

[Erik Schnetter]

On Aug 28, 2007, at 07:08:13, Yosef Zlochower wrote:

> Hi,
>
>   What is the best procedure for obtaining the timestep for the  
> finest  grid
> currently used in a Carpet driven simulation?

Yosef,

I would probably write a small piece of C++ code which looks into  
Carpet's semi-public data structures.  Semi-public because this code  
works only if Carpet is active and you wouldn't want to have too much  
of such code around.

I can also introduce an aliased function or a GH extension where you  
can access this information.  But see below for some important  
points; let's make sure that you really get the information you need  
before I extend Carpet's API.

> Would something like the code below work for arbitrary numbers
>  of processors, or would some processors find different values of dt
>  because they do not own a segment of the finest grid?

No, all processors iterate over the same levels in the same way.  In  
fact, they currently all have to own the same number of components  
per level.  (This number is usually 1 if there are many processors.)

> Can ResetDt be
>  scheduled in the same timebin as FindDt (I am thinking about the  
> case where
>  a new refinement level may be added, or perhaps taken away)?

Yes.  I would schedule this in postregridinitial and in postregrid.   
Then this routine is only executed when the refinement structure  
changes.  Be aware that the thorn Time can also change dt, so it may  
be best if this routine only returns the time refinement factor, not  
the time step itself.  Otherwise it also needs to be scheduled after  
Time may have changed something.

> Would
>  the correct value of dt be used in MainCode as scheduled?

I think the code you present below is correct.

Pay attention to the relative order of this global mode routine to  
other level mode routines.  If you access grid functions which are  
evolved in time, e.g. via interpolation or reduction operations, then  
you need to ensure that your global mode routine is called after  
these other routines.  This is not guaranteed by global mode --  
global mode routines are executed first, not last, so that grid  
function evolution routines can access the result of global mode  
routines.

You can schedule your routine in level mode instead of global mode,  
and do nothing unless the refinement level is the finest refinement  
level.  This lets you choose the order with AFTER and BEFORE.  AFTER  
and BEFORE do not work well with global mode routines, because global  
mode routines are not always called, and then these specifications  
are ignored.  This also set automatically the time delta correctly:

USES INCLUDE HEADER: carpet.h

#include <carpet.h>

extern "C"
int
number_of_refinement_levels ()
{
   return Carpet::reflevels;
}

extern "C"
int
current_refinement_level ()
{
   return Carpet::reflevel;
}

if (current_refinement_level() == number_of_refinement_levels() - 1)

>  Is there a much cleaner way to do this?

There probably is.  You are interested in evolving grid arrays in  
time; one problem with this is that MoL wants to evolve them once for  
every refinement level, not only on the finest level.  I would begin  
by updating MoL.  MoL would then evolve grid arrays together with the  
finest level of grid functions, and would then automatically use the  
correct time step size there.  Then you only need to provide the RHS.



Instead of scheduling routines, you can also write a short piece of C+ 
+ code which access the maximum time refinement factor:

USES INCLUDE HEADER: carpet.h

#include <carpet.h>

extern "C"
int
max_time_refinement_factor ()
{
   return Carpet::timereffacts.at (reflevels - 1);
}

CCTK_REAL dt = CCTK_DELTA_TIME / max_time_refinement_factor();

-erik

> schedule.ccl:
>
> schedule ResetDT AT EVOL BEFORE FindDt
> {
>   LANG: C
>   OPTIONS: GLOBAL
> } "set dt to some big value"
>
> schedule FindDt AT EVOL BEFORE MainCode
> {
>   LANG: C
>   OPTIONS: GLOBAL LOOP-LOCAL
> } "find dt"
>
> schedule MainCode AT EVOL AFTER FindDT
> {
>    LANG: C
>    OPTIONS: GLOBAL
>  } "use the smallest value of dt"
>
> source.c:
>
> static CCTK_REAL smallest_dt
>
> void ResetDT(CCTK_ARGUMENTS)
> {
>   smallest_dt = BIG_NUMBER;
> }
>
> void FindDt(CCTK_ARGUMENTS)
> {
>   DECLARE_CCTK_ARGUMENTS
>   if (fabs(smallest_dt)> fabs(CCTK_DELTA_TIME))
>   {
>      smallest_dt = CCTK_DELTA_TIME;
>    }
> }
>
> void MainCode(CCTK_ARGUMENTS)
> {
>   ...
>   fdot = INTERPOLATE_SOME_VALUE_FROM_GFS;
>
>   /* f is a cactus grid array */
>   f = f_p + smallest_dt * fdot;
> }
>
> Thanks
>
> Yosef
> _______________________________________________
> developers mailing list
> developers at lists.carpetcode.org
> http://lists.carpetcode.org/listinfo/developers
>


-- 
Erik Schnetter <schnetter at cct.lsu.edu>

My email is as private as my paper mail.  I therefore support encrypting
and signing email messages.  Get my PGP key from www.keyserver.net.



-------------- next part --------------
A non-text attachment was scrubbed...
Name: PGP.sig
Type: application/pgp-signature
Size: 186 bytes
Desc: This is a digitally signed message part
Url : /archives/developers/attachments/20070828/ddeebc0e/attachment.pgp