Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/PDE/CompFlow/Problem/GaussHumpCompflow.cpp
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 defines 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 : :
15 : : #include "GaussHumpCompflow.hpp"
16 : : #include "FieldOutput.hpp"
17 : :
18 : : using inciter::CompFlowProblemGaussHump;
19 : :
20 : : tk::InitializeFn::result_type
21 : 650155 : CompFlowProblemGaussHump::initialize( ncomp_t ncomp,
22 : : const std::vector< EOS >& mat_blk,
23 : : tk::real x,
24 : : tk::real y,
25 : : tk::real,
26 : : tk::real t )
27 : : // *****************************************************************************
28 : : //! Evaluate analytical solution at (x,y,z,t) for all components
29 : : //! \param[in] ncomp Number of scalar components in this PDE system
30 : : //! \param[in] x X coordinate where to evaluate the solution
31 : : //! \param[in] y Y coordinate where to evaluate the solution
32 : : //! \param[in] t Time where to evaluate the solution
33 : : //! \return Values of all components evaluated at (x)
34 : : //! \note The function signature must follow tk::InitializeFn
35 : : // *****************************************************************************
36 : : {
37 : : Assert( ncomp == 5, "Number of scalar components must be 5" );
38 : :
39 : 650155 : const auto vel = prescribedVelocity( ncomp, x, y, 0.0 );
40 : :
41 : : tk::real r, p, u, v, w, rE;
42 : :
43 : : // center of the hump
44 : 650155 : auto x0 = 0.25 + vel[0][0]*t;
45 : 650155 : auto y0 = 0.25 + vel[0][1]*t;
46 : :
47 : : // density
48 : 650155 : r = 1.0 + exp( -((x-x0)*(x-x0)
49 : 650155 : + (y-y0)*(y-y0))/(2.0 * 0.005) );
50 : : // pressure
51 : 650155 : p = 1.0;
52 : : // velocity
53 : 650155 : u = 1;
54 : 650155 : v = 1;
55 : 650155 : w = 0;
56 : : // total specific energy
57 [ + - ]: 650155 : rE = mat_blk[0].compute< EOS::totalenergy >( r, u, v, w, p );
58 : :
59 [ + - ][ + - ]: 1300310 : return {{ r, r*u, r*v, r*w, rE }};
[ - - ]
60 : : }
61 : :
62 : : tk::InitializeFn::result_type
63 : 131286 : CompFlowProblemGaussHump::analyticSolution( ncomp_t ncomp,
64 : : const std::vector< EOS >& mat_blk,
65 : : tk::real x,
66 : : tk::real y,
67 : : tk::real,
68 : : tk::real t )
69 : : // *****************************************************************************
70 : : //! Evaluate analytical solution at (x,y,z,t) for all components
71 : : //! \param[in] ncomp Number of scalar components in this PDE system
72 : : //! \param[in] x X coordinate where to evaluate the solution
73 : : //! \param[in] y Y coordinate where to evaluate the solution
74 : : //! \param[in] t Time where to evaluate the solution
75 : : //! \return Values of all components evaluated at (x)
76 : : //! \note The function signature must follow tk::InitializeFn
77 : : // *****************************************************************************
78 : : {
79 : : Assert( ncomp == 5, "Number of scalar components must be 5" );
80 : :
81 : 131286 : const auto vel = prescribedVelocity( ncomp, x, y, 0.0 );
82 : :
83 : : // center of the hump
84 : 131286 : auto x0 = 0.25 + vel[0][0]*t;
85 : 131286 : auto y0 = 0.25 + vel[0][1]*t;
86 : :
87 : : // density
88 : 131286 : auto r = 1.0 + exp( -((x-x0)*(x-x0) + (y-y0)*(y-y0))/(2.0 * 0.005) );
89 : : // pressure
90 : 131286 : auto p = 1.0;
91 : : // velocity
92 : 131286 : auto u = 1.0;
93 : 131286 : auto v = 1.0;
94 : 131286 : auto w = 0.0;
95 : : // total specific energy
96 [ + - ]: 131286 : auto E = mat_blk[0].compute< EOS::totalenergy >( r, u, v, w, p ) / r;
97 : :
98 [ + - ][ + - ]: 262572 : return {{ r, u, v, w, E, p }};
[ - - ]
99 : : }
100 : :
101 : : std::vector< std::string >
102 [ + - ]: 138 : CompFlowProblemGaussHump::analyticFieldNames( ncomp_t ) const
103 : : // *****************************************************************************
104 : : // Return analytic field names to be output to file
105 : : //! \return Vector of strings labelling fields output in file
106 : : // *****************************************************************************
107 : : {
108 : : std::vector< std::string > n;
109 [ + - ]: 138 : n.push_back( "density_analytical" );
110 [ + - ]: 138 : n.push_back( "x-velocity_analytical" );
111 [ + - ]: 138 : n.push_back( "y-velocity_analytical" );
112 [ + - ]: 138 : n.push_back( "z-velocity_analytical" );
113 [ + - ]: 138 : n.push_back( "specific_total_energy_analytical" );
114 [ + - ]: 138 : n.push_back( "pressure_analytical" );
115 : :
116 : 138 : return n;
117 : : }
118 : :
119 : : std::vector< std::string >
120 : 5 : CompFlowProblemGaussHump::names( ncomp_t ) const
121 : : // *****************************************************************************
122 : : // Return names of integral variables to be output to diagnostics file
123 : : //! \return Vector of strings labelling integral variables output
124 : : // *****************************************************************************
125 : : {
126 [ + - ][ + - ]: 30 : return { "r", "ru", "rv", "rw", "re" };
[ + - ][ + - ]
[ + - ][ - + ]
[ + + ][ - + ]
[ - - ][ - - ]
127 : : }
128 : :
129 : : std::vector< std::array< tk::real, 3 > >
130 : 781441 : CompFlowProblemGaussHump::prescribedVelocity( ncomp_t ncomp, tk::real,
131 : : tk::real, tk::real )
132 : : // *****************************************************************************
133 : : //! Assign prescribed velocity at a point
134 : : //! \param[in] ncomp Number of components in this transport equation
135 : : //! \return Velocity assigned to all vertices of a tetrehedron, size:
136 : : //! ncomp * ndim = [ncomp][3]
137 : : // *****************************************************************************
138 : : {
139 : 781441 : std::vector< std::array< tk::real, 3 > > vel( ncomp );
140 : :
141 [ + + ]: 4688646 : for (ncomp_t c=0; c<ncomp; ++c)
142 : 3907205 : vel[c] = {{ 1, 1, 0.0 }};
143 : :
144 : 781441 : return vel;
145 : : }
|