1 // Copyright (c) 2017-2021, Lawrence Livermore National Security, LLC and 2 // other Axom Project Developers. See the top-level LICENSE file for details. 3 // 4 // SPDX-License-Identifier: (BSD-3-Clause) 5 6 7 // Hydro driver 8 // Fri Mar 27 14:14:28 PDT 2015 9 10 // HYDRO C: Physics-wise, this hydro is still pretty dumb: multiple 11 // materials, gamma-law gas. The mesh is an arbitrary collection of 12 // polygons. Architecturally, we use as few temporaries as possible. 13 // No sweet returning of objects when you can instead pass in an array 14 // and modify it in place. No STL containers. Explicit memory 15 // allocation of all the memory we ever need at problem start, and 16 // explicit deletes when the hydro object is destroyed. 17 18 #include "State.hpp" 19 20 #include "TinyHydroTypes.hpp" 21 22 namespace tinyHydro { 23 24 class Hydro 25 { 26 public: 27 28 public: 29 Hydro(State * s); // using a pointer here with PyBindGen prevents 30 // State being garbage collected by Python 31 ~Hydro(); 32 33 // problem control and advance methods 34 void step(double dt); 35 double newDT(); 36 void steps(int numSteps); 37 void advance(double stopTime); 38 void initialize(void); // call once before taking steps 39 void setBC(const char * boundary, double xVel, double yVel); 40 41 // compute time derivatives from first arg, store rates in 2nd arg 42 void calcDerivs(const State & state, State & dState, double dt); 43 44 // physics methods 45 void calcQ(const State & state); 46 void calcForce(const State & state); 47 void calcDivU(const NodalVectorField& u); 48 void applyVelocityBC(NodalVectorField& u); 49 double totalEnergy(const State & s) const; 50 double totalEnergy(void) const; 51 void calcMaxCs(void); 52 53 // update all the algebraically determinable stuff, like density 54 void updateConstitutives(State & state); 55 56 // setter methods (used from Python) setState(const State s)57 void setState(const State s) {state = s; } 58 59 // getter methods getState(void)60 State * getState(void) {return &state; } getQ(int i)61 double getQ(int i) {return Q[i]; } getForce(int i)62 VectorXY getForce(int i) {return force[i]; } 63 int numBCnodes(int bc); 64 int bcNode(int bc, int node); 65 fuzzyEqual(double a,double b,double tol=1.0e-15)66 bool fuzzyEqual(double a, double b, double tol = 1.0e-15) 67 { 68 return (fabs(a - b) < tol); 69 } 70 71 public: 72 PolygonMeshXY& mesh; 73 74 State& state; 75 State halfStep; // re-used temporary for RK2 76 State dState; // re-used temporary for RK2 77 78 ZonalScalarField Q; // re-used for Q 79 ZonalScalarField divu; // for dt, Q, and some PdV work schemes 80 81 NodalVectorField force; // re-used for force calc 82 NodalVectorField halfStepVelocity; // re-used for work calc 83 84 BoundaryEdgeSet boundaryEdgeSet; // SLAM -- new set for domain boundaries (of size 4) 85 BoundaryEdgeVectorField bcVelocity; // velocity BC values // SLAM -- map to velocity boundary conditions for each edge 86 NodeSubset bcNodes[NUM_DOMAIN_BOUNDARIES]; // SLAM -- indirection sets of nodes on each boundary edge 87 88 89 ZonalScalarField* initialMass; // per material initial mass of zones 90 91 ZonalScalarField totalMass; // overall mass in zone 92 NodalScalarField nodeMass; // also keep node masses 93 94 ZonalScalarField* pressure; // per material zone pressure 95 96 ZonalScalarField totalPressure; // sum of pressures in a zone 97 ZonalScalarField maxCs; // 98 99 100 int cycle; 101 int dtZone; 102 double time; 103 double Cq, Cl; // Q coefficients 104 double cfl; 105 int printCycle; 106 107 private: 108 // Unused variables... 109 // int numBCs; // <UNUSED>? 110 111 #ifdef TINY_HYDRO_USE_CORNER_DATA 112 VectorXY * cornerforce; // re-used for force calc // <UNUSED>? 113 #endif 114 }; 115 116 } // end namespace tinyHydro 117