Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/PDE/Transport/Problem/ShearDiff.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 transport equations
9 : : \details This file defines a Problem policy class for the transport
10 : : equations, defined in PDE/Transport/CGTransport.h implementing
11 : : node-centered continuous Galerkin (CG) and PDE/Transport/DGTransport.h
12 : : implementing cell-centered discontinuous Galerkin (DG) discretizations.
13 : : See PDE/Transport/Problem.h for general requirements on Problem policy
14 : : classes for cg::Transport and dg::Transport.
15 : : */
16 : : // *****************************************************************************
17 : :
18 : : #include "ShearDiff.hpp"
19 : : #include "Inciter/InputDeck/InputDeck.hpp"
20 : :
21 : : namespace inciter {
22 : :
23 : : extern ctr::InputDeck g_inputdeck;
24 : :
25 : : } // ::inciter
26 : :
27 : : using inciter::TransportProblemShearDiff;
28 : :
29 : : std::vector< tk::real >
30 : 722103 : TransportProblemShearDiff::initialize( ncomp_t system, ncomp_t ncomp,
31 : : tk::real x, tk::real y, tk::real z, tk::real t )
32 : : // *****************************************************************************
33 : : // Evaluate analytical solution at (x,y,z,t) for all components
34 : : //! \param[in] system Equation system index
35 : : //! \param[in] ncomp Number of components in this transport equation system
36 : : //! \param[in] x X coordinate where to evaluate the solution
37 : : //! \param[in] y Y coordinate where to evaluate the solution
38 : : //! \param[in] z Z coordinate where to evaluate the solution
39 : : //! \param[in] t Time where to evaluate the solution
40 : : //! \return Values of all components evaluated at (x,y,t)
41 : : // *****************************************************************************
42 : : {
43 : : using tag::param;
44 : :
45 : 722103 : const auto& u0 = g_inputdeck.get< param, eq, tag::u0 >()[system];
46 : 722103 : const auto& d = g_inputdeck.get< param, eq, tag::diffusivity >()[system];
47 : 722103 : const auto& l = g_inputdeck.get< param, eq, tag::lambda >()[system];
48 : :
49 [ + - ]: 722103 : std::vector< tk::real > r( ncomp );
50 [ + + ]: 1513289 : for (ncomp_t c=0; c<ncomp; ++c) {
51 : 791186 : const auto li = 2*c;
52 : 791186 : const auto di = 3*c;
53 : 791186 : const auto phi3s = (l[li+0]*l[li+0]*d[di+1]/d[di+0] +
54 : 791186 : l[li+1]*l[li+1]*d[di+2]/d[di+0]) / 12.0;
55 : 791186 : r[c] =
56 : 791186 : 1.0 / ( 8.0 * std::pow(M_PI,3.0/2.0) *
57 : 791186 : std::sqrt(d[di+0]*d[di+1]*d[di+2]) *
58 : 1582372 : std::pow(t,3.0/2.0) * std::sqrt(1.0+phi3s*t*t) ) *
59 : 791186 : exp( -std::pow( x - u0[c]*t -
60 : 791186 : 0.5*(l[li+0]*y + l[li+1]*z)*t, 2.0 ) /
61 : 791186 : ( 4.0 * d[di+0] * t * (1.0 + phi3s*t*t) )
62 : 791186 : -y*y / ( 4.0 * d[di+1] * t )
63 : 791186 : -z*z / ( 4.0 * d[di+2] * t ) );
64 : : }
65 : :
66 : 722103 : return r;
67 : : }
68 : :
69 : : void
70 : 13 : TransportProblemShearDiff::errchk( ncomp_t system, ncomp_t ncomp ) const
71 : : // *****************************************************************************
72 : : // Do error checking on PDE parameters
73 : : //! \param[in] system Equation system index, i.e., which transport equation
74 : : //! system we operate on among the systems of PDEs
75 : : //! \param[in] ncomp Number of components in this transport equation
76 : : // *****************************************************************************
77 : : {
78 : : using tag::param;
79 : :
80 : 13 : const auto& u0 = g_inputdeck.get< param, eq, tag::u0 >()[system];
81 [ - + ][ - - ]: 13 : ErrChk( ncomp == u0.size(),
[ - - ][ - - ]
82 : : "Wrong number of advection-diffusion PDE parameters 'u0'" );
83 : :
84 : 13 : const auto& lambda = g_inputdeck.get< param, eq, tag::lambda >()[system];
85 [ - + ][ - - ]: 13 : ErrChk( 2*ncomp == lambda.size(),
[ - - ][ - - ]
86 : : "Wrong number of advection-diffusion PDE parameters 'lambda'" );
87 : :
88 : 13 : const auto& d = g_inputdeck.get< param, eq, tag::diffusivity >()[system];
89 [ - + ][ - - ]: 13 : ErrChk( 3*ncomp == d.size(),
[ - - ][ - - ]
90 : : "Wrong number of advection-diffusion PDE parameters 'diffusivity'" );
91 : 13 : }
92 : :
93 : : std::vector< std::array< tk::real, 3 > >
94 : 7878780 : TransportProblemShearDiff::prescribedVelocity( ncomp_t system, ncomp_t ncomp,
95 : : tk::real, tk::real y, tk::real z,
96 : : tk::real )
97 : : // *****************************************************************************
98 : : // Assign prescribed shear velocity at a point
99 : : //! \param[in] system Equation system index, i.e., which transport equation
100 : : //! system we operate on among the systems of PDEs
101 : : //! \param[in] ncomp Number of components in this transport equation
102 : : //! \param[in] y y coordinate at which to assign velocity
103 : : //! \param[in] z Z coordinate at which to assign velocity
104 : : //! \return Velocity assigned to all vertices of a tetrehedron, size:
105 : : //! ncomp * ndim = [ncomp][3]
106 : : // *****************************************************************************
107 : : {
108 : : using tag::param;
109 : :
110 : 7878780 : const auto& u0 = g_inputdeck.get< param, eq, tag::u0 >()[ system ];
111 : 7878780 : const auto& l = g_inputdeck.get< param, eq, tag::lambda >()[ system ];
112 : :
113 [ + - ]: 7878780 : std::vector< std::array< tk::real, 3 > > vel( ncomp );
114 [ + + ]: 16883100 : for (ncomp_t c=0; c<ncomp; ++c)
115 : 9004320 : vel[c] = {{ u0[c] + l[2*c+0]*y + l[2*c+1]*z, 0.0, 0.0 }};
116 : :
117 : 7878780 : return vel;
118 : : }
|
