Blockstructured Adaptive Mesh Refinement in object-oriented C++

---

       

---

GridFunction1.c

Go to the documentation of this file.

#ifndef _included_GridFunction1_c
#define _included_GridFunction1_c

/*****************************************************************************/
/**** Constructors ****/
/*****************************************************************************/

/*****************************************************************************/
/**** Simple Constructors - the GF is explicitly setup using GF_Set*() ****/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[], 
                                              GridHierarchy& gh)
        : GridFunctionVoid(gh.dagh_type(),1,name,gh,0,1,
                           gh.default_alignment(1),
                           gh.boundarywidth(),
                           gh.externalghostwidth(),
                           DAGHComm,
                           gh.boundarytype(),
                           gh.adaptboundarytype(),
                           DAGHNoExternalGhost),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) { 

  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;

  const int pnum = comm_service::proc_num();
  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }
  
  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/*** Constructor for backward compatibility ***/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int s_sten,
                                              GridHierarchy& gh,
                                              const int cflag, 
                                              const int extghflag)
        : GridFunctionVoid(gh.dagh_type(),1,name,t_sten,s_sten,
                           gh,0,1,DAGH_All,gh.boundarywidth(),
                           gh.externalghostwidth(),
                           cflag,gh.boundarytype(),
                           gh.adaptboundarytype(),extghflag), 
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) { 

  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
  
  const int pnum = comm_service::proc_num();
  register int p = 0;

  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int s_sten,
                                              const int cfac,
                                              GridHierarchy& gh,
                                              const int cflag, 
                                              const int extghflag)
        : GridFunctionVoid(gh.dagh_type(),1,name,t_sten,s_sten,
                           gh,0,cfac,DAGH_All,gh.boundarywidth(),
                           gh.externalghostwidth(),
                           cflag,gh.boundarytype(),
                           gh.adaptboundarytype(),extghflag), 
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) { 

  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
  
  const int pnum = comm_service::proc_num();
  register int p = 0;

  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/*** Constructor to set GF type (Cell/NonCell Centerted). ***/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int s_sten,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,0,1,DAGH_All,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0),
          iofunc(0) {
 
  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
  
  const int pnum = comm_service::proc_num();
  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
   
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int* s_sten,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,0,1,DAGH_All,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) {
 
  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);

  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
  
  const int pnum = comm_service::proc_num();
  register int p = 0;

  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }
  
  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/*** Constructor for GF's with gfRank < dagh.rank ***/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int s_sten,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int loff,
                                              const int cfac,
                                              const int align,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,loff,cfac,align,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0),
          iofunc(0) { 

  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
  
  const int pnum = comm_service::proc_num();

  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int* s_sten,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int loff,
                                              const int cfac,
                                              const int align,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,loff,cfac,align,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0),
          iofunc(0) {
 
  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;
   
  const int pnum = comm_service::proc_num();
  register int p = 0;

  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;

  GF_Compose();
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/*** Constructor for GF's at a given time and level ***/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int s_sten,
                                              const int time, 
                                              const int level,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int loff,
                                              const int cfac,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,loff,cfac,DAGH_All,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0),
          iofunc(0) {
 
  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;   

  const int pnum = comm_service::proc_num();
  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose(time,level);
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[],
                                              const int t_sten, 
                                              const int* s_sten,
                                              const int time, 
                                              const int level,
                                              GridHierarchy& gh, 
                                              const int type,
                                              const int loff,
                                              const int cfac,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(((type == DAGHNull) ? gh.dagh_type() : type),
                           1,name,t_sten,s_sten,gh,loff,cfac,DAGH_All,
                           gh.boundarywidth(),gh.externalghostwidth(),
                           cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) { 

  gfid = gh.DAGH_AddGridFunction((GridFunctionVoid *) this);
  
  length = new int[gh.totallevels()];
  for (register int l=0; l<gh.totallevels(); l++)
    length[l] = 0;   
  
  const int pnum = comm_service::proc_num();
  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose(time,level);
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
  
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/*** Constructor from another GF ***/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::GridFunction(1)(const char name[], 
                         const GridFunction(1)<DAGH_GFType>& gf,
                                              const int time, 
                                              const int level,
                                              const int cflag, 
                                              const int bflag,
                                              const int adptbflag,
                                              const int extghflag)
        : GridFunctionVoid(name,1,gf,cflag,bflag,adptbflag,extghflag),
          gdb(0), bvalue(0),
          ifunc(0), ufunc(0), bfunc(0), abfunc(0), pfunc(0), rfunc(0), 
          iofunc(0) {
 
  length = new int[dagh.totallevels()];
  for (register int l=0; l<dagh.totallevels(); l++)
    length[l] = 0;   
  
  for (register int i=0;i<2*time_sten_rad+1;i++) 
    time_alias[i] = gf.time_alias[i];

  gfid = dagh.DAGH_AddGridFunction((GridFunctionVoid *) this);

  const int pnum = comm_service::proc_num(); 
  register int p = 0;
  
  /* set up storage for ghost comm */
  if (comm()) {
    ghost_send_info = new GF_Interaction** [pnum];
    ghost_recv_info = new GF_Interaction** [pnum];
    ghost_recv_server = new GridTableGhostRcv** [pnum];
    for (p=0;p<pnum;p++) {
      ghost_send_info[p] = (GF_Interaction**) 0;
      ghost_recv_info[p] = (GF_Interaction**) 0;
      ghost_recv_server[p] = (GridTableGhostRcv**) 0;
    }
  }

  /* setup storage for data comm */
  data_recv_server = new GridTableDataRcv** [pnum];
  for (p=0;p<pnum;p++) data_recv_server[p] = (GridTableDataRcv**) 0;
  
  GF_Compose(gf,time,level);
  
  if (dagh.chkpt_restart() && checkpoint()) 
    GF_CheckpointRestart();
   
  myargc = 0;
  myargs = new char[1];
  myargs[1] = 0;
}

/*****************************************************************************/
/**** Destructor ****/
/*****************************************************************************/
template <class DAGH_GFType>
GridFunction(1)<DAGH_GFType>::~GridFunction(1)() {
  GF_DeleteGDBStorage();
  
  if (comm()) {
    GF_DeleteGhostCommInfo();
    if (ghost_send_info) delete [] ghost_send_info;
    ghost_send_info = (GF_Interaction***) 0;
    if (ghost_recv_cnt) delete [] ghost_recv_cnt;
    ghost_recv_cnt = 0;
    if (ghost_recv_info) delete [] ghost_recv_info;
    ghost_recv_info = (GF_Interaction***) 0;
    if (ghost_recv_server) delete [] ghost_recv_server;
    ghost_recv_server = (GridTableGhostRcv***) 0;
  }
  GF_DeleteGhostBndryInfo();
  
  GF_DeleteDataCommInfo();
  if (data_recv_server) delete [] data_recv_server;
  data_recv_server = (GridTableDataRcv***) 0;
  
  if (gt) delete gt;
  
  dagh.DAGH_DelGridFunction(gfid);
  
  delete myargs;
  delete [] length;
}
   
/*****************************************************************************/
/**** DeleteInfo ***/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_DeleteGDBStorage() {
  const int timesteps = 2*time_sten_rad+1;
  if (gdb) {
    for (register int t = 0 ; t < timesteps ; t++) 
      GF_DeleteGDBStorage(t);
    delete [] gdb; 
    gdb = (GridDataBlock(1)<DAGH_GFType> ****) 0;
  }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_DeleteGDBStorage(const int t) {
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();

  if (gdb) if (gdb[t]) {
    for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
      GF_DeleteGDBStorage(t, l);
    delete [] gdb[t]; 
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) 0;
  }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_DeleteGDBStorage(const int t, const int l) {
  assert (l < dagh.totallevels());
  
  if (gdb) if (gdb[t]) 
    if (gdb[t][l]) {
      for (register int i = 0 ; i < length[l] ; i++) 
        if (gdb[t][l][i]) { 
          delete gdb[t][l][i]; 
          gdb[t][l][i] = (GridDataBlock(1)<DAGH_GFType> *) 0; 
        }
      delete [] gdb[t][l]; 
      gdb[t][l] = (GridDataBlock(1)<DAGH_GFType> **) 0;
    }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CreateGDBStorage(const int t) {
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();

  assert (t!=time_sten_rad);
  if (gdb) if (!gdb[t]) {
    gdb[t] = new GridDataBlock(1)<DAGH_GFType> **[levels];
    for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      gdb[t][l] = (GridDataBlock(1)<DAGH_GFType> **) NULL;
      GF_CreateGDBStorage(t, l);
    }
  }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CreateGDBStorage(const int t, const int l) {
  assert (l < dagh.totallevels());
  
  if (gdb) if (gdb[t]) 
    if (!gdb[t][l]) {
      GridBoxList* lgbl = dagh.lgbl(l);
      assert (length[l] == Length(lgbl));
      GridBox** mgb = new GridBox*[length[l]];
      register int m=0; 
      for (register GridBox* tmpgb=lgbl->first();tmpgb;tmpgb=lgbl->next()) 
        mgb[m] = tmpgb; m++;
      
      gdb[t][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];
      const int mindex = dagh.GlobalMaxIndex(l)+1;
      for (register int i=0;i<length[l];i++) {
        if (gdb[time_sten_rad][l][i]) {      
          gdb[t][l][i] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt,
                                                               interactions,
                                                               i, *mgb[i],
                                                               t, l, bwidth, 
                                                               bndrybboxlist[l],
                                                               prolongbboxlist[l],
                                                               extghostwidth,
                                                               overlap,
                                                               space_sten_rad,
                                                               alignment, mindex,
                                                               comm(t,l), 
                                                               gdb[time_sten_rad][l][i]);
             
        }
        else
          gdb[t][l][i] = (GridDataBlock(1)<DAGH_GFType>*) NULL;
      } 
      delete [] mgb;
    }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_DeleteGhostCommInfo() {
  const int pnum = comm_service::proc_num();

  if (comm()) {
    const int levels = dagh.totallevels();
    const int cnt = DAGHMaxAxis * DAGHMaxDirs * levels;
    
    register int i = 0, j = 0;
    
    for (i=pnum-1;i>=0;i--) {
      if (ghost_send_info[i]) {
        for (j=cnt-1;j>=0;j--) {
          if (ghost_send_info[i] && ghost_send_info[i][j]) {
            delete ghost_send_info[i][j];
            ghost_send_info[i][j] = (GF_Interaction*) 0;
          }
        }
        delete [] ghost_send_info[i];
        ghost_send_info[i] = (GF_Interaction**) 0;
      }
      if (ghost_recv_info[i] || ghost_recv_server[i]) {
        for (j=cnt-1;j>=0;j--) {
          if (ghost_recv_info[i] && ghost_recv_info[i][j]) {
            delete ghost_recv_info[i][j];
            ghost_recv_info[i][j] = (GF_Interaction*) 0;
          }
          if (ghost_recv_server[i] && ghost_recv_server[i][j]) {
            delete ghost_recv_server[i][j];
            ghost_recv_server[i][j] = (GridTableGhostRcv*) 0;
          }
        }
        if (ghost_recv_info[i]) delete [] ghost_recv_info[i];
        ghost_recv_info[i] = (GF_Interaction**) 0;
        if (ghost_recv_server[i]) delete [] ghost_recv_server[i];
        ghost_recv_server[i] = (GridTableGhostRcv**) 0;
      }
    }
  }
}
 
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_DeleteDataCommInfo() {
  const int pnum = comm_service::proc_num();
  
  const int timesteps = 2*time_sten_rad+1;
  for (register int i=pnum-1;i>=0;i--) {
    if (data_recv_server[i]) {
      for (register int j=timesteps-1;j>=0;j--) {
        if (data_recv_server[i][j]) {
          delete data_recv_server[i][j];
          data_recv_server[i][j] = (GridTableDataRcv*) 0;
        }
      }
      delete [] data_recv_server[i];
      data_recv_server[i] = (GridTableDataRcv**) 0;
    }
  }
}

/*****************************************************************************/
/**** Compose() ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Compose() {
  AllocError::SetTexts("GridFunction(1)","GF_Compose()");      
  
  if (comm_flag >= DAGHTemplateComm) {
    template_flag = comm_flag - DAGHTemplateComm;
#ifdef DEBUG_PRINT
    assert(dagh.gflist[template_flag]);
#endif
    comm_flag = dagh.gflist[template_flag]->comm_type();
  }

  if (composed()) return;
  /* set up the interactions */
  interactions.compute_interactions(GhostInteraction::face_flag(comm_flag),
                                    GhostInteraction::corner_flag(comm_flag),
                                    GhostInteraction::edge_flag(comm_flag),
                                    GhostInteraction::simple_flag(comm_flag));

  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();

  register int t = 0, l = 0, i = 0;  
  const int tsteps = 2*time_sten_rad+1;
  
  gt = new GridTable(gfid,levels,tsteps);
  
  GridBoxList** lgbl = dagh.lgbl();
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
    length[l] = Length(lgbl[l]);

#ifdef DEBUG_PRINT_GF_GT
  ( comm_service::log() << "\n************* Grid Table *************\n"
                        << "------ [" << gfname << "] ------" << "\n"
                        << *gt
                        << "\n************* ***************** *************\n"
                          ).flush();
#endif

  /* If memory is already allocated, delete it.... */
  GF_DeleteGDBStorage();

  /* Allocate memory */
  gdb = new GridDataBlock(1)<DAGH_GFType> ***[tsteps];
  for (t=0;t<tsteps;t++) {
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
  }
  for (t=0;t<tsteps;t++) {
    if (!gdb[time_alias[t]]) {
      gdb[time_alias[t]] = new GridDataBlock(1)<DAGH_GFType> **[levels];
      for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
        gdb[time_alias[t]][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];
        for (i=0;i<length[l];i++) 
          gdb[time_alias[t]][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
      }
    }
  }
   
  GridFunction(1)<DAGH_GFType>* gftemplate = 0;
  if (template_flag != DAGHNull) {
    gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
    if (gftemplate && !gftemplate->composed()) gftemplate = 0;
  }
  const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;

  GF_GatherGhostBndryInfo();

  register int gbi = 0;
  GridBox* gb = 0;
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    if (!lgbl[l]) continue;
    const int mindex = dagh.GlobalMaxIndex(l)+1;
    for (gb = lgbl[l]->first(), gbi = 0; gb; gb = lgbl[l]->next(), gbi++) {
      const int ct = dagh.getCurrentTime(l);
      int tmplcti = DAGHNull;
      if (template_flag != DAGHNull && 
          gftemplate && 
          gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
          && gftemplate->gdb[tmplcti][l])
        gdbtmpl = gftemplate->gdb[tmplcti][l][gbi];
      else gdbtmpl = 0;
      for (int tc=0;tc<tsteps;tc++) { 
        t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
        if (gdb[t]) if (gdb[t][l]) {
          gdb[t][l][gbi] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt, 
                                                                 interactions,
                                                                 gbi, *gb,
                                                                 t, l, bwidth,
                                                                 bndrybboxlist[l],
                                                                 prolongbboxlist[l],
                                                                 extghostwidth,
                                                                 overlap, 
                                                                 space_sten_rad, 
                                                                 alignment, mindex,
                                                                 comm(t,l), 
                                                                 gdbtmpl);
          if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][gbi];
        }
      }
    }
  }

  if (comm()) {
    GF_GatherGhostCommInfo();
    GF_SetUpGhostCommServers();
  }
  compose_flag = DAGHTrue;
}

/*****************************************************************************/
/**** Compose(const int T, const int L) ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Compose(const int T, const int L) {
  AllocError::SetTexts("GridFunction(1)","GF_Compose(const int T, const int L)");      

  if (comm_flag >= DAGHTemplateComm) {
    template_flag = comm_flag - DAGHTemplateComm;
#ifdef DEBUG_PRINT
    assert(dagh.gflist[template_flag]);
#endif
    comm_flag = dagh.gflist[template_flag]->comm_type();
  }

  if (composed()) return;
  /* set up the interactions */
  interactions.compute_interactions(GhostInteraction::face_flag(comm_flag),
                                    GhostInteraction::corner_flag(comm_flag),
                                    GhostInteraction::edge_flag(comm_flag),
                                    GhostInteraction::simple_flag(comm_flag));

  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();

  register int t = 0, l = 0, i = 0;

  const int tsteps = 2*time_sten_rad+1;
  const int fl = dagh.finelevel();

  if (L == DAGHAllLevels) {
    gt = new GridTable(gfid,levels,tsteps);
  }
  else if (L <= fl) {
    gt = new GridTable(gfid,levels,tsteps);
  }
  else {
    // right now I will allow only one extra level
#ifdef DEBUG_PRINT
    assert((L-fl) == 1); 
#endif
    gt = new GridTable(gfid,levels,tsteps);
  }

#ifdef DEBUG_PRINT_GF_GT
  ( comm_service::log() << "\n************* Grid Table *************\n"
                        << "------ [" << gfname << "] ------" << "\n"
                        << *gt
                        << "\n************* ***************** *************\n"
                         ).flush();
#endif

  GridBoxList** lgbl = dagh.lgbl();
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
    if (L == DAGHAllLevels || L == l) 
      length[l] = Length(lgbl[l]);

  /* If memory is already allocated, delete it.... */
  GF_DeleteGDBStorage();

  /* Allocate memory */
  gdb = new GridDataBlock(1)<DAGH_GFType> ***[tsteps];
  for (t=0;t<tsteps;t++) {
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
  }
  for (t=0;t<tsteps;t++) {
    for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      if (T != DAGHAllTimes && dagh_timeindex(T,l) != time_alias[t]) continue;
      if (!gdb[time_alias[t]]) {
        gdb[time_alias[t]] = new GridDataBlock(1)<DAGH_GFType> **[levels];
        for (register int ll=Max(crslev,crslev+loffset); ll<Min(levels,levels+loffset); ll++)
          gdb[time_alias[t]][ll] = (GridDataBlock(1)<DAGH_GFType> **) NULL;
      }
      if ((L == DAGHAllLevels || L == l) && !gdb[time_alias[t]][l]) {
        gdb[time_alias[t]][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];

        for (i=0;i<length[l];i++)
          gdb[time_alias[t]][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
      }
    }
  }

  GridFunction(1)<DAGH_GFType>* gftemplate = 0;
  if (template_flag != DAGHNull) {
    gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
    if (gftemplate && !gftemplate->composed()) gftemplate = 0;
  }
  const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;

  GF_GatherGhostBndryInfo();

  register int gbi = 0;
  GridBox* gb = 0;
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    if (!lgbl[l]) continue;
    if (L != DAGHAllLevels && L != l) continue;
    const int mindex = dagh.GlobalMaxIndex(l)+1;
    const int ct = dagh.getCurrentTime(l);
    int tmplcti = DAGHNull;
    for (gb = lgbl[l]->first(), gbi = 0; gb; gb = lgbl[l]->next(), gbi++) {
      if (template_flag != DAGHNull && 
          gftemplate && 
          gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
          && gftemplate->gdb[tmplcti][l])
        gdbtmpl = gftemplate->gdb[tmplcti][l][gbi];
      else gdbtmpl = 0;
      for (int tc=0;tc<tsteps;tc++) { 
        t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
        if (gdb[t]) if (gdb[t][l]) {
          gdb[t][l][gbi] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt, 
                                                                 interactions,
                                                                 gbi, *gb,
                                                                 t, l, bwidth, 
                                                                 bndrybboxlist[l],
                                                                 prolongbboxlist[l],
                                                                 extghostwidth,
                                                                 overlap, 
                                                                 space_sten_rad,  
                                                                 alignment, mindex,
                                                                 comm(t,l), 
                                                                 gdbtmpl);
          if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][gbi];
        }
      }
    }
  }

  if (comm()) {
    GF_GatherGhostCommInfo();
    GF_SetUpGhostCommServers();
  }

  compose_flag = DAGHTrue;
}

/*****************************************************************************/
/*** Compose(const GridFunction(1)<DAGH_GFType>&const int T, const int L) ***/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Compose(
                                        const GridFunction(1)<DAGH_GFType>& gf,
                                        const int T, const int L) {
  AllocError::SetTexts("GridFunction(1)",
                       "GF_Compose(GridFunction(1)&gf,const int T, const int L)");      

  if (comm_flag >= DAGHTemplateComm) {
    template_flag = comm_flag - DAGHTemplateComm;
#ifdef DEBUG_PRINT
    assert(dagh.gflist[template_flag]);
#endif
    comm_flag = dagh.gflist[template_flag]->comm_type();
  }

  if (composed()) return;
  /* set up the interactions */
  interactions.compute_interactions(GhostInteraction::face_flag(comm_flag),
                                    GhostInteraction::corner_flag(comm_flag),
                                    GhostInteraction::edge_flag(comm_flag),
                                    GhostInteraction::simple_flag(comm_flag));

  const int me = comm_service::proc_me();
  const int pnum = comm_service::proc_num();
  
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();

  register int t = 0, l = 0, i = 0;
  const int tsteps = 2*time_sten_rad+1;
  const int fl = dagh.finelevel();

  if (L == DAGHAllLevels) {
    gt = new GridTable(gfid,levels,tsteps);
   }
  else if (L <= fl) {
    gt = new GridTable(gfid,levels,tsteps);
  }
  else {
    // right now I will allow only one extra level
#ifdef DEBUG_PRINT
    assert((L-fl) == 1); 
#endif
    gt = new GridTable(gfid,levels,tsteps);
  }

  length = gf.length;

#ifdef DEBUG_PRINT_GF_GT
  ( comm_service::log() << "\n************* Grid Table *************\n"
                        << "------ [" << gfname << "] ------" << "\n"
                        << *gt
                        << "\n************* ***************** *************\n"
                         ).flush();
#endif

  /* If memory is already allocated, delete it.... */
  GF_DeleteGDBStorage();

  /* Allocate memory */
  gdb = new GridDataBlock(1)<DAGH_GFType> ***[tsteps];
  for (t=0;t<tsteps;t++) {
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
  }
  for (t=0;t<tsteps;t++) {
    for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      if (T != DAGHAllTimes && dagh_timeindex(T,l) != time_alias[t]) continue;
      if (!gdb[time_alias[t]]) {
        gdb[time_alias[t]] = new GridDataBlock(1)<DAGH_GFType> **[levels];
        for (register int ll=Max(crslev,crslev+loffset); ll<Min(levels,levels+loffset); ll++)
          gdb[time_alias[t]][ll] = (GridDataBlock(1)<DAGH_GFType> **) NULL;
      }
      if ((L == DAGHAllLevels || L == l) && !gdb[time_alias[t]][l]) {
        gdb[time_alias[t]][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];

        for (i=0;i<length[l];i++)
          gdb[time_alias[t]][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
      }
    }
  }

  GridFunction(1)<DAGH_GFType>* gftemplate = 0;
  if (template_flag != DAGHNull) {
    gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
    if (gftemplate && !gftemplate->composed()) gftemplate = 0;
  }
  const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;

  GridBoxList** lgbl = dagh.lgbl();
  GridBox* gb = 0;
  GDB_Interaction* pi;     
  int npi;

  GF_GatherGhostBndryInfo();

  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    if (L != DAGHAllLevels && L != l) continue;
    if (!lgbl[l]) continue;
    const int mindex = dagh.GlobalMaxIndex(l)+1;
    const int cti = dagh_timeindex(dagh.getCurrentTime(l),l);
    const int ct = dagh.getCurrentTime(l);
    int tmplcti = DAGHNull;
    gb = lgbl[l]->first();
    for (i=0;i<length[l];i++) {
      if (template_flag != DAGHNull && 
          gftemplate && 
          gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
          && gftemplate->gdb[tmplcti][l])
        gdbtmpl = gftemplate->gdb[tmplcti][l][i];
      else gdbtmpl = 0;
      for (int tc=0;tc<tsteps;tc++) { 
        t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
        const int ll = (l < fl) ? l : fl;
        if (gdb[t]) if (gdb[t][l]) 
          if (gf.gdb[t]) if (gf.gdb[t][ll]) if (gf.gdb[t][ll][i]) { 
            gdb[t][l][i] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt, 
                                                                 interactions,
                                                                 i, *gb,
                                                                 t, l, bwidth, 
                                                                 bndrybboxlist[l],
                                                                 prolongbboxlist[l],
                                                                 extghostwidth,
                                                                 overlap, 
                                                                 space_sten_rad, 
                                                                 alignment, mindex,
                                                                 comm(t,l), gdbtmpl);

            if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][i];

            // Initialize data from the other GF ....
            if (DAGHInitGFOnCreation == DAGHTrue) { 
              if (l-1==ll && pfunc && prolong() && 
                  (npi=gdb[t][l][i]->gdb_parent_index)>0 && 
                  (pi=gdb[t][l][i]->gdb_parent_info)) { /* prolong */
                for (register int j=0; j<npi; j++)
                  if (gf.gdb[t][ll][pi[j].idx]) 
                    (*pfunc)(FORTRAN_ARGS(1,gf.gdb[t][ll][pi[j].idx]->data),
                             FORTRAN_ARGS(1,gdb[t][l][i]->data),
                             BOUNDING_BOX(pi[j].bbox),myargs,&myargc);
              }
              else if (l == ll)  /* copy... */
                gdb[t][l][i]->data.copy(gf.gdb[t][ll][i]->data);
              
            }
          } // end inititalize from other GF
      }
      gb = lgbl[l]->next();
    }
  }

  if (comm()) {
    GF_GatherGhostCommInfo();
    GF_SetUpGhostCommServers();
  }
  
  compose_flag = DAGHTrue;
}

/*****************************************************************************/
/**** .. & Recompose  ****/ 
/* Needs consistent data. Ensure calling Syncs() before recomposing !            */
/* Reads and writes also values of ghostcells. During redistribution the values */
/* of ghostcells might be copied into interior cells. This is why they have to be */
/* consistent.  */
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Recompose(int* reuse_possible,
                                                    GridBoxList** ollist,
                                                    GridBoxList** rlist, 
                                                    GridBoxList** slist,
                                                    GridBoxList** olist)

{
#ifdef DEBUG_PRINT_GF_RG
  ( comm_service::log() << "\n************* GF_Recompose() [" << gfid << "] *************\n").flush();
#endif
   AllocError::SetTexts("GridFunction(1)", "GF_Recompose(...)");    
   const int me = comm_service::proc_me();
   const int pnum = comm_service::proc_num();

   register int t=0,i=0,oi=0,p=0,m=0;
   int l;
   register GridBox* tmpgb = 0;

   GridDataBlock(1)<DAGH_GFType> ****oldgdb = gdb; gdb = 0;
   int* oldlength = length;

   gt->resettable();

   const int crslev = dagh.coarselevel();
   const int levels = dagh.totallevels();
   const int timesteps = 2*time_sten_rad+1;
   length = new int[levels];
   GridBoxList** lgbl = dagh.lgbl();
   GridBox*** mgb = new GridBox**[levels];
   for (l=0; l<levels; l++) {
     if (l>=Max(crslev,crslev+loffset) && l<Min(levels,levels+loffset))
       length[l] = Length(lgbl[l]);
     else
       length[l] = 0;
     mgb[l] = new GridBox*[length[l]];
     if (length[l]>0)
       for (tmpgb=lgbl[l]->first(),m=0;tmpgb;tmpgb=lgbl[l]->next()) {
         mgb[l][m] = tmpgb; m++;
       }
   }
   int* reuse = new int[levels];
   for (l=0; l<levels; l++) reuse[l] = 0;
   for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
     for (t=0;t<timesteps;t++) {
       if (!oldgdb[t]) break;
       if (!oldgdb[t][l]) break;
     }
     if (t==timesteps) reuse[l] = reuse_possible[l];
   }

   gdb = new GridDataBlock(1)<DAGH_GFType> ***[timesteps];
   for (t=0;t<timesteps;t++) {
     gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
     if (oldgdb[t]) {
       gdb[t] = new GridDataBlock(1)<DAGH_GFType> **[levels];
       for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
         gdb[t][l] = (GridDataBlock(1)<DAGH_GFType> **) NULL;
         if (oldgdb[t][l]) {
           if (reuse[l])
             gdb[t][l] = oldgdb[t][l];
           else {
             gdb[t][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];
             for (i=0;i<length[l];i++)
               gdb[t][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
           }
         }
       }
     }
   }

   /* destroy old storage for data comm */
   GF_DeleteDataCommInfo();     

   /* destroy old storage for ghost comm */
   GF_DeleteGhostCommInfo();
   GF_DeleteGhostBndryInfo();
   
   GF_GatherGhostBndryInfo();

   register int idx = 0;
   register int oldowner = DAGHNoBody;

   /* Comm: Post Rcvs */
   if (comm_service::dce() && pnum > 1) {
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n*************     Recv List     *************\n" ).flush();
#endif
     if (rlist) {

       idx = 0;
       oldowner = DAGHNoBody;
     
       if (template_flag == DAGHNull) {
    
         if (rcvsize) delete [] rcvsize;     
         rcvsize = new unsigned[timesteps*pnum];     
         for (p=0;p<timesteps*pnum;p++) 
           rcvsize[p] = 0;

         const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;
         
         for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
           if (!rlist[l] || !lgbl[l]) continue;
           if (rlist[l]->isempty()) continue;
#ifdef DEBUG_PRINT_GF_RG
           ( comm_service::log() << "   on l=" << l << "\n" << *rlist[l] ).flush();
#endif
           const int mindex = dagh.GlobalMaxIndex(l)+1;
           for (tmpgb=rlist[l]->first();tmpgb;
                tmpgb=rlist[l]->next()) {
             oldowner = tmpgb->gbOwner();
             idx = lgbl[l]->index(tmpgb->gbBBox());
             gdbtmpl = 0;
             for (int tc=0;tc<timesteps;tc++) { 
               t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
               if (!gdb[t]) continue; 
               if (!gdb[t][l]) continue;
               if (!gdb[t][l][idx]) {
                 gdb[t][l][idx] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt,
                                                                        interactions,
                                                                        idx, *mgb[l][idx],
                                                                        t, l, bwidth, 
                                                                        bndrybboxlist[l],
                                                                        prolongbboxlist[l],
                                                                        extghostwidth,
                                                                        overlap,
                                                                        space_sten_rad,
                                                                        alignment, mindex,
                                                                        comm(t,l), 
                                                                        gdbtmpl);
                 if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][idx];
               } 
               if (l<=lmax_recompose) {
                 BBox gb = usedbbox(*tmpgb,l,alignment);
                 rcvsize[oldowner+t*pnum] += gdhdr::gdbsize(sizeof(DAGH_GFType)*gb.size());

#ifdef DEBUG_PRINT_GF_RG
                 ( comm_service::log() << "GF_Recompose::RecvInfo "
                   << "[From:" << oldowner << "]"
                   << "[BB:" << gb << "]"
                   << "[Si:" << tmpgb->gbIndex() << "]"
                   << "[Ri:" << idx << "]"
                   << "[Sl:" << l << "]"
                   << endl ).flush();
#endif
               }
             }
           }
         }

       } else {

         GridFunction(1)<DAGH_GFType>* gftemplate = 
           (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];

         rcvsize = gftemplate->rcvsize;
       }

       for (t=0;t<timesteps;t++) {
         if (!gdb[t]) continue;
         for (p=0;p<pnum;p++) {
           if (me == p) continue;
           if (rcvsize[p+t*pnum] > 0) {
             if (!data_recv_server[p]) {
               data_recv_server[p] = new GridTableDataRcv *[timesteps];
               for (i=0; i<timesteps; i++) data_recv_server[p][i] = 0; 
             }
             data_recv_server[p][t] = 
               new GridTableDataRcv(*gt, (DAGHDataTag | t), rcvsize[p+t*pnum], p);
             data_recv_server[p][t]->postrcv();
           }
         }
       }
       
       if (template_flag != DAGHNull) {
         rcvsize = 0;
       }
     }
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n************* ***************** *************\n" ).flush();
#endif
   } 
   
   /* Comm: Do Sends */
   if (comm_service::dce() && pnum > 1) {
     int newowner = DAGHNoBody;
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n*************     Send List     *************\n" ).flush();
#endif
     if (slist) {

       idx = 0;

       int mlindex = 0;
       for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
         if (slist[l]) if (!slist[l]->isempty())
           mlindex += slist[l]->number();
       
       if (template_flag == DAGHNull) {
         
         if (sndcnt) delete [] sndcnt;
         sndcnt = new short[timesteps*pnum];
         for (p=0;p<timesteps*pnum;p++) sndcnt[p] = 0;

         if (sndsize) delete [] sndsize;
         if (sndbbox) delete [] sndbbox;
         if (sndindex) delete [] sndindex;
         if (rcvindex) delete [] rcvindex;
         if (sndlevel) delete [] sndlevel;

         sndsize = new unsigned[timesteps*pnum*mlindex];
         sndbbox = new BBox[timesteps*pnum*mlindex];
         sndindex = new short[timesteps*pnum*mlindex];
         rcvindex = new short[timesteps*pnum*mlindex];
         sndlevel = new short[timesteps*pnum*mlindex];

         register int tmpid = 0;       
         for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
           if (!slist[l] || !ollist[l]) continue;
           if (slist[l]->isempty()) continue;
#ifdef DEBUG_PRINT_GF_RG
           ( comm_service::log() << "   on l=" << l << "\n" << *slist[l] ).flush();
#endif
           for(tmpgb=slist[l]->first();tmpgb;
               tmpgb=slist[l]->next()) {
             GridBox oldgb = ollist[l]->find(tmpgb->gbBBox());
             idx = oldgb.gbIndex();
             newowner = tmpgb->gbOwner();
             if (l>lmax_recompose) break;
             BBox gb = usedbbox(*tmpgb,l,alignment);
             for (t=0;t<timesteps;t++) {
               if (!gdb[t]) continue;
               if (!gdb[t][l]) continue;

               tmpid = t*pnum*mlindex + newowner*mlindex + sndcnt[newowner+t*pnum];
               sndbbox[tmpid] = gb;
               sndsize[tmpid] = sizeof(DAGH_GFType)*sndbbox[tmpid].size();
               sndindex[tmpid] = idx;
               rcvindex[tmpid] = tmpgb->gbIndex();
               sndlevel[tmpid] = l;
                   
#ifdef DEBUG_PRINT_GF_RG
               ( comm_service::log() << "GF_Recompose::SendInfo [" 
                 << sndcnt[newowner+t*pnum] << "]"
                 << "[To:" << newowner << "]"
                 << "[BB:" << sndbbox[tmpid] << "]"
                 << "[Sz:" << sndsize[tmpid] << "]"
                 << "[Si:" << sndindex[tmpid] << "]"
                 << "[Ri:" << rcvindex[tmpid] << "]"
                 << "[Sl:" << sndlevel[tmpid] << "]"
                 << endl ).flush();
#endif
               sndcnt[newowner+t*pnum]++;                
             }
           }
         }

       } else {

         GridFunction(1)<DAGH_GFType>* gftemplate = 
           (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];

         sndsize = gftemplate->sndsize;
         sndbbox = gftemplate->sndbbox;
         sndindex = gftemplate->sndindex;
         rcvindex = gftemplate->rcvindex;
         sndlevel = gftemplate->sndlevel;
         sndcnt = gftemplate->sndcnt;
       }

       register int s = 0;
       register int off = 0;
       for (t=0;t<timesteps;t++) {
         if (!oldgdb[t] || !gdb[t]) continue;
         for (p=0;p<pnum;p++) {
           if (p == me) continue;
               
#ifdef DEBUG_PRINT_GF_RG
           ( comm_service::log() << "sndcnt[p]:" << sndcnt[p+t*pnum] << endl ).flush();
#endif   
           
           if (sndcnt[p+t*pnum] == 0) continue;
           off = t*pnum*mlindex + p*mlindex;
           GridDataBucket<DAGH_GFType> *gdbkt = 
             new GridDataBucket<DAGH_GFType>(sndcnt[p+t*pnum],(sndsize+off),DAGHPacked);
           s = 0;
           while (s<sndcnt[p+t*pnum]) {
#ifdef DEBUG_PRINT
             assert(oldgdb[t][sndlevel[off+s]][sndindex[off+s]]);
#endif
             const int time_value = dagh_timevalue(time_sten_rad,sndlevel[off+s]);
             oldgdb[t][sndlevel[off+s]][sndindex[off+s]]->gdbWriteData(time_value,sndbbox[off+s],
                                                                       *gdbkt,s);
             (gdbkt->head(s))->index = rcvindex[off+s];
             (gdbkt->head(s))->level = sndlevel[off+s];
#ifdef DEBUG_PRINT_GF_RG
             ( comm_service::log() << "GF_Recompose::Sending"
               << *gdbkt->head(s) << " called BBox: " << sndbbox[off+s]
               << endl ).flush();
#endif
             s++;
           }
#ifdef DEBUG_PRINT_GF_RG
           ( comm_service::log() << "send" << endl ).flush();
#endif
           gt->send((DAGHDataTag | t),gdbkt,p);
         }
       }
     
       if (template_flag != DAGHNull) {       
         sndsize = 0;
         sndbbox = 0;
         sndindex = 0;
         rcvindex = 0;
         sndlevel = 0;
         sndcnt = 0;
       }
     } 
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n************* ***************** *************\n" ).flush();
#endif
   }

   /* Allocate remaining GDB's */
   {
     GridFunction(1)<DAGH_GFType>* gftemplate = 0;
     if (template_flag != DAGHNull) {
       gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
     }
     const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;
     
     for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
       if (!lgbl[l] || reuse[l]) continue;
       if (lgbl[l]->isempty()) continue;
       const int mindex = dagh.GlobalMaxIndex(l)+1;
       const int ct = dagh.getCurrentTime(l);
       int tmplcti = DAGHNull;       
       for (i = 0 ; i < length[l] ; i++) {
         if (template_flag != DAGHNull &&
             gftemplate && 
             gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
             && gftemplate->gdb[tmplcti][l])
           gdbtmpl = gftemplate->gdb[tmplcti][l][i];
         else gdbtmpl = 0;
         for (int tc=0;tc<timesteps;tc++) { 
           t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
           if (!gdb[t]) continue;
           if (!gdb[t][l]) continue;
           if (!gdb[t][l][i]) {
#ifdef DEBUG_PRINT_GF_RG
             ( comm_service::log() << "Creating for " << mgb[l][i]->gbBBox() 
               << " on l=" << l << " for l=" << l-loffset   ).flush();
#endif
             gdb[t][l][i] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt,
                                                                  interactions,
                                                                  i, *mgb[l][i],
                                                                  t, l, bwidth, 
                                                                  bndrybboxlist[l],
                                                                  prolongbboxlist[l],
                                                                  extghostwidth,
                                                                  overlap,
                                                                  space_sten_rad,
                                                                  alignment, mindex,
                                                                  comm(t,l), 
                                                                  gdbtmpl);
             assert (!gdb[t][l][i]->boundingbox().empty());
#ifdef DEBUG_PRINT_GF_RG
             ( comm_service::log() << "  Done " << gdb[t][l][i]->boundingbox() << "\n" ).flush();
#endif
             if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][i];
           }
         }
       }
     }
   }
   
   // Setup new serves
   if (comm()) {
     GF_GatherGhostCommInfo();
     GF_SetUpGhostCommServers();
   }

   /* All Data: interpolate and/or copy */
   {
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n************* Overlap List *************\n" ).flush();
#endif
     for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
       if (!lgbl[l] || reuse[l]) continue;
       if (lgbl[l]->isempty()) continue;
       if (l>lmax_recompose) break;
#ifdef DEBUG_PRINT_GF_RG
       ( comm_service::log() << "   on l=" << l << "\n" << *olist[l] ).flush();
#endif
       GDB_Interaction* pi;
       int npi;

       for (t=0;t<timesteps;t++) {
         if (!gdb[t]) continue;
         if (!gdb[t][l]) continue;
         
         for (i=0; i<length[l]; i++) if (l>Max(crslev,crslev+loffset) && pfunc && prolong())
           if (gdb[t][l][i] && gdb[t][l-1]) 
             if ((npi=gdb[t][l][i]->gdb_parent_index)>0 && 
                 (pi=gdb[t][l][i]->gdb_parent_info)) {
               for (register int j=0; j<npi; j++)
                 if (gdb[t][l-1][pi[j].idx]) {
#ifdef DEBUG_PRINT_GF_RG
                   ( comm_service::log() << "Prolonging to " << gdb[t][l][i]->boundingbox() 
                     << " from " << gdb[t][l-1][pi[j].idx]->boundingbox() << " using " 
                     << pi[j].bbox << "\n" ).flush();
#endif
                   if (t==dagh_timeindex(dagh.getCurrentTime(l-1),l-1) &&
                       t==dagh_timeindex(dagh.getCurrentTime(l),l))
                     (*pfunc)(FORTRAN_ARGS(1,gdb[t][l-1][pi[j].idx]->data),
                              FORTRAN_ARGS(1,gdb[t][l][i]->data),
                              BOUNDING_BOX(pi[j].bbox),
                              myargs,&myargc);
                   else
                     gdb[t][l][i]->data.equals(DAGH_GFType(0.0),pi[j].bbox);    
                 }             
             }

         /* Local copies - use only internal cells */
         for (tmpgb=olist[l]->first();
              tmpgb;tmpgb=olist[l]->next()) {
#ifdef DEBUG_PRINT
           assert (tmpgb->gbOwner() == me);
#endif
           oi = tmpgb->gbIndex();
           i = lgbl[l]->index(tmpgb->gbBBox());
           if (!gdb[t][l][i]) continue;
#ifdef DEBUG_PRINT
           assert(oldgdb[t][l][oi]);
#endif

#ifdef DEBUG_PRINT_GF_RG
           ( comm_service::log() << "Copying from " << oldgdb[t][l][oi]->boundingbox() 
             << " to " << gdb[t][l][i]->boundingbox() << " using " 
             << oldgdb[t][l][oi]->interiorbox() << "\n" ).flush();
#endif
           if (gdb[t][l][i]->boundingbox().intersects(
               oldgdb[t][l][oi]->boundingbox())) 
             gdb[t][l][i]->data.copy(oldgdb[t][l][oi]->data,
                                     oldgdb[t][l][oi]->interiorbox());
         }
        
         /* Receives*/
         if (comm_service::dce() && pnum > 1) {
           for (p=0; p<pnum; p++) {
             if (data_recv_server[p] && data_recv_server[p][t]) {
               if (!data_recv_server[p][t]->received())
                 comm_service::serve(*data_recv_server[p][t]->req());
               delete data_recv_server[p][t];
               data_recv_server[p][t] = 0;
             }
           }
           GF_ReadData(t,l);
         }
        
         /* Now synchronize all grids and set ghost cells */
         int tnum = dagh_timevalue(t,l);
         GF_BndryUpdate(tnum,l);
       }  
     }
#ifdef DEBUG_PRINT_GF_RG
     ( comm_service::log() << "\n************* ***************** *************\n" ).flush();
#endif
     /* Empty GridTable after ALL levels got their data. */
     for (t=0;t<timesteps;t++) {
       List<GridDataBucketVoid*>& lgdbkt = gt->data(t);
       GridDataBucketVoid** gdbkt = 0;
       DAGHListLoop(lgdbkt, gdbkt, GridDataBucketVoid*) {
         GridDataBucket<DAGH_GFType>* rcvbkt =
           (GridDataBucket<DAGH_GFType> *) *gdbkt;
         delete rcvbkt; rcvbkt = 0;
       } DAGHEndLoop
       lgdbkt.empty();
     }
     for (p=0;p<pnum;p++) {
       if (data_recv_server[p]) {
         delete [] data_recv_server[p];
         data_recv_server[p] = 0;
       }
     }
   }   
   
#ifdef DEBUG_PRINT_GF_MEMORY
   int newmem = 0; int newgrids=0;
   int oldmem = 0; int oldgrids=0;
   for (t=0;t<timesteps;t++) 
     if (gdb[t])   
       for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
         if (gdb[t][l]) 
           for (i=0;i<length[l];i++) 
             if (gdb[t][l][i]) { 
               newmem += gdb[t][l][i]->MemoryUsage();
               newgrids++;
             }
   if (oldgdb) {
     for (t=0;t<timesteps;t++) 
       if (oldgdb[t])   
         for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
           if (oldgdb[t][l]) 
             for (i=0;i<oldlength[l];i++) 
               if (oldgdb[t][l][i]) { 
                 oldmem += oldgdb[t][l][i]->MemoryUsage();
                 oldgrids++;
               }
   }
   ( comm_service::log() << "\n************* GF_Recompose() [" << gfid << "] *************\n"
     << "Memory GridDataBlocks: New=" << newmem 
     << " Bytes, Remaining Old=" << oldmem << " Bytes, Whole=" << newmem+oldmem 
     << "\nGrids: New=" << newgrids << " Remaining Old=" << oldgrids 
     << "\n************* ******************* *************\n"
     ).flush();
#endif

   /* Delete Old Data */
   if (oldgdb) {
     for (t=0;t<timesteps;t++) {
       if (oldgdb[t])   {
         for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
           if (oldgdb[t][l] && !reuse[l]) {
             for (i=0;i<oldlength[l];i++) {
               if (oldgdb[t][l][i]){ 
                 delete oldgdb[t][l][i];
               }
             }
             delete [] oldgdb[t][l];
           }
         }
         delete [] oldgdb[t];
       } 
     } 
     delete [] oldgdb;
   }
   delete [] oldlength;
   if (mgb) {
     for (l=0; l<levels; l++) 
       if (mgb[l]) delete [] mgb[l];
     delete [] mgb;
   }
   delete [] reuse;
}

/*****************************************************************************/
/**** .Recompose  using a checkpoint file...****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CheckpointRecompose() {

  AllocError::SetTexts("GridFunction(1)", "GF_CheckpointRecompose()");    
  
  register int t=0,l=0,i=0;
  
  const int levels = dagh.totallevels();
  const int timesteps = 2*time_sten_rad+1;
  const int crslev = dagh.coarselevel();
  
  GridDataBlock(1)<DAGH_GFType> ****oldgdb = gdb; gdb = 0;
  int* oldlength = new int[levels];
  for (i=0; i<levels; i++)
    oldlength[i] = length[i];

  gt->resettable();
   
  GridBoxList** lgbl = dagh.lgbl();
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
    length[l] = Length(lgbl[l]);
     
  gdb = new GridDataBlock(1)<DAGH_GFType> ***[timesteps];
  for (t=0;t<timesteps;t++) {
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
    if (oldgdb[t]) {
      gdb[t] = new GridDataBlock(1)<DAGH_GFType> **[levels];
      for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
        gdb[t][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];
        for (i=0;i<length[l];i++)
          gdb[t][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
      }
    }
  }

  /* destroy old storage for ghost comm */
  GF_DeleteGhostCommInfo();
  GF_DeleteGhostBndryInfo();
   
  /* destroy old storage for data comm */
  GF_DeleteDataCommInfo();

  GF_GatherGhostBndryInfo();

  /* Delete Old Data */
  if (oldgdb) {
    for (t=0;t<timesteps;t++) {
      if (oldgdb[t])   {
        for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
          if (oldgdb[t][l]) {
            for (i=0;i<oldlength[l];i++) 
              if (oldgdb[t][l][i]) delete oldgdb[t][l][i];
            
            delete [] oldgdb[t][l];
          }
        }
        delete [] oldgdb[t];
      } 
    } 
    delete [] oldgdb;
  }
  delete [] oldlength;

  /* Allocate GDB's */
  {
    GridFunction(1)<DAGH_GFType>* gftemplate = 0;
    if (template_flag != DAGHNull) {
      gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
    }
    const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;

    register GridBox* tmpgb = 0;
    register int idx = 0;
    for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      if (!lgbl[l]) continue;
      const int mindex = dagh.GlobalMaxIndex(l)+1;
      const int ct = dagh.getCurrentTime(l);
      int tmplcti = DAGHNull;

      for (tmpgb=lgbl[l]->first(),idx=0;tmpgb;tmpgb=lgbl[l]->next(),idx++) {   
        if (template_flag != DAGHNull &&
            gftemplate && 
            gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
            && gftemplate->gdb[tmplcti][l])
          gdbtmpl = gftemplate->gdb[tmplcti][l][idx];
        else gdbtmpl = 0;
        for (int tc=0;tc<timesteps;tc++) { 
          t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
          if (!gdb[t]) continue;
          if (!gdb[t][l]) continue;

          if (!gdb[t][l][idx]) {
            gdb[t][l][idx] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt,
                                                                   interactions,
                                                                   idx, *tmpgb,
                                                                   t, l, bwidth, 
                                                                   bndrybboxlist[l],
                                                                   prolongbboxlist[l],
                                                                   extghostwidth,
                                                                   overlap,
                                                                   space_sten_rad,
                                                                   alignment, mindex,
                                                                   comm(t,l), 
                                                                   gdbtmpl);
            if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][idx];
          }
        }
      }
    }
  }
  
  /* Setup new serves */
  if (comm()) {
    GF_GatherGhostCommInfo();
    GF_SetUpGhostCommServers();
  }
  
  /* Read data from checkpoint file */
  if (checkpoint()) GF_CheckpointRestart();
}


/*****************************************************************************/
/**** Recompose from memory ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CheckpointRecompose(strstream& ifs) {

  AllocError::SetTexts("GridFunction(1)",  "GF_CheckpointRecompose(strstream& ifs)");    
  
  register int t=0,l=0,i=0;

  const int levels = dagh.totallevels();
  const int timesteps = 2*time_sten_rad+1;  
  const int crslev = dagh.coarselevel();
  
  GridDataBlock(1)<DAGH_GFType> ****oldgdb = gdb; gdb = 0;
  int* oldlength = new int[levels];
  for (i=0; i<levels; i++)
    oldlength[i] = length[i];
  
  gt->resettable();
   
  GridBoxList** lgbl = dagh.lgbl();
  for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) 
    length[l] = Length(lgbl[l]);
  
  gdb = new GridDataBlock(1)<DAGH_GFType> ***[timesteps];
  for (t=0;t<timesteps;t++) {
    gdb[t] = (GridDataBlock(1)<DAGH_GFType> ***) NULL;
    if (oldgdb[t]) {
      gdb[t] = new GridDataBlock(1)<DAGH_GFType> **[levels];
      for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
        gdb[t][l] = new GridDataBlock(1)<DAGH_GFType> *[length[l]];
        for (i=0;i<length[l];i++)
          gdb[t][l][i] = (GridDataBlock(1)<DAGH_GFType> *) NULL;
      }
    }
  }

  /* destroy old storage for ghost comm */
  GF_DeleteGhostCommInfo();
  GF_DeleteGhostBndryInfo();
  
  /* destroy old storage for data comm */
  GF_DeleteDataCommInfo();
  
  GF_GatherGhostBndryInfo();

  /* Delete Old Data */
  if (oldgdb) {
    for (t=0;t<timesteps;t++) {
      if (oldgdb[t])   {
        for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
          if (oldgdb[t][l]) {
            for (i=0;i<oldlength[l];i++) 
              if (oldgdb[t][l][i]) delete oldgdb[t][l][i];

            delete [] oldgdb[t][l];
          }
        }
        delete [] oldgdb[t];
      } 
    } 
    delete [] oldgdb;
  }

  /* Allocate GDB's */
  {
    GridFunction(1)<DAGH_GFType>* gftemplate = 0;
    if (template_flag != DAGHNull) {
      gftemplate = (GridFunction(1)<DAGH_GFType>*) dagh.gflist[template_flag];
    }
    const GridDataBlock(1)<DAGH_GFType>* gdbtmpl = 0;
    
    register GridBox* tmpgb = 0;
    register int idx = 0;
    for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      if (!lgbl[l]) continue;
      const int mindex = dagh.GlobalMaxIndex(l)+1;
      const int ct = dagh.getCurrentTime(l);
      int tmplcti = DAGHNull;
      
      for (tmpgb=lgbl[l]->first(),idx=0;tmpgb;tmpgb=lgbl[l]->next(),idx++) {   
        if (template_flag != DAGHNull &&
            gftemplate && 
            gftemplate->gdb[(tmplcti=gftemplate->dagh_timeindex(ct,l))] 
            && gftemplate->gdb[tmplcti][l])
          gdbtmpl = gftemplate->gdb[tmplcti][l][idx];
        else gdbtmpl = 0;
        for (int tc=0;tc<timesteps;tc++) { 
          t = (tc<=time_sten_rad ? tc+time_sten_rad : tc-time_sten_rad-1);
          if (!gdb[t]) continue;
          if (!gdb[t][l]) continue;

          if (!gdb[t][l][idx]) {
            gdb[t][l][idx] = new GridDataBlock(1)<DAGH_GFType>(*this, dagh, *gt,
                                                                   interactions,
                                                                   idx, *tmpgb,
                                                                   t, l, bwidth, 
                                                                   bndrybboxlist[l],
                                                                   prolongbboxlist[l],
                                                                   extghostwidth,
                                                                   overlap,
                                                                   space_sten_rad,
                                                                   alignment, mindex,
                                                                   comm(t,l), 
                                                                   gdbtmpl);
            if (gdbtmpl == 0 || comm(t,l)) gdbtmpl = gdb[t][l][idx];
          }
        }
      }
    }
  }

  /* Setup new serves */
  if (comm()) {
    GF_GatherGhostCommInfo();
    GF_SetUpGhostCommServers();
  }
  
  /* Read data from checkpoint file */
  if (checkpoint()) GF_CheckpointRestart(ifs);
}

/*****************************************************************************/
/* Set Up Ghost Comm Servers */
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_SetUpGhostCommServers() {
  if (!comm()) return;

  const int pnum = comm_service::proc_num();
  const int me = comm_service::proc_me();  
  const int levels = dagh.totallevels();
  const int cnt = levels * DAGHMaxAxis * DAGHMaxDirs ;
  
  for (register int i=0;i<pnum;i++) {
    if (i == me || !ghost_recv_info[i]) continue;
    if (!ghost_recv_server[i]) {
      ghost_recv_server[i] = new GridTableGhostRcv *[cnt];
      for (register int c=0; c<cnt; c++) ghost_recv_server[i][c] = 0;
    }
    for (register int j=0;j<cnt;j++) {
      if (!ghost_recv_info[i][j]) continue;
      const GF_Interaction* gri = ghost_recv_info[i][j];
      ghost_recv_server[i][j] =
        new GridTableGhostRcv(*gt, (DAGHGhostTag | j), gri->tsize, i);
    }
  }
}

/*****************************************************************************/
/**** Gather Ghost Communication Info ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_GatherGhostCommInfo() {
  if (!comm()) return;
  
  const int me = comm_service::proc_me();
  
  register int t=0, l=0, c=0, i=0, p=0, m=0, k=0;
  for (t=0;t<=2*time_sten_rad;t++) if (gdb[t]) break;
  
  const int pnum = comm_service::proc_num();
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();
  const int cnt1 = levels * DAGHMaxAxis * DAGHMaxDirs ;
  const int cnt2 = DAGHMaxAxis * DAGHMaxDirs ;

  if (ghost_recv_cnt) delete [] ghost_recv_cnt;
  ghost_recv_cnt = new short[cnt1];
  for (k=0;k<cnt1;k++) ghost_recv_cnt[k] = 0;
  
  List<unsigned> rsize;
  List<unsigned> wsize;
  unsigned r=0, w=0;

  GridDataBlock(1)<DAGH_GFType> *g = 0;
  GDB_Interaction *gdbi;
  int ngdbi;
  
  for (p=0;p<pnum;p++) {
    for (l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      for (c=0;c<cnt2;c++) {
        const int idx1 = l*cnt2 + c;
        for (i=0; i<length[l]; i++) {
          if (gdb[t]) if (gdb[t][l]) if (gdb[t][l][i]) {
            g = gdb[t][l][i];
            /* Write Info */
            if (g->gdb_write_info[p] && g->gdb_write_index[p])
              if ((ngdbi=g->gdb_write_index[p][c])>0 && 
                  (gdbi=g->gdb_write_info[p][c]))
                for (m=0; m<ngdbi; m++) {
                  /* Allocate storage */
                  if (!ghost_send_info[p]) {
                    ghost_send_info[p] = new GF_Interaction* [cnt1];
                    for (k=0;k<cnt1;k++) ghost_send_info[p][k] = 0;
                  }
                  if (!ghost_send_info[p][idx1]) 
                    ghost_send_info[p][idx1] = new GF_Interaction;

#ifdef DEBUG_PRINT_GF_COMM
                  ( comm_service::log() << "sendbox " << idx1 << " " 
                    << gdbi[m].bbox << "  " << gdbi[m].idx << endl ).flush();
#endif
                  w = sizeof(DAGH_GFType)*gdbi[m].bbox.size();
                  ghost_send_info[p][idx1]->tsize += gdhdr::gdbsize(w);
#ifdef DEBUG_PRINT
                  assert( w != 0 );
#endif
                  wsize.add(w);
                }
              
            /* Read Info */
            if (g->gdb_read_info[p] && g->gdb_read_index[p])
              if ((ngdbi=g->gdb_read_index[p][c])>0 && 
                  (gdbi=g->gdb_read_info[p][c]))
                for (m=0; m<ngdbi; m++) {
                  /* Allocate storage */
                  if (!ghost_recv_info[p]) {
                    ghost_recv_info[p] = new GF_Interaction* [cnt1];
                    for (k=0;k<cnt1;k++) ghost_recv_info[p][k] = 0;
                  }
                  if (!ghost_recv_info[p][idx1]) 
                    ghost_recv_info[p][idx1] = new GF_Interaction;

#ifdef DEBUG_PRINT_GF_COMM
                  ( comm_service::log() << "receivebox " << idx1 << " " 
                    << gdbi[m].bbox << "  " << gdbi[m].idx << endl ).flush();
#endif
  
                  r = sizeof(DAGH_GFType)*gdbi[m].bbox.size();
                  ghost_recv_info[p][idx1]->tsize += gdhdr::gdbsize(r);
#ifdef DEBUG_PRINT
                  assert( r != 0 );
#endif
                  rsize.add(r); 
                }
          }
        }
        if (ghost_send_info[p] && ghost_send_info[p][idx1])
          wsize.array(ghost_send_info[p][idx1]->size,ghost_send_info[p][idx1]->cnt);
        if (ghost_recv_info[p] && ghost_recv_info[p][idx1]) {
          rsize.array(ghost_recv_info[p][idx1]->size,ghost_recv_info[p][idx1]->cnt);
          if (p != me) ghost_recv_cnt[idx1] += ghost_recv_info[p][idx1]->cnt;
        }
        wsize.empty(); rsize.empty();
      }
    }
#ifdef DEBUG_PRINT_GF_COMM
    ( comm_service::log() << "[GF_Interaction:[" << p << "]" << endl ).flush();
    for (k=0;k<cnt1;k++) {
      if (ghost_recv_info[p] && ghost_recv_info[p][k])
        ( comm_service::log() << "\t [Rcv[" << k << "]" 
                              << *ghost_recv_info[p][k] << "]"
                              << endl ).flush();
    }
    ( comm_service::log() << endl ).flush();
    for (k=0;k<cnt1;k++) {
      if (ghost_send_info[p] && ghost_send_info[p][k])
        ( comm_service::log() << "\t [Snd[" << k << "]" 
                              << *ghost_send_info[p][k] << "]"
                              << endl ).flush();
    }
    ( comm_service::log() << "]" << endl ).flush();
#endif 
  }
}

/*****************************************************************************/
/**** Swap storage for time levels  ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_SwapTimeLevels(const int l, 
                                                     const int t1, 
                                                     const int t2) {
  const int ti1 = dagh_timeindex(t1,l); 
  const int ti2 = dagh_timeindex(t2,l); 

  if (ti1 == ti2) return;

#ifdef DEBUG_PRINT
  assert(gdb[ti1] && gdb[ti2]);
#endif

  GridDataBlock(1)<DAGH_GFType>** gdb1 = gdb[ti1][l];
  GridDataBlock(1)<DAGH_GFType>** gdb2 = gdb[ti2][l];

  if (gdb1) {
    for (register int i=0;i<length[l];i++)
      if (gdb1[i]) gdb1[i]->timenum = ti2;
  }
  if (gdb2) {
    for (register int i=0;i<length[l];i++)
      if (gdb2[i]) gdb2[i]->timenum = ti1;
  }
  gdb[ti1][l] = gdb2;
  gdb[ti2][l] = gdb1;
}

/* Move next -> current && current -> previous && previous -> next */
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CycleTimeLevels(const int l) {
  if (time_sten_rad == 0)  return;

  const int cti = dagh_timeindex(dagh.getCurrentTime(l),l);
  const int pti = dagh_timeindex(dagh.getPreviousTime(l),l);
  const int nti = dagh_timeindex(dagh.getNextTime(l),l);

#ifdef DEBUG_PRINT
  assert(gdb[cti] && gdb[pti] && gdb[nti]);
#endif

  GridDataBlock(1)<DAGH_GFType>** cgdb = gdb[cti][l];
  GridDataBlock(1)<DAGH_GFType>** pgdb = gdb[pti][l];
  GridDataBlock(1)<DAGH_GFType>** ngdb = gdb[nti][l];

  register int i = 0;

  if (cgdb && (cti != pti)) {
    for (i=0;i<length[l];i++)
      if (cgdb[i]) cgdb[i]->timenum = pti;
  }
  if (pgdb && (pti != nti)) {
    for (i=0;i<length[l];i++)
      if (pgdb[i]) pgdb[i]->timenum = nti;
  }
  if (ngdb && (nti != cti)) {
    for (i=0;i<length[l];i++)
      if (ngdb[i]) ngdb[i]->timenum = cti;
  }
  if (cti != nti) gdb[cti][l] = ngdb;
  if (pti != cti) gdb[pti][l] = cgdb;
  if (nti != pti) gdb[nti][l] = pgdb;
}

/*****************************************************************************/
/**** BBox Queries ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::intbboxlist(BBoxList &bbl, const int l) {
  if (!bbl.isempty()) bbl.empty();
  register int t = 0, i = 0;
  for (t=0;t<=2*time_sten_rad;t++) if (gdb[t]) break;
  for (i=0; i<length[l]; i++)
    if (gdb[t][l][i]) bbl.add(gdb[t][l][i]->interiorbox());
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::databboxlist(BBoxList &bbl, const int l) {
  if (!bbl.isempty()) bbl.empty();
  register int t = 0, i = 0;
  for (t=0;t<=2*time_sten_rad;t++) if (gdb[t]) break;
  for (i=0; i<length[l]; i++)
    if (gdb[t][l][i]) bbl.add(gdb[t][l][i]->databox());
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::boundingbboxlist(BBoxList &bbl, const int l) { 
  if (!bbl.isempty()) bbl.empty();
  register int t = 0, i = 0;
  for (t=0;t<=2*time_sten_rad;t++) if (gdb[t]) break;
  for (i=0; i<length[l]; i++)
    if (gdb[t][l][i]) bbl.add(gdb[t][l][i]->boundingbox());
}


/*****************************************************************************/
/**** Checkpoint/Restart ****/
/*****************************************************************************/
template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Checkpoint(ofstream& ofs) {
#ifdef DEBUG_PRINT
  assert (dagh.chkpt_restart() && checkpoint());
#endif

  const int crslev = dagh.coarselevel();
  const short levels = dagh.totallevels();
  const short timesteps = 2*time_sten_rad + 1;
  ofs.write((char*)&levels,sizeof(short));
  int maxlength = 0;
  for (register int i=0; i<levels; i++) 
    if (maxlength<length[i]) maxlength=length[i];
  ofs.write((char*)&maxlength,sizeof(int));
  ofs.write((char*)&timesteps,sizeof(short));

  const int dirnum = levels*maxlength*timesteps;
  streampos* chkptdir = new streampos[dirnum];
  streampos chkptdirpos = ofs.tellp();
  ofs.write((char*)chkptdir,sizeof(streampos)*dirnum);

  for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    for (register int t=0; t<timesteps; t++) {
      for (register int c=0; c<maxlength; c++) {
        chkptdir[c*levels*timesteps+l*timesteps+t] = DAGHNull;
        if (gdb[t]) if (gdb[t][l]) if (c<length[l] && gdb[t][l][c]) {
          chkptdir[c*levels*timesteps+l*timesteps+t] = ofs.tellp();
          ofs << *gdb[t][l][c];
        }
      }
    }
  }
  streampos curpos = ofs.tellp();
  ofs.seekp(chkptdirpos);
  ofs.write((char*)chkptdir,sizeof(streampos)*dirnum);
  ofs.seekp(curpos);
  if (chkptdir) delete [] chkptdir;
}

template <class DAGH_GFType>
int GridFunction(1)<DAGH_GFType>::GF_Checkpoint_StrStream_Memory() {
#ifdef DEBUG_PRINT
  assert (checkpoint());
#endif
  const int crslev = dagh.coarselevel();
  const short levels = dagh.totallevels();
  const short timesteps = 2*time_sten_rad + 1;
  int maxlength = 0;
  for (register int i=0; i<levels; i++) 
    if (maxlength<length[i]) maxlength=length[i];
  const int dirnum = levels*maxlength*timesteps;

  int mem = 0;
  mem += 2*sizeof(short)+sizeof(int);
  mem += sizeof(streampos)*dirnum;
  for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    int tc = dagh_timeindex(dagh.getCurrentTime(l),l);
    for (register int t=0; t<timesteps; t++) 
      for (register int c=0; c<maxlength; c++) 
        if (tc==t && gdb[t]) if (gdb[t][l]) if (c<length[l] && gdb[t][l][c]) 
          mem += gdb[t][l][c]->databox().size()*sizeof(DAGH_GFType)+
                 sizeof(BBox)+2*sizeof(short);    
  }
  return mem;
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_Checkpoint(strstream& ofs) {
#ifdef DEBUG_PRINT
  assert (dagh.chkpt_restart() && checkpoint());
#endif
  
  const int crslev = dagh.coarselevel();
  const short levels = dagh.totallevels();
  const short timesteps = 2*time_sten_rad + 1;
  ofs.write((char*)&levels,sizeof(short));
  int maxlength = 0;
  for (register int i=0; i<levels; i++) 
    if (maxlength<length[i]) maxlength=length[i];
  ofs.write((char*)&maxlength,sizeof(int));
  ofs.write((char*)&timesteps,sizeof(short));
  
  const int dirnum = levels*maxlength*timesteps;
  streampos* chkptdir = new streampos[dirnum];
  streampos chkptdirpos = ofs.tellp();
  ofs.write((char*)chkptdir,sizeof(streampos)*dirnum);

  for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    int tc = dagh_timeindex(dagh.getCurrentTime(l),l);
    for (register int t=0; t<timesteps; t++) {
      for (register int c=0; c<maxlength; c++) {
        chkptdir[c*levels*timesteps+l*timesteps+t] = DAGHNull;
        if (tc==t && gdb[t]) if (gdb[t][l]) if (c<length[l] && gdb[t][l][c]) {
          chkptdir[c*levels*timesteps+l*timesteps+t] = ofs.tellp();
          ofs << *gdb[t][l][c];
        }
      }
    }
  }
  streampos curpos = ofs.tellp();

  /* Set get-pointer on current end to allow seekp(curpos)     */
  /* Absolutely necessary for newer strstreams implementations */
  ofs.seekg((streampos)0,ios::end);

  ofs.seekp(chkptdirpos);
  ofs.write((char*)chkptdir,sizeof(streampos)*dirnum);   
  ofs.seekp(curpos);
  if (chkptdir) delete [] chkptdir;
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CheckpointRestart(int proc) {
#ifdef DEBUG_PRINT
  assert (dagh.chkpt_restart() && checkpoint());
#endif

  const int num = comm_service::proc_num();
  const int me = comm_service::proc_me();
  const int oldpnum = dagh.chkpt_pnum();
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();
  const int timesteps = 2*time_sten_rad + 1;
  
  if (oldpnum == num) { /* Old dist == new dist */
    register streampos pos;
    
    ifstream ifs;
    int eflag = dagh.DAGH_GetGFChkptIStream(me,gfname,gfid,ifs);
    if (eflag == DAGHFalse) return;
    
    short hislevels = 0, histsteps = 0; 
    int hismaxlength = 0;
    ifs.read((char*)&hislevels,sizeof(short));
#ifdef DEBUG_PRINT
    assert (hislevels == levels);
#endif
    ifs.read((char*)&hismaxlength,sizeof(int));
    ifs.read((char*)&histsteps,sizeof(short));
    
    const int dirnum = hislevels*hismaxlength*histsteps;
    streampos* chkptdir = new streampos[dirnum];
    ifs.read((char*)chkptdir,sizeof(streampos)*dirnum);    
    for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      for (register int t=0; t<timesteps; t++) {
        for (register int c=0; c<hismaxlength; c++) {
          if (gdb[t]) if (gdb[t][l]) if (c<length[l] && gdb[t][l][c] && t < histsteps &&
             (pos=chkptdir[c*hislevels*histsteps+l*histsteps+t])!=(streampos) DAGHNull) {
            ifs.seekg(pos);
            ifs >> *gdb[t][l][c];
          } 
        }
      }
    }
    if (chkptdir) delete [] chkptdir;
    dagh.DAGH_CloseChkptIStream(ifs);
  }
  else { /* A new dist. */
    register streampos pos;

    if (proc<0) proc=me;
    ifstream ifs;
    int eflag = dagh.DAGH_GetGFChkptIStream(proc,gfname,gfid,ifs);
    if (eflag == DAGHFalse) return;
    
    short hislevels = 0, histsteps = 0; 
    int hismaxlength = 0;
    ifs.read((char*)&hislevels,sizeof(short));
#ifdef DEBUG_PRINT
    assert (hislevels == levels);
#endif
    ifs.read((char*)&hismaxlength,sizeof(int));
    ifs.read((char*)&histsteps,sizeof(short));
    
    const int dirnum = hislevels*hismaxlength*histsteps;
    streampos* chkptdir = new streampos[dirnum];
    ifs.read((char*)chkptdir,sizeof(streampos)*dirnum);    
    for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
      for (register int t=0; t<timesteps; t++) {
        for (register int hisc=0; hisc<hismaxlength; hisc++) {
          if (gdb[t]) if (gdb[t][l]) if (t<histsteps &&
             (pos=chkptdir[hisc*hislevels*histsteps+l*histsteps+t])!=(streampos) DAGHNull) {
            ifs.seekg(pos);

            short histimenum, hislevelnum;
            ifs.read((char*)&histimenum,sizeof(short));
            ifs.read((char*)&hislevelnum,sizeof(short));
//          ( comm_service::log() << "Read in: " << hisc << "  " << " on l=" << hislevelnum 
//            << " for t=" << histimenum << "\n" ).flush();
   
#ifdef DEBUG_PRINT
            assert (histimenum==t && hislevelnum == l);
#endif
            GridData(1)<DAGH_GFType> gdbhelp;       
            ifs >> gdbhelp;
            for (register int c=0; c<length[l]; c++) 
              if (gdb[t][l][c]) {
                BBox bb = gdbhelp.bbox();
                bb.grow(-extghostwidth);                
                bb.shrinkupperbydim(overlap);
                bb.shrinkbydim(space_sten_rad);   
                gdb[t][l][c]->data.copy(gdbhelp,bb);
              }
          } 
        }
        /* Now synchronize all grids and set ghost cells */
        if (oldpnum-1 == proc) {
          int tnum = dagh_timevalue(t,l);
          GF_BndryUpdate(tnum,l);
        }
      }
    }
    if (chkptdir) delete [] chkptdir;
    dagh.DAGH_CloseChkptIStream(ifs);
  }
}

template <class DAGH_GFType>
void GridFunction(1)<DAGH_GFType>::GF_CheckpointRestart(strstream& ifs) {
#ifdef DEBUG_PRINT
  assert (dagh.chkpt_restart() && checkpoint());
#endif

  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();
  const int timesteps = 2*time_sten_rad + 1;
  
  register streampos pos;
   
  int eflag = dagh.DAGH_GetGFChkptStrStream(gfname,gfid,ifs);
  if (eflag == DAGHFalse) return;
   
  short hislevels = 0, histsteps = 0; 
  int hismaxlength = 0;
  ifs.read((char*)&hislevels,sizeof(short));
#ifdef DEBUG_PRINT
  assert (hislevels == levels);
#endif
  ifs.read((char*)&hismaxlength,sizeof(int));
  ifs.read((char*)&histsteps,sizeof(short));
  
  const int dirnum = hislevels*hismaxlength*histsteps;
  streampos* chkptdir = new streampos[dirnum];
  ifs.read((char*)chkptdir,sizeof(streampos)*dirnum);    
  for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
    int tc = dagh_timeindex(dagh.getCurrentTime(l),l);
    for (register int t=0; t<timesteps; t++) {
      for (register int c=0; c<hismaxlength; c++) {
        if (tc==t && gdb[t]) if (gdb[t][l]) if (c<length[l] && gdb[t][l][c])
          if (t < histsteps &&
              (pos=chkptdir[c*hislevels*histsteps+l*histsteps+t])!=(streampos) DAGHNull) {
            ifs.seekg(pos);
            ifs >> *gdb[t][l][c];
          }     
          else 
            gdb[t][l][c]->griddata() = 0;
      }
    }
  }
  
  if (chkptdir) delete [] chkptdir;
}

/*****************************************************************************/
/**** Some debugging ****/
/*****************************************************************************/
template <class DAGH_GFType>
ostream &GridFunction(1)<DAGH_GFType>::GF_DebugPrintIntData(ostream &os,
                                                            const int time, 
                                                            const int level) {
  os << *this;
  os << "(" << time << ") ";
  os << "(" << level << ") ";
  os << "\n";
  register const int t = dagh_timeindex(time,level);
  register const int l = level;
  for (register int i=0; i<length[l]; i++) { 
    if (gdb[t][l][i]) {
      os << "Interior: " << gdb[t][l][i]->boundingbox() << "\n";
      GridData(1)<DAGH_GFType> tmpgd(gdb[t][l][i]->boundingbox());
      tmpgd = gdb[t][l][i]->griddata();
      os << tmpgd;
      os << "\n";
    }
  }
  return os;
}

template <class DAGH_GFType>
ostream &GridFunction(1)<DAGH_GFType>::GF_DebugPrintData(ostream &os, 
                                                         const int time, 
                                                         const int level) {
  os << *this;
  os << "(" << time << ") ";
  os << "(" << level << ") ";
  os << "\n";
  register const int t = dagh_timeindex(time,level);
  register const int l = level;
  for (register int i=0; i<length[l]; i++)
    if (gdb[t][l][i]) os << gdb[t][l][i]->griddata();
  return os;
}

template <class DAGH_GFType>
ostream &GridFunction(1)<DAGH_GFType>::GF_DebugPrintDataBlk(ostream &os,
                                                            const int time, 
                                                            const int level) {
  os << *this;
  os << "(" << time << ") ";
  os << "(" << level << ") ";
  os << "\n";
  register const int t = dagh_timeindex(time,level);
  register const int l = level;
  for (register int i=0; i<length[l]; i++)
    if (gdb[t][l][i]) os << *gdb[t][l][i] << "\n";
  return os;
}

/*****************************************************************************/
/**** Overloaded stdout ****/
/*****************************************************************************/

template <class DAGH_GFType>
ostream&  operator << (ostream& os, const GridFunction(1)<DAGH_GFType>& gf) {

  if (&gf == (GridFunction(1)<DAGH_GFType> *) NULL) return os;
  
  os << "GridFunction: " << gf.GF_Name();
  
  return os;
}


template <class DAGH_GFType>
ostream &GridFunction(1)<DAGH_GFType>::GF_DebugPrintTheData(ostream &os) {
    
  const int crslev = dagh.coarselevel();
  const int levels = dagh.totallevels();
  const int timesteps = 2*time_sten_rad+1;

  os << "******* THE DATA *********" << endl;
  if (gdb) {
    for (register int t=0;t<timesteps;t++) {
      if (gdb[t])   {
        for (register int l=Max(crslev,crslev+loffset); l<Min(levels,levels+loffset); l++) {
          if (gdb[t][l]) {
            for (register int i=0;i<length[l];i++) {
              if (gdb[t][l][i]){ 
                os << "gdb[ " << t << "][" 
                   << l << "][" << i << "]:" << endl;
                os << *gdb[t][l][i] << endl;
              }
            }
          }
        }
      } 
    } 
  }
  os << "******* END THE DATA *********" << endl;
  return os;
}

#include "GridFunctionComm1.h"
#include "GridFunctionFuncs1.h"
#include "GridFunctionOps1.h"
#include "GridFunctionOpsRed1.h"
#include "GridFunctionBndry1.h"
#include "GridFunctionIO1.h"

#endif

---

Quickstart     Users Guide     Programmers Reference     Installation      Examples     Download

---

AMROC Main      Home      Contact

last update: 06/01/04