29 #include "grins_config.h"
38 #include "libmesh/getpot.h"
39 #include "libmesh/quadrature.h"
40 #include "libmesh/boundary_info.h"
41 #include "libmesh/fem_system.h"
42 #include "libmesh/fe_interface.h"
46 template<
typename StressStrainLaw>
48 bool is_compressible )
50 _stress_strain_law(input),
52 _is_compressible(is_compressible)
55 if( !input.have_variable(
"Physics/"+physics_name+
"/h0") )
57 std::cerr <<
"Error: Must specify initial thickness for "+physics_name << std::endl
58 <<
" Input the option Physics/"+physics_name+
"/h0" << std::endl;
63 (
_h0, input,
"Physics/"+physics_name+
"/h0",
_h0 );
72 template<
typename StressStrainLaw>
78 template<
typename StressStrainLaw>
88 _lambda_sq_var = system->add_variable(
"lambda_sq", GRINSEnums::FIRST, GRINSEnums::LAGRANGE);
94 template<
typename StressStrainLaw>
100 if( input.have_variable(section) )
102 unsigned int n_vars = input.vector_variable_size(section);
104 for(
unsigned int v = 0; v < n_vars; v++ )
106 std::string name = input(section,
"DIE!",v);
108 if( name == std::string(
"stress") )
112 _stress_indices.resize(7);
128 else if( name == std::string(
"stress_zz" ) )
133 else if( name == std::string(
"strain") )
136 _strain_indices.resize(3);
146 std::cerr <<
"Error: Invalue output_vars value for "+
elastic_membrane << std::endl
147 <<
" Found " << name << std::endl
148 <<
" Acceptable values are: stress" << std::endl
149 <<
" strain" << std::endl;
158 template<
typename StressStrainLaw>
163 const unsigned int n_u_dofs = context.get_dof_indices(_disp_vars.u_var()).size();
165 const std::vector<libMesh::Real> &JxW =
166 this->get_fe(context)->get_JxW();
169 libMesh::DenseSubVector<libMesh::Number> &Fu = context.get_elem_residual(_disp_vars.u_var());
170 libMesh::DenseSubVector<libMesh::Number> &Fv = context.get_elem_residual(_disp_vars.v_var());
171 libMesh::DenseSubVector<libMesh::Number> &Fw = context.get_elem_residual(_disp_vars.w_var());
173 libMesh::DenseSubMatrix<libMesh::Number>& Kuu = context.get_elem_jacobian(_disp_vars.u_var(),_disp_vars.u_var());
174 libMesh::DenseSubMatrix<libMesh::Number>& Kuv = context.get_elem_jacobian(_disp_vars.u_var(),_disp_vars.v_var());
175 libMesh::DenseSubMatrix<libMesh::Number>& Kuw = context.get_elem_jacobian(_disp_vars.u_var(),_disp_vars.w_var());
177 libMesh::DenseSubMatrix<libMesh::Number>& Kvu = context.get_elem_jacobian(_disp_vars.v_var(),_disp_vars.u_var());
178 libMesh::DenseSubMatrix<libMesh::Number>& Kvv = context.get_elem_jacobian(_disp_vars.v_var(),_disp_vars.v_var());
179 libMesh::DenseSubMatrix<libMesh::Number>& Kvw = context.get_elem_jacobian(_disp_vars.v_var(),_disp_vars.w_var());
181 libMesh::DenseSubMatrix<libMesh::Number>& Kwu = context.get_elem_jacobian(_disp_vars.w_var(),_disp_vars.u_var());
182 libMesh::DenseSubMatrix<libMesh::Number>& Kwv = context.get_elem_jacobian(_disp_vars.w_var(),_disp_vars.v_var());
183 libMesh::DenseSubMatrix<libMesh::Number>& Kww = context.get_elem_jacobian(_disp_vars.w_var(),_disp_vars.w_var());
185 unsigned int n_qpoints = context.get_element_qrule().n_points();
188 const std::vector<std::vector<libMesh::Real> >& dphi_dxi =
189 this->get_fe(context)->get_dphidxi();
191 const std::vector<std::vector<libMesh::Real> >& dphi_deta =
192 this->get_fe(context)->get_dphideta();
194 const libMesh::DenseSubVector<libMesh::Number>& u_coeffs = context.get_elem_solution( _disp_vars.u_var() );
195 const libMesh::DenseSubVector<libMesh::Number>& v_coeffs = context.get_elem_solution( _disp_vars.v_var() );
196 const libMesh::DenseSubVector<libMesh::Number>& w_coeffs = context.get_elem_solution( _disp_vars.w_var() );
199 const std::vector<libMesh::RealGradient>& dxdxi = this->get_fe(context)->get_dxyzdxi();
200 const std::vector<libMesh::RealGradient>& dxdeta = this->get_fe(context)->get_dxyzdeta();
202 const unsigned int dim = 2;
204 for (
unsigned int qp=0; qp != n_qpoints; qp++)
207 libMesh::Gradient grad_u, grad_v, grad_w;
208 for(
unsigned int d = 0; d < n_u_dofs; d++ )
210 libMesh::RealGradient u_gradphi( dphi_dxi[d][qp], dphi_deta[d][qp] );
211 grad_u += u_coeffs(d)*u_gradphi;
212 grad_v += v_coeffs(d)*u_gradphi;
213 grad_w += w_coeffs(d)*u_gradphi;
216 libMesh::RealGradient grad_x( dxdxi[qp](0), dxdeta[qp](0) );
217 libMesh::RealGradient grad_y( dxdxi[qp](1), dxdeta[qp](1) );
218 libMesh::RealGradient grad_z( dxdxi[qp](2), dxdeta[qp](2) );
222 libMesh::Real lambda_sq = 0;
224 this->compute_metric_tensors( qp, *(this->get_fe(context)), context,
225 grad_u, grad_v, grad_w,
226 a_cov, a_contra, A_cov, A_contra,
232 _stress_strain_law.compute_stress_and_elasticity(dim,a_contra,a_cov,A_contra,A_cov,tau,C);
234 libMesh::Real jac = JxW[qp];
236 for (
unsigned int i=0; i != n_u_dofs; i++)
238 libMesh::RealGradient u_gradphi( dphi_dxi[i][qp], dphi_deta[i][qp] );
240 for(
unsigned int alpha = 0; alpha < dim; alpha++ )
242 for(
unsigned int beta = 0; beta < dim; beta++ )
244 Fu(i) -= 0.5*tau(alpha,beta)*_h0*( (grad_x(beta) + grad_u(beta))*u_gradphi(alpha) +
245 (grad_x(alpha) + grad_u(alpha))*u_gradphi(beta) )*jac;
247 Fv(i) -= 0.5*tau(alpha,beta)*_h0*( (grad_y(beta) + grad_v(beta))*u_gradphi(alpha) +
248 (grad_y(alpha) + grad_v(alpha))*u_gradphi(beta) )*jac;
250 Fw(i) -= 0.5*tau(alpha,beta)*_h0*( (grad_z(beta) + grad_w(beta))*u_gradphi(alpha) +
251 (grad_z(alpha) + grad_w(alpha))*u_gradphi(beta) )*jac;
256 if( compute_jacobian )
258 for (
unsigned int i=0; i != n_u_dofs; i++)
260 libMesh::RealGradient u_gradphi_i( dphi_dxi[i][qp], dphi_deta[i][qp] );
262 for (
unsigned int j=0; j != n_u_dofs; j++)
264 libMesh::RealGradient u_gradphi_j( dphi_dxi[j][qp], dphi_deta[j][qp] );
266 for(
unsigned int alpha = 0; alpha < dim; alpha++ )
268 for(
unsigned int beta = 0; beta < dim; beta++ )
270 const libMesh::Real diag_term = 0.5*_h0*jac*tau(alpha,beta)*( u_gradphi_j(beta)*u_gradphi_i(alpha) +
271 u_gradphi_j(alpha)*u_gradphi_i(beta) );
272 Kuu(i,j) -= diag_term;
274 Kvv(i,j) -= diag_term;
276 Kww(i,j) -= diag_term;
278 for(
unsigned int lambda = 0; lambda < dim; lambda++ )
280 for(
unsigned int mu = 0; mu < dim; mu++ )
282 const libMesh::Real dgamma_du = 0.5*( u_gradphi_j(lambda)*(grad_x(mu)+grad_u(mu)) +
283 (grad_x(lambda)+grad_u(lambda))*u_gradphi_j(mu) );
285 const libMesh::Real dgamma_dv = 0.5*( u_gradphi_j(lambda)*(grad_y(mu)+grad_v(mu)) +
286 (grad_y(lambda)+grad_v(lambda))*u_gradphi_j(mu) );
288 const libMesh::Real dgamma_dw = 0.5*( u_gradphi_j(lambda)*(grad_z(mu)+grad_w(mu)) +
289 (grad_z(lambda)+grad_w(lambda))*u_gradphi_j(mu) );
291 const libMesh::Real C1 = 0.5*_h0*jac*C(alpha,beta,lambda,mu);
293 const libMesh::Real x_term = C1*( (grad_x(beta)+grad_u(beta))*u_gradphi_i(alpha) +
294 (grad_x(alpha)+grad_u(alpha))*u_gradphi_i(beta) );
296 const libMesh::Real y_term = C1*( (grad_y(beta)+grad_v(beta))*u_gradphi_i(alpha) +
297 (grad_y(alpha)+grad_v(alpha))*u_gradphi_i(beta) );
299 const libMesh::Real z_term = C1*( (grad_z(beta)+grad_w(beta))*u_gradphi_i(alpha) +
300 (grad_z(alpha)+grad_w(alpha))*u_gradphi_i(beta) );
302 Kuu(i,j) -= x_term*dgamma_du;
304 Kuv(i,j) -= x_term*dgamma_dv;
306 Kuw(i,j) -= x_term*dgamma_dw;
308 Kvu(i,j) -= y_term*dgamma_du;
310 Kvv(i,j) -= y_term*dgamma_dv;
312 Kvw(i,j) -= y_term*dgamma_dw;
314 Kwu(i,j) -= z_term*dgamma_du;
316 Kwv(i,j) -= z_term*dgamma_dv;
318 Kww(i,j) -= z_term*dgamma_dw;
332 template<
typename StressStrainLaw>
338 if( _is_compressible )
340 unsigned int n_qpoints = context.get_element_qrule().n_points();
342 const unsigned int n_u_dofs = context.get_dof_indices(_disp_vars.u_var()).size();
344 const std::vector<libMesh::Real> &JxW = context.get_element_fe(this->_lambda_sq_var)->get_JxW();
346 libMesh::DenseSubVector<libMesh::Number>& Fl = context.get_elem_residual(this->_lambda_sq_var);
348 const std::vector<std::vector<libMesh::Real> >& phi =
349 context.get_element_fe(this->_lambda_sq_var)->get_phi();
351 const unsigned int n_lambda_sq_dofs = context.get_dof_indices(this->_lambda_sq_var).size();
353 const libMesh::DenseSubVector<libMesh::Number>& u_coeffs = context.get_elem_solution( _disp_vars.u_var() );
354 const libMesh::DenseSubVector<libMesh::Number>& v_coeffs = context.get_elem_solution( _disp_vars.v_var() );
355 const libMesh::DenseSubVector<libMesh::Number>& w_coeffs = context.get_elem_solution( _disp_vars.w_var() );
358 const std::vector<std::vector<libMesh::Real> >& dphi_dxi =
359 this->get_fe(context)->get_dphidxi();
361 const std::vector<std::vector<libMesh::Real> >& dphi_deta =
362 this->get_fe(context)->get_dphideta();
364 for (
unsigned int qp=0; qp != n_qpoints; qp++)
366 libMesh::Real jac = JxW[qp];
368 libMesh::Gradient grad_u, grad_v, grad_w;
369 for(
unsigned int d = 0; d < n_u_dofs; d++ )
371 libMesh::RealGradient u_gradphi( dphi_dxi[d][qp], dphi_deta[d][qp] );
372 grad_u += u_coeffs(d)*u_gradphi;
373 grad_v += v_coeffs(d)*u_gradphi;
374 grad_w += w_coeffs(d)*u_gradphi;
378 libMesh::Real lambda_sq = 0;
380 this->compute_metric_tensors( qp, *(this->get_fe(context)), context,
381 grad_u, grad_v, grad_w,
382 a_cov, a_contra, A_cov, A_contra,
385 libMesh::Real stress_33 = _stress_strain_law.compute_33_stress( a_contra, a_cov, A_contra, A_cov );
387 for (
unsigned int i=0; i != n_lambda_sq_dofs; i++)
389 Fl(i) += stress_33*phi[i][qp]*jac;
392 if( compute_jacobian )
394 libmesh_not_implemented();
402 template<
typename StressStrainLaw>
422 template<
typename StressStrainLaw>
427 libmesh_not_implemented();
431 template<
typename StressStrainLaw>
434 const libMesh::Point& point,
435 libMesh::Real& value )
437 bool is_stress =
false;
438 if( !_stress_indices.empty() )
439 is_stress= ( _stress_indices[0] == quantity_index ||
440 _stress_indices[1] == quantity_index ||
441 _stress_indices[2] == quantity_index ||
442 _stress_indices[3] == quantity_index ||
443 _stress_indices[4] == quantity_index ||
444 _stress_indices[5] == quantity_index ||
445 _stress_indices[6] == quantity_index ||
446 _stress_zz_index == quantity_index );
448 bool is_strain =
false;
449 if( !_strain_indices.empty() )
450 is_strain = ( _strain_indices[0] == quantity_index ||
451 _strain_indices[1] == quantity_index ||
452 _strain_indices[2] == quantity_index );
454 if( is_stress || is_strain )
456 const unsigned int n_u_dofs = context.get_dof_indices(_disp_vars.u_var()).size();
458 const libMesh::DenseSubVector<libMesh::Number>& u_coeffs = context.get_elem_solution( _disp_vars.u_var() );
459 const libMesh::DenseSubVector<libMesh::Number>& v_coeffs = context.get_elem_solution( _disp_vars.v_var() );
460 const libMesh::DenseSubVector<libMesh::Number>& w_coeffs = context.get_elem_solution( _disp_vars.w_var() );
463 libMesh::AutoPtr<libMesh::FEGenericBase<libMesh::Real> > fe_new =
464 this->build_new_fe( context.get_elem(), this->get_fe(context),
467 const std::vector<std::vector<libMesh::Real> >& dphi_dxi =
468 fe_new->get_dphidxi();
470 const std::vector<std::vector<libMesh::Real> >& dphi_deta =
471 fe_new->get_dphideta();
474 const std::vector<libMesh::RealGradient>& dxdxi = fe_new->get_dxyzdxi();
475 const std::vector<libMesh::RealGradient>& dxdeta = fe_new->get_dxyzdeta();
478 libMesh::Gradient grad_u, grad_v, grad_w;
479 for(
unsigned int d = 0; d < n_u_dofs; d++ )
481 libMesh::RealGradient u_gradphi( dphi_dxi[d][0], dphi_deta[d][0] );
482 grad_u += u_coeffs(d)*u_gradphi;
483 grad_v += v_coeffs(d)*u_gradphi;
484 grad_w += w_coeffs(d)*u_gradphi;
487 libMesh::RealGradient grad_x( dxdxi[0](0), dxdeta[0](0) );
488 libMesh::RealGradient grad_y( dxdxi[0](1), dxdeta[0](1) );
489 libMesh::RealGradient grad_z( dxdxi[0](2), dxdeta[0](2) );
492 libMesh::Real lambda_sq = 0;
495 this->compute_metric_tensors(0, *fe_new, context, grad_u, grad_v, grad_w,
496 a_cov, a_contra, A_cov, A_contra, lambda_sq );
501 if( _strain_indices[0] == quantity_index )
503 value = 0.5*(A_cov(0,0) - a_cov(0,0));
505 else if( _strain_indices[1] == quantity_index )
507 value = 0.5*(A_cov(0,1) - a_cov(0,1));
509 else if( _strain_indices[2] == quantity_index )
511 value = 0.5*(A_cov(1,1) - a_cov(1,1));
523 libMesh::Real det_a = a_cov(0,0)*a_cov(1,1) - a_cov(0,1)*a_cov(1,0);
524 libMesh::Real det_A = A_cov(0,0)*A_cov(1,1) - A_cov(0,1)*A_cov(1,0);
526 libMesh::Real I3 = lambda_sq*det_A/det_a;
529 _stress_strain_law.compute_stress(2,a_contra,a_cov,A_contra,A_cov,tau);
531 if( _stress_indices[0] == quantity_index )
534 value = tau(0,0)/std::sqrt(I3);
536 else if( _stress_indices[1] == quantity_index )
539 value = tau(0,1)/std::sqrt(I3);
541 else if( _stress_indices[2] == quantity_index )
544 value = tau(1,1)/std::sqrt(I3);
546 else if( _stress_indices[3] == quantity_index )
548 value = 0.5*(tau(0,0) + tau(1,1)) + std::sqrt(0.25*(tau(0,0)-tau(1,1))*(tau(0,0)-tau(1,1))
549 + tau(0,1)*tau(0,1) );
551 else if( _stress_indices[4] == quantity_index ||
552 _stress_indices[5] == quantity_index ||
553 _stress_indices[6] == quantity_index )
557 tau(2,2) = _stress_strain_law.compute_33_stress( a_contra, a_cov, A_contra, A_cov );
560 libMesh::Real stress_I1 = tau(0,0) + tau(1,1) + tau(2,2);
561 libMesh::Real stress_I2 = 0.5*(stress_I1*stress_I1 - (tau(0,0)*tau(0,0) + tau(1,1)*tau(1,1)
562 + tau(2,2)*tau(2,2) + tau(0,1)*tau(0,1)
563 + tau(1,0)*tau(1,0)) );
565 libMesh::Real stress_I3 = tau(2,2)*( tau(0,0)*tau(1,1) - tau(1,0)*tau(0,1) );
571 libMesh::Real C1 = (stress_I1*stress_I1 - 3*stress_I2);
574 libMesh::Real C2 = (2*stress_I1*stress_I1*stress_I1 - 9*stress_I1*stress_I2 + 27*stress_I3)/54;
576 libMesh::Real theta = std::acos( C2/(2*std::sqrt(C1*C1*C1)) )/3.0;
578 if( _stress_indices[4] == quantity_index )
580 value = (stress_I1 + 2.0*std::sqrt(C1)*std::cos(theta))/3.0;
583 if( _stress_indices[5] == quantity_index )
585 value = (stress_I1 + 2.0*std::sqrt(C1)*std::cos(theta+
Constants::two_pi/3.0))/3.0;
588 if( _stress_indices[6] == quantity_index )
590 value = (stress_I1 + 2.0*std::sqrt(C1)*std::cos(theta+2.0*
Constants::two_pi/3.0))/3.0;
593 else if( _stress_zz_index == quantity_index )
595 value = _stress_strain_law.compute_33_stress( a_contra, a_cov, A_contra, A_cov );
609 template<
typename StressStrainLaw>
611 const libMesh::FEGenericBase<libMesh::Real>* fe,
612 const libMesh::Point p )
615 FEType fe_type = fe->get_fe_type();
620 libmesh_assert(&elem || fe_type.family == SCALAR);
622 unsigned int elem_dim = &elem ? elem.dim() : 0;
624 AutoPtr<FEGenericBase<libMesh::Real> >
625 fe_new(FEGenericBase<libMesh::Real>::build(elem_dim, fe_type));
629 Point master_point = &elem ?
630 FEInterface::inverse_map(elem_dim, fe_type, &elem, p) :
633 std::vector<Point> coor(1, master_point);
636 fe_new->reinit (&elem, &coor);
641 template<
typename StressStrainLaw>
643 const libMesh::FEBase& elem,
645 const libMesh::Gradient& grad_u,
646 const libMesh::Gradient& grad_v,
647 const libMesh::Gradient& grad_w,
652 libMesh::Real& lambda_sq )
654 const std::vector<libMesh::RealGradient>& dxdxi = elem.get_dxyzdxi();
655 const std::vector<libMesh::RealGradient>& dxdeta = elem.get_dxyzdeta();
657 const std::vector<libMesh::Real>& dxidx = elem.get_dxidx();
658 const std::vector<libMesh::Real>& dxidy = elem.get_dxidy();
659 const std::vector<libMesh::Real>& dxidz = elem.get_dxidz();
661 const std::vector<libMesh::Real>& detadx = elem.get_detadx();
662 const std::vector<libMesh::Real>& detady = elem.get_detady();
663 const std::vector<libMesh::Real>& detadz = elem.get_detadz();
665 libMesh::RealGradient dxi( dxidx[qp], dxidy[qp], dxidz[qp] );
666 libMesh::RealGradient deta( detadx[qp], detady[qp], detadz[qp] );
668 libMesh::RealGradient dudxi( grad_u(0), grad_v(0), grad_w(0) );
669 libMesh::RealGradient dudeta( grad_u(1), grad_v(1), grad_w(1) );
673 a_cov(0,0) = dxdxi[qp]*dxdxi[qp];
674 a_cov(0,1) = dxdxi[qp]*dxdeta[qp];
675 a_cov(1,0) = dxdeta[qp]*dxdxi[qp];
676 a_cov(1,1) = dxdeta[qp]*dxdeta[qp];
678 libMesh::Real det_a = a_cov(0,0)*a_cov(1,1) - a_cov(0,1)*a_cov(1,0);
682 A_cov(0,0) = (dxdxi[qp] + dudxi)*(dxdxi[qp] + dudxi);
683 A_cov(0,1) = (dxdxi[qp] + dudxi)*(dxdeta[qp] + dudeta);
684 A_cov(1,0) = (dxdeta[qp] + dudeta)*(dxdxi[qp] + dudxi);
685 A_cov(1,1) = (dxdeta[qp] + dudeta)*(dxdeta[qp] + dudeta);
689 a_contra(0,0) = dxi*dxi;
690 a_contra(0,1) = dxi*deta;
691 a_contra(1,0) = deta*dxi;
692 a_contra(1,1) = deta*deta;
695 libMesh::Real det_A = A_cov(0,0)*A_cov(1,1) - A_cov(0,1)*A_cov(1,0);
698 A_contra(0,0) = A_cov(1,1)/det_A;
699 A_contra(0,1) = -A_cov(0,1)/det_A;
700 A_contra(1,0) = -A_cov(1,0)/det_A;
701 A_contra(1,1) = A_cov(0,0)/det_A;
708 if( _is_compressible )
710 lambda_sq = context.interior_value(this->_lambda_sq_var, qp);
715 lambda_sq = det_a/det_A;
718 A_cov(2,2) = lambda_sq;
719 A_contra(2,2) = 1.0/lambda_sq;
virtual void set_parameter(libMesh::Number ¶m_variable, const GetPot &input, const std::string ¶m_name, libMesh::Number param_default)
Each subclass can simultaneously read a parameter value from.
unsigned int register_quantity(std::string name)
Register quantity to be postprocessed.
GRINS::ICHandlingBase * _ic_handler
GRINS::BCHandlingBase * _bc_handler
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...
Base class for reading and handling initial conditions for physics classes.
void compute_metric_tensors(unsigned int qp, const libMesh::FEBase &elem, const AssemblyContext &context, const libMesh::Gradient &grad_u, const libMesh::Gradient &grad_v, const libMesh::Gradient &grad_w, libMesh::TensorValue< libMesh::Real > &a_cov, libMesh::TensorValue< libMesh::Real > &a_contra, libMesh::TensorValue< libMesh::Real > &A_cov, libMesh::TensorValue< libMesh::Real > &A_contra, libMesh::Real &lambda_sq)
virtual void element_constraint(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Constraint part(s) of physics for element interiors.
const PhysicsName elastic_membrane
virtual void init_variables(libMesh::FEMSystem *system)
Initialize variables for this physics.
virtual void element_time_derivative(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Time dependent part(s) of physics for element interiors.
virtual void register_postprocessing_vars(const GetPot &input, PostProcessedQuantities< libMesh::Real > &postprocessing)
Register postprocessing variables for ElasticMembrane.
virtual void compute_postprocessed_quantity(unsigned int quantity_index, const AssemblyContext &context, const libMesh::Point &point, libMesh::Real &value)
Compute the registered postprocessed quantities.
virtual void init_variables(libMesh::FEMSystem *system)
Initialize variables for this physics.
const libMesh::Real two_pi
virtual ~ElasticMembrane()
libMesh::AutoPtr< libMesh::FEGenericBase< libMesh::Real > > build_new_fe(const libMesh::Elem &elem, const libMesh::FEGenericBase< libMesh::Real > *fe, const libMesh::Point p)
virtual void side_time_derivative(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Time dependent part(s) of physics for boundaries of elements on the domain boundary.