00001
00002
00003
00004
00005
00006 #ifndef AMROC_PROBLEM_H
00007 #define AMROC_PROBLEM_H
00008
00009 #include "eulerrhok2.h"
00010
00011 #define NEQUATIONS 14
00012 #define NEQUSED 13
00013 #define NFIXUP 12
00014 #define NAUX 4
00015
00016 #include "ClpProblem.h"
00017
00018 #define OWN_FLAGGING
00019 #define OWN_GFMAMRSOLVER
00020 #include "ClpStdGFMProblem.h"
00021 #include "F77Interfaces/F77UpdateLevelTransfer.h"
00022 #include "SpecialOutput/AMRGridAccumulation.h"
00023
00024 class FlaggingSpecific :
00025 public AMRFlagging<VectorType,FixupType,FlagType,DIM> {
00026 typedef AMRFlagging<VectorType,FixupType,FlagType,DIM> base;
00027 public:
00028 FlaggingSpecific(solver_type& solver) : base(solver) {
00029 base::AddCriterion(new ScaledGradient<VectorType,FlagType,DIM>());
00030 base::AddCriterion(new LimiterType<VectorType,FlagType,DIM>());
00031 base::AddCriterion(new AbsoluteError<VectorType,FixupType,FlagType,DIM>(solver));
00032 base::AddCriterion(new RelativeError<VectorType,FixupType,FlagType,DIM>(solver));
00033 base::AddCriterion(new F77ScaledGradient<VectorType,FlagType,DIM>(f_flg));
00034 base::AddCriterion(new F77LimiterType<VectorType,FlagType,DIM>(f_flg));
00035 base::AddCriterion(new F77AbsoluteError<VectorType,FixupType,FlagType,DIM>(solver,f_flg));
00036 base::AddCriterion(new F77RelativeError<VectorType,FixupType,FlagType,DIM>(solver,f_flg));
00037 }
00038 ~FlaggingSpecific() { DeleteAllCriterions(); }
00039 };
00040
00041 class SolverSpecific :
00042 public AMRGFMSolver<VectorType,FixupType,FlagType,DIM> {
00043 typedef VectorType::InternalDataType DataType;
00044 typedef AMRGFMSolver<VectorType,FixupType,FlagType,DIM> base;
00045 typedef AMRGridAccumulation<DataType,DIM> accumulation_type;
00046 typedef F77GFMFileOutput<VectorType,FixupType,FlagType,DIM> output_type;
00047 public:
00048 SolverSpecific(IntegratorSpecific& integ,
00049 base::initial_condition_type& init,
00050 base::boundary_conditions_type& bc) : base(integ, init, bc) {
00051 SetLevelTransfer(new F77UpdateLevelTransfer<VectorType,DIM>(f_prolong, f_restrict, f_tupdate));
00052 SetFileOutput(new F77GFMFileOutput<VectorType,FixupType,FlagType,DIM>(*this,f_out));
00053 SetFixup(new FixupSpecific(integ));
00054 SetFlagging(new FlaggingSpecific(*this));
00055 AddGFM(new GhostFluidMethod<VectorType,DIM>(
00056 new F77GFMBoundary<VectorType,DIM>(f_ibndrfl,f_itrans),
00057 new F77GFMLevelSet<DataType,DIM>(f_lset)));
00058 _Accumulation = new accumulation_type();
00059 MaxType = 1;
00060 }
00061
00062 ~SolverSpecific() {
00063 delete _LevelTransfer;
00064 delete _Flagging;
00065 delete _Fixup;
00066 delete _FileOutput;
00067 delete _Accumulation;
00068 }
00069
00070 virtual void register_at(ControlDevice& Ctrl) { base::register_at(Ctrl); }
00071 virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
00072 base::register_at(Ctrl, prefix);
00073 RegisterAt(base::LocCtrl,"MaxType",MaxType);
00074 _Accumulation->register_at(base::LocCtrl, prefix);
00075 }
00076
00077 virtual void SetupData() {
00078 base::SetupData();
00079 _Accumulation->SetupData(base::PGH(), base::NGhosts());
00080 }
00081
00082 virtual void Initialize_(const double& dt_start) {
00083 base::Initialize_(dt_start);
00084 _Accumulation->Initialize();
00085 }
00086
00087 virtual void Output() {
00088 int Time = CurrentTime(base::GH(),0);
00089 if (Time == base::LastOutputTime()) return;
00090 _Accumulation->GridCombine(VizServer);
00091 int me = MY_PROC;
00092 if (me == VizServer) {
00093 int Level = _Accumulation->AccLevel();
00094 base::FileOutput_().WriteOut(_Accumulation->AccGrid(),_Accumulation->AccName(),
00095 Time,Level,base::t[0]);
00096 }
00097 base::Output();
00098 }
00099
00100 virtual void Checkpointing_(const char* CheckpointFile) {
00101 base::Checkpointing_(CheckpointFile);
00102 _Accumulation->Checkpointing();
00103 }
00104
00105 virtual bool Restart_(const char* CheckpointFile) {
00106 if (base::Restart_(CheckpointFile)) {
00107 _Accumulation->Restart();
00108 return true;
00109 }
00110 else
00111 return false;
00112 }
00113
00114 virtual void AdvanceLevel(const int Level, int RegridEvery, bool RegridDone,
00115 bool ShadowAllowed, bool DoFixup,
00116 bool RecomposeBaseLev, bool RecomposeHighLev) {
00117
00118 base::AdvanceLevel(Level,RegridEvery,RegridDone,ShadowAllowed,
00119 DoFixup,RecomposeBaseLev,RecomposeHighLev);
00120
00121 int Time = CurrentTime(base::GH(),Level);
00122 if (MaxType) {
00123 AdaptiveBoundaryUpdate(base::U(),Time,Level);
00124 Sync(base::U(),Time,Level);
00125 ExternalBoundaryUpdate(base::U(),Time,Level);
00126
00127 int irho = 1;
00128 int iu = NEQUSED-DIM-2;
00129 int iv = NEQUSED-DIM-1;
00130 int TStep = TimeStep(base::U(),Level);
00131 ((output_type*) base::_FileOutput)->Transform(base::U(), base::Work(), Time+TStep,
00132 Level, irho, base::t[Level]);
00133 forall (base::U(),Time,Level,c)
00134 Coords ss(base::U()(Time,Level,c).bbox().stepsize());
00135 BBox bbi(base::U().interiorbbox(Time,Level,c));
00136 BeginFastIndex2(U, base::U()(Time,Level,c).bbox(),
00137 base::U()(Time,Level,c).data(), VectorType);
00138 BeginFastIndex2(Work, Work()(Time,Level,c).bbox(),
00139 Work()(Time,Level,c).data(), DataType);
00140 BeginFastIndex2(WorkN, Work()(Time+TStep,Level,c).bbox(),
00141 Work()(Time+TStep,Level,c).data(), DataType);
00142 DataType uy, vx;
00143 DCoords dx = base::GH().worldStep(ss);
00144 for_2 (i, j, bbi, ss)
00145 vx = (FastIndex2(U,i+ss(0),j)(iv)/FastIndex2(WorkN,i+ss(0),j)-
00146 FastIndex2(U,i-ss(0),j)(iv)/FastIndex2(WorkN,i-ss(0),j))/(2.0*dx(0));
00147 uy = (FastIndex2(U,i,j+ss(1))(iu)/FastIndex2(WorkN,i,j+ss(1))-
00148 FastIndex2(U,i,j-ss(1))(iu)/FastIndex2(WorkN,i,j-ss(1)))/(2.0*dx(1));
00149 FastIndex2(Work,i,j) = std::fabs(vx-uy);
00150 end_for
00151 EndFastIndex2(U);
00152 EndFastIndex2(Work);
00153 EndFastIndex2(WorkN);
00154 end_forall
00155 }
00156 else {
00157 int ipress = base::Dim()+4;
00158 ((output_type*) base::_FileOutput)->Transform(base::U(), base::Work(), Time, Level,
00159 ipress, base::t[Level]);
00160 }
00161
00162 _Accumulation->Accumulate(base::Work(), Time, Level, base::t[Level]);
00163 }
00164
00165 protected:
00166 int MaxType;
00167 accumulation_type* _Accumulation;
00168 };
00169
00170 #endif
00171
00172