36 #include "libmesh/quadrature.h"
70 #ifdef GRINS_USE_GRVY_TIMERS
71 this->_timer->BeginTimer(
"SpalartAllmarasSPGSMStabilization::element_time_derivative");
76 libMesh::Elem &elem_pointer = context.get_elem();
79 const unsigned int n_nu_dofs = context.get_dof_indices(this->_turbulence_vars.nu()).size();
82 const std::vector<libMesh::Real> &JxW =
83 context.get_element_fe(this->_turbulence_vars.nu())->get_JxW();
87 const std::vector<std::vector<libMesh::RealGradient> >& nu_gradphi =
88 context.get_element_fe(this->_turbulence_vars.nu())->get_dphi();
96 libMesh::DenseSubVector<libMesh::Number> &Fnu = context.get_elem_residual(this->_turbulence_vars.nu());
98 libMesh::FEBase* fe = context.get_element_fe(this->_turbulence_vars.nu());
100 unsigned int n_qpoints = context.get_element_qrule().n_points();
103 libMesh::UniquePtr< libMesh::DenseVector<libMesh::Real> > distance_qp;
106 distance_qp = this->distance_function->interpolate(&elem_pointer, context.get_element_qrule().get_points());
108 for (
unsigned int qp=0; qp != n_qpoints; qp++)
110 libMesh::Gradient grad_nu;
111 grad_nu = context.interior_gradient(this->_turbulence_vars.nu(), qp);
113 libMesh::Real jac = JxW[qp];
116 libMesh::Real _mu_qp = this->_mu(context, qp);
124 u = context.interior_value(this->_flow_vars.u(), qp);
125 v = context.interior_value(this->_flow_vars.v(), qp);
127 libMesh::NumberVectorValue U(u,v);
129 U(2) = context.interior_value(this->_flow_vars.w(), qp);
133 libMesh::RealGradient g = this->_stab_helper.compute_g( fe, context, qp );
134 libMesh::RealTensor G = this->_stab_helper.compute_G( fe, context, qp );
136 libMesh::Real tau_spalart = this->_stab_helper.compute_tau_spalart( context, qp, g, G, this->_rho, U, _mu_qp, this->_is_steady );
138 libMesh::Number RM_spalart = this->_stab_helper.compute_res_spalart_steady( context, qp, this->_rho, _mu_qp, (*distance_qp)(qp), this->_infinite_distance );
140 for (
unsigned int i=0; i != n_nu_dofs; i++)
142 Fnu(i) += jac*( -tau_spalart*RM_spalart*this->_rho*(U*nu_gradphi[i][qp]) );
145 if( compute_jacobian )
147 libmesh_not_implemented();
152 #ifdef GRINS_USE_GRVY_TIMERS
153 this->_timer->EndTimer(
"SpalartAllmarasSPGSMStabilization::element_time_derivative");
164 #ifdef GRINS_USE_GRVY_TIMERS
165 this->_timer->BeginTimer(
"SpalartAllmarasSPGSMStabilization::mass_residual");
170 libMesh::Elem &elem_pointer = context.get_elem();
173 const unsigned int n_nu_dofs = context.get_dof_indices(this->_turbulence_vars.nu()).size();
176 const std::vector<libMesh::Real> &JxW =
177 context.get_element_fe(this->_turbulence_vars.nu())->get_JxW();
180 const std::vector<std::vector<libMesh::RealGradient> >& nu_gradphi =
181 context.get_element_fe(this->_turbulence_vars.nu())->get_dphi();
183 libMesh::DenseSubVector<libMesh::Number> &Fnu = context.get_elem_residual(this->_turbulence_vars.nu());
185 libMesh::FEBase* fe = context.get_element_fe(this->_turbulence_vars.nu());
187 unsigned int n_qpoints = context.get_element_qrule().n_points();
190 libMesh::UniquePtr< libMesh::DenseVector<libMesh::Real> > distance_qp;
193 distance_qp = this->distance_function->interpolate(&elem_pointer, context.get_element_qrule().get_points());
195 for (
unsigned int qp=0; qp != n_qpoints; qp++)
197 libMesh::RealGradient g = this->_stab_helper.compute_g( fe, context, qp );
198 libMesh::RealTensor G = this->_stab_helper.compute_G( fe, context, qp );
200 libMesh::RealGradient U( context.fixed_interior_value( this->_flow_vars.u(), qp ),
201 context.fixed_interior_value( this->_flow_vars.v(), qp ) );
203 libMesh::Real _mu_qp = this->_mu(context, qp);
205 if( this->_dim == 3 )
207 U(2) = context.fixed_interior_value( this->_flow_vars.w(), qp );
210 libMesh::Real tau_spalart = this->_stab_helper.compute_tau_spalart( context, qp, g, G, this->_rho, U, _mu_qp, this->_is_steady );
212 libMesh::Real RM_spalart = this->_stab_helper.compute_res_spalart_transient( context, qp, this->_rho );
214 for (
unsigned int i=0; i != n_nu_dofs; i++)
216 Fnu(i) += -JxW[qp]*tau_spalart*RM_spalart*this->_rho*(U*nu_gradphi[i][qp]);
219 if( compute_jacobian )
221 libmesh_not_implemented();
226 #ifdef GRINS_USE_GRVY_TIMERS
227 this->_timer->EndTimer(
"SpalartAllmarasSPGSMStabilization::mass_residual");
INSTANTIATE_TURBULENCE_MODELS_SUBCLASS(SpalartAllmarasSPGSMStabilization)
SpalartAllmarasSPGSMStabilization()
virtual void init_variables(libMesh::FEMSystem *system)
Initialize variables for this physics.
virtual void init_variables(libMesh::FEMSystem *system)
Initialize variables for this physics.
virtual void register_parameter(const std::string ¶m_name, libMesh::ParameterMultiAccessor< libMesh::Number > ¶m_pointer) const
Each subclass will register its copy of an independent.
virtual void element_time_derivative(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Time dependent part(s) of physics for element interiors.
virtual void register_parameter(const std::string ¶m_name, libMesh::ParameterMultiAccessor< libMesh::Number > ¶m_pointer) const
Each subclass will register its copy of an independent.
virtual void mass_residual(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Mass matrix part(s) for element interiors. All boundary terms lie within the time_derivative part...