clawpack/applications/euler_chem/1d/ReflecDet/src/Problem.h

// -*- C++ -*-

// Copyright (C) 2002 Ralf Deiterding
// Brandenburgische Universitaet Cottbus

#ifndef AMROC_PROBLEM_H
#define AMROC_PROBLEM_H

#include "eulerrhok1.h"

#define NEQUATIONS 15
#define NEQUSED    12
#define NFIXUP     11
#define NAUX        0

#include "ClpProblem.h"

#define OWN_FLAGGING
#define OWN_AMRSOLVER
#include "ClpStdProblem.h"
#include "F77Interfaces/F77UpdateLevelTransfer.h"
 
#define MAXREGIONS   5

template <class VectorType, class FlagType, int dim>
class PermanentFlagging : public StdCriterion<VectorType,FlagType,dim> {
  typedef typename VectorType::InternalDataType DataType;
  typedef StdCriterion<VectorType,FlagType,dim> base;
public:
  typedef typename base::vec_grid_fct_type vec_grid_fct_type;  
  typedef typename base::grid_fct_type grid_fct_type;
  typedef typename base::flag_fct_type flag_fct_type;

  PermanentFlagging() {
    base::_IsUsed = true;
  }

  virtual void register_at(ControlDevice& Ctrl) { register_at(Ctrl, ""); }
  virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {     
    base::LocCtrl = Ctrl.getSubDevice(prefix+"PermanentFlagging");
    char Name[16];
    for (int i=0; i<MAXREGIONS; i++) {  
      x1[i] = -1.0; x2[i] = -1.0; lev[i] = -1;
      std::sprintf(Name,"x1(%d)",i+1);
      RegisterAt(base::LocCtrl,Name,x1[i]);
      std::sprintf(Name,"x2(%d)",i+1);
      RegisterAt(base::LocCtrl,Name,x2[i]);
      std::sprintf(Name,"l(%d)",i+1);
      RegisterAt(base::LocCtrl,Name,lev[i]);
    }    
  }

  virtual bool SetFlags(vec_grid_fct_type& u, grid_fct_type& work, flag_fct_type& flags, 
                        const int& cnt, const int& Time, const int& Level, const double& t, 
                        const FlagType& FlagValue) {
    forall(flags, Time, Level, c)
      BeginFastIndex1(err_flag, flags(Time,Level,c).bbox(), 
                      flags(Time,Level,c).data(), FlagType);
      BBox OpBox = flags.interiorbbox(Time,Level,c);
      Coords& OpBox_stepsize = OpBox.stepsize();
      Coords cell=OpBox.lower();
      for_1 (n, OpBox, OpBox_stepsize)
        cell(0) = n; 
        DCoords xc = base::GH().worldCoords(cell, OpBox_stepsize);
        for (int i=0; i<MAXREGIONS; i++)  
          if (xc(0)>=x1[i] && xc(0)<=x2[i] && Level<=lev[i]) 
            FastIndex1(err_flag,n) = BadPoint; 
      end_for
      EndFastIndex1(err_flag);
    end_forall
    return true;
  } 

  virtual void OutputName(char* name, int cnt) { std::sprintf(name,"permanent_%d",cnt+1); }

private:
  double x1[MAXREGIONS];
  double x2[MAXREGIONS];
  int lev[MAXREGIONS];
};

class FlaggingSpecific : 
  public AMRFlagging<VectorType,FixupType,FlagType,DIM> {
  typedef AMRFlagging<VectorType,FixupType,FlagType,DIM> base;
public:
  FlaggingSpecific(base::solver_type& solver) : base(solver) {
      base::AddCriterion(new F77ScaledGradient<VectorType,FlagType,DIM>(f_flg));
      base::AddCriterion(new F77LimiterType<VectorType,FlagType,DIM>(f_flg));
      base::AddCriterion(new F77AbsoluteError<VectorType,FixupType,FlagType,DIM>(solver,f_flg));
      base::AddCriterion(new F77RelativeError<VectorType,FixupType,FlagType,DIM>(solver,f_flg));
      base::AddCriterion(new PermanentFlagging<VectorType,FlagType,DIM>());
    }

  ~FlaggingSpecific() { DeleteAllCriterions(); }
};

class SolverSpecific : 
  public AMRSolver<VectorType,FixupType,FlagType,DIM> {
  typedef AMRSolver<VectorType,FixupType,FlagType,DIM> base;
public:
  SolverSpecific(IntegratorSpecific& integ, 
                 base::initial_condition_type& init,
                 base::boundary_conditions_type& bc) :
    AMRSolver<VectorType,FixupType,FlagType,DIM>(integ, init, bc) {
      SetLevelTransfer(new F77UpdateLevelTransfer<VectorType,DIM>(f_prolong, f_restrict, f_tupdate));
    SetFileOutput(new F77FileOutput<VectorType,DIM>(f_out)); 
    SetFixup(new FixupSpecific(integ));
    SetFlagging(new FlaggingSpecific(*this)); 
  }  

  ~SolverSpecific() {
    delete _LevelTransfer;
    delete _Flagging;
    delete _Fixup;
    delete _FileOutput;