1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163 | // *****************************************************************************
/*!
\file src/PDE/MultiSpecies/MiscMultiSpeciesFns.hpp
\copyright 2012-2015 J. Bakosi,
2016-2018 Los Alamos National Security, LLC.,
2019-2021 Triad National Security, LLC.
All rights reserved. See the LICENSE file for details.
\brief Misc multi-species system functions
\details This file defines functions that required for multi-species
compressible fluid dynamics.
*/
// *****************************************************************************
#include <iostream>
#include "MiscMultiSpeciesFns.hpp"
#include "Inciter/InputDeck/InputDeck.hpp"
#include "Integrate/Basis.hpp"
#include "MultiSpecies/MultiSpeciesIndexing.hpp"
#include "EoS/GetMatProp.hpp"
namespace inciter {
extern ctr::InputDeck g_inputdeck;
void initializeSpeciesEoS( std::vector< EOS >& mat_blk )
// *****************************************************************************
// Initialize the species block with configured EOS
//! \param[in,out] mat_blk Material block that gets initialized
// *****************************************************************************
{
// EoS initialization
// ---------------------------------------------------------------------------
auto nspec = g_inputdeck.get< tag::multispecies, tag::nspec >();
const auto& matprop = g_inputdeck.get< tag::material >();
const auto& matidxmap = g_inputdeck.get< tag::matidxmap >();
// assume only one type of species
auto mateos = matprop[matidxmap.get< tag::eosidx >()[0]].get<tag::eos>();
for (std::size_t k=0; k<nspec; ++k) {
mat_blk.emplace_back(mateos, EqType::multispecies, k);
}
}
tk::real
timeStepSizeMultiSpecies(
const std::vector< EOS >& mat_blk,
const std::vector< int >& esuf,
const tk::Fields& geoFace,
const tk::Fields& geoElem,
const std::size_t nelem,
std::size_t nspec,
const tk::Fields& U,
const tk::Fields& /*P*/ )
// *****************************************************************************
// Time step restriction for multi species cell-centered schemes
//! \param[in] mat_blk EOS species block
//! \param[in] esuf Elements surrounding elements array
//! \param[in] geoFace Face geometry array
//! \param[in] geoElem Element geometry array
//! \param[in] nelem Number of elements
//! \param[in] nspec Number of speciess in this PDE system
//! \param[in] U High-order solution vector
// //! \param[in] P High-order vector of primitives
//! \return Maximum allowable time step based on cfl criterion
// *****************************************************************************
{
const auto ndof = g_inputdeck.get< tag::ndof >();
const auto rdof = g_inputdeck.get< tag::rdof >();
std::size_t ncomp = U.nprop()/rdof;
tk::real u, v, w, a, vn, dSV_l, dSV_r;<--- The scope of the variable 'u' can be reduced. [+]The scope of the variable 'u' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level. <--- The scope of the variable 'v' can be reduced. [+]The scope of the variable 'v' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level. <--- The scope of the variable 'w' can be reduced. [+]The scope of the variable 'w' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level. <--- The scope of the variable 'a' can be reduced. [+]The scope of the variable 'a' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level. <--- The scope of the variable 'vn' can be reduced. [+]The scope of the variable 'vn' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level. <--- The scope of the variable 'dSV_l' can be reduced. [+]The scope of the variable 'dSV_l' can be reduced. Warning: Be careful when fixing this message, especially when there are inner loops. Here is an example where cppcheck will write that the scope for 'i' can be reduced:
void f(int x)
{
int i = 0;
if (x) {
// it's safe to move 'int i = 0;' here
for (int n = 0; n < 10; ++n) {
// it is possible but not safe to move 'int i = 0;' here
do_something(&i);
}
}
}
When you see this message it is always safe to reduce the variable scope 1 level.
std::vector< tk::real > delt(U.nunk(), 0.0);
std::vector< tk::real > ugp(ncomp, 0.0);<--- Variable 'ugp' is assigned a value that is never used.
// compute maximum characteristic speed at all internal element faces
for (std::size_t f=0; f<esuf.size()/2; ++f)
{
std::size_t el = static_cast< std::size_t >(esuf[2*f]);
auto er = esuf[2*f+1];
std::array< tk::real, 3 > fn;
fn[0] = geoFace(f,1);
fn[1] = geoFace(f,2);
fn[2] = geoFace(f,3);
// left element
// Compute the basis function for the left element
std::vector< tk::real > B_l(rdof, 0.0);
B_l[0] = 1.0;
// get conserved quantities
ugp = eval_state(ncomp, rdof, ndof, el, U, B_l);
// mixture density
tk::real rhob(0.0);
for (std::size_t k=0; k<nspec; ++k)
rhob += ugp[multispecies::densityIdx(nspec, k)];
// advection velocity
u = ugp[multispecies::momentumIdx(nspec, 0)]/rhob;
v = ugp[multispecies::momentumIdx(nspec, 1)]/rhob;
w = ugp[multispecies::momentumIdx(nspec, 2)]/rhob;
vn = u*geoFace(f,1) + v*geoFace(f,2) + w*geoFace(f,3);
// acoustic speed
auto p = mat_blk[0].compute< EOS::pressure >( rhob, u, v, w,
ugp[multispecies::energyIdx(nspec, 0)] );
a = mat_blk[0].compute< EOS::soundspeed >( rhob, p );
dSV_l = geoFace(f,0) * (std::fabs(vn) + a);
// right element
if (er > -1) {
std::size_t eR = static_cast< std::size_t >( er );
// Compute the basis function for the right element
std::vector< tk::real > B_r(rdof, 0.0);
B_r[0] = 1.0;
// get conserved quantities
ugp = eval_state( ncomp, rdof, ndof, eR, U, B_r);
// mixture density
rhob = 0.0;
for (std::size_t k=0; k<nspec; ++k)
rhob += ugp[multispecies::densityIdx(nspec, k)];
// advection velocity
u = ugp[multispecies::momentumIdx(nspec, 0)]/rhob;
v = ugp[multispecies::momentumIdx(nspec, 1)]/rhob;
w = ugp[multispecies::momentumIdx(nspec, 2)]/rhob;
vn = u*geoFace(f,1) + v*geoFace(f,2) + w*geoFace(f,3);
// acoustic speed
p = mat_blk[0].compute< EOS::pressure >( rhob, u, v, w,
ugp[multispecies::energyIdx(nspec, 0)] );
a = mat_blk[0].compute< EOS::soundspeed >( rhob, p );
dSV_r = geoFace(f,0) * (std::fabs(vn) + a);
delt[eR] += std::max( dSV_l, dSV_r );
} else {
dSV_r = dSV_l;
}
delt[el] += std::max( dSV_l, dSV_r );
}
tk::real mindt = std::numeric_limits< tk::real >::max();
// compute allowable dt
for (std::size_t e=0; e<nelem; ++e)
{
mindt = std::min( mindt, geoElem(e,0)/delt[e] );
}
return mindt;
}
} //inciter::
|