www.mooseframework.org
CNSFVSlopeLimitingOneD.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
9 
10 template <>
11 InputParameters
13 {
14  InputParameters params = validParams<SlopeLimitingBase>();
15 
16  params.addClassDescription("A use object that serves as base class for slope limiting to get the "
17  "limited slopes of cell average variables in 1-D.");
18 
19  params.addRequiredCoupledVar("rho", "Density at P0 (constant monomial)");
20 
21  params.addRequiredCoupledVar("rhou", "X-momentum at P0 (constant monomial)");
22 
23  params.addRequiredCoupledVar("rhoe", "Total energy at P0 (constant monomial)");
24 
25  params.addRequiredParam<UserObjectName>("fluid_properties",
26  "Name for fluid properties user object");
27 
28  MooseEnum scheme("none minmod mc superbee", "none");
29 
30  params.addParam<MooseEnum>("scheme", scheme, "Slope limiting scheme");
31 
32  return params;
33 }
34 
35 CNSFVSlopeLimitingOneD::CNSFVSlopeLimitingOneD(const InputParameters & parameters)
36  : SlopeLimitingBase(parameters),
37  _rho(getVar("rho", 0)),
38  _rhou(getVar("rhou", 0)),
39  _rhoe(getVar("rhoe", 0)),
40  _fp(getUserObject<SinglePhaseFluidProperties>("fluid_properties")),
41  _scheme(getParam<MooseEnum>("scheme"))
42 {
43 }
44 
45 std::vector<RealGradient>
47 {
48  const Elem * elem = _current_elem;
49 
51  unsigned int iv;
52 
54  unsigned int is;
55 
57  unsigned int in;
58 
60  unsigned int nside = elem->n_sides();
61 
63  unsigned int nsten = nside + 1;
64 
66  unsigned int nvars = 5;
67 
69  std::vector<RealGradient> ugrad(nvars, RealGradient(0.0, 0.0, 0.0));
70 
72  std::vector<Real> xc(nsten, 0.);
73 
75  xc[0] = elem->centroid()(0);
76 
78  std::vector<std::vector<Real>> ucell(nsten, std::vector<Real>(nvars, 0.0));
79 
82 
87 
92 
97 
98  std::vector<std::vector<Real>> sigma(nsten, std::vector<Real>(nvars, 0.0));
99 
104 
105  Real rho = _rho->getElementalValue(elem);
106  Real rhou = _rhou->getElementalValue(elem);
107  Real rhoe = _rhoe->getElementalValue(elem);
108 
113 
114  ucell[0][1] = rhou / rho;
115 
116  Real v = 1. / rho;
117  Real e = rhoe / rho - 0.5 * ucell[0][1] * ucell[0][1];
118 
119  ucell[0][0] = _fp.pressure(v, e);
120  ucell[0][4] = _fp.temperature(v, e);
121 
123 
124  unsigned int bflag = 0;
125 
127 
128  for (is = 0; is < nside; is++)
129  {
130  in = is + 1;
131 
133  if (elem->neighbor(is) != NULL)
134  {
135  const Elem * neig = elem->neighbor(is);
136 
137  xc[in] = neig->centroid()(0);
138 
143 
144  rho = _rho->getElementalValue(neig);
145  rhou = _rhou->getElementalValue(neig);
146  rhoe = _rhoe->getElementalValue(neig);
147 
152 
153  ucell[in][1] = rhou / rho;
154 
155  v = 1. / rho;
156  e = rhoe / rho - 0.5 * ucell[in][1] * ucell[in][1];
157 
158  ucell[in][0] = _fp.pressure(v, e);
159  ucell[in][4] = _fp.temperature(v, e);
160 
162 
163  for (iv = 0; iv < nvars; iv++)
164  sigma[in][iv] = (ucell[0][iv] - ucell[in][iv]) / (xc[0] - xc[in]);
165  }
171  else
172  {
173  bflag = in;
174  }
175  }
176 
179 
180  switch (_scheme)
181  {
183  case 0:
184  break;
185 
187  case 1:
188 
189  if (bflag == 0)
190  {
191  for (iv = 0; iv < nvars; iv++)
192  {
193  if ((sigma[1][iv] * sigma[2][iv]) > 0.)
194  {
195  if (std::abs(sigma[1][iv]) < std::abs(sigma[2][iv]))
196  ugrad[iv](0) = sigma[1][iv];
197  else
198  ugrad[iv](0) = sigma[2][iv];
199  }
200  }
201  }
202  break;
203 
205  case 2:
206 
207  if (bflag == 0)
208  {
209  for (iv = 0; iv < nvars; iv++)
210  sigma[0][iv] = (ucell[1][iv] - ucell[2][iv]) / (xc[1] - xc[2]);
211 
212  for (iv = 0; iv < nvars; iv++)
213  {
214  if (sigma[0][iv] > 0. && sigma[1][iv] > 0. && sigma[2][iv] > 0.)
215  ugrad[iv](0) = std::min(sigma[0][iv], 2. * std::min(sigma[1][iv], sigma[2][iv]));
216  else if (sigma[0][iv] < 0. && sigma[1][iv] < 0. && sigma[2][iv] < 0.)
217  ugrad[iv](0) = std::max(sigma[0][iv], 2. * std::max(sigma[1][iv], sigma[2][iv]));
218  }
219  }
220  break;
221 
223  case 3:
224 
225  if (bflag == 0)
226  {
227  for (iv = 0; iv < nvars; iv++)
228  {
229  Real sigma1 = 0.;
230  Real sigma2 = 0.;
231 
233  if (sigma[2][iv] > 0. && sigma[1][iv] > 0.)
234  sigma1 = std::min(sigma[2][iv], 2. * sigma[1][iv]);
235  else if (sigma[2][iv] < 0. && sigma[1][iv] < 0.)
236  sigma1 = std::max(sigma[2][iv], 2. * sigma[1][iv]);
237 
239  if (sigma[2][iv] > 0. && sigma[1][iv] > 0.)
240  sigma2 = std::min(2. * sigma[2][iv], sigma[1][iv]);
241  else if (sigma[2][iv] < 0. && sigma[1][iv] < 0.)
242  sigma2 = std::max(2. * sigma[2][iv], sigma[1][iv]);
243 
245  if (sigma1 > 0. && sigma2 > 0.)
246  ugrad[iv](0) = std::max(sigma1, sigma2);
247  else if (sigma1 < 0. && sigma2 < 0.)
248  ugrad[iv](0) = std::min(sigma1, sigma2);
249  }
250  }
251  break;
252 
253  default:
254  mooseError("Unknown 1D TVD slope limiter scheme");
255  break;
256  }
257 
258  return ugrad;
259 }
const SinglePhaseFluidProperties & _fp
fluid properties user object
MooseEnum _scheme
One-D slope limiting scheme.
virtual Real temperature(Real v, Real u) const =0
Temperature as a function of specific internal energy and specific volume.
MooseVariable * _rho
the input density
virtual std::vector< RealGradient > limitElementSlope() const
compute the limited slope of the cell
MooseVariable * _rhou
the input x-momentum
CNSFVSlopeLimitingOneD(const InputParameters &parameters)
InputParameters validParams< CNSFVSlopeLimitingOneD >()
Common class for single phase fluid properties.
Base class for slope limiting to limit the slopes of cell average variables.
InputParameters validParams< SlopeLimitingBase >()
MooseVariable * _rhoe
the input total energy
virtual Real pressure(Real v, Real u) const =0
Pressure as a function of specific internal energy and specific volume.