Blockstructured Adaptive Mesh Refinement in object-oriented C++

---

       

---

AMRFlaggingBase.h

Go to the documentation of this file.

#ifndef AMROC_AMRFLAGGING_BASE_H
#define AMROC_AMRFLAGGING_BASE_H

#include "Flagging.h"
#include "AMRFlaggingBaseInterface.h"

#ifndef AMRFlaggingBaseName
#define AMRFlaggingBase(dim)      name2(AMRFlaggingBase,dim)
#define AMRFlaggingBaseName
#endif

template <class VectorType, class FixupType, class FlagType>
class AMRFlaggingBase(DIM) : public Flagging(DIM)<VectorType,FixupType,FlagType>,
  public AMRFlaggingBaseInterface(DIM)<VectorType,FlagType> {
  typedef typename VectorType::InternalDataType DataType;
  typedef Solver(DIM)<VectorType,FixupType,FlagType> solver_type;    
  typedef GridData(DIM)<FlagType> flag_data_type;
  typedef GridData(DIM)<VectorType>   vec_grid_data_type;
 
public:
  AMRFlaggingBase(DIM)(solver_type& solver) :
    Flagging(DIM)<VectorType,FixupType,FlagType>(solver), _PlotGrad(0), _PlotErrEst(0), 
    _BaseFlag(0) {
    Tolerance   = 0.0;
    ToleranceSp = 0.0;
    _NFlags = 2*NEquations();
  }

  virtual void register_at(ControlDevice& Ctrl) {}
  virtual void register_at(ControlDevice& Ctrl,const string& prefix) {
    LocCtrl = Ctrl.getSubDevice(prefix+"AMRFlagging");
    char ToleranceName[16];
    for (int d=0; d<NEquations(); d++) {
      sprintf(ToleranceName,"TolGrad(%d)",d+1);
      RegisterAt(LocCtrl,ToleranceName,ToleranceSp(d));
      sprintf(ToleranceName,"TolErrEst(%d)",d+1);
      RegisterAt(LocCtrl,ToleranceName,Tolerance(d));
      RegisterAt(LocCtrl,"PlotGrad",_PlotGrad);
      RegisterAt(LocCtrl,"PlotErrEst",_PlotErrEst);
    }    
  }
  virtual void update() {
    SetErrorEstimation(Tolerance>0.0 || ErrorEstimation());
  }

  virtual void SetFlags(const int Time, const int Level, double t, double dt) {
    FlagCellsByDifference(Time, Level);
    if (ErrorEstimation() && Time>0) {
      //---------------------------------------------------------------
      // Take a step on the U(). U()(Time,Level) has to remain 
      // unchanged in a fractional step method!
      //---------------------------------------------------------------
      END_WATCH_FLAGGING
        Solver_().IntegrateLevel(U(), Time, Level, t, dt, false, 0.0, 1); 
      START_WATCH
      FlagCellsByError(Time, Level);
    }
  }   

  virtual int NFlags() const { return _NFlags; }
protected:  
  virtual void FlagCellsByDifference(const int Time, const int Level) { 

    if (Time < 0) return;
    int me = MY_PROC(GH); 
    int TStep = TimeStep(U(),Level);
    for (int k=0; k<NEquations(); k++) {
      if (ToleranceSp(k)<=0.0) continue;
      forall(U(),Time,Level,c)
        equals_from(Work()(Time,Level,c), U()(Time,Level,c), k);
        Work()(Time+TStep,Level,c) = 0.0;
      end_forall      
      FlagByDifference(Time, Level, Work(), ToleranceSp(k), k+1+BaseFlag());

      if (PlotGrad() && Solver_().LastOutputTime()==Time) {
        char name[DAGHBktGFNameWidth+16];
        sprintf(name,"grad_%d_%d",k+1,Time); 
        if (Level==0 && me==VizServer) 
          cout << " *** Writing " << name << endl;   
        Write(Work(),Time+TStep,Level,DAGH_Double,name); 
#ifndef DAGH_NO_MPI
        MPI_Barrier(comm_service::comm()); 
#endif
      }
    }
  }

  virtual void FlagCellsByError(const int Time, const int Level) { 
    
    if (Time <= 0) return;
    int me = MY_PROC(GH); 
    int TStep = TimeStep(U(),Level);
    int ShTStep = TimeStep(Ush(),Level);
    int myargc = sizeof(int); int myargs = NEquations();
    forall (U(),Time,Level,c)      
      //------------------------------------------------------------------------
      // First update Ush() from U()
      //------------------------------------------------------------------------
      f_restrict_amr(FA(DIM,U()(Time+TStep,Level,c)),
                     FA(DIM,Ush()(Time+ShTStep,Level,c)),
                     FBA(U().interiorbbox(Time+TStep,Level,c)), 
                     (char *)&myargs, &myargc); 
    end_forall

    //------------------------------------------------------------------------
    // Estimate error by Richardson-Extrapolation
    //------------------------------------------------------------------------
    DataType Order = pow(2.0,Integrator_().NMethodOrder()+1.0) - 2.0;
    for (int k=0; k<NEquations(); k++) {
      if (Tolerance(k)<=0.0) continue;
      forall(Ush(),Time,Level,c)
        equals_from(Worksh()(Time+ShTStep,Level,c), Ush()(Time+ShTStep,Level,c), k);
        minus_from(Worksh()(Time+ShTStep,Level,c), Ush()(Time,Level,c), k);
        Worksh()(Time+ShTStep,Level,c).divide(Order);
      end_forall      
      FlagByErrorEstimation(Time, Level, Worksh(), Tolerance(k), k+1+NEquations()+BaseFlag());

      if (PlotErrEst() && Solver_().LastOutputTime()==Time) {
        char name[DAGHBktGFNameWidth+16];
        sprintf(name,"errest_%d_%d",k+1,Time); 
        if (Level==0 && me==VizServer) 
          cout << " *** Writing " << name << endl;   
        Write(Worksh(),Time+ShTStep,Level,DAGH_Double,name); 
#ifndef DAGH_NO_MPI
        MPI_Barrier(comm_service::comm()); 
#endif
      }
    }
  } 

  inline const int& PlotGrad() const { return _PlotGrad; }
  inline const int& PlotErrEst() const { return _PlotErrEst; }
  inline const int& BaseFlag() const { return _BaseFlag; }
  void SetBaseFlag(const int bf) { _BaseFlag = bf; }

protected:
  VectorType Tolerance, ToleranceSp;
  int _PlotGrad, _PlotErrEst;
  int _BaseFlag, _NFlags;
};

#endif

---

Quickstart     Users Guide     Programmers Reference     Installation      Examples     Download

---

AMROC Main      Home      Contact

last update: 06/01/04