AMROC - Source code documentation

00001 #ifndef AMROC_CLP_INTEGRATOR2_H
00002 #define AMROC_CLP_INTEGRATOR2_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(2) : public ClpIntegratorBase(2)<VectorType,AuxVectorType> {
00030   typedef typename VectorType::InternalDataType DataType;
00031   typedef ClpIntegratorBase(2)<VectorType,AuxVectorType> base_type;
00032   typedef GridData(2)<VectorType> vec_grid_data_type;
00033   typedef GridFunction(2)<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(2)(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(2)","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     DataType* LeftStateSl  = new DataType[maxm*NEqUsed()];
00057     DataType* RightStateSl = new DataType[maxm*NEqUsed()];
00058     DataType* auxlSl = new DataType[maxm*NAux()];
00059     DataType* auxrSl = new DataType[maxm*NAux()];
00060 
00061     AllocError::SetTexts("ClpIntegrator(2)","calculate_Riemann_Problem(): fluctuations");
00062     DataType* apdqSl = new DataType[maxm*NEqUsed()];
00063     DataType* amdqSl = new DataType[maxm*NEqUsed()];
00064 
00065     copySlice(LeftStateSl, LeftState, NEqUsed(), 1);
00066     copySlice(RightStateSl, RightState, NEqUsed(), 0);
00067         
00068     copyAuxSlice(auxlSl, auxl, NAux(), 1);
00069     copyAuxSlice(auxrSl, auxr, NAux(), 0);
00070 
00071     (*normal_flux_func)(direction, mx, NEqUsed(), NWaves(), NGhosts(), mx,
00072                         LeftStateSl, RightStateSl,
00073                         NAux(), auxlSl, auxrSl, wave, s, amdqSl, apdqSl);
00074 
00075     addSlices(NewState, apdqSl, amdqSl, NEqUsed(), 1);
00076     
00077     delete [] LeftStateSl; delete [] RightStateSl;
00078     delete [] auxlSl; delete [] auxrSl;
00079     delete [] apdqSl; delete [] amdqSl;
00080     delete [] wave;
00081     delete [] s;
00082   }
00083 
00084   virtual double ComputeGrid(vec_grid_data_type& StateVec, 
00085                              vec_grid_data_type* Flux[], DataType* aux,
00086                              double dt, DCoords& dx) {
00087 
00088     AllocError::SetTexts("ClpIntegrator(2)","calculate_time_step(): work");
00089     Coords ex = StateVec.extents();
00090     int mx = ex(0)-2*NGhosts();
00091     int my = ex(1)-2*NGhosts();      
00092     int maxm = Max(mx,my);
00093     int msl = 0;
00094 #ifdef EXTENDED_INTEGRATOR
00095     msl = 4;
00096 #endif  
00097 
00098     int mwork = (maxm + 2*NGhosts())*((12+msl+NWaves())*NEqUsed() + 
00099                                       NWaves() + 3*NAux() + 2);      
00100     DataType* work = new DataType[mwork];
00101     
00102     int faddm = 0;
00103     int faddp = faddm + (maxm+2*NGhosts())*NEqUsed();
00104     int gaddm = faddp + (maxm+2*NGhosts())*NEqUsed();
00105     int gaddp = gaddm + 2*(maxm+2*NGhosts())*NEqUsed();
00106     int q1d =   gaddp + 2*(maxm+2*NGhosts())*NEqUsed();
00107     int dtdx1d = q1d + (maxm+2*NGhosts())*NEqUsed();
00108     int dtdy1d = dtdx1d + (maxm+2*NGhosts());
00109     int aux1 =   dtdy1d + (maxm+2*NGhosts());
00110     int aux2 = aux1 + (maxm+2*NGhosts())*NAux();
00111     int aux3 = aux2 + (maxm+2*NGhosts())*NAux();
00112     int next = aux3 + (maxm+2*NGhosts())*NAux();
00113     int mwork1 = mwork - next;
00114     double cfl = 0.0;
00115 
00116     // The ghostcells should remain untouched during integration!
00117     f_step(maxm,mx,my,NEquations(),NEqUsed(),
00118            NAux(),NWaves(),NGhosts(),mx,my,
00119            StateVec.data(),
00120            aux,dx(0),dx(1),dt,method,mthlim,cfl, 
00121            Flux[0]->data(),Flux[1]->data(),
00122            Flux[2]->data(),Flux[3]->data(), 
00123            &(work[faddm]),&(work[faddp]),
00124            &(work[gaddm]),&(work[gaddp]),
00125            &(work[q1d]),
00126            &(work[dtdx1d]),&(work[dtdy1d]),
00127            &(work[aux1]),&(work[aux2]),&(work[aux3]), 
00128            &(work[next]),mwork1,
00129            normal_flux_func,transverse_flux_func);
00130 
00131     delete [] work;
00132     return cfl;
00133   }
00134 
00135 private:
00136   void copySlice(DataType* target, const ld_vec_grid_data_type &source, 
00137                  const int meqn, const int move) {
00138     const BBox& bbox = source.bbox();
00139     BeginFastIndex1(src, bbox, source.data(), const VectorType);
00140     for (int m=0; m<meqn; m++) {
00141       if (move) target++;
00142       BeginFastIndex1(tgt, bbox, target, DataType);
00143       for_1 (k, bbox, bbox.stepsize()) 
00144         FastIndex1(tgt, k) = FastIndex1(src, k)(m);
00145       end_for
00146       EndFastIndex1(tgt);
00147       if (!move) target++;
00148       target = &(target[bbox.size()]);
00149     }
00150     EndFastIndex1(src);
00151   }
00152 
00153   void copyAuxSlice(DataType* target, const ld_aux_grid_data_type &source, 
00154                     const int meqn, const int move) {
00155     const BBox& bbox = source.bbox();
00156     BeginFastIndex1(src, bbox, source.data(), const AuxVectorType);
00157     for (int m=0; m<meqn; m++) {
00158       if (move) target++;
00159       BeginFastIndex1(tgt, bbox, target, DataType);
00160       for_1 (k, bbox, bbox.stepsize()) 
00161         FastIndex1(tgt, k) = FastIndex1(src, k)(m);
00162       end_for
00163       EndFastIndex1(tgt);
00164       if (!move) target++;
00165       target = &(target[bbox.size()]);
00166     }
00167     EndFastIndex1(src);
00168   }
00169 
00170   void addSlices(ld_vec_grid_data_type &target,
00171                  DataType* source1, DataType* source2, 
00172                  const int meqn, const int move) {
00173     const BBox& bbox = target.bbox();
00174     BeginFastIndex1(tgt, bbox, target.data(), VectorType);
00175     for (int m=0; m<meqn; m++) {
00176       if (move) { source1++; source2++; } 
00177       BeginFastIndex1(src1, bbox, source1, const DataType);
00178       BeginFastIndex1(src2, bbox, source2, const DataType);
00179       for_1 (k, bbox, bbox.stepsize()) 
00180         FastIndex1(tgt, k)(m) = FastIndex1(src1, k) + FastIndex1(src2, k);
00181       end_for
00182       EndFastIndex1(src1);
00183       EndFastIndex1(src2);
00184       if (!move) { source1++; source2++; } 
00185       source1 = &(source1[bbox.size()]);
00186       source2 = &(source2[bbox.size()]);
00187     }
00188     EndFastIndex1(tgt);
00189   }
00190 
00191   transverse_func_type transverse_flux_func;
00192 };
00193 
00194 
00195 
00196 #endif
ClpIntegrator2.h

last update: 06/01/04
