amroc/amr/SpecialOutput/AMRGridAccumulation.h

Go to the documentation of this file.

// -*- C++ -*-

#ifndef AMROC_AMR_GRIDACCUMULATION_H
#define AMROC_AMR_GRIDACCUMULATION_H

#include <string>
#include <ios>
#include <fstream>
#include <cstdio>
#include <cmath>
#include <unistd.h>
#include <sys/stat.h>
#include "SpecialOutput/AMRAccumulation.h"

template <class DataType, int dim>
class AMRGridAccumulation : public AMRAccumulation<DataType,dim> {
  typedef AMRAccumulation<DataType,dim> base;
public:
  typedef typename base::grid_fct_type grid_fct_type;
  typedef typename base::grid_data_type grid_data_type;

  typedef void (*veloc_func_type) ( const DOUBLE v[] );

  AMRGridAccumulation() : base(), f_veloc(0) { InitVar();}  
  AMRGridAccumulation(veloc_func_type veloc) : base(), f_veloc(veloc) { InitVar();}  

  void InitVar() {  
    _AccGrid = (grid_data_type *) 0;
    for (int d=0; d<dim; d++) {
      accoffset[d] = 0.0;
      accvelocity[d] = 0.0;
      accminshape[d] = 0;
      accmaxshape[d] = 0;
    }
    CheckpointSave = 1;
    char Name[16];
    std::sprintf(Name,"acc");
    _CheckpointName = Name;
  }

  virtual ~AMRGridAccumulation() {
    if (_AccGrid) delete _AccGrid;
  }

  virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
    base::register_at(Ctrl, prefix);
    char VariableName[32];
    for (int d=0; d<dim; d++) {
      std::sprintf(VariableName,"CellMin(%d)",d+1);
      RegisterAt(base::LocCtrl,VariableName,accminshape[d]);
      std::sprintf(VariableName,"CellMax(%d)",d+1);
      RegisterAt(base::LocCtrl,VariableName,accmaxshape[d]);
      std::sprintf(VariableName,"Offset(%d)",d+1);
      RegisterAt(base::LocCtrl,VariableName,accoffset[d]);
      std::sprintf(VariableName,"Velocity(%d)",d+1);
      RegisterAt(base::LocCtrl,VariableName,accvelocity[d]);
    }
    RegisterAt(base::LocCtrl,"CheckpointSave",CheckpointSave);
    RegisterAt(base::LocCtrl,"CheckpointName",_CheckpointName);
  }
  virtual void register_at(ControlDevice& Ctrl) { register_at(Ctrl, ""); }

  virtual void finish() { 
    if (_AccGrid) delete _AccGrid;
    _AccGrid = (grid_data_type *) 0;    
  }

  virtual void SetupData(GridHierarchy* gh, const int& ghosts) {
    base::_Hierarchy = gh;
    base::_Ghosts = ghosts;

    if (base::_AccLevel > TotalLevels(base::GH())-1) base::_AccLevel = TotalLevels(base::GH())-1;
    Coords accl = Coords(dim,accminshape)*base::GH().wholebbox().stepsize();
    Coords accu = Coords(dim,accmaxshape)*base::GH().wholebbox().stepsize();
    _AccBBox = BBox(accl,accu,base::GH().wholebbox().stepsize());
    if (base::_AccLevel<0) _AccBBox.setempty();
    if (_AccBBox.empty()) return;
    _AccBBox.refine(base::GH().refinedby(base::_AccLevel));
    _AccGrid = new grid_data_type(_AccBBox);
    if (f_veloc) {
      DataType veloc[dim]; 
      f_veloc(veloc);
      for (int d=0; d<dim; d++) 
        accvelocity[d] = veloc[d];
    }
  }

  virtual void Initialize() {
    base::GridInitializationOp(AccGrid());    
  }

  virtual void Restart() {
    int me = MY_PROC; 
    char fname[32];
    std::sprintf(fname,"%s.%d.cp",CheckpointName(),me);
    std::ifstream infile;
    infile.open(fname, std::ios::in);
    infile >> *_AccGrid;   
    infile.close();
  }

  virtual void Checkpointing() {
    int me = MY_PROC; 
    char fname[256];
    std::sprintf(fname,"%s.%d.cp",CheckpointName(),me);
    int len = std::strlen(fname);
    for (int i=0; i<len; i++) 
      if (fname[i]=='/') {
        fname[i]='\0';
        mkdir(fname, S_IRWXU | S_IRGRP | S_IXGRP | S_IROTH | S_IXOTH);
        fname[i]='/';
      }
    if (CheckpointSave) {
      char ofname[256];
      int last = -1;
      for (int i=0; i<len-1; i++) 
        if (fname[i]=='/') 
          last = i;
      if (last>=0) std::strncpy(ofname,fname,last+1);
      ofname[++last] = '.';
      std::strncpy(ofname+last+1,fname+last,len-last);
      ofname[len+1] = '\0';
      rename(fname,ofname);
    }
    std::ofstream outfile;
    outfile.open(fname, std::ios::out);
    outfile << *_AccGrid;   
    outfile.close();
  }

  virtual void Accumulate(grid_fct_type& u, const int& Time, 
                          const int& Level, const double& t) {
    if (AccBBox().empty()) return;
    Coords offset(dim,0);
    for (int d=0; d<dim; d++) 
      offset(d) = int((accoffset[d]+accvelocity[d]*t/DeltaX(base::GH(),d,base::AccLevel())))*
        AccBBox().stepsize(d);

    forall (u,Time,Level,c)
      base::AccumulateGrid(u(Time,Level,c),u.interiorbbox(Time,Level,c),
                           AccGrid(),AccBBox(),offset);
    end_forall
  }    

  void GridCombine(const int TargetNode) {
    if (AccBBox().empty()) return;
#ifdef DAGH_NO_MPI
#else
    if (comm_service::dce() && comm_service::proc_num() > 1) {  
      int num = comm_service::proc_num();
      MPI_Status status;
      int R;
      int size = sizeof(DataType)*AccBBox().size();
      int tag = 10001;
      if (MY_PROC!=TargetNode) {
        R = MPI_Send((void *)AccGrid().data(), size, MPI_BYTE, TargetNode, tag, comm_service::comm());
        if ( MPI_SUCCESS != R ) 
          comm_service::error_die( "GridCombine", "MPI_Send", R );
      }
      else {
        grid_data_type hGrid(AccBBox());
        for (int proc=0; proc<num; proc++) {
          if (proc==TargetNode) continue;
          R = MPI_Recv((void *)hGrid.data(), size, MPI_BYTE, proc, tag, comm_service::comm(), &status);
          if ( MPI_SUCCESS != R ) 
            comm_service::error_die( "GridCombine", "MPI_Recv", R );
          base::GridAccumulationOp(AccGrid(),hGrid);
        }
      }
    }
#endif
  }

  inline grid_data_type& AccGrid() { return *_AccGrid;}
  inline const BBox& AccBBox() const { return _AccBBox;}
  inline BBox AccBBox() { return _AccBBox;}
  inline const char* CheckpointName() { return _CheckpointName.c_str(); }

  inline void SetFunc(veloc_func_type veloc) { f_veloc = veloc; }
  veloc_func_type GetFunc() const { return f_veloc; }


protected:
  veloc_func_type f_veloc;
  grid_data_type *_AccGrid;
  BBox _AccBBox;
  std::string _CheckpointName;
  DataType accoffset[dim];
  DataType accvelocity[dim];
  int accminshape[dim];  
  int accmaxshape[dim];  
  int CheckpointSave;  
}; 
 

#endif

---