00001
00002
00003 #ifndef AMROC_PROBLEM_H
00004 #define AMROC_PROBLEM_H
00005
00006 #define DIM 2
00007
00008 #define OWN_AMRSOLVER
00009 #define OWN_FLAGGING
00010
00011 #include "WENOProblem.h"
00012 #include "WENOStdProblem.h"
00013 #include "AMRInterpolation.h"
00014 #include "WENOStatistics.h"
00015
00016 #ifdef OWN_FLAGGING
00017 #define NZONES 3
00018
00019 class FlaggingSpecific :
00020 public AMRFlagging<VectorType,FixupType,FlagType,DIM> {
00021 typedef AMRFlagging<VectorType,FixupType,FlagType,DIM> base;
00022 public:
00023 FlaggingSpecific(base::solver_type& solver) : base(solver) {
00024 base::AddCriterion(new ByValue<VectorType,FlagType,DIM>());
00025 base::AddCriterion(new ScaledGradient<VectorType,FlagType,DIM>());
00026 base::AddCriterion(new LimiterType<VectorType,FlagType,DIM>());
00027 base::AddCriterion(new AbsoluteError<VectorType,FixupType,FlagType,DIM>(solver));
00028 base::AddCriterion(new RelativeError<VectorType,FixupType,FlagType,DIM>(solver));
00029 #ifdef f_flgout
00030 base::AddCriterion(new F77ByValue<VectorType,FlagType,DIM>(f_flgout));
00031 base::AddCriterion(new F77ScaledGradient<VectorType,FlagType,DIM>(f_flgout));
00032 base::AddCriterion(new F77LimiterType<VectorType,FlagType,DIM>(f_flgout));
00033 base::AddCriterion(new F77AbsoluteError<VectorType,FixupType,FlagType,DIM>(solver,f_flgout));
00034 base::AddCriterion(new F77RelativeError<VectorType,FixupType,FlagType,DIM>(solver,f_flgout));
00035 #endif
00036 }
00037
00038 virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
00039 base::register_at(Ctrl, prefix);
00040 char VariableName[32];
00041 for (int iz=0; iz < NZONES ; iz++)
00042 for (int d=0; d<DIM; d++) {
00043 sprintf(VariableName,"DCellsMin%d(%d)",iz+1,d+1);
00044 RegisterAt(base::LocCtrl,VariableName,dcmin[iz][d]);
00045 sprintf(VariableName,"DCellsMax%d(%d)",iz+1,d+1);
00046 RegisterAt(base::LocCtrl,VariableName,dcmax[iz][d]);
00047 }
00048 }
00049 virtual void register_at(ControlDevice& Ctrl) {
00050 register_at(Ctrl, "");
00051 }
00052
00053 virtual void SetFlags(const int Time, const int Level, double t, double dt) {
00054 base::SetFlags(Time, Level, t, dt);
00055 int max_level = MaxLevel(base::GH());
00056 for ( int iz = 0 ; iz < NZONES ; iz ++ ) {
00057 Coords lb(base::Dim(), dcmin[iz]), ub(base::Dim(), dcmax[iz]),
00058 s(base::Dim(),1);
00059 BBox bbox(lb*StepSize(base::GH(),0), ub*StepSize(base::GH(),0),
00060 s*StepSize(base::GH(),0));
00061 bbox.refine(StepSize(base::GH(),0)/StepSize(base::GH(),Level));
00062
00063 forall (base::Flags(),Time,Level,c)
00064 if ( Level > max_level - (2+iz) )
00065 base::Flags()(Time,Level,c).equals(GoodPoint, bbox);
00066 end_forall
00067 }
00068 }
00069
00070 ~FlaggingSpecific() { DeleteAllCriterions(); }
00071 private:
00072 int dcmin[NZONES][DIM], dcmax[NZONES][DIM];
00073 };
00074 #endif
00075
00076 #ifdef OWN_AMRSOLVER
00077
00078 #define f_initialize_geometry FORTRAN_NAME_(init_geom,INIT_GEOM)
00079 extern "C" {
00080 void f_initialize_geometry (const INTEGER&,const DOUBLE*);
00081 }
00082
00083 class SolverSpecific :
00084 public AMRSolver<VectorType,FixupType,FlagType,DIM> {
00085 typedef AMRSolver<VectorType,FixupType,FlagType,DIM> base;
00086 typedef AMRInterpolation<VectorType,DIM> interpolation_type;
00087 typedef F77FileOutput<VectorType,DIM> output_type;
00088 typedef WENOStatistics<VectorType,interpolation_type,output_type,DIM> stat_type;
00089 public:
00090 SolverSpecific(IntegratorSpecific& integ,
00091 base::initial_condition_type& init,
00092 base::boundary_conditions_type& bc) :
00093 AMRSolver<VectorType,FixupType,FlagType,DIM>(integ, init, bc) {
00094 SetLevelTransfer(new F77LevelTransfer<VectorType,DIM>(f_prolong, f_restrict));
00095 #ifdef f_flgout
00096 SetFileOutput(new F77FileOutput<VectorType,DIM>(f_flgout));
00097 #else
00098 SetFileOutput(new FileOutput<VectorType,DIM>());
00099 #endif
00100 SetFixup(new FixupSpecific());
00101 SetFlagging(new FlaggingSpecific(*this));
00102 _Interpolation = new interpolation_type();
00103 _Stats = new stat_type(_Interpolation, (output_type*)_FileOutput);
00104 }
00105
00106 ~SolverSpecific() {
00107 delete _LevelTransfer;
00108 delete _Flagging;
00109 delete _Fixup;
00110 delete _FileOutput;
00111 delete _Interpolation;
00112 delete _Stats;
00113 }
00114
00115 virtual void register_at(ControlDevice& Ctrl) { base::register_at(Ctrl); }
00116 virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
00117 base::register_at(Ctrl, prefix);
00118 if (_Stats) _Stats->register_at(base::LocCtrl, "");
00119 }
00120
00121 virtual void SetupData() {
00122 base::SetupData();
00123 _Interpolation->SetupData(base::PGH(), base::NGhosts());
00124 _Stats->Setup(base::PGH(), base::NGhosts(), base::shape, base::geom);
00125 f_initialize_geometry (DIM,geom);
00126 }
00127
00128 virtual void Initialize_(const double& dt_start) {
00129 base::Initialize_(dt_start);
00130 _Stats->Evaluate(base::t, base::U(), base::Work());
00131 }
00132
00133 virtual void Advance(double& t, double& dt) {
00134 base::Advance(t, dt);
00135 _Stats->Evaluate(base::t, base::U(), base::Work());
00136 }
00137
00138 protected:
00139 interpolation_type* _Interpolation;
00140 stat_type* _Stats;
00141 };