amroc/lbm/LBMLevelTransfer.h

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 // Copyright (C) 2013 Ralf Deiterding
00004 // German Aerospace Center (DLR)
00005 
00006 #ifndef LBM_LEVELTRANSFER_H
00007 #define LBM_LEVELTRANSFER_H
00008 
00017 #include "F77Interfaces/F77LevelTransfer.h"  
00018 #include "LBMFixupOps.h"
00019 
00020 #if DIM == 1 
00021 
00022 #define f_restrict FORTRAN_NAME(restrict1, RESTRICT1)
00023 #define f_prolong FORTRAN_NAME(prolong1, PROLONG1)
00024 
00025 #elif DIM == 2 
00026 
00027 #define f_restrict FORTRAN_NAME(restrict2, RESTRICT2)
00028 #define f_prolong FORTRAN_NAME(prolong2, PROLONG2)
00029 
00030 #elif DIM == 3 
00031 
00032 #define f_restrict FORTRAN_NAME(restrict3, RESTRICT3)
00033 #define f_prolong FORTRAN_NAME(prolong3, PROLONG3)
00034 
00035 #endif
00036 
00037 extern "C" {
00038   void f_restrict();
00039   void f_prolong();
00040 }
00041 
00049 template <class LBMType, int dim>
00050 class LBMF77LevelTransfer : public F77LevelTransfer<typename LBMType::MicroType,dim>,
00051                             public LBMFixupOps<LBMType,dim> {
00052   typedef F77LevelTransfer<typename LBMType::MicroType,dim> base;
00053   typedef LBMFixupOps<LBMType,dim> lbmfixupops;
00054   typedef typename LBMType::MicroType MicroType;
00055   typedef typename MicroType::InternalDataType DataType;
00056 
00057 public:
00058   typedef typename base::vec_grid_data_type vec_grid_data_type;
00059   typedef typename base::generic_func_type generic_func_type;
00060 
00061   LBMF77LevelTransfer() : base(), _TimeInterpolate(0), _RescaleNonEq(0) {}
00062   LBMF77LevelTransfer(LBMType &lbm, generic_func_type prolong, generic_func_type restrct) : 
00063     base(prolong,restrct), _lbm(lbm), _TimeInterpolate(0), _RescaleNonEq(0) {
00064     base::SetAdaptiveBoundaryType(DAGHAdaptBoundaryUserDef);
00065     for (int d=0; d<2*dim; d++) 
00066       indices[d] = (int *) 0;
00067   }
00068   
00069   virtual ~LBMF77LevelTransfer() {
00070     for (int d=0; d<2*dim; d++) 
00071       if (indices[d]) delete [] indices[d];
00072   }
00073 
00074   virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
00075     base::LocCtrl = Ctrl.getSubDevice(prefix+"LevelTransfer");
00076     RegisterAt(base::LocCtrl,"TimeInterpolate",_TimeInterpolate);
00077     RegisterAt(base::LocCtrl,"RescaleNonEq",_RescaleNonEq);
00078   }
00079 
00080   virtual void SetupData(GridHierarchy* gh, const int& ghosts) {
00081     base::SetupData(gh, ghosts);
00082     for (int d=0; d<2*dim; d++) {
00083       indices[d] = new int[LBM().NMicroVar()];
00084       Nindices[d] = LBM().IncomingIndices(d,indices[d]);
00085       for (register int n=Nindices[d]; n<LBM().NMicroVar(); n++) 
00086         indices[d][n]=-1;
00087     }
00088   }
00089 
00090   virtual void SetAdaptBndry(vec_grid_data_type &target, vec_grid_data_type &target_help_previous, 
00091                              vec_grid_data_type &target_help_next, vec_grid_data_type &source_next, const double &frac,
00092                              vec_grid_data_type &source_previous, const double &oneminusfrac, 
00093                              const int &target_level, const BBox &bb, const int side) {
00094     int source_level = target_level-1;
00095     if (frac==0.0) base::Prolong(source_previous, source_level, target, target_level, bb);
00096     else if (frac==1.0) base::Prolong(source_next, source_level, target, target_level, bb);
00097     else {
00098       if (_TimeInterpolate) {
00099         vec_grid_data_type work(bb);
00100         base::Prolong(source_previous, source_level, target_help_previous, target_level, bb);
00101         base::Prolong(source_next, source_level, target_help_next, target_level, bb);
00102         work.lin_interp(target_help_previous,oneminusfrac,target_help_next,frac,bb);
00103         if (_RescaleNonEq) {
00104           DataType fac=0.;
00105           if (DeltaT(base::GH(),target_level)>0.)
00106             fac=LBM().Omega(DeltaT(base::GH(),source_level))/
00107               (RefineFactor(base::GH(),source_level)*LBM().Omega(DeltaT(base::GH(),target_level)));
00108           lbmfixupops::equilibrium(LBM(), work, bb, indices[side], Nindices[side]);
00109           lbmfixupops::scale_indices(target, work, bb, indices[side], Nindices[side], fac);
00110         }
00111         else 
00112           lbmfixupops::copy_indices(target, work, bb, indices[side], Nindices[side]);
00113       }   
00114     }
00115   }
00116 
00117   inline LBMType& LBM() { return _lbm; }
00118   inline const LBMType& LBM() const { return _lbm; }
00119 
00120 protected:
00121   LBMType &_lbm;
00122   int _TimeInterpolate, _RescaleNonEq;
00123   int Nindices[2*dim];
00124   int *indices[2*dim];
00125 };
00126 
00127 #endif