Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Compound Members   File Members  

AMRNavierStokes.H

Go to the documentation of this file.
00001 /* _______              __
00002   / ___/ /  ___  __ _  / /  ___
00003  / /__/ _ \/ _ \/  ' \/ _ \/ _ \
00004  \___/_//_/\___/_/_/_/_.__/\___/ 
00005 */
00006 
00007 // AMRNavierStokes.H
00008 
00009 //
00010 // This software is copyright (C) by the Lawrence Berkeley
00011 // National Laboratory.  Permission is granted to reproduce
00012 // this software for non-commercial purposes provided that
00013 // this notice is left intact.
00014 // 
00015 // It is acknowledged that the U.S. Government has rights to
00016 // this software under Contract DE-AC03-765F00098 between
00017 // the U.S.  Department of Energy and the University of
00018 // California.
00019 //
00020 // This software is provided as a professional and academic
00021 // contribution for joint exchange. Thus it is experimental,
00022 // is provided ``as is'', with no warranties of any kind
00023 // whatsoever, no support, no promise of updates, or printed
00024 // documentation. By using this software, you acknowledge
00025 // that the Lawrence Berkeley National Laboratory and
00026 // Regents of the University of California shall have no
00027 // liability with respect to the infringement of other
00028 // copyrights by any part of this software.
00029 //
00030 // Dan Martin, Fri, Jan 14, 2000
00031 
00032 #ifndef _AMRNavierStokes_H_
00033 #define _AMRNavierStokes_H_
00034 
00035 #include "AMRLevel.H"
00036 #include "LevelFluxRegister.H"
00037 #include "CoarseAverage.H"
00038 #include "PiecewiseLinearFillPatch.H"
00039 #include "CCProjector.H"
00040 #include "LevelHelmholtzSolver.H"
00041 #include "PoissonOp.H"
00042 #include "VelocityPoissonOp.H"
00043 #include "PhysBCUtil.H"
00045 enum viscousSolverTypes {
00047   backwardEuler = 0,
00049   CrankNicolson,
00051   TGA, 
00053   NUM_SOLVER_TYPES};
00054 
00055 
00056 
00058 enum scalTypes {
00060   xMarker = 0,
00062   yMarker,
00063 #if CH_SPACEDIM >= 3
00064 
00065   zMarker,
00066 #endif
00067 
00068   NUM_SCAL_TYPES};
00069 
00070 
00071 
00072 
00074 class AMRNavierStokes: public AMRLevel
00075 {
00076 public:
00077   
00079   AMRNavierStokes();
00080   
00082   AMRNavierStokes(AMRLevel* a_coarser_level_ptr, 
00083            const Box& a_prob_domain,
00084            int a_level, int a_ref_ratio);
00085 
00087   AMRNavierStokes(AMRLevel* a_coarser_level_ptr, 
00088            const ProblemDomain& a_prob_domain,
00089            int a_level, int a_ref_ratio);
00090 
00092   virtual ~AMRNavierStokes();
00093 
00094 
00096   virtual AMRLevel* makeAMRLevel(AMRLevel* a_coarser_level_ptr,
00097                                  const Box& a_problem_domain,
00098                                  int a_level, int a_ref_ratio) const;
00099 
00101   virtual AMRLevel* makeAMRLevel(AMRLevel* a_coarser_level_ptr,
00102                                  const ProblemDomain& a_problem_domain,
00103                                  int a_level, int a_ref_ratio) const;
00104 
00106   virtual void define(AMRLevel* a_coarse_level_ptr, 
00107                       const Box& a_problem_domain,
00108                       int a_level, int a_ref_ratio);
00109 
00111   virtual void define(AMRLevel* a_coarse_level_ptr, 
00112                       const ProblemDomain& a_problem_domain,
00113                       int a_level, int a_ref_ratio);
00114 
00116   virtual Real advance();
00117 
00119   virtual void postTimeStep();
00120 
00122   virtual void tagCells(IntVectSet& a_tags);
00123 
00125   virtual void tagCellsInit(IntVectSet& a_tags);
00126 
00128   virtual void regrid(const Vector<Box>& a_new_grids);
00129 
00130 
00131 
00133   virtual void postRegrid(int a_lBase);
00134 
00135 
00137   virtual void initialGrid(const Vector<Box>& a_new_grids);
00138 
00140   virtual void initialData();
00141 
00143   virtual void postInitialize();
00144 
00145 #ifdef HDF5
00146 
00147   virtual void writeCheckpointHeader(HDF5Handle& a_handle) const;
00148 
00150   virtual void writeCheckpointLevel(HDF5Handle& a_handle) const;
00151 
00153   virtual void readCheckpointHeader(HDF5Handle& a_handle);
00154 
00156   virtual void readCheckpointLevel(HDF5Handle& a_handle);
00157 
00159   virtual void writePlotHeader(HDF5Handle& a_handle) const;
00160 
00162   virtual void writePlotLevel(HDF5Handle& a_handle) const;
00163 #endif
00164 
00166   virtual Real computeDt();
00167 
00169   virtual Real computeInitialDt();
00170 
00172   virtual void CFL(Real a_cfl);
00173 
00175   virtual void refinementThreshold(Real a_refine_threshold);
00176 
00178   bool finestLevel() const;
00179   
00181   bool isEmpty() const;
00182   
00183   // compute vorticity
00184   void computeVorticity(LevelData<FArrayBox>& a_vorticity) const;
00185 
00186   // compute kinetic energy (only on interior cells)
00187   void computeKineticEnergy(LevelData<FArrayBox>& a_energy) const;
00188 
00190 
00192   Real Dx() const;
00193   
00195   Real Nu() const;
00196 
00198   AMRNavierStokes* finerNSPtr() const;
00199 
00201   AMRNavierStokes* crseNSPtr() const;
00202 
00204   LevelData<FArrayBox>& newVel();
00205 
00207   const LevelData<FArrayBox>& newVel() const;
00208 
00210   LevelData<FArrayBox>& oldVel();
00211 
00213   const LevelData<FArrayBox>& oldVel() const;
00214 
00216   void velocity(LevelData<FArrayBox>& a_vel, Real a_time) const;
00217 
00219   LevelData<FArrayBox>& newLambda();
00220 
00222   const LevelData<FArrayBox>& newLambda() const;
00223 
00225   LevelData<FArrayBox>& oldLambda();
00226 
00228   const LevelData<FArrayBox>& oldLambda() const;
00229 
00231   void lambda(LevelData<FArrayBox>& a_lambda, Real a_time) const;
00232   
00234   LevelData<FArrayBox>& newScal(const int a_comp);
00235 
00237   const LevelData<FArrayBox>& newScal(const int a_comp) const;
00238 
00240   LevelData<FArrayBox>& oldScal(const int a_comp);
00241 
00243   const LevelData<FArrayBox>& oldScal(const int a_comp) const;
00244 
00246   void scalars(LevelData<FArrayBox>& a_scalars, const Real a_time, 
00247                const int a_comp) const;
00248 
00249   
00251   void setPhysBC(const PhysBCUtil& a_BC);
00252 
00254   PhysBCUtil* getPhysBCPtr() const;
00255 
00256 
00257 protected:
00259   DisjointBoxLayout loadBalance(const Vector<Box>& a_grids);
00260 
00262   void levelSetup (const DisjointBoxLayout& level_domain);
00263 
00265   void readParameters();
00266 
00268   void finestLevel(bool a_finest_level);
00269 
00271   void swapOldAndNewStates();
00272 
00274   void resetStates(const Real a_time);
00275 
00277   void computeAdvectionVelocities(LevelData<FluxBox>& a_adv_vel);
00278 
00279 
00280   /*@ManMemo: do predictor -- returns approx. to (U dot del)U.
00281     Also updates flux registers with momentum fluxes and
00282     part of viscous fluxes which contains the old-time velocity.
00283   */
00284   void predictVelocities(LevelData<FArrayBox>& a_uDelU,
00285                          LevelData<FluxBox>& a_advVel);
00286 
00287   /*@ManMemo: This function computes uStar -- the approximation
00288     of the new-time velocity which will be projected.  The
00289     advection term (u-dot-del-U) is passed in in a_uStar, which 
00290     is modified in place to contain uStar. Note that both viscous
00291     and inviscid cases, this returns the approximation to 
00292     u^{n+1} + dt*gradPi */
00293   void computeUStar(LevelData<FArrayBox>& a_uStar);
00294 
00295 
00296 
00298   void advectScalar(LevelData<FArrayBox>& a_new_scal,
00299                     LevelData<FArrayBox>& a_old_scal,
00300                     LevelData<FluxBox>& a_adv_vel,
00301                     LevelFluxRegister* a_crseFluxRegPtr,
00302                     LevelFluxRegister& a_fineFluxReg,
00303                     Real a_dt);
00304 
00306   void advectDiffuseScalar(LevelData<FArrayBox>& a_new_scal,
00307                            LevelData<FArrayBox>& a_old_scal,
00308                            LevelData<FluxBox>& a_adv_vel,
00309                            const Real& a_diffusive_coeff,
00310                            const LevelData<FArrayBox>* a_crse_scal_old,
00311                            const LevelData<FArrayBox>* a_crse_scal_new,
00312                            const Real a_old_crse_time,
00313                            const Real a_new_crse_time,
00314                            LevelFluxRegister* a_crseFluxRegPtr,
00315                            LevelFluxRegister& a_fineFluxReg,
00316                            DomainGhostBC& a_scalPhysBC,
00317                            Real a_dt);
00318 
00320   void initializeGlobalPressure();
00321 
00323   void initializeLevelPressure(Real a_currentTime, Real a_dtInit);
00324 
00326   void smoothVelocityField(int lbase);
00327 
00329   void computeLapVel(LevelData<FArrayBox>& a_lapVel,
00330                      LevelData<FArrayBox>& a_vel,
00331                      const LevelData<FArrayBox>* a_crseVelPtr);
00332 
00334   void computeLapScal(LevelData<FArrayBox>& a_lapScal,
00335                       LevelData<FArrayBox>& a_scal,
00336                       DomainGhostBC& a_physBC,
00337                       const LevelData<FArrayBox>* a_crseScalPtr);
00338 
00340   void fillVelocity(LevelData<FArrayBox>& a_vel, Real a_time);
00341 
00343   void fillVelocity(LevelData<FArrayBox>& a_vel, Real a_time,  
00344                     int vel_comp, int dest_comp, int num_comp);
00345 
00348   void fillLambda(LevelData<FArrayBox>& a_lambda, Real a_time) const;
00349 
00351   void fillScalars(LevelData<FArrayBox>& a_scal, Real a_time,
00352                    const int a_comp) const;
00353 
00354 
00356   virtual int numPlotComps() const;
00357 
00359   virtual void getPlotData(LevelData<FArrayBox>& a_plot_data) const;
00360 
00362   void  defineRegridAMRSolver(AMRSolver& a_solver, 
00363                               const Vector<DisjointBoxLayout>& a_grids,
00364                               const Vector<ProblemDomain>& a_domains,
00365                               const Vector<Real>& a_amrDx,
00366                               const Vector<int>& a_refRatios,
00367                               const int& a_lBase);
00368   
00369 
00370 
00371 
00372 
00373 protected:
00375   LevelData<FArrayBox>* m_vel_old_ptr;
00376 
00378   LevelData<FArrayBox>* m_vel_new_ptr;
00379 
00380 #ifdef DEBUG
00381 
00382   LevelData<FArrayBox>* m_vel_save_ptr;  
00383 
00384   Real m_saved_time;
00385 
00386 #endif
00387 
00389   LevelData<FArrayBox>* m_lambda_old_ptr;
00390   
00392   LevelData<FArrayBox>* m_lambda_new_ptr;
00393 
00394 #ifdef DEBUG
00395 
00396   LevelData<FArrayBox>* m_lambda_save_ptr;  
00397 
00398 
00399 #endif
00400 
00402   Vector<LevelData<FArrayBox>*> m_scal_old;
00403   
00405   Vector<LevelData<FArrayBox>*> m_scal_new;
00406 
00407 #ifdef DEBUG
00408 
00409   Vector<LevelData<FArrayBox>*> m_scal_save;
00410 #endif
00411 
00413   static const int s_num_vel_comps = SpaceDim;
00414 
00416   static int s_num_scal_comps;
00417 
00419   static const int s_max_scal_comps=SpaceDim;
00420 
00422   static const char* s_vel_names[s_num_vel_comps];
00423 
00425   static const char* s_scal_names[s_max_scal_comps];
00426 
00428   static Vector<Real> s_scal_coeffs;
00429 
00431   static Vector<Real> s_domLength;
00432   
00434   static bool s_ppInit;
00435 
00437   static Real s_init_shrink;
00438 
00440   static Real s_max_dt;
00441 
00442   static Real s_prescribedDt;
00443 
00445   static bool s_project_initial_vel;
00446 
00448   static bool s_initialize_pressures;
00449 
00451   static bool s_smooth_after_regrid;
00452 
00454   static bool s_dump_smoothing_plots;
00455 
00457   static Real s_regrid_smoothing_coeff;
00458 
00460   static int s_num_init_passes;
00461 
00463   static bool s_reflux_momentum;
00464 
00466   static bool s_set_bogus_values;
00467 
00469   static Real s_bogus_value;
00470 
00472   static bool s_tag_vorticity;
00473 
00475   static bool s_tag_theta;
00476 
00478   static Real s_vort_factor;
00479 
00481   static int s_tags_grow;
00482 
00484   static bool s_implicit_reflux;
00485 
00487   static bool s_applyFreestreamCorrection;
00488 
00490   static bool s_reflux_scal;
00491 
00493   static bool s_implicit_scal_reflux;
00494 
00496   static int s_viscous_solver_type;
00497 
00499   static Real s_viscous_solver_tol;
00500 
00502   static Real s_nu;
00503 
00505   static bool s_boussinesq;
00506 
00508   static Real s_theta_coeff;
00509 
00511   static Real s_rayleighNo;
00512 
00514   static Real s_prandtlNo;
00515 
00517   static Real s_wallTemp;
00518 
00520   static bool s_specifyInitialGrids;
00521 
00523   static string s_initialGridFile;
00524 
00526   static bool s_init_vel_from_vorticity;
00527   
00529   static Real s_backgroundVel;
00530   
00531 
00533   Vector<Tuple<Real, SpaceDim> > m_dataPoints;
00534 
00536   Real m_dx;
00537 
00539   Real m_cfl;
00540 
00542   Real m_dt_save;
00543 
00545   int m_level_steps;
00546 
00548   CoarseAverage m_coarse_average;
00549 
00551   CoarseAverage m_coarse_average_lambda;
00552 
00553 #if 0
00554 
00555   CoarseAverage m_coarse_average_scal;
00556 #endif
00557 
00559   CCProjector m_projection;
00560   
00562   LevelFluxRegister m_flux_register;
00563 
00565   LevelFluxRegister m_lambda_flux_reg;
00566 
00568   Vector<LevelFluxRegister*> m_scal_fluxreg_ptrs;
00569 
00571   QuadCFInterp m_velCFInterp; 
00572 
00574   VelocityPoissonOp m_velocityPoissonOp;
00575 
00577   PoissonOp m_scalarsPoissonOp;
00578 
00580   LevelHelmholtzSolver m_viscous_solver;
00581 
00583   Real m_refine_threshold;
00584 
00586   bool m_finest_level;
00587 
00589   bool m_is_empty;
00590 
00592   bool m_regrid_smoothing_done;
00593 
00595   PhysBCUtil* m_physBCPtr;
00596 
00597 
00598 };
00599 
00600 #endif
00601 
00602 

Generated on Wed Apr 30 18:00:39 2003 for Chombo&INS by doxygen1.2.16