Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/PDE/CompFlow/Problem/NLEnergyGrowth.hpp
4 : : \copyright 2012-2015 J. Bakosi,
5 : : 2016-2018 Los Alamos National Security, LLC.,
6 : : 2019-2021 Triad National Security, LLC.
7 : : All rights reserved. See the LICENSE file for details.
8 : : \brief Problem configuration for the compressible flow equations
9 : : \details This file declares a problem policy class for the compressible flow
10 : : equations, defined in PDE/CompFlow/CompFlow.h. See PDE/CompFlow/Problem.h
11 : : for general requirements on Problem policy classes for CompFlow.
12 : : */
13 : : // *****************************************************************************
14 : : #ifndef CompFlowProblemNLEnergyGrowth_h
15 : : #define CompFlowProblemNLEnergyGrowth_h
16 : :
17 : : #include <string>
18 : : #include <unordered_set>
19 : :
20 : : #include "Types.hpp"
21 : : #include "Fields.hpp"
22 : : #include "FunctionPrototypes.hpp"
23 : : #include "SystemComponents.hpp"
24 : : #include "Inciter/Options/Problem.hpp"
25 : : #include "Inciter/InputDeck/InputDeck.hpp"
26 : : #include "EoS/EoS.hpp"
27 : :
28 : : namespace inciter {
29 : :
30 : : extern ctr::InputDeck g_inputdeck;
31 : :
32 : : //! CompFlow system of PDEs problem: nonlinear energy growth (NLEG)
33 : : //! \see Waltz, et. al, "Manufactured solutions for the three-dimensional Euler
34 : : //! equations with relevance to Inertial Confinement Fusion", Journal of
35 : : //! Computational Physics 267 (2014) 196-209.
36 : : class CompFlowProblemNLEnergyGrowth {
37 : :
38 : : private:
39 : : using ncomp_t = tk::ctr::ncomp_t;
40 : : using eq = tag::compflow;
41 : :
42 : : //! Compute internal energy parameter
43 : : static tk::real hx( tk::real bx, tk::real by, tk::real bz,
44 : : tk::real x, tk::real y, tk::real z );
45 : :
46 : : //! Compute a power of the internal energy
47 : : static tk::real ec( tk::real ce, tk::real kappa, tk::real t, tk::real h,
48 : : tk::real p );
49 : :
50 : : public:
51 : : //! Initialize numerical solution
52 : : static tk::InitializeFn::result_type
53 : : initialize( ncomp_t system, ncomp_t, tk::real x, tk::real y,
54 : : tk::real z, tk::real t );
55 : :
56 : : //! Evaluate analytical solution at (x,y,z,t) for all components
57 : : static tk::InitializeFn::result_type
58 : : analyticSolution( ncomp_t system, ncomp_t, tk::real x, tk::real y,
59 : : tk::real z, tk::real t );
60 : :
61 : : //! Compute and return source term for NLEG manufactured solution
62 : : //! \param[in] system Equation system index, i.e., which compressible
63 : : //! flow equation system we operate on among the systems of PDEs
64 : : //! \param[in] x X coordinate where to evaluate the solution
65 : : //! \param[in] y Y coordinate where to evaluate the solution
66 : : //! \param[in] z Z coordinate where to evaluate the solution
67 : : //! \param[in] t Physical time at which to evaluate the source
68 : : //! \param[in,out] r Density source
69 : : //! \param[in,out] ru X momentum source
70 : : //! \param[in,out] rv Y momentum source
71 : : //! \param[in,out] rw Z momentum source
72 : : //! \param[in,out] re Specific total energy source
73 : : //! \note The function signature must follow tk::SrcFn
74 : : static tk::CompFlowSrcFn::result_type
75 : 3162360 : src( ncomp_t system, tk::real x, tk::real y, tk::real z,
76 : : tk::real t, tk::real& r, tk::real& ru, tk::real& rv, tk::real& rw,
77 : : tk::real& re )
78 : : {
79 : : using tag::param; using std::sin; using std::cos;
80 : : // manufactured solution parameters
81 : 3162360 : const auto a = g_inputdeck.get< param, eq, tag::alpha >()[system];
82 : 3162360 : const auto bx = g_inputdeck.get< param, eq, tag::betax >()[system];
83 : 3162360 : const auto by = g_inputdeck.get< param, eq, tag::betay >()[system];
84 : 3162360 : const auto bz = g_inputdeck.get< param, eq, tag::betaz >()[system];
85 : 3162360 : const auto ce = g_inputdeck.get< param, eq, tag::ce >()[system];
86 : 3162360 : const auto kappa = g_inputdeck.get< param, eq, tag::kappa >()[system];
87 : 3162360 : const auto r0 = g_inputdeck.get< param, eq, tag::r0 >()[system];
88 : : // ratio of specific heats
89 : : const auto g = gamma< tag::compflow >(system);
90 : : // spatial component of density field
91 : 3162360 : const auto gx = 1.0 - x*x - y*y - z*z;
92 : : // derivative of spatial component of density field
93 : 3162360 : const std::array< tk::real, 3 > dg{{ -2.0*x, -2.0*y, -2.0*z }};
94 : : // spatial component of energy field
95 : 3162360 : const auto h = hx( bx, by, bz, x, y, z );
96 : : // derivative of spatial component of energy field
97 : : std::array< tk::real, 3 >
98 : 3162360 : dh{{ -bx*M_PI*sin(bx*M_PI*x)*cos(by*M_PI*y)*cos(bz*M_PI*z),
99 : 3162360 : -by*M_PI*cos(bx*M_PI*x)*sin(by*M_PI*y)*cos(bz*M_PI*z),
100 : 3162360 : -bz*M_PI*cos(bx*M_PI*x)*cos(by*M_PI*y)*sin(bz*M_PI*z) }};
101 : : // temporal function f and its derivative
102 : 3162360 : const auto ft = std::exp(-a*t);
103 : 3162360 : const auto dfdt = -a*ft;
104 : : // density and its derivatives
105 : 3162360 : const auto rho = r0 + ft*gx;
106 : 3162360 : const std::array< tk::real, 3 > drdx{{ ft*dg[0], ft*dg[1], ft*dg[2] }};
107 : 3162360 : const auto drdt = gx*dfdt;
108 : : // internal energy and its derivatives
109 : 3162360 : const auto ie = ec( ce, kappa, t, h, -1.0/3.0 );
110 : : const std::array< tk::real, 3 > dedx{{
111 : 3162360 : 2.0*std::pow(ie,4.0)*kappa*h*dh[0]*t,
112 : 3162360 : 2.0*std::pow(ie,4.0)*kappa*h*dh[1]*t,
113 : 3162360 : 2.0*std::pow(ie,4.0)*kappa*h*dh[2]*t }};
114 : 3162360 : const auto dedt = kappa*h*h*std::pow(ie,4.0);
115 : : // density source
116 : 3162360 : r = drdt;
117 : : // momentum source
118 : 3162360 : ru = (g-1.0)*(rho*dedx[0] + ie*drdx[0]);
119 : 3162360 : rv = (g-1.0)*(rho*dedx[1] + ie*drdx[1]);
120 : 3162360 : rw = (g-1.0)*(rho*dedx[2] + ie*drdx[2]);
121 : : // energy source
122 : 3162360 : re = rho*dedt + ie*drdt;
123 : 3162360 : }
124 : :
125 : : //! Return analytic field names to be output to file
126 : : std::vector< std::string > analyticFieldNames( ncomp_t ) const;
127 : :
128 : : //! Return names of integral variables to be output to diagnostics file
129 : : std::vector< std::string > names( ncomp_t ) const;
130 : :
131 : : //! Return problem type
132 : : static ctr::ProblemType type() noexcept
133 : : { return ctr::ProblemType::NL_ENERGY_GROWTH; }
134 : : };
135 : :
136 : : } // inciter::
137 : :
138 : : #endif // CompFlowProblemNLEnergyGrowth_h
|