www.mooseframework.org
Public Member Functions | Protected Attributes | List of all members
CNSFVHLLCInflowOutflowBoundaryFlux Class Reference

A user object that computes inflow/outflow boundary flux using the HLLC approximate Riemann solver. More...

#include <CNSFVHLLCInflowOutflowBoundaryFlux.h>

Inheritance diagram for CNSFVHLLCInflowOutflowBoundaryFlux:
[legend]

Public Member Functions

 CNSFVHLLCInflowOutflowBoundaryFlux (const InputParameters &parameters)
 
virtual ~CNSFVHLLCInflowOutflowBoundaryFlux ()
 
virtual void calcFlux (unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, std::vector< Real > &flux) const
 Solve the Riemann problem on the boundary face. More...
 
virtual void calcJacobian (unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, DenseMatrix< Real > &jac1) const
 Compute the Jacobian matrix on the boundary face. More...
 
virtual void execute ()
 
virtual void initialize ()
 
virtual void finalize ()
 
virtual const std::vector< Real > & getFlux (unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, THREAD_ID tid) const
 Get the boundary flux vector. More...
 
virtual const DenseMatrix< Real > & getJacobian (unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, THREAD_ID tid) const
 Get the boundary Jacobian matrix. More...
 

Protected Attributes

const SinglePhaseFluidProperties_fp
 
Real _inf_rho
 
Real _inf_uadv
 
Real _inf_vadv
 
Real _inf_wadv
 
Real _inf_pres
 
unsigned int _cached_side_id
 
dof_id_type _cached_elem_id
 
std::vector< std::vector< Real > > _flux
 Threaded storage for fluxes. More...
 
std::vector< DenseMatrix< Real > > _jac1
 Threaded storage for jacobians. More...
 

Detailed Description

A user object that computes inflow/outflow boundary flux using the HLLC approximate Riemann solver.

Definition at line 24 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Constructor & Destructor Documentation

CNSFVHLLCInflowOutflowBoundaryFlux::CNSFVHLLCInflowOutflowBoundaryFlux ( const InputParameters &  parameters)

Definition at line 36 of file CNSFVHLLCInflowOutflowBoundaryFlux.C.

38  : BoundaryFluxBase(parameters),
39  _fp(getUserObject<SinglePhaseFluidProperties>("fluid_properties")),
40  _inf_rho(getParam<Real>("infinity_density")),
41  _inf_uadv(getParam<Real>("infinity_x_velocity")),
42  _inf_vadv(getParam<Real>("infinity_y_velocity")),
43  _inf_wadv(getParam<Real>("infinity_z_velocity")),
44  _inf_pres(getParam<Real>("infinity_pressure"))
45 {
46 }
BoundaryFluxBase(const InputParameters &parameters)
CNSFVHLLCInflowOutflowBoundaryFlux::~CNSFVHLLCInflowOutflowBoundaryFlux ( )
virtual

Definition at line 48 of file CNSFVHLLCInflowOutflowBoundaryFlux.C.

48 {}

Member Function Documentation

void CNSFVHLLCInflowOutflowBoundaryFlux::calcFlux ( unsigned int  iside,
dof_id_type  ielem,
const std::vector< Real > &  uvec1,
const RealVectorValue &  dwave,
std::vector< Real > &  flux 
) const
virtual

Solve the Riemann problem on the boundary face.

Parameters
[in]isidelocal index of current side
[in]ielemglobal index of the current element
[in]uvec1vector of variables on the host side
[in]dwavevector of unit normal
[out]fluxflux vector for conservation equations

pass the inputs to local

assign the proper size for the flux vector

compute the local Mach number with state variables on the left

subsonic inflow case – fix pressure (but also OK to extrapolate pressure) – fix temperature – fix velocity – fix mass fractions

state variables required at ghost cell <– inifinity state

get the enthalpy

get the so-called Roe averaged variables

averaged speed of sound

compute the eigenvalues at the left, right, and Roe states

get the S_L and S_R defined by Eq. (13).

compute the S_M defined by Eq. (12).

compute the Omega_l, Omega_r, and p^*

compute the U_l^, U_r^

compute the fluxes according to the wave speed

subsonic outflow case – fix pressure – extrapolate temperature – extrapolate velocities – extrapolate mass fractions

supersonic inflow case – fix pressure – fix temperature – fix velocity – fix mass fractions

supersonic outflow case – extrapolate pressure – extrapolate temperature – extrapolate velocity – extrapolate mass fractions

Implements BoundaryFluxBase.

Definition at line 51 of file CNSFVHLLCInflowOutflowBoundaryFlux.C.

56 {
57  Real eps = 1e-6;
58 
60 
61  Real rho1 = uvec1[0];
62  Real rhou1 = uvec1[1];
63  Real rhov1 = uvec1[2];
64  Real rhow1 = uvec1[3];
65  Real rhoe1 = uvec1[4];
66 
67  Real nx = dwave(0);
68  Real ny = dwave(1);
69  Real nz = dwave(2);
70 
72 
73  flux.resize(5);
74 
76 
77  Real rhom1 = 1. / rho1;
78  Real uadv1 = rhou1 * rhom1;
79  Real vadv1 = rhov1 * rhom1;
80  Real wadv1 = rhow1 * rhom1;
81  Real vdov1 = uadv1 * uadv1 + vadv1 * vadv1 + wadv1 * wadv1;
82  Real eint1 = rhoe1 * rhom1 - 0.5 * vdov1;
83  Real pres1 = _fp.pressure(rhom1, eint1);
84  Real csou1 = _fp.c(rhom1, eint1);
85  Real gamma = _fp.gamma(rhom1, eint1);
86  Real mach1 = std::sqrt(vdov1) / csou1;
87 
88  if (mach1 > -1. && mach1 < 0.)
89  {
95 
97 
98  Real rho2 = _inf_rho;
99  Real uadv2 = _inf_uadv;
100  Real vadv2 = _inf_vadv;
101  Real wadv2 = _inf_wadv;
102  Real pres2 = _inf_pres;
103  Real rhom2 = 1. / rho2;
104  Real rhou2 = rho2 * uadv2;
105  Real rhov2 = rho2 * vadv2;
106  Real rhow2 = rho2 * wadv2;
107  Real vdov2 = uadv2 * uadv2 + vadv2 * vadv2 + wadv2 * wadv2;
108  Real eint2 = _fp.e(pres2, rho2);
109  Real rhoe2 = rho2 * eint2 + 0.5 * rho2 * vdov2;
110  Real csou2 = _fp.c(rhom2, eint2);
111 
113 
114  Real enth1 = (rhoe1 + pres1) * rhom1;
115  Real enth2 = (rhoe2 + pres2) * rhom2;
116 
118 
119  Real rhsca = std::sqrt(rho2 / rho1);
120  Real rmden = 1. / (rhsca + 1.);
121 
122  // Real rhoav = rhsca * rho1;
123  Real uaver = (rhsca * uadv2 + uadv1) * rmden;
124  Real vaver = (rhsca * vadv2 + vadv1) * rmden;
125  Real waver = (rhsca * wadv2 + wadv1) * rmden;
126  Real entav = (rhsca * enth2 + enth1) * rmden;
127 
129 
130  Real qave5 = 0.5 * (uaver * uaver + vaver * vaver + waver * waver);
131  Real cssa2 = std::max(eps, (gamma - 1.) * (entav - qave5));
132  Real cssav = std::sqrt(cssa2);
133 
135 
136  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
137  Real vdon2 = uadv2 * nx + vadv2 * ny + wadv2 * nz;
138  Real vnave = uaver * nx + vaver * ny + waver * nz;
139 
141 
142  Real s1 = std::min(vdon1 - csou1, vnave - cssav);
143  Real s2 = std::max(vdon2 + csou2, vnave + cssav);
144 
145  Real dsv1 = s1 - vdon1;
146  Real dsv2 = s2 - vdon2;
147 
149 
150  Real sm =
151  (rho2 * vdon2 * dsv2 - rho1 * vdon1 * dsv1 + pres1 - pres2) / (rho2 * dsv2 - rho1 * dsv1);
152 
154 
155  Real omeg1 = 1. / (s1 - sm);
156  Real omeg2 = 1. / (s2 - sm);
157  Real prsta = rho1 * dsv1 * (sm - vdon1) + pres1;
158 
159  Real prst1 = prsta - pres1;
160  Real prst2 = prsta - pres2;
161 
163 
164  Real rhol = omeg1 * dsv1 * rho1;
165  Real rhoul = omeg1 * (dsv1 * rhou1 + prst1 * nx);
166  Real rhovl = omeg1 * (dsv1 * rhov1 + prst1 * ny);
167  Real rhowl = omeg1 * (dsv1 * rhow1 + prst1 * nz);
168  Real rhoel = omeg1 * (dsv1 * rhoe1 - pres1 * vdon1 + prsta * sm);
169 
170  Real rhor = omeg2 * dsv2 * rho2;
171  Real rhour = omeg2 * (dsv2 * rhou2 + prst2 * nx);
172  Real rhovr = omeg2 * (dsv2 * rhov2 + prst2 * ny);
173  Real rhowr = omeg2 * (dsv2 * rhow2 + prst2 * nz);
174  Real rhoer = omeg2 * (dsv2 * rhoe2 - pres2 * vdon2 + prsta * sm);
175 
177 
178  if (s1 > 0.)
179  {
180  flux[0] = vdon1 * rho1;
181  flux[1] = vdon1 * rhou1 + pres1 * nx;
182  flux[2] = vdon1 * rhov1 + pres1 * ny;
183  flux[3] = vdon1 * rhow1 + pres1 * nz;
184  flux[4] = vdon1 * (rhoe1 + pres1);
185  }
186  else if (s1 <= 0. && sm > 0.)
187  {
188  flux[0] = sm * rhol;
189  flux[1] = sm * rhoul + prsta * nx;
190  flux[2] = sm * rhovl + prsta * ny;
191  flux[3] = sm * rhowl + prsta * nz;
192  flux[4] = sm * (rhoel + prsta);
193  }
194  else if (sm <= 0. && s2 >= 0.)
195  {
196  flux[0] = sm * rhor;
197  flux[1] = sm * rhour + prsta * nx;
198  flux[2] = sm * rhovr + prsta * ny;
199  flux[3] = sm * rhowr + prsta * nz;
200  flux[4] = sm * (rhoer + prsta);
201  }
202  else if (s2 < 0.)
203  {
204  flux[0] = vdon2 * rho2;
205  flux[1] = vdon2 * rhou2 + pres2 * nx;
206  flux[2] = vdon2 * rhov2 + pres2 * ny;
207  flux[3] = vdon2 * rhow2 + pres2 * nz;
208  flux[4] = vdon2 * (rhoe2 + pres2);
209  }
210  else
211  {
212  mooseError("Weird wave speed occured in",
213  "\n",
214  name(),
215  ": ",
216  __FUNCTION__,
217  "\n",
218  "iside = ",
219  iside,
220  "\n",
221  "ielem = ",
222  ielem,
223  "\n",
224  "rho1 = ",
225  rho1,
226  "\n",
227  "rhou1 = ",
228  rhou1,
229  "\n",
230  "rhov1 = ",
231  rhov1,
232  "\n",
233  "rhow1 = ",
234  rhow1,
235  "\n",
236  "rhoe1 = ",
237  rhoe1,
238  "\n",
239  "pres1 = ",
240  pres1,
241  "\n",
242  "enth1 = ",
243  enth1,
244  "\n",
245  "csou1 = ",
246  csou1,
247  "\n",
248  "rho2 = ",
249  rho2,
250  "\n",
251  "rhou2 = ",
252  rhou2,
253  "\n",
254  "rhov2 = ",
255  rhov2,
256  "\n",
257  "rhoe2 = ",
258  rhoe2,
259  "\n",
260  "pres2 = ",
261  pres2,
262  "\n",
263  "enth2 = ",
264  enth2,
265  "\n",
266  "csou2 = ",
267  csou2,
268  "\n",
269  "vdon1 = ",
270  vdon1,
271  "\n",
272  "vdon2 = ",
273  vdon2,
274  "\n",
275  "vnave = ",
276  vnave,
277  "\n",
278  "cssav = ",
279  cssav,
280  "\n",
281  "s1 = ",
282  s1,
283  "\n",
284  "s2 = ",
285  s2,
286  "\n",
287  "sm = ",
288  sm,
289  "\n",
290  "omeg1 = ",
291  omeg1,
292  "\n",
293  "omeg2 = ",
294  omeg2,
295  "\n",
296  "prsta = ",
297  prsta,
298  "\n",
299  "Please check before continuing!\n");
300  }
301  }
302  else if (mach1 >= 0. && mach1 < 1.)
303  {
309 
310  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
311 
312  flux[0] = vdon1 * rho1;
313  flux[1] = vdon1 * rhou1 + _inf_pres * nx;
314  flux[2] = vdon1 * rhov1 + _inf_pres * ny;
315  flux[3] = vdon1 * rhow1 + _inf_pres * nz;
316  flux[4] = vdon1 * (rhoe1 + _inf_pres);
317  }
318  else if (mach1 <= -1.)
319  {
325 
326  Real rho2 = _inf_rho;
327  Real uadv2 = _inf_uadv;
328  Real vadv2 = _inf_vadv;
329  Real wadv2 = _inf_wadv;
330  Real pres2 = _inf_pres;
331  Real rhou2 = rho2 * uadv2;
332  Real rhov2 = rho2 * vadv2;
333  Real rhow2 = rho2 * wadv2;
334  Real vdov2 = uadv2 * uadv2 + vadv2 * vadv2 + wadv2 * wadv2;
335  Real eint2 = _fp.e(pres2, rho2);
336  Real rhoe2 = rho2 * eint2 + 0.5 * rho2 * vdov2;
337  Real vdon2 = uadv2 * nx + vadv2 * ny + wadv2 * nz;
338 
339  flux[0] = vdon2 * rho2;
340  flux[1] = vdon2 * rhou2 + pres2 * nx;
341  flux[2] = vdon2 * rhov2 + pres2 * ny;
342  flux[3] = vdon2 * rhow2 + pres2 * nz;
343  flux[4] = vdon2 * (rhoe2 + pres2);
344  }
345  else if (mach1 >= 1.)
346  {
352 
353  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
354 
355  flux[0] = vdon1 * rho1;
356  flux[1] = vdon1 * rhou1 + pres1 * nx;
357  flux[2] = vdon1 * rhov1 + pres1 * ny;
358  flux[3] = vdon1 * rhow1 + pres1 * nz;
359  flux[4] = vdon1 * (rhoe1 + pres1);
360  }
361  else
362  mooseError("Something is wrong in ",
363  name(),
364  ": ",
365  __FUNCTION__,
366  "\n",
367  "ielem = ",
368  ielem,
369  "\n",
370  "iside = ",
371  iside,
372  "\n",
373  "mach1 = ",
374  mach1,
375  "\n");
376 }
virtual Real c(Real v, Real u) const =0
Sound speed.
virtual Real gamma(Real v, Real u) const
Compute the ratio of specific heats.
virtual Real pressure(Real v, Real u) const =0
Pressure as a function of specific internal energy and specific volume.
virtual Real e(Real pressure, Real rho) const =0
Computes internal energy from pressure and density.
void CNSFVHLLCInflowOutflowBoundaryFlux::calcJacobian ( unsigned int  iside,
dof_id_type  ielem,
const std::vector< Real > &  uvec1,
const RealVectorValue &  dwave,
DenseMatrix< Real > &  jac1 
) const
virtual

Compute the Jacobian matrix on the boundary face.

Parameters
[in]isidelocal index of current side
[in]ielemglobal index of the current element
[in]uvec1vector of variables on the host side
[in]dwavevector of unit normal
[out]jac1Jacobian matrix contribution

pass the inputs to local

assign the proper size for the Jacobian matrix

compute the local Mach number with state variables on the left

subsonic inflow case – fix pressure (but also OK to extrapolate pressure) – fix temperature – fix velocity – fix mass fractions

state variables required at ghost cell <– inifinity state

get the enthalpy

get the so-called Roe averaged variables

averaged speed of sound

compute the eigenvalues at the left, right, and Roe states

get the S_L and S_R defined by Eq. (13).

compute the S_M defined by Eq. (12).

compute the flux Jacobians according to the wave speed

get the Jacobian matrix on the left only

compute the Omega_l and p^* by Eq. (10). and (11).

compute the U_l^ by Eq. (7).

compute rho^

compute a(S_M)/a(U_l) by Eq. (42).

compute a(p^*)/a(U_l) by Eq. (44).

compute a(rho_l^)/a(U_l) by Eq. (46).

compute a(rhou_l^)/a(U_l) by Eq. (47).

compute a(rhov_l^)/a(U_l) by Eq. (49).

compute a(rhow_l^)/a(U_l) by Eq. (51).

compute a(rhoe_l^)/a(U_l) by Eq. (53). note that enth1 = (rhoe1 + pres1)*rhom1 and rhoepl = rhoel + prsta

compute the HLLC Jacobians a(F_l^)/a(U_l) by Eq. (40).

compute the Omega_r, and p^* by Eq. (10). and (11).

compute the U_r^ by Eq. (8).

compute rho^

compute a(S_M)/a(U_l) by Eq. (42).

compute a(p^*)/a(U_l) by Eq. (44).

compute a(rho_r^)/a(U_l) by Eq. (46a*).

compute a(rhou_r^)/a(U_l) by Eq. (48*).

compute a(rhov_r^)/a(U_l) by Eq. (50*).

compute a(rhow_r^)/a(U_l) by Eq. (52*).

compute a(rhoe_r^)/a(U_l) by Eq. (54*).

compute the HLLC Jacobians a(F_r^)/a(U_l) by Eq. (41*).

do nothing

compute the Omega_l and p^* by Eq. (10). and (11).

subsonic outflow case – fix pressure – extrapolate temperature – extrapolate velocities – extrapolate mass fractions

supersonic inflow case – fix pressure – fix temperature – fix velocity – fix mass fractions

do nothing

supersonic outflow case – extrapolate pressure – extrapolate temperature – extrapolate velocity – extrapolate mass fractions

Implements BoundaryFluxBase.

Definition at line 379 of file CNSFVHLLCInflowOutflowBoundaryFlux.C.

384 {
385  Real eps = 1e-6;
386 
388 
389  Real rho1 = uvec1[0];
390  Real rhou1 = uvec1[1];
391  Real rhov1 = uvec1[2];
392  Real rhow1 = uvec1[3];
393  Real rhoe1 = uvec1[4];
394 
395  Real nx = dwave(0);
396  Real ny = dwave(1);
397  Real nz = dwave(2);
398 
400 
401  jac1.resize(5, 5);
402 
404 
405  Real rhom1 = 1. / rho1;
406  Real uadv1 = rhou1 * rhom1;
407  Real vadv1 = rhov1 * rhom1;
408  Real wadv1 = rhow1 * rhom1;
409  Real vdov1 = uadv1 * uadv1 + vadv1 * vadv1 + wadv1 * wadv1;
410  Real eint1 = rhoe1 * rhom1 - 0.5 * vdov1;
411  Real pres1 = _fp.pressure(rhom1, eint1);
412  Real csou1 = _fp.c(rhom1, eint1);
413  Real gamma = _fp.gamma(rhom1, eint1);
414  Real gamm1 = gamma - 1.;
415  Real gamm2 = 2. - gamma;
416  Real rq051 = 0.5 * gamm1 * vdov1;
417  Real mach1 = std::sqrt(vdov1) / csou1;
418 
419  if (mach1 > -1. && mach1 < 0.)
420  {
426 
428 
429  Real rho2 = _inf_rho;
430  Real uadv2 = _inf_uadv;
431  Real vadv2 = _inf_vadv;
432  Real wadv2 = _inf_wadv;
433  Real pres2 = _inf_pres;
434  Real rhom2 = 1. / rho2;
435  Real rhou2 = rho2 * uadv2;
436  Real rhov2 = rho2 * vadv2;
437  Real rhow2 = rho2 * wadv2;
438  Real vdov2 = uadv2 * uadv2 + vadv2 * vadv2 + wadv2 * wadv2;
439  Real eint2 = _fp.e(pres2, rho2);
440  Real rhoe2 = rho2 * eint2 + 0.5 * rho2 * vdov2;
441  Real csou2 = _fp.c(rhom2, eint2);
442 
444 
445  Real enth1 = (rhoe1 + pres1) * rhom1;
446  Real enth2 = (rhoe2 + pres2) * rhom2;
447 
449 
450  Real rhsca = std::sqrt(rho2 / rho1);
451  Real rmden = 1. / (rhsca + 1.);
452 
453  // Real rhoav = rhsca * rho1;
454  Real uaver = (rhsca * uadv2 + uadv1) * rmden;
455  Real vaver = (rhsca * vadv2 + vadv1) * rmden;
456  Real waver = (rhsca * wadv2 + wadv1) * rmden;
457  Real entav = (rhsca * enth2 + enth1) * rmden;
458 
460 
461  Real qave5 = 0.5 * (uaver * uaver + vaver * vaver + waver * waver);
462  Real cssa2 = std::max(eps, (gamma - 1.) * (entav - qave5));
463  Real cssav = std::sqrt(cssa2);
464 
466 
467  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
468  Real vdon2 = uadv2 * nx + vadv2 * ny + wadv2 * nz;
469  Real vnave = uaver * nx + vaver * ny + waver * nz;
470 
472 
473  Real s1 = std::min(vdon1 - csou1, vnave - cssav);
474  Real s2 = std::max(vdon2 + csou2, vnave + cssav);
475 
476  Real dsv1 = s1 - vdon1;
477  Real dsv2 = s2 - vdon2;
478 
480 
481  Real sm =
482  (rho2 * vdon2 * dsv2 - rho1 * vdon1 * dsv1 + pres1 - pres2) / (rho2 * dsv2 - rho1 * dsv1);
483 
485 
486  if (s1 > 0.)
487  {
489 
490  jac1(0, 0) = 0.;
491  jac1(0, 1) = nx;
492  jac1(0, 2) = ny;
493  jac1(0, 3) = nz;
494  jac1(0, 4) = 0.;
495 
496  jac1(1, 0) = rq051 * nx - uadv1 * vdon1;
497  jac1(1, 1) = gamm2 * nx * uadv1 + vdon1;
498  jac1(1, 2) = ny * uadv1 - vadv1 * gamm1 * nx;
499  jac1(1, 3) = nz * uadv1 - wadv1 * gamm1 * nx;
500  jac1(1, 4) = gamm1 * nx;
501 
502  jac1(2, 0) = rq051 * ny - vadv1 * vdon1;
503  jac1(2, 1) = nx * vadv1 - uadv1 * gamm1 * ny;
504  jac1(2, 2) = gamm2 * ny * vadv1 + vdon1;
505  jac1(2, 3) = nz * vadv1 - wadv1 * gamm1 * ny;
506  jac1(2, 4) = gamm1 * ny;
507 
508  jac1(3, 0) = rq051 * nz - wadv1 * vdon1;
509  jac1(3, 1) = nx * wadv1 - uadv1 * gamm1 * nz;
510  jac1(3, 2) = ny * wadv1 - vadv1 * gamm1 * nz;
511  jac1(3, 3) = gamm2 * nz * wadv1 + vdon1;
512  jac1(3, 4) = gamm1 * nz;
513 
514  jac1(4, 0) = (rq051 - enth1) * vdon1;
515  jac1(4, 1) = nx * enth1 - gamm1 * uadv1 * vdon1;
516  jac1(4, 2) = ny * enth1 - gamm1 * vadv1 * vdon1;
517  jac1(4, 3) = nz * enth1 - gamm1 * wadv1 * vdon1;
518  jac1(4, 4) = gamma * vdon1;
519  }
520  else if (s1 <= 0. && sm > 0.)
521  {
523 
524  Real omeg1 = 1. / (s1 - sm);
525  Real prsta = rho1 * dsv1 * (sm - vdon1) + pres1;
526 
527  Real prst1 = prsta - pres1;
528 
530 
531  Real rhol = omeg1 * dsv1 * rho1;
532  Real rhoul = omeg1 * (dsv1 * rhou1 + prst1 * nx);
533  Real rhovl = omeg1 * (dsv1 * rhov1 + prst1 * ny);
534  Real rhowl = omeg1 * (dsv1 * rhow1 + prst1 * nz);
535  Real rhoel = omeg1 * (dsv1 * rhoe1 - pres1 * vdon1 + prsta * sm);
536 
537  Real rhoepl = rhoel + prsta;
538 
540 
541  Real rhotd = std::max(eps, rho2 * (s2 - vdon2) - rho1 * (s1 - vdon1));
542  Real rhotm = 1. / rhotd;
543 
545 
546  Real sm_rho1 = rhotm * (-vdon1 * vdon1 + rq051 + sm * s1);
547  Real sm_rhou1 = rhotm * (nx * (2. * vdon1 - s1 - sm) - gamm1 * uadv1);
548  Real sm_rhov1 = rhotm * (ny * (2. * vdon1 - s1 - sm) - gamm1 * vadv1);
549  Real sm_rhow1 = rhotm * (nz * (2. * vdon1 - s1 - sm) - gamm1 * wadv1);
550  Real sm_rhoe1 = rhotm * gamm1;
551 
553 
554  Real ps_rho1 = rho2 * dsv2 * sm_rho1;
555  Real ps_rhou1 = rho2 * dsv2 * sm_rhou1;
556  Real ps_rhov1 = rho2 * dsv2 * sm_rhov1;
557  Real ps_rhow1 = rho2 * dsv2 * sm_rhow1;
558  Real ps_rhoe1 = rho2 * dsv2 * sm_rhoe1;
559 
561 
562  Real rhol_rho1 = omeg1 * (s1 + rhol * sm_rho1);
563  Real rhol_rhou1 = omeg1 * (-nx + rhol * sm_rhou1);
564  Real rhol_rhov1 = omeg1 * (-ny + rhol * sm_rhov1);
565  Real rhol_rhow1 = omeg1 * (-nz + rhol * sm_rhow1);
566  Real rhol_rhoe1 = omeg1 * (rhol * sm_rhoe1);
567 
569 
570  Real rhoul_rho1 = omeg1 * (uadv1 * vdon1 - nx * rq051 + nx * ps_rho1 + rhoul * sm_rho1);
571  Real rhoul_rhou1 = omeg1 * (dsv1 - gamm2 * nx * uadv1 + nx * ps_rhou1 + rhoul * sm_rhou1);
572  Real rhoul_rhov1 =
573  omeg1 * (-uadv1 * ny + gamm1 * nx * vadv1 + nx * ps_rhov1 + rhoul * sm_rhov1);
574  Real rhoul_rhow1 =
575  omeg1 * (-uadv1 * nz + gamm1 * nx * wadv1 + nx * ps_rhow1 + rhoul * sm_rhow1);
576  Real rhoul_rhoe1 = omeg1 * (-gamm1 * nx + nx * ps_rhoe1 + rhoul * sm_rhoe1);
577 
579 
580  Real rhovl_rho1 = omeg1 * (vadv1 * vdon1 - ny * rq051 + ny * ps_rho1 + rhovl * sm_rho1);
581  Real rhovl_rhou1 =
582  omeg1 * (-vadv1 * nx + gamm1 * ny * uadv1 + ny * ps_rhou1 + rhovl * sm_rhou1);
583  Real rhovl_rhov1 = omeg1 * (dsv1 - gamm2 * ny * vadv1 + ny * ps_rhov1 + rhovl * sm_rhov1);
584  Real rhovl_rhow1 =
585  omeg1 * (-vadv1 * nz + gamm1 * ny * wadv1 + ny * ps_rhow1 + rhovl * sm_rhow1);
586  Real rhovl_rhoe1 = omeg1 * (-gamm1 * ny + ny * ps_rhoe1 + rhovl * sm_rhoe1);
587 
589 
590  Real rhowl_rho1 = omeg1 * (wadv1 * vdon1 - nz * rq051 + nz * ps_rho1 + rhowl * sm_rho1);
591  Real rhowl_rhou1 =
592  omeg1 * (-wadv1 * nx + gamm1 * nz * uadv1 + nz * ps_rhou1 + rhowl * sm_rhou1);
593  Real rhowl_rhov1 =
594  omeg1 * (-wadv1 * ny + gamm1 * nz * vadv1 + nz * ps_rhov1 + rhowl * sm_rhov1);
595  Real rhowl_rhow1 = omeg1 * (dsv1 - gamm2 * nz * wadv1 + nz * ps_rhow1 + rhowl * sm_rhow1);
596  Real rhowl_rhoe1 = omeg1 * (-gamm1 * nz + nz * ps_rhoe1 + rhowl * sm_rhoe1);
597 
601 
602  Real rhoel_rho1 = omeg1 * (vdon1 * enth1 - vdon1 * rq051 + sm * ps_rho1 + rhoepl * sm_rho1);
603  Real rhoel_rhou1 =
604  omeg1 * (-nx * enth1 + gamm1 * vdon1 * uadv1 + sm * ps_rhou1 + rhoepl * sm_rhou1);
605  Real rhoel_rhov1 =
606  omeg1 * (-ny * enth1 + gamm1 * vdon1 * vadv1 + sm * ps_rhov1 + rhoepl * sm_rhov1);
607  Real rhoel_rhow1 =
608  omeg1 * (-nz * enth1 + gamm1 * vdon1 * wadv1 + sm * ps_rhow1 + rhoepl * sm_rhow1);
609  Real rhoel_rhoe1 = omeg1 * (s1 - vdon1 * gamma + sm * ps_rhoe1 + rhoepl * sm_rhoe1);
610 
612 
613  jac1(0, 0) = sm * rhol_rho1 + rhol * sm_rho1;
614  jac1(0, 1) = sm * rhol_rhou1 + rhol * sm_rhou1;
615  jac1(0, 2) = sm * rhol_rhov1 + rhol * sm_rhov1;
616  jac1(0, 3) = sm * rhol_rhow1 + rhol * sm_rhow1;
617  jac1(0, 4) = sm * rhol_rhoe1 + rhol * sm_rhoe1;
618 
619  jac1(1, 0) = sm * rhoul_rho1 + rhoul * sm_rho1 + nx * ps_rho1;
620  jac1(1, 1) = sm * rhoul_rhou1 + rhoul * sm_rhou1 + nx * ps_rhou1;
621  jac1(1, 2) = sm * rhoul_rhov1 + rhoul * sm_rhov1 + nx * ps_rhov1;
622  jac1(1, 3) = sm * rhoul_rhow1 + rhoul * sm_rhow1 + nx * ps_rhow1;
623  jac1(1, 4) = sm * rhoul_rhoe1 + rhoul * sm_rhoe1 + nx * ps_rhoe1;
624 
625  jac1(2, 0) = sm * rhovl_rho1 + rhovl * sm_rho1 + ny * ps_rho1;
626  jac1(2, 1) = sm * rhovl_rhou1 + rhovl * sm_rhou1 + ny * ps_rhou1;
627  jac1(2, 2) = sm * rhovl_rhov1 + rhovl * sm_rhov1 + ny * ps_rhov1;
628  jac1(2, 3) = sm * rhovl_rhow1 + rhovl * sm_rhow1 + ny * ps_rhow1;
629  jac1(2, 4) = sm * rhovl_rhoe1 + rhovl * sm_rhoe1 + ny * ps_rhoe1;
630 
631  jac1(3, 0) = sm * rhowl_rho1 + rhowl * sm_rho1 + nz * ps_rho1;
632  jac1(3, 1) = sm * rhowl_rhou1 + rhowl * sm_rhou1 + nz * ps_rhou1;
633  jac1(3, 2) = sm * rhowl_rhov1 + rhowl * sm_rhov1 + nz * ps_rhov1;
634  jac1(3, 3) = sm * rhowl_rhow1 + rhowl * sm_rhow1 + nz * ps_rhow1;
635  jac1(3, 4) = sm * rhowl_rhoe1 + rhowl * sm_rhoe1 + nz * ps_rhoe1;
636 
637  jac1(4, 0) = sm * (rhoel_rho1 + ps_rho1) + rhoepl * sm_rho1;
638  jac1(4, 1) = sm * (rhoel_rhou1 + ps_rhou1) + rhoepl * sm_rhou1;
639  jac1(4, 2) = sm * (rhoel_rhov1 + ps_rhov1) + rhoepl * sm_rhov1;
640  jac1(4, 3) = sm * (rhoel_rhow1 + ps_rhow1) + rhoepl * sm_rhow1;
641  jac1(4, 4) = sm * (rhoel_rhoe1 + ps_rhoe1) + rhoepl * sm_rhoe1;
642  }
643  else if (sm <= 0. && s2 >= 0.)
644  {
646 
647  Real omeg2 = 1. / (s2 - sm);
648  Real prsta = rho2 * dsv2 * (sm - vdon2) + pres2;
649 
650  Real prst2 = prsta - pres2;
651 
653 
654  Real rhor = omeg2 * dsv2 * rho2;
655  Real rhour = omeg2 * (dsv2 * rhou2 + prst2 * nx);
656  Real rhovr = omeg2 * (dsv2 * rhov2 + prst2 * ny);
657  Real rhowr = omeg2 * (dsv2 * rhow2 + prst2 * nz);
658  Real rhoer = omeg2 * (dsv2 * rhoe2 - pres2 * vdon2 + prsta * sm);
659 
660  Real rhoepr = rhoer + prsta;
661 
663 
664  Real rhotd = std::max(eps, rho2 * (s2 - vdon2) - rho1 * (s1 - vdon1));
665  Real rhotm = 1. / rhotd;
666 
668 
669  Real sm_rho1 = rhotm * (-vdon1 * vdon1 + rq051 + sm * s1);
670  Real sm_rhou1 = rhotm * (nx * (2. * vdon1 - s1 - sm) - gamm1 * uadv1);
671  Real sm_rhov1 = rhotm * (ny * (2. * vdon1 - s1 - sm) - gamm1 * vadv1);
672  Real sm_rhow1 = rhotm * (nz * (2. * vdon1 - s1 - sm) - gamm1 * wadv1);
673  Real sm_rhoe1 = rhotm * gamm1;
674 
676 
677  Real ps_rho1 = rho2 * dsv2 * sm_rho1;
678  Real ps_rhou1 = rho2 * dsv2 * sm_rhou1;
679  Real ps_rhov1 = rho2 * dsv2 * sm_rhov1;
680  Real ps_rhow1 = rho2 * dsv2 * sm_rhow1;
681  Real ps_rhoe1 = rho2 * dsv2 * sm_rhoe1;
682 
684 
685  Real rhor_rho1 = omeg2 * rhor * sm_rho1;
686  Real rhor_rhou1 = omeg2 * rhor * sm_rhou1;
687  Real rhor_rhov1 = omeg2 * rhor * sm_rhov1;
688  Real rhor_rhow1 = omeg2 * rhor * sm_rhow1;
689  Real rhor_rhoe1 = omeg2 * rhor * sm_rhoe1;
690 
692 
693  Real rhour_rho1 = omeg2 * (nx * ps_rho1 + rhour * sm_rho1);
694  Real rhour_rhou1 = omeg2 * (nx * ps_rhou1 + rhour * sm_rhou1);
695  Real rhour_rhov1 = omeg2 * (nx * ps_rhov1 + rhour * sm_rhov1);
696  Real rhour_rhow1 = omeg2 * (nx * ps_rhow1 + rhour * sm_rhow1);
697  Real rhour_rhoe1 = omeg2 * (nx * ps_rhoe1 + rhour * sm_rhoe1);
698 
700 
701  Real rhovr_rho1 = omeg2 * (ny * ps_rho1 + rhovr * sm_rho1);
702  Real rhovr_rhou1 = omeg2 * (ny * ps_rhou1 + rhovr * sm_rhou1);
703  Real rhovr_rhov1 = omeg2 * (ny * ps_rhov1 + rhovr * sm_rhov1);
704  Real rhovr_rhow1 = omeg2 * (ny * ps_rhow1 + rhovr * sm_rhow1);
705  Real rhovr_rhoe1 = omeg2 * (ny * ps_rhoe1 + rhovr * sm_rhoe1);
706 
708 
709  Real rhowr_rho1 = omeg2 * (nz * ps_rho1 + rhowr * sm_rho1);
710  Real rhowr_rhou1 = omeg2 * (nz * ps_rhou1 + rhowr * sm_rhou1);
711  Real rhowr_rhov1 = omeg2 * (nz * ps_rhov1 + rhowr * sm_rhov1);
712  Real rhowr_rhow1 = omeg2 * (nz * ps_rhow1 + rhowr * sm_rhow1);
713  Real rhowr_rhoe1 = omeg2 * (nz * ps_rhoe1 + rhowr * sm_rhoe1);
714 
716 
717  Real rhoer_rho1 = omeg2 * (sm * ps_rho1 + rhoepr * sm_rho1);
718  Real rhoer_rhou1 = omeg2 * (sm * ps_rhou1 + rhoepr * sm_rhou1);
719  Real rhoer_rhov1 = omeg2 * (sm * ps_rhov1 + rhoepr * sm_rhov1);
720  Real rhoer_rhow1 = omeg2 * (sm * ps_rhow1 + rhoepr * sm_rhow1);
721  Real rhoer_rhoe1 = omeg2 * (sm * ps_rhoe1 + rhoepr * sm_rhoe1);
722 
724 
725  jac1(0, 0) = sm * rhor_rho1 + rhor * sm_rho1;
726  jac1(0, 1) = sm * rhor_rhou1 + rhor * sm_rhou1;
727  jac1(0, 2) = sm * rhor_rhov1 + rhor * sm_rhov1;
728  jac1(0, 3) = sm * rhor_rhow1 + rhor * sm_rhow1;
729  jac1(0, 4) = sm * rhor_rhoe1 + rhor * sm_rhoe1;
730 
731  jac1(1, 0) = sm * rhour_rho1 + rhour * sm_rho1 + nx * ps_rho1;
732  jac1(1, 1) = sm * rhour_rhou1 + rhour * sm_rhou1 + nx * ps_rhou1;
733  jac1(1, 2) = sm * rhour_rhov1 + rhour * sm_rhov1 + nx * ps_rhov1;
734  jac1(1, 3) = sm * rhour_rhow1 + rhour * sm_rhow1 + nx * ps_rhow1;
735  jac1(1, 4) = sm * rhour_rhoe1 + rhour * sm_rhoe1 + nx * ps_rhoe1;
736 
737  jac1(2, 0) = sm * rhovr_rho1 + rhovr * sm_rho1 + ny * ps_rho1;
738  jac1(2, 1) = sm * rhovr_rhou1 + rhovr * sm_rhou1 + ny * ps_rhou1;
739  jac1(2, 2) = sm * rhovr_rhov1 + rhovr * sm_rhov1 + ny * ps_rhov1;
740  jac1(2, 3) = sm * rhovr_rhow1 + rhovr * sm_rhow1 + ny * ps_rhow1;
741  jac1(2, 4) = sm * rhovr_rhoe1 + rhovr * sm_rhoe1 + ny * ps_rhoe1;
742 
743  jac1(3, 0) = sm * rhowr_rho1 + rhowr * sm_rho1 + nz * ps_rho1;
744  jac1(3, 1) = sm * rhowr_rhou1 + rhowr * sm_rhou1 + nz * ps_rhou1;
745  jac1(3, 2) = sm * rhowr_rhov1 + rhowr * sm_rhov1 + nz * ps_rhov1;
746  jac1(3, 3) = sm * rhowr_rhow1 + rhowr * sm_rhow1 + nz * ps_rhow1;
747  jac1(3, 4) = sm * rhowr_rhoe1 + rhowr * sm_rhoe1 + nz * ps_rhoe1;
748 
749  jac1(4, 0) = sm * (rhoer_rho1 + ps_rho1) + rhoepr * sm_rho1;
750  jac1(4, 1) = sm * (rhoer_rhou1 + ps_rhou1) + rhoepr * sm_rhou1;
751  jac1(4, 2) = sm * (rhoer_rhov1 + ps_rhov1) + rhoepr * sm_rhov1;
752  jac1(4, 3) = sm * (rhoer_rhow1 + ps_rhow1) + rhoepr * sm_rhow1;
753  jac1(4, 4) = sm * (rhoer_rhoe1 + ps_rhoe1) + rhoepr * sm_rhoe1;
754  }
755  else if (s2 < 0.)
756  {
758  }
759  else
760  {
762 
763  Real omeg1 = 1. / (s1 - sm);
764  Real omeg2 = 1. / (s2 - sm);
765  Real prsta = rho1 * dsv1 * (sm - vdon1) + pres1;
766 
767  mooseError("Weird wave speed occured in ",
768  name(),
769  ": ",
770  __FUNCTION__,
771  "\n",
772  "iside = ",
773  iside,
774  "\n",
775  "ielem = ",
776  ielem,
777  "\n",
778  "rho1 = ",
779  rho1,
780  "\n",
781  "rhou1 = ",
782  rhou1,
783  "\n",
784  "rhov1 = ",
785  rhov1,
786  "\n",
787  "rhow1 = ",
788  rhow1,
789  "\n",
790  "rhoe1 = ",
791  rhoe1,
792  "\n",
793  "pres1 = ",
794  pres1,
795  "\n",
796  "enth1 = ",
797  enth1,
798  "\n",
799  "csou1 = ",
800  csou1,
801  "\n",
802  "rho2 = ",
803  rho2,
804  "\n",
805  "rhou2 = ",
806  rhou2,
807  "\n",
808  "rhov2 = ",
809  rhov2,
810  "\n",
811  "rhoe2 = ",
812  rhoe2,
813  "\n",
814  "pres2 = ",
815  pres2,
816  "\n",
817  "enth2 = ",
818  enth2,
819  "\n",
820  "csou2 = ",
821  csou2,
822  "\n",
823  "vdon1 = ",
824  vdon1,
825  "\n",
826  "vdon2 = ",
827  vdon2,
828  "\n",
829  "vnave = ",
830  vnave,
831  "\n",
832  "cssav = ",
833  cssav,
834  "\n",
835  "s1 = ",
836  s1,
837  "\n",
838  "s2 = ",
839  s2,
840  "\n",
841  "sm = ",
842  sm,
843  "\n",
844  "omeg1 = ",
845  omeg1,
846  "\n",
847  "omeg2 = ",
848  omeg2,
849  "\n",
850  "prsta = ",
851  prsta,
852  "\n",
853  "Please check before continuing!\n");
854  }
855  }
856  else if (mach1 >= 0. && mach1 < 1.)
857  {
863 
864  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
865  Real enth1 = (rhoe1 + _inf_pres) * rhom1;
866 
867  jac1(0, 0) = 0.;
868  jac1(0, 1) = nx;
869  jac1(0, 2) = ny;
870  jac1(0, 3) = nz;
871  jac1(0, 4) = 0.;
872 
873  jac1(1, 0) = -vdon1 * uadv1;
874  jac1(1, 1) = nx * uadv1 + vdon1;
875  jac1(1, 2) = ny * uadv1;
876  jac1(1, 3) = nz * uadv1;
877  jac1(1, 4) = 0.;
878 
879  jac1(2, 0) = -vdon1 * vadv1;
880  jac1(2, 1) = nx * vadv1;
881  jac1(2, 2) = ny * vadv1 + vdon1;
882  jac1(2, 3) = nz * vadv1;
883  jac1(2, 4) = 0.;
884 
885  jac1(3, 0) = -vdon1 * wadv1;
886  jac1(3, 1) = nx * wadv1;
887  jac1(3, 2) = ny * wadv1;
888  jac1(3, 3) = nz * wadv1 + vdon1;
889  jac1(3, 4) = 0.;
890 
891  jac1(4, 0) = -vdon1 * enth1;
892  jac1(4, 1) = nx * enth1;
893  jac1(4, 2) = ny * enth1;
894  jac1(4, 3) = nz * enth1;
895  jac1(4, 4) = vdon1;
896  }
897  else if (mach1 <= -1.)
898  {
904 
906  }
907  else if (mach1 >= 1.)
908  {
914 
915  Real vdon1 = uadv1 * nx + vadv1 * ny + wadv1 * nz;
916  Real enth1 = (rhoe1 + pres1) * rhom1;
917 
918  jac1(0, 0) = 0.;
919  jac1(0, 1) = nx;
920  jac1(0, 2) = ny;
921  jac1(0, 3) = nz;
922  jac1(0, 4) = 0.;
923 
924  jac1(1, 0) = rq051 * nx - uadv1 * vdon1;
925  jac1(1, 1) = gamm2 * nx * uadv1 + vdon1;
926  jac1(1, 2) = ny * uadv1 - vadv1 * gamm1 * nx;
927  jac1(1, 3) = nz * uadv1 - wadv1 * gamm1 * nx;
928  jac1(1, 4) = gamm1 * nx;
929 
930  jac1(2, 0) = rq051 * ny - vadv1 * vdon1;
931  jac1(2, 1) = nx * vadv1 - uadv1 * gamm1 * ny;
932  jac1(2, 2) = gamm2 * ny * vadv1 + vdon1;
933  jac1(2, 3) = nz * vadv1 - wadv1 * gamm1 * ny;
934  jac1(2, 4) = gamm1 * ny;
935 
936  jac1(3, 0) = rq051 * nz - wadv1 * vdon1;
937  jac1(3, 1) = nx * wadv1 - uadv1 * gamm1 * nz;
938  jac1(3, 2) = ny * wadv1 - vadv1 * gamm1 * nz;
939  jac1(3, 3) = gamm2 * nz * wadv1 + vdon1;
940  jac1(3, 4) = gamm1 * nz;
941 
942  jac1(4, 0) = (rq051 - enth1) * vdon1;
943  jac1(4, 1) = nx * enth1 - gamm1 * uadv1 * vdon1;
944  jac1(4, 2) = ny * enth1 - gamm1 * vadv1 * vdon1;
945  jac1(4, 3) = nz * enth1 - gamm1 * wadv1 * vdon1;
946  jac1(4, 4) = gamma * vdon1;
947  }
948  else
949  mooseError("Something is wrong in ",
950  name(),
951  ": ",
952  __FUNCTION__,
953  "\n",
954  "ielem = ",
955  ielem,
956  "\n",
957  "iside = ",
958  iside,
959  "\n",
960  "mach1 = ",
961  mach1,
962  "\n");
963 }
virtual Real c(Real v, Real u) const =0
Sound speed.
virtual Real gamma(Real v, Real u) const
Compute the ratio of specific heats.
virtual Real pressure(Real v, Real u) const =0
Pressure as a function of specific internal energy and specific volume.
virtual Real e(Real pressure, Real rho) const =0
Computes internal energy from pressure and density.
void BoundaryFluxBase::execute ( )
virtualinherited

Definition at line 36 of file BoundaryFluxBase.C.

37 {
38 }
void BoundaryFluxBase::finalize ( )
virtualinherited

Definition at line 41 of file BoundaryFluxBase.C.

42 {
43 }
const std::vector< Real > & BoundaryFluxBase::getFlux ( unsigned int  iside,
dof_id_type  ielem,
const std::vector< Real > &  uvec1,
const RealVectorValue &  dwave,
THREAD_ID  tid 
) const
virtualinherited

Get the boundary flux vector.

Parameters
[in]isidelocal index of current side
[in]ielemglobal index of the current element
[in]uvec1vector of variables on the host side
[in]dwavevector of unit normal

Definition at line 46 of file BoundaryFluxBase.C.

Referenced by CNSFVBC::computeQpResidual(), and AEFVBC::computeQpResidual().

51 {
52  Threads::spin_mutex::scoped_lock lock(_mutex);
53  if (_cached_elem_id != ielem || _cached_side_id != iside)
54  {
55  _cached_elem_id = ielem;
56  _cached_side_id = iside;
57 
58  calcFlux(iside, ielem, uvec1, dwave, _flux[tid]);
59  }
60  return _flux[tid];
61 }
unsigned int _cached_side_id
static Threads::spin_mutex _mutex
std::vector< std::vector< Real > > _flux
Threaded storage for fluxes.
dof_id_type _cached_elem_id
virtual void calcFlux(unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, std::vector< Real > &flux) const =0
Solve the Riemann problem on the boundary face.
const DenseMatrix< Real > & BoundaryFluxBase::getJacobian ( unsigned int  iside,
dof_id_type  ielem,
const std::vector< Real > &  uvec1,
const RealVectorValue &  dwave,
THREAD_ID  tid 
) const
virtualinherited

Get the boundary Jacobian matrix.

Parameters
[in]isidelocal index of current side
[in]ielemglobal index of the current element
[in]uvec1vector of variables on the host side
[in]dwavevector of unit normal

Definition at line 64 of file BoundaryFluxBase.C.

Referenced by CNSFVBC::computeQpJacobian(), AEFVBC::computeQpJacobian(), and CNSFVBC::computeQpOffDiagJacobian().

69 {
70  Threads::spin_mutex::scoped_lock lock(_mutex);
71  if (_cached_elem_id != ielem || _cached_side_id != iside)
72  {
73  _cached_elem_id = ielem;
74  _cached_side_id = iside;
75 
76  calcJacobian(iside, ielem, uvec1, dwave, _jac1[tid]);
77  }
78  return _jac1[tid];
79 }
unsigned int _cached_side_id
static Threads::spin_mutex _mutex
std::vector< DenseMatrix< Real > > _jac1
Threaded storage for jacobians.
virtual void calcJacobian(unsigned int iside, dof_id_type ielem, const std::vector< Real > &uvec1, const RealVectorValue &dwave, DenseMatrix< Real > &jac1) const =0
Compute the Jacobian matrix on the boundary face.
dof_id_type _cached_elem_id
void BoundaryFluxBase::initialize ( )
virtualinherited

Definition at line 29 of file BoundaryFluxBase.C.

30 {
31  _cached_elem_id = 0;
32  _cached_side_id = libMesh::invalid_uint;
33 }
unsigned int _cached_side_id
dof_id_type _cached_elem_id

Member Data Documentation

dof_id_type BoundaryFluxBase::_cached_elem_id
mutableprotectedinherited
unsigned int BoundaryFluxBase::_cached_side_id
mutableprotectedinherited
std::vector<std::vector<Real> > BoundaryFluxBase::_flux
mutableprotectedinherited

Threaded storage for fluxes.

Definition at line 98 of file BoundaryFluxBase.h.

Referenced by BoundaryFluxBase::BoundaryFluxBase(), and BoundaryFluxBase::getFlux().

const SinglePhaseFluidProperties& CNSFVHLLCInflowOutflowBoundaryFlux::_fp
protected

Definition at line 43 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

Real CNSFVHLLCInflowOutflowBoundaryFlux::_inf_pres
protected

Definition at line 49 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

Real CNSFVHLLCInflowOutflowBoundaryFlux::_inf_rho
protected

Definition at line 45 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

Real CNSFVHLLCInflowOutflowBoundaryFlux::_inf_uadv
protected

Definition at line 46 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

Real CNSFVHLLCInflowOutflowBoundaryFlux::_inf_vadv
protected

Definition at line 47 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

Real CNSFVHLLCInflowOutflowBoundaryFlux::_inf_wadv
protected

Definition at line 48 of file CNSFVHLLCInflowOutflowBoundaryFlux.h.

Referenced by calcFlux(), and calcJacobian().

std::vector<DenseMatrix<Real> > BoundaryFluxBase::_jac1
mutableprotectedinherited

Threaded storage for jacobians.

Definition at line 101 of file BoundaryFluxBase.h.

Referenced by BoundaryFluxBase::BoundaryFluxBase(), and BoundaryFluxBase::getJacobian().


The documentation for this class was generated from the following files: