00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
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
00184 void computeVorticity(LevelData<FArrayBox>& a_vorticity) const;
00185
00186
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
00281
00282
00283
00284 void predictVelocities(LevelData<FArrayBox>& a_uDelU,
00285 LevelData<FluxBox>& a_advVel);
00286
00287
00288
00289
00290
00291
00292
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
00532
00534 static bool s_write_divergence;
00535
00537 static bool s_write_lambda;
00538
00540 static bool s_write_time_derivatives;
00541
00543 static bool s_write_vorticity;
00544
00546 static bool s_write_scalars;
00547
00549 static bool s_write_dScalar_dt;
00550
00552 static bool s_write_strains;
00553
00555 static bool s_write_grad_eLambda;
00556
00558 static bool s_write_error;
00559
00561
00562
00564 static bool s_write_proc_ids;
00565
00567 Vector<Tuple<Real, SpaceDim> > m_dataPoints;
00568
00570 Real m_dx;
00571
00573 Real m_cfl;
00574
00576 Real m_dt_save;
00577
00579 int m_level_steps;
00580
00582 CoarseAverage m_coarse_average;
00583
00585 CoarseAverage m_coarse_average_lambda;
00586
00587 #if 0
00588
00589 CoarseAverage m_coarse_average_scal;
00590 #endif
00591
00593 CCProjector m_projection;
00594
00596 LevelFluxRegister m_flux_register;
00597
00599 LevelFluxRegister m_lambda_flux_reg;
00600
00602 Vector<LevelFluxRegister*> m_scal_fluxreg_ptrs;
00603
00605 QuadCFInterp m_velCFInterp;
00606
00608 VelocityPoissonOp m_velocityPoissonOp;
00609
00611 PoissonOp m_scalarsPoissonOp;
00612
00614 LevelHelmholtzSolver m_viscous_solver;
00615
00617 Real m_refine_threshold;
00618
00620 bool m_finest_level;
00621
00623 bool m_is_empty;
00624
00626 bool m_regrid_smoothing_done;
00627
00629 PhysBCUtil* m_physBCPtr;
00630
00631
00632
00633 };
00634
00635 #endif
00636
00637