Blockstructured Adaptive Mesh Refinement in object-oriented C++
Blockstructured Adaptive Mesh Refinement in object-oriented C++
00001 #ifndef AMROC_CLP_INTEGRATOR3_H 00002 #define AMROC_CLP_INTEGRATOR3_H 00003 00011 #include "ClpIntegratorBase.h" 00012 00013 #ifndef ClpIntegratorName 00014 #define ClpIntegrator(dim) name2(ClpIntegrator,dim) 00015 #define ClpIntegratorName 00016 #endif 00017 00028 template <class VectorType, class AuxVectorType> 00029 class ClpIntegrator(3) : public ClpIntegratorBase(3)<VectorType,AuxVectorType> { 00030 typedef typename VectorType::InternalDataType DataType; 00031 typedef ClpIntegratorBase(3)<VectorType,AuxVectorType> base_type; 00032 typedef GridData(3)<VectorType> vec_grid_data_type; 00033 typedef GridFunction(3)<VectorType> vec_grid_fct_type; 00034 typedef GridData(DIM_1)<VectorType> ld_vec_grid_data_type; 00035 typedef GridData(DIM_1)<AuxVectorType> ld_aux_grid_data_type; 00036 public: 00037 ClpIntegrator(3)(const int equ, const int wv, const int gh) : 00038 base_type(equ, wv, gh) { 00039 transverse_flux_func = f_transverseflux; 00040 } 00041 00042 //------------------------------------------------------------------ 00043 // Solve Riemann problems. Copy left and right data into appropriate 00044 // format for normal_flux_func(). Note that left state has to shifted 00045 // one element. This is done in copySlice(). 00046 //------------------------------------------------------------------ 00047 virtual void ComputeRP(ld_vec_grid_data_type& LeftState, ld_aux_grid_data_type& auxl, 00048 ld_vec_grid_data_type& RightState, ld_aux_grid_data_type& auxr, 00049 ld_vec_grid_data_type& NewState, const int& direction) { 00050 00051 AllocError::SetTexts("ClpIntegrator(3)","ComputeRP(): work"); 00052 int maxm = LeftState.extents(0) + 1; 00053 int mx = maxm-2*NGhosts(); 00054 DataType* wave = new DataType[maxm*NEqUsed()*NWaves()]; 00055 DataType* s = new DataType[maxm*NWaves()]; 00056 00057 AllocError::SetTexts("ClpIntegrator(3)","calculate_Riemann_Problem(): slices"); 00058 DataType* LeftStateSl = new DataType[maxm*NEqUsed()]; 00059 DataType* RightStateSl = new DataType[maxm*NEqUsed()]; 00060 DataType* apdqSl = new DataType[maxm*NEqUsed()]; 00061 DataType* amdqSl = new DataType[maxm*NEqUsed()]; 00062 DataType* auxlSl = new DataType[maxm*3*NAux()]; 00063 DataType* auxrSl = new DataType[maxm*3*NAux()]; 00064 00065 BBox bb(LeftState.bbox()); 00066 BBox bbSl(bb); 00067 AlignBBoxToSlice(bbSl); 00068 for (register int y=bb.lower(1); y<=bb.upper(1); y+=bb.stepsize(1)) { 00069 00070 copySlice(LeftStateSl, bbSl, LeftState, y, NEqUsed(), 1); 00071 copySlice(RightStateSl, bbSl, RightState, y, NEqUsed(), 0); 00072 00073 // auxiliary data along neighboring slices generally aren't needed in 00074 // the rpn3 routine. 00075 copyAuxSlice(&(auxlSl[maxm*NAux()]), bbSl, auxl, y, NAux(), 1); 00076 copyAuxSlice(&(auxrSl[maxm*NAux()]), bbSl, auxr, y, NAux(), 0); 00077 00078 (*normal_flux_func)(direction, mx, NEqUsed(), NWaves(), NGhosts(), mx, 00079 LeftStateSl, RightStateSl, NAux(), 00080 (const DataType *)auxlSl, 00081 (const DataType *)auxrSl, 00082 wave, s, amdqSl, apdqSl); 00083 00084 addSlices(NewState, apdqSl, amdqSl, bbSl, y, NEqUsed(), 1); 00085 } 00086 00087 delete [] LeftStateSl; delete [] RightStateSl; 00088 delete [] auxlSl; delete [] auxrSl; 00089 delete [] apdqSl; delete [] amdqSl; 00090 delete [] wave; 00091 delete [] s; 00092 } 00093 00094 virtual double ComputeGrid(vec_grid_data_type& StateVec, 00095 vec_grid_data_type* Flux[], DataType* aux, 00096 double dt, DCoords& dx) { 00097 00098 AllocError::SetTexts("ClpIntegrator(2)","calculate_time_step(): work"); 00099 Coords ex = StateVec.extents(); 00100 int mx = ex(0)-2*NGhosts(); 00101 int my = ex(1)-2*NGhosts(); 00102 int mz = ex(2)-2*NGhosts(); 00103 int maxm = Max(Max(mx,my),mz); 00104 int msl = 0; 00105 #ifdef EXTENDED_INTEGRATOR 00106 msl = 4; 00107 #endif 00108 int mwork = (maxm + 2*NGhosts())*((54+msl+NWaves())*NEqUsed() + 00109 NWaves() + 9*NAux() + 3); 00110 DataType* work = new DataType[mwork]; 00111 00112 int faddm = 0; 00113 int faddp = faddm + (maxm+2*NGhosts())*NEqUsed(); 00114 int gaddm = faddp + (maxm+2*NGhosts())*NEqUsed(); 00115 int gaddp = gaddm + 6*(maxm+2*NGhosts())*NEqUsed(); 00116 int haddm = gaddp + 6*(maxm+2*NGhosts())*NEqUsed(); 00117 int haddp = haddm + 6*(maxm+2*NGhosts())*NEqUsed(); 00118 int q1d = haddp + 6*(maxm+2*NGhosts())*NEqUsed(); 00119 int dtdx1d = q1d + (maxm+2*NGhosts())*NEqUsed(); 00120 int dtdy1d = dtdx1d + (maxm+2*NGhosts()); 00121 int dtdz1d = dtdy1d + (maxm+2*NGhosts()); 00122 int aux1 = dtdz1d + (maxm+2*NGhosts()); 00123 int aux2 = aux1 + 3*(maxm+2*NGhosts())*NAux(); 00124 int aux3 = aux2 + 3*(maxm+2*NGhosts())*NAux(); 00125 int next = aux3 + 3*(maxm+2*NGhosts())*NAux(); 00126 int mwork1 = mwork - next; 00127 double cfl = 0.0; 00128 00129 // The ghostcells must remain untouched during integration! 00130 f_step(maxm,mx,my,mz,NEquations(),NEqUsed(), 00131 NAux(),NWaves(),NGhosts(),mx,my,mz, 00132 StateVec.data(), 00133 aux,dx(0),dx(1),dx(2),dt,method,mthlim,cfl, 00134 Flux[0]->data(),Flux[1]->data(), 00135 Flux[2]->data(),Flux[3]->data(), 00136 Flux[4]->data(),Flux[5]->data(), 00137 &(work[faddm]),&(work[faddp]), 00138 &(work[gaddm]),&(work[gaddp]), 00139 &(work[haddm]),&(work[haddp]), 00140 &(work[q1d]), 00141 &(work[dtdx1d]),&(work[dtdy1d]),&(work[dtdz1d]), 00142 &(work[aux1]),&(work[aux2]),&(work[aux3]), 00143 &(work[next]),mwork1, 00144 normal_flux_func,transverse_flux_func); 00145 00146 delete [] work; 00147 return cfl; 00148 } 00149 00150 private: 00151 void AlignBBoxToSlice(BBox &bb) { 00152 gdbAlignBBox(1, bb, DAGH_X); 00153 bb.rank = 1; 00154 bb.lower().rank = 1; 00155 bb.upper().rank = 1; 00156 bb.stepsize().rank = 1; 00157 } 00158 00159 void copySlice(DataType* target, const BBox &targetbbox, 00160 const ld_vec_grid_data_type &source, 00161 const int y, const int meqn, const int move) { 00162 BBox bbox = targetbbox * source.bbox(); 00163 BeginFastIndex2(src, source.bbox(), source.data(), const VectorType); 00164 for (int m=0; m<meqn; m++) { 00165 if (move) target++; 00166 BeginFastIndex1(tgt, targetbbox, target, DataType); 00167 for_1 (k, bbox, bbox.stepsize()) 00168 FastIndex1(tgt, k) = FastIndex2(src, k, y)(m); 00169 end_for 00170 EndFastIndex1(tgt); 00171 if (!move) target++; 00172 target = &(target[targetbbox.size()]); 00173 } 00174 EndFastIndex2(src); 00175 } 00176 00177 void copyAuxSlice(DataType* target, const BBox &targetbbox, 00178 const ld_aux_grid_data_type &source, 00179 const int y, const int meqn, const int move) { 00180 BBox bbox = targetbbox * source.bbox(); 00181 BeginFastIndex2(src, source.bbox(), source.data(), const AuxVectorType); 00182 for (int m=0; m<meqn; m++) { 00183 if (move) target++; 00184 BeginFastIndex1(tgt, targetbbox, target, DataType); 00185 for_1 (k, bbox, bbox.stepsize()) 00186 FastIndex1(tgt, k) = FastIndex2(src, k, y)(m); 00187 end_for 00188 EndFastIndex1(tgt); 00189 if (!move) target++; 00190 target = &(target[targetbbox.size()]); 00191 } 00192 EndFastIndex2(src); 00193 } 00194 00195 void addSlices(ld_vec_grid_data_type &target, 00196 DataType* source1, DataType* source2, const BBox &sourcebbox, 00197 const int y, const int meqn, const int move) { 00198 BBox bbox = target.bbox() * sourcebbox; 00199 BeginFastIndex2(tgt, target.bbox(), target.data(), VectorType); 00200 for (int m=0; m<meqn; m++) { 00201 if (move) { source1++; source2++; } 00202 BeginFastIndex1(src1, sourcebbox, source1, const DataType); 00203 BeginFastIndex1(src2, sourcebbox, source2, const DataType); 00204 for_1 (k, bbox, bbox.stepsize()) 00205 FastIndex2(tgt, k, y)(m) = FastIndex1(src1, k) + FastIndex1(src2, k); 00206 end_for 00207 EndFastIndex1(src1); 00208 EndFastIndex1(src2); 00209 if (!move) { source1++; source2++; } 00210 source1 = &(source1[sourcebbox.size()]); 00211 source2 = &(source2[sourcebbox.size()]); 00212 } 00213 EndFastIndex2(tgt); 00214 } 00215 00216 transverse_func_type transverse_flux_func; 00217 }; 00218 00219 00220 #endif
Quickstart Users Guide Programmers Reference Installation Examples Download
AMROC Main Home Contactlast update: 06/01/04