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

#include <InertialForceBeam.h>

Inheritance diagram for InertialForceBeam:
[legend]

Public Member Functions

 InertialForceBeam (const InputParameters &parameters)
 
virtual void computeResidual () override
 
virtual void computeJacobian () override
 
virtual void computeOffDiagJacobian (MooseVariableFEBase &jvar) override
 

Protected Member Functions

virtual Real computeQpResidual () override
 

Private Attributes

const MaterialProperty< Real > & _density
 Density of the beam. More...
 
unsigned int _nrot
 Number of coupled rotational variables. More...
 
unsigned int _ndisp
 Number of coupled displacement variables. More...
 
std::vector< unsigned int > _rot_num
 Variable numbers corresponding to rotational variables. More...
 
std::vector< unsigned int > _disp_num
 Variable numbers corresponding to displacement variables. More...
 
std::vector< unsigned int > _vel_num
 Variable numbers corresponding to velocity aux variables. More...
 
std::vector< unsigned int > _accel_num
 Variable numbers corresponding to acceleraion aux variables. More...
 
std::vector< unsigned int > _rot_vel_num
 Variable numbers corresponding to rotational velocity aux variables. More...
 
std::vector< unsigned int > _rot_accel_num
 Variable numbers corresponding to rotational acceleration aux variables. More...
 
const VariableValue & _area
 Coupled variable for beam cross-sectional area. More...
 
const VariableValue & _Ay
 Coupled variable for first moment of area of beam in y direction, i.e., integral of y*dA over the cross-section. More...
 
const VariableValue & _Az
 Coupled variable for first moment of area of beam in z direction, i.e., integral of z*dA over the cross-section. More...
 
const VariableValue & _Ix
 Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA over the cross-section. More...
 
const VariableValue & _Iy
 Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the cross-section. More...
 
const VariableValue & _Iz
 Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the cross-section. More...
 
const Real _beta
 Newmark time integration parameter. More...
 
const Real _gamma
 Newmark time integraion parameter. More...
 
const MaterialProperty< Real > & _eta
 Mass proportional Rayleigh damping parameter. More...
 
const Real _alpha
 HHT time integration parameter. More...
 
const MaterialProperty< RankTwoTensor > & _original_local_config
 Rotational transformation from global to initial beam local coordinate system. More...
 
const MaterialProperty< Real > & _original_length
 Initial length of beam. More...
 
const unsigned int _component
 Direction along which residual is calculated. More...
 
RealVectorValue _vel_old_0
 Old translational and rotational velocities at the two nodes of the beam in the global coordinate system. More...
 
RealVectorValue _vel_old_1
 
RealVectorValue _rot_vel_old_0
 
RealVectorValue _rot_vel_old_1
 
RealVectorValue _vel_0
 Current translational and rotational velocities at the two nodes of the beam in the global coordinate system. More...
 
RealVectorValue _vel_1
 
RealVectorValue _rot_vel_0
 
RealVectorValue _rot_vel_1
 
RealVectorValue _accel_0
 Current translational and rotational accelerations at the two nodes of the beam in the global coordinate system. More...
 
RealVectorValue _accel_1
 
RealVectorValue _rot_accel_0
 
RealVectorValue _rot_accel_1
 
RealVectorValue _local_vel_old_0
 Old translational and rotational velocities at the two nodes of the beam in the initial beam local coordinate system. More...
 
RealVectorValue _local_vel_old_1
 
RealVectorValue _local_rot_vel_old_0
 
RealVectorValue _local_rot_vel_old_1
 
RealVectorValue _local_vel_0
 Current translational and rotational velocities at the two nodes of the beam in the initial beam local coordinate system. More...
 
RealVectorValue _local_vel_1
 
RealVectorValue _local_rot_vel_0
 
RealVectorValue _local_rot_vel_1
 
RealVectorValue _local_accel_0
 Current translational and rotational accelerations at the two nodes of the beam in the initial beam local coordinate system. More...
 
RealVectorValue _local_accel_1
 
RealVectorValue _local_rot_accel_0
 
RealVectorValue _local_rot_accel_1
 
std::vector< RealVectorValue > _local_force
 Forces and moments at the two end nodes of the beam in the initial beam local configuration. More...
 
std::vector< RealVectorValue > _local_moment
 
RealVectorValue _global_force_0
 Forces and moments at the two end nodes of the beam in the global coordinate system. More...
 
RealVectorValue _global_force_1
 
RealVectorValue _global_moment_0
 
RealVectorValue _global_moment_1
 
const VariableValue * _du_dot_du
 Coupled variable for du_dot_du calculated by time integrator. More...
 
const VariableValue * _du_dotdot_du
 Coupled variable for du_dotdot_du calculated by time integrator. More...
 

Detailed Description

Definition at line 20 of file InertialForceBeam.h.

Constructor & Destructor Documentation

◆ InertialForceBeam()

InertialForceBeam::InertialForceBeam ( const InputParameters &  parameters)

Definition at line 71 of file InertialForceBeam.C.

72  : TimeKernel(parameters),
73  _density(getMaterialProperty<Real>("density")),
74  _nrot(coupledComponents("rotations")),
75  _ndisp(coupledComponents("displacements")),
76  _rot_num(_nrot),
82  _area(coupledValue("area")),
83  _Ay(coupledValue("Ay")),
84  _Az(coupledValue("Az")),
85  _Ix(isParamValid("Ix") ? coupledValue("Ix") : _zero),
86  _Iy(coupledValue("Iy")),
87  _Iz(coupledValue("Iz")),
88  _beta(isParamValid("beta") ? getParam<Real>("beta") : 0.1),
89  _gamma(isParamValid("gamma") ? getParam<Real>("gamma") : 0.1),
90  _eta(getMaterialProperty<Real>("eta")),
91  _alpha(getParam<Real>("alpha")),
92  _original_local_config(getMaterialPropertyByName<RankTwoTensor>("initial_rotation")),
93  _original_length(getMaterialPropertyByName<Real>("original_length")),
94  _component(getParam<unsigned int>("component")),
95  _local_force(2),
96  _local_moment(2)
97 {
98  // Checking for consistency between the length of the provided rotations and displacements vector
99  if (_ndisp != _nrot)
100  mooseError("InertialForceBeam: The number of variables supplied in 'displacements' and "
101  "'rotations' must match.");
102 
103  if (isParamValid("beta") && isParamValid("gamma") && isParamValid("velocities") &&
104  isParamValid("accelerations") && isParamValid("rotational_velocities") &&
105  isParamValid("rotational_accelerations"))
106  {
107  if ((coupledComponents("velocities") != _ndisp) ||
108  (coupledComponents("accelerations") != _ndisp) ||
109  (coupledComponents("rotational_velocities") != _ndisp) ||
110  (coupledComponents("rotational_accelerations") != _ndisp))
111  mooseError(
112  "InertialForceBeam: The number of variables supplied in 'velocities', "
113  "'accelerations', 'rotational_velocities' and 'rotational_accelerations' must match "
114  "the number of displacement variables.");
115 
116  // fetch coupled velocities and accelerations
117  for (unsigned int i = 0; i < _ndisp; ++i)
118  {
119  MooseVariable * vel_variable = getVar("velocities", i);
120  _vel_num[i] = vel_variable->number();
121 
122  MooseVariable * accel_variable = getVar("accelerations", i);
123  _accel_num[i] = accel_variable->number();
124 
125  MooseVariable * rot_vel_variable = getVar("rotational_velocities", i);
126  _rot_vel_num[i] = rot_vel_variable->number();
127 
128  MooseVariable * rot_accel_variable = getVar("rotational_accelerations", i);
129  _rot_accel_num[i] = rot_accel_variable->number();
130  }
131  }
132  else if (!isParamValid("beta") && !isParamValid("gamma") && !isParamValid("velocities") &&
133  !isParamValid("accelerations") && !isParamValid("rotational_velocities") &&
134  !isParamValid("rotational_accelerations"))
135  {
136  _du_dot_du = &coupledDotDu("displacements", 0);
137  _du_dotdot_du = &coupledDotDotDu("displacements", 0);
138  }
139  else
140  mooseError("InertialForceBeam: Either all or none of `beta`, `gamma`, `velocities`, "
141  "`accelerations`, `rotational_velocities` and `rotational_accelerations` should be "
142  "provided as input.");
143 
144  // fetch coupled displacements and rotations
145  for (unsigned int i = 0; i < _ndisp; ++i)
146  {
147  MooseVariable * disp_variable = getVar("displacements", i);
148  _disp_num[i] = disp_variable->number();
149 
150  MooseVariable * rot_variable = getVar("rotations", i);
151  _rot_num[i] = rot_variable->number();
152  }
153 }
const VariableValue * _du_dot_du
Coupled variable for du_dot_du calculated by time integrator.
std::vector< RealVectorValue > _local_moment
unsigned int _nrot
Number of coupled rotational variables.
const Real _gamma
Newmark time integraion parameter.
std::vector< unsigned int > _accel_num
Variable numbers corresponding to acceleraion aux variables.
const Real _beta
Newmark time integration parameter.
const MaterialProperty< RankTwoTensor > & _original_local_config
Rotational transformation from global to initial beam local coordinate system.
unsigned int _ndisp
Number of coupled displacement variables.
const Real _alpha
HHT time integration parameter.
const VariableValue & _Iz
Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the...
const unsigned int _component
Direction along which residual is calculated.
std::vector< RealVectorValue > _local_force
Forces and moments at the two end nodes of the beam in the initial beam local configuration.
std::vector< unsigned int > _vel_num
Variable numbers corresponding to velocity aux variables.
std::vector< unsigned int > _rot_accel_num
Variable numbers corresponding to rotational acceleration aux variables.
std::vector< unsigned int > _rot_vel_num
Variable numbers corresponding to rotational velocity aux variables.
const VariableValue & _Ix
Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA ove...
const VariableValue & _Ay
Coupled variable for first moment of area of beam in y direction, i.e., integral of y*dA over the cro...
const VariableValue & _area
Coupled variable for beam cross-sectional area.
const MaterialProperty< Real > & _original_length
Initial length of beam.
const VariableValue & _Iy
Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the ...
const MaterialProperty< Real > & _density
Density of the beam.
const MaterialProperty< Real > & _eta
Mass proportional Rayleigh damping parameter.
std::vector< unsigned int > _disp_num
Variable numbers corresponding to displacement variables.
const VariableValue & _Az
Coupled variable for first moment of area of beam in z direction, i.e., integral of z*dA over the cro...
std::vector< unsigned int > _rot_num
Variable numbers corresponding to rotational variables.
const VariableValue * _du_dotdot_du
Coupled variable for du_dotdot_du calculated by time integrator.

Member Function Documentation

◆ computeJacobian()

void InertialForceBeam::computeJacobian ( )
overridevirtual

Definition at line 419 of file InertialForceBeam.C.

Referenced by computeOffDiagJacobian().

420 {
421  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
422  _local_ke.resize(ke.m(), ke.n());
423  _local_ke.zero();
424 
425  mooseAssert(_beta > 0.0, "InertialForceBeam: Beta parameter should be positive.");
426 
427  Real factor = 0.0;
428  if (isParamValid("beta"))
429  factor = 1.0 / (_beta * _dt * _dt) + _eta[0] * (1.0 + _alpha) * _gamma / _beta / _dt;
430  else
431  factor = (*_du_dotdot_du)[0] + _eta[0] * (1.0 + _alpha) * (*_du_dot_du)[0];
432 
433  for (unsigned int i = 0; i < _test.size(); ++i)
434  {
435  for (unsigned int j = 0; j < _phi.size(); ++j)
436  {
437  if (_component < 3)
438  _local_ke(i, j) = (i == j ? 1.0 / 3.0 : 1.0 / 6.0) * _density[0] * _area[0] *
439  _original_length[0] * factor;
440  else if (_component > 2)
441  {
442  RankTwoTensor I;
443  if (isParamValid("Ix"))
444  I(0, 0) = _Ix[0];
445  else
446  I(0, 0) = _Iy[0] + _Iz[0];
447  I(1, 1) = _Iz[0];
448  I(2, 2) = _Iy[0];
449 
450  // conversion from local config to global coordinate system
451  RankTwoTensor Ig = _original_local_config[0].transpose() * I * _original_local_config[0];
452 
453  _local_ke(i, j) = (i == j ? 1.0 / 3.0 : 1.0 / 6.0) * _density[0] *
454  Ig(_component - 3, _component - 3) * _original_length[0] * factor;
455  }
456  }
457  }
458 
459  ke += _local_ke;
460 
461  if (_has_diag_save_in)
462  {
463  unsigned int rows = ke.m();
464  DenseVector<Number> diag(rows);
465  for (unsigned int i = 0; i < rows; ++i)
466  diag(i) = _local_ke(i, i);
467 
468  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
469  for (unsigned int i = 0; i < _diag_save_in.size(); ++i)
470  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
471  }
472 }
const VariableValue * _du_dot_du
Coupled variable for du_dot_du calculated by time integrator.
const Real _gamma
Newmark time integraion parameter.
const Real _beta
Newmark time integration parameter.
const MaterialProperty< RankTwoTensor > & _original_local_config
Rotational transformation from global to initial beam local coordinate system.
const Real _alpha
HHT time integration parameter.
const VariableValue & _Iz
Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the...
const unsigned int _component
Direction along which residual is calculated.
const VariableValue & _Ix
Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA ove...
const VariableValue & _area
Coupled variable for beam cross-sectional area.
const MaterialProperty< Real > & _original_length
Initial length of beam.
const VariableValue & _Iy
Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the ...
const MaterialProperty< Real > & _density
Density of the beam.
const MaterialProperty< Real > & _eta
Mass proportional Rayleigh damping parameter.

◆ computeOffDiagJacobian()

void InertialForceBeam::computeOffDiagJacobian ( MooseVariableFEBase &  jvar)
overridevirtual

Definition at line 475 of file InertialForceBeam.C.

476 {
477  mooseAssert(_beta > 0.0, "InertialForceBeam: Beta parameter should be positive.");
478 
479  Real factor = 0.0;
480  if (isParamValid("beta"))
481  factor = 1.0 / (_beta * _dt * _dt) + _eta[0] * (1.0 + _alpha) * _gamma / _beta / _dt;
482  else
483  factor = (*_du_dotdot_du)[0] + _eta[0] * (1.0 + _alpha) * (*_du_dot_du)[0];
484 
485  size_t jvar_num = jvar.number();
486  if (jvar_num == _var.number())
487  computeJacobian();
488  else
489  {
490  unsigned int coupled_component = 0;
491  bool disp_coupled = false;
492  bool rot_coupled = false;
493 
494  for (unsigned int i = 0; i < _ndisp; ++i)
495  {
496  if (jvar_num == _disp_num[i] && _component > 2)
497  {
498  coupled_component = i;
499  disp_coupled = true;
500  break;
501  }
502  }
503 
504  for (unsigned int i = 0; i < _nrot; ++i)
505  {
506  if (jvar_num == _rot_num[i])
507  {
508  coupled_component = i + 3;
509  rot_coupled = true;
510  break;
511  }
512  }
513 
514  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_num);
515 
516  if (disp_coupled || rot_coupled)
517  {
518  for (unsigned int i = 0; i < _test.size(); ++i)
519  {
520  for (unsigned int j = 0; j < _phi.size(); ++j)
521  {
522  if (_component < 3 && coupled_component > 2)
523  {
524  RankTwoTensor A;
525  A(0, 1) = _Az[0];
526  A(0, 2) = -_Ay[0];
527  A(1, 0) = -_Az[0];
528  A(2, 0) = _Ay[0];
529 
530  // conversion from local config to global coordinate system
531  const RankTwoTensor Ag =
532  _original_local_config[0].transpose() * A * _original_local_config[0];
533 
534  ke(i, j) += (i == j ? 1.0 / 3.0 : 1.0 / 6.0) * _density[0] *
535  Ag(_component, coupled_component - 3) * _original_length[0] * factor;
536  }
537  else if (_component > 2 && coupled_component < 3)
538  {
539  RankTwoTensor A;
540  A(0, 1) = -_Az[0];
541  A(0, 2) = _Ay[0];
542  A(1, 0) = _Az[0];
543  A(2, 0) = -_Ay[0];
544 
545  // conversion from local config to global coordinate system
546  const RankTwoTensor Ag =
547  _original_local_config[0].transpose() * A * _original_local_config[0];
548 
549  ke(i, j) += (i == j ? 1.0 / 3.0 : 1.0 / 6.0) * _density[0] *
550  Ag(_component - 3, coupled_component) * _original_length[0] * factor;
551  }
552  else if (_component > 2 && coupled_component > 2)
553  {
554  RankTwoTensor I;
555  if (isParamValid("Ix"))
556  I(0, 0) = _Ix[0];
557  else
558  I(0, 0) = _Iy[0] + _Iz[0];
559  I(1, 1) = _Iz[0];
560  I(2, 2) = _Iy[0];
561 
562  // conversion from local config to global coordinate system
563  const RankTwoTensor Ig =
564  _original_local_config[0].transpose() * I * _original_local_config[0];
565 
566  ke(i, j) += (i == j ? 1.0 / 3.0 : 1.0 / 6.0) * _density[0] *
567  Ig(_component - 3, coupled_component - 3) * _original_length[0] * factor;
568  }
569  }
570  }
571  }
572  }
573 }
virtual void computeJacobian() override
const VariableValue * _du_dot_du
Coupled variable for du_dot_du calculated by time integrator.
unsigned int _nrot
Number of coupled rotational variables.
const Real _gamma
Newmark time integraion parameter.
const Real _beta
Newmark time integration parameter.
const MaterialProperty< RankTwoTensor > & _original_local_config
Rotational transformation from global to initial beam local coordinate system.
unsigned int _ndisp
Number of coupled displacement variables.
const Real _alpha
HHT time integration parameter.
const VariableValue & _Iz
Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the...
const unsigned int _component
Direction along which residual is calculated.
const VariableValue & _Ix
Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA ove...
const VariableValue & _Ay
Coupled variable for first moment of area of beam in y direction, i.e., integral of y*dA over the cro...
const MaterialProperty< Real > & _original_length
Initial length of beam.
const VariableValue & _Iy
Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the ...
const MaterialProperty< Real > & _density
Density of the beam.
const MaterialProperty< Real > & _eta
Mass proportional Rayleigh damping parameter.
std::vector< unsigned int > _disp_num
Variable numbers corresponding to displacement variables.
const VariableValue & _Az
Coupled variable for first moment of area of beam in z direction, i.e., integral of z*dA over the cro...
std::vector< unsigned int > _rot_num
Variable numbers corresponding to rotational variables.

◆ computeQpResidual()

virtual Real InertialForceBeam::computeQpResidual ( )
inlineoverrideprotectedvirtual

Definition at line 33 of file InertialForceBeam.h.

33 { return 0.0; };

◆ computeResidual()

void InertialForceBeam::computeResidual ( )
overridevirtual

Definition at line 156 of file InertialForceBeam.C.

157 {
158  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
159  mooseAssert(re.size() == 2, "Beam element only has two nodes.");
160  _local_re.resize(re.size());
161  _local_re.zero();
162 
163  if (_dt != 0.0)
164  {
165  // fetch the two end nodes for _current_elem
166  std::vector<Node *> node;
167  for (unsigned int i = 0; i < 2; ++i)
168  node.push_back(_current_elem->get_node(i));
169 
170  // Fetch the solution for the two end nodes at time t
171  NonlinearSystemBase & nonlinear_sys = _fe_problem.getNonlinearSystemBase();
172 
173  if (isParamValid("beta"))
174  {
175  const NumericVector<Number> & sol = *nonlinear_sys.currentSolution();
176  const NumericVector<Number> & sol_old = nonlinear_sys.solutionOld();
177 
178  AuxiliarySystem & aux = _fe_problem.getAuxiliarySystem();
179  const NumericVector<Number> & aux_sol_old = aux.solutionOld();
180 
181  mooseAssert(_beta > 0.0, "InertialForceBeam: Beta parameter should be positive.");
182  mooseAssert(_eta[0] >= 0.0, "InertialForceBeam: Eta parameter should be non-negative.");
183 
184  for (unsigned int i = 0; i < _ndisp; ++i)
185  {
186  // obtain delta displacement
187  unsigned int dof_index_0 = node[0]->dof_number(nonlinear_sys.number(), _disp_num[i], 0);
188  unsigned int dof_index_1 = node[1]->dof_number(nonlinear_sys.number(), _disp_num[i], 0);
189  const Real disp_0 = sol(dof_index_0) - sol_old(dof_index_0);
190  const Real disp_1 = sol(dof_index_1) - sol_old(dof_index_1);
191 
192  dof_index_0 = node[0]->dof_number(nonlinear_sys.number(), _rot_num[i], 0);
193  dof_index_1 = node[1]->dof_number(nonlinear_sys.number(), _rot_num[i], 0);
194  const Real rot_0 = sol(dof_index_0) - sol_old(dof_index_0);
195  const Real rot_1 = sol(dof_index_1) - sol_old(dof_index_1);
196 
197  // obtain new translational and rotational velocities and accelerations using newmark-beta
198  // time integration
199  _vel_old_0(i) = aux_sol_old(node[0]->dof_number(aux.number(), _vel_num[i], 0));
200  _vel_old_1(i) = aux_sol_old(node[1]->dof_number(aux.number(), _vel_num[i], 0));
201  const Real accel_old_0 = aux_sol_old(node[0]->dof_number(aux.number(), _accel_num[i], 0));
202  const Real accel_old_1 = aux_sol_old(node[1]->dof_number(aux.number(), _accel_num[i], 0));
203 
204  _rot_vel_old_0(i) = aux_sol_old(node[0]->dof_number(aux.number(), _rot_vel_num[i], 0));
205  _rot_vel_old_1(i) = aux_sol_old(node[1]->dof_number(aux.number(), _rot_vel_num[i], 0));
206  const Real rot_accel_old_0 =
207  aux_sol_old(node[0]->dof_number(aux.number(), _rot_accel_num[i], 0));
208  const Real rot_accel_old_1 =
209  aux_sol_old(node[1]->dof_number(aux.number(), _rot_accel_num[i], 0));
210 
211  _accel_0(i) =
212  1. / _beta * (disp_0 / (_dt * _dt) - _vel_old_0(i) / _dt - accel_old_0 * (0.5 - _beta));
213  _accel_1(i) =
214  1. / _beta * (disp_1 / (_dt * _dt) - _vel_old_1(i) / _dt - accel_old_1 * (0.5 - _beta));
215  _rot_accel_0(i) =
216  1. / _beta *
217  (rot_0 / (_dt * _dt) - _rot_vel_old_0(i) / _dt - rot_accel_old_0 * (0.5 - _beta));
218  _rot_accel_1(i) =
219  1. / _beta *
220  (rot_1 / (_dt * _dt) - _rot_vel_old_1(i) / _dt - rot_accel_old_1 * (0.5 - _beta));
221 
222  _vel_0(i) = _vel_old_0(i) + (_dt * (1 - _gamma)) * accel_old_0 + _gamma * _dt * _accel_0(i);
223  _vel_1(i) = _vel_old_1(i) + (_dt * (1 - _gamma)) * accel_old_1 + _gamma * _dt * _accel_1(i);
224  _rot_vel_0(i) = _rot_vel_old_0(i) + (_dt * (1 - _gamma)) * rot_accel_old_0 +
225  _gamma * _dt * _rot_accel_0(i);
226  _rot_vel_1(i) = _rot_vel_old_1(i) + (_dt * (1 - _gamma)) * rot_accel_old_1 +
227  _gamma * _dt * _rot_accel_1(i);
228  }
229  }
230  else
231  {
232  if (!nonlinear_sys.solutionUDot())
233  mooseError("InertialForceBeam: Time derivative of solution (`u_dot`) is not stored. Please "
234  "set uDotRequested() to true in FEProblemBase before requesting `u_dot`.");
235 
236  if (!nonlinear_sys.solutionUDotOld())
237  mooseError("InertialForceBeam: Old time derivative of solution (`u_dot_old`) is not "
238  "stored. Please set uDotOldRequested() to true in FEProblemBase before "
239  "requesting `u_dot_old`.");
240 
241  if (!nonlinear_sys.solutionUDotDot())
242  mooseError("InertialForceBeam: Second time derivative of solution (`u_dotdot`) is not "
243  "stored. Please set uDotDotRequested() to true in FEProblemBase before "
244  "requesting `u_dotdot`.");
245 
246  const NumericVector<Number> & vel = *nonlinear_sys.solutionUDot();
247  const NumericVector<Number> & vel_old = *nonlinear_sys.solutionUDotOld();
248  const NumericVector<Number> & accel = *nonlinear_sys.solutionUDotDot();
249 
250  for (unsigned int i = 0; i < _ndisp; ++i)
251  {
252  // translational velocities and accelerations
253  unsigned int dof_index_0 = node[0]->dof_number(nonlinear_sys.number(), _disp_num[i], 0);
254  unsigned int dof_index_1 = node[1]->dof_number(nonlinear_sys.number(), _disp_num[i], 0);
255  _vel_0(i) = vel(dof_index_0);
256  _vel_1(i) = vel(dof_index_1);
257  _vel_old_0(i) = vel_old(dof_index_0);
258  _vel_old_1(i) = vel_old(dof_index_1);
259  _accel_0(i) = accel(dof_index_0);
260  _accel_1(i) = accel(dof_index_1);
261 
262  // rotational velocities and accelerations
263  dof_index_0 = node[0]->dof_number(nonlinear_sys.number(), _rot_num[i], 0);
264  dof_index_1 = node[1]->dof_number(nonlinear_sys.number(), _rot_num[i], 0);
265  _rot_vel_0(i) = vel(dof_index_0);
266  _rot_vel_1(i) = vel(dof_index_1);
267  _rot_vel_old_0(i) = vel_old(dof_index_0);
268  _rot_vel_old_1(i) = vel_old(dof_index_1);
269  _rot_accel_0(i) = accel(dof_index_0);
270  _rot_accel_1(i) = accel(dof_index_1);
271  }
272  }
273 
274  // transform translational and rotational velocities and accelerations to the initial local
275  // configuration of the beam
282 
289 
290  // local residual
291  for (unsigned int i = 0; i < _ndisp; ++i)
292  {
293  if (_component < 3)
294  {
295  _local_force[0](i) = _density[0] * _area[0] * _original_length[0] / 3.0 *
296  (_local_accel_0(i) + _local_accel_1(i) / 2.0 +
297  _eta[0] * (1.0 + _alpha) * (_local_vel_0(i) + _local_vel_1(i) / 2.0) -
298  _alpha * _eta[0] * (_local_vel_old_0(i) + _local_vel_old_1(i) / 2.0));
299  _local_force[1](i) = _density[0] * _area[0] * _original_length[0] / 3.0 *
300  (_local_accel_1(i) + _local_accel_0(i) / 2.0 +
301  _eta[0] * (1.0 + _alpha) * (_local_vel_1(i) + _local_vel_0(i) / 2.0) -
302  _alpha * _eta[0] * (_local_vel_old_1(i) + _local_vel_old_0(i) / 2.0));
303  }
304 
305  if (_component > 2)
306  {
307  Real I = _Iy[0] + _Iz[0];
308  if (isParamValid("Ix") && (i == 0))
309  I = _Ix[0];
310  if (i == 1)
311  I = _Iz[0];
312  else if (i == 2)
313  I = _Iy[0];
314 
315  _local_moment[0](i) =
316  _density[0] * I * _original_length[0] / 3.0 *
317  (_local_rot_accel_0(i) + _local_rot_accel_1(i) / 2.0 +
318  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_0(i) + _local_rot_vel_1(i) / 2.0) -
319  _alpha * _eta[0] * (_local_rot_vel_old_0(i) + _local_rot_vel_old_1(i) / 2.0));
320  _local_moment[1](i) =
321  _density[0] * I * _original_length[0] / 3.0 *
322  (_local_rot_accel_1(i) + _local_rot_accel_0(i) / 2.0 +
323  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_1(i) + _local_rot_vel_0(i) / 2.0) -
324  _alpha * _eta[0] * (_local_rot_vel_old_1(i) + _local_rot_vel_old_0(i) / 2.0));
325  }
326  }
327 
328  // If Ay or Az are non-zero, contribution of rotational accelerations to translational forces
329  // and vice versa have to be added
330  if (_component < 3)
331  {
332  _local_force[0](0) +=
333  _density[0] * _original_length[0] / 3.0 *
334  (_Az[0] * (_local_rot_accel_0(1) + _local_rot_accel_1(1) / 2.0 +
335  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_0(1) + _local_rot_vel_1(1) / 2.0) -
336  _alpha * _eta[0] * (_local_rot_vel_old_0(1) + _local_rot_vel_old_1(1) / 2.0)) -
337  _Ay[0] * (_local_rot_accel_0(2) + _local_rot_accel_1(2) / 2.0 +
338  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_0(2) + _local_rot_vel_1(2) / 2.0) -
339  _alpha * _eta[0] * (_local_rot_vel_old_0(2) + _local_rot_vel_old_1(2) / 2.0)));
340  _local_force[1](0) +=
341  _density[0] * _original_length[0] / 3.0 *
342  (_Az[0] * (_local_rot_accel_1(1) + _local_rot_accel_0(1) / 2.0 +
343  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_1(1) + _local_rot_vel_0(1) / 2.0) -
344  _alpha * _eta[0] * (_local_rot_vel_old_1(1) + _local_rot_vel_old_0(1) / 2.0)) -
345  _Ay[0] * (_local_rot_accel_1(2) + _local_rot_accel_0(2) / 2.0 +
346  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_1(2) + _local_rot_vel_0(2) / 2.0) -
347  _alpha * _eta[0] * (_local_rot_vel_old_1(2) + _local_rot_vel_old_0(2) / 2.0)));
348 
349  _local_force[0](1) +=
350  -_density[0] * _original_length[0] / 3.0 * _Az[0] *
351  (_local_rot_accel_0(0) + _local_rot_accel_1(0) / 2.0 +
352  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_0(0) + _local_rot_vel_1(0) / 2.0) -
353  _alpha * _eta[0] * (_local_rot_vel_old_0(0) + _local_rot_vel_old_1(0) / 2.0));
354  _local_force[1](1) +=
355  -_density[0] * _original_length[0] / 3.0 * _Az[0] *
356  (_local_rot_accel_1(0) + _local_rot_accel_0(0) / 2.0 +
357  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_1(0) + _local_rot_vel_0(0) / 2.0) -
358  _alpha * _eta[0] * (_local_rot_vel_old_1(0) + _local_rot_vel_old_0(0) / 2.0));
359 
360  _local_force[0](2) +=
361  _density[0] * _original_length[0] / 3.0 * _Ay[0] *
362  (_local_rot_accel_0(0) + _local_rot_accel_1(0) / 2.0 +
363  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_0(0) + _local_rot_vel_1(0) / 2.0) -
364  _alpha * _eta[0] * (_local_rot_vel_old_0(0) + _local_rot_vel_old_1(0) / 2.0));
365  _local_force[1](2) +=
366  _density[0] * _original_length[0] / 3.0 * _Ay[0] *
367  (_local_rot_accel_1(0) + _local_rot_accel_0(0) / 2.0 +
368  _eta[0] * (1.0 + _alpha) * (_local_rot_vel_1(0) + _local_rot_vel_0(0) / 2.0) -
369  _alpha * _eta[0] * (_local_rot_vel_old_1(0) + _local_rot_vel_old_0(0) / 2.0));
370  }
371  else
372  {
373  _local_moment[0](0) += _density[0] * _original_length[0] / 3.0 *
374  (-_Az[0] * (_local_accel_0(1) + _local_accel_1(1) / 2.0) +
375  _Ay[0] * (_local_accel_0(1) + _local_accel_1(1) / 2.0));
376  _local_moment[1](0) += _density[0] * _original_length[0] / 3.0 *
377  (-_Az[0] * (_local_accel_1(1) + _local_accel_0(1) / 2.0) +
378  _Ay[0] * (_local_accel_1(1) + _local_accel_0(1) / 2.0));
379 
380  _local_moment[0](1) += _density[0] * _original_length[0] / 3.0 * _Az[0] *
381  (_local_accel_0(0) + _local_accel_1(0) / 2.0);
382  _local_moment[1](1) += _density[0] * _original_length[0] / 3.0 * _Az[0] *
383  (_local_accel_1(0) + _local_accel_0(0) / 2.0);
384 
385  _local_moment[0](2) += -_density[0] * _original_length[0] / 3.0 * _Ay[0] *
386  (_local_accel_0(0) + _local_accel_1(0) / 2.0);
387  _local_moment[1](2) += -_density[0] * _original_length[0] / 3.0 * _Ay[0] *
388  (_local_accel_1(0) + _local_accel_0(0) / 2.0);
389  }
390 
391  // Global force and moments
392  if (_component < 3)
393  {
396  _local_re(0) = _global_force_0(_component);
397  _local_re(1) = _global_force_1(_component);
398  }
399  else
400  {
403  _local_re(0) = _global_moment_0(_component - 3);
404  _local_re(1) = _global_moment_1(_component - 3);
405  }
406  }
407 
408  re += _local_re;
409 
410  if (_has_save_in)
411  {
412  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
413  for (unsigned int i = 0; i < _save_in.size(); ++i)
414  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
415  }
416 }
RealVectorValue _local_vel_old_0
Old translational and rotational velocities at the two nodes of the beam in the initial beam local co...
RealVectorValue _accel_1
RealVectorValue _rot_vel_0
RealVectorValue _local_rot_vel_old_1
RealVectorValue _local_rot_vel_1
std::vector< RealVectorValue > _local_moment
RealVectorValue _local_vel_0
Current translational and rotational velocities at the two nodes of the beam in the initial beam loca...
RealVectorValue _rot_accel_1
const Real _gamma
Newmark time integraion parameter.
std::vector< unsigned int > _accel_num
Variable numbers corresponding to acceleraion aux variables.
RealVectorValue _vel_0
Current translational and rotational velocities at the two nodes of the beam in the global coordinate...
const Real _beta
Newmark time integration parameter.
const MaterialProperty< RankTwoTensor > & _original_local_config
Rotational transformation from global to initial beam local coordinate system.
RealVectorValue _local_rot_vel_0
RealVectorValue _local_accel_1
unsigned int _ndisp
Number of coupled displacement variables.
const Real _alpha
HHT time integration parameter.
const VariableValue & _Iz
Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the...
RealVectorValue _local_rot_accel_0
RealVectorValue _local_rot_vel_old_0
const unsigned int _component
Direction along which residual is calculated.
std::vector< RealVectorValue > _local_force
Forces and moments at the two end nodes of the beam in the initial beam local configuration.
std::vector< unsigned int > _vel_num
Variable numbers corresponding to velocity aux variables.
std::vector< unsigned int > _rot_accel_num
Variable numbers corresponding to rotational acceleration aux variables.
std::vector< unsigned int > _rot_vel_num
Variable numbers corresponding to rotational velocity aux variables.
RealVectorValue _rot_vel_1
const VariableValue & _Ix
Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA ove...
const VariableValue & _Ay
Coupled variable for first moment of area of beam in y direction, i.e., integral of y*dA over the cro...
const VariableValue & _area
Coupled variable for beam cross-sectional area.
RealVectorValue _vel_old_1
const MaterialProperty< Real > & _original_length
Initial length of beam.
RealVectorValue _rot_vel_old_0
RealVectorValue _global_moment_1
const VariableValue & _Iy
Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the ...
RealVectorValue _local_vel_1
RealVectorValue _rot_vel_old_1
RealVectorValue _vel_old_0
Old translational and rotational velocities at the two nodes of the beam in the global coordinate sys...
RealVectorValue _global_force_0
Forces and moments at the two end nodes of the beam in the global coordinate system.
const MaterialProperty< Real > & _density
Density of the beam.
RealVectorValue _accel_0
Current translational and rotational accelerations at the two nodes of the beam in the global coordin...
RealVectorValue _global_moment_0
const MaterialProperty< Real > & _eta
Mass proportional Rayleigh damping parameter.
RealVectorValue _local_rot_accel_1
std::vector< unsigned int > _disp_num
Variable numbers corresponding to displacement variables.
const VariableValue & _Az
Coupled variable for first moment of area of beam in z direction, i.e., integral of z*dA over the cro...
RealVectorValue _local_vel_old_1
std::vector< unsigned int > _rot_num
Variable numbers corresponding to rotational variables.
RealVectorValue _vel_1
RealVectorValue _local_accel_0
Current translational and rotational accelerations at the two nodes of the beam in the initial beam l...
RealVectorValue _global_force_1
RealVectorValue _rot_accel_0

Member Data Documentation

◆ _accel_0

RealVectorValue InertialForceBeam::_accel_0
private

Current translational and rotational accelerations at the two nodes of the beam in the global coordinate system.

Definition at line 136 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _accel_1

RealVectorValue InertialForceBeam::_accel_1
private

Definition at line 136 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _accel_num

std::vector<unsigned int> InertialForceBeam::_accel_num
private

Variable numbers corresponding to acceleraion aux variables.

Definition at line 55 of file InertialForceBeam.h.

Referenced by computeResidual(), and InertialForceBeam().

◆ _alpha

const Real InertialForceBeam::_alpha
private

HHT time integration parameter.

Definition at line 106 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _area

const VariableValue& InertialForceBeam::_area
private

Coupled variable for beam cross-sectional area.

Definition at line 64 of file InertialForceBeam.h.

Referenced by computeJacobian(), and computeResidual().

◆ _Ay

const VariableValue& InertialForceBeam::_Ay
private

Coupled variable for first moment of area of beam in y direction, i.e., integral of y*dA over the cross-section.

Definition at line 70 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), and computeResidual().

◆ _Az

const VariableValue& InertialForceBeam::_Az
private

Coupled variable for first moment of area of beam in z direction, i.e., integral of z*dA over the cross-section.

Definition at line 76 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), and computeResidual().

◆ _beta

const Real InertialForceBeam::_beta
private

Newmark time integration parameter.

Definition at line 97 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _component

const unsigned int InertialForceBeam::_component
private

Direction along which residual is calculated.

Definition at line 118 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _density

const MaterialProperty<Real>& InertialForceBeam::_density
private

Density of the beam.

Definition at line 33 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _disp_num

std::vector<unsigned int> InertialForceBeam::_disp_num
private

Variable numbers corresponding to displacement variables.

Definition at line 49 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), computeResidual(), and InertialForceBeam().

◆ _du_dot_du

const VariableValue* InertialForceBeam::_du_dot_du
private

Coupled variable for du_dot_du calculated by time integrator.

Definition at line 171 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and InertialForceBeam().

◆ _du_dotdot_du

const VariableValue* InertialForceBeam::_du_dotdot_du
private

Coupled variable for du_dotdot_du calculated by time integrator.

Definition at line 176 of file InertialForceBeam.h.

Referenced by InertialForceBeam().

◆ _eta

const MaterialProperty<Real>& InertialForceBeam::_eta
private

Mass proportional Rayleigh damping parameter.

Definition at line 103 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _gamma

const Real InertialForceBeam::_gamma
private

Newmark time integraion parameter.

Definition at line 100 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _global_force_0

RealVectorValue InertialForceBeam::_global_force_0
private

Forces and moments at the two end nodes of the beam in the global coordinate system.

Definition at line 166 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _global_force_1

RealVectorValue InertialForceBeam::_global_force_1
private

Definition at line 166 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _global_moment_0

RealVectorValue InertialForceBeam::_global_moment_0
private

Definition at line 166 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _global_moment_1

RealVectorValue InertialForceBeam::_global_moment_1
private

Definition at line 166 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _Ix

const VariableValue& InertialForceBeam::_Ix
private

Coupled variable for second moment of area of beam in x direction, i.e., integral of (y^2+z^2)*dA over the cross-section.

Definition at line 82 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _Iy

const VariableValue& InertialForceBeam::_Iy
private

Coupled variable for second moment of area of beam in y direction, i.e., integral of y^2*dA over the cross-section.

Definition at line 88 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _Iz

const VariableValue& InertialForceBeam::_Iz
private

Coupled variable for second momemnt of area of beam in z direction, i.e., integral of z^2*dA over the cross-section.

Definition at line 94 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _local_accel_0

RealVectorValue InertialForceBeam::_local_accel_0
private

Current translational and rotational accelerations at the two nodes of the beam in the initial beam local coordinate system.

Definition at line 154 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_accel_1

RealVectorValue InertialForceBeam::_local_accel_1
private

Definition at line 154 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_force

std::vector<RealVectorValue> InertialForceBeam::_local_force
private

Forces and moments at the two end nodes of the beam in the initial beam local configuration.

Definition at line 160 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_moment

std::vector<RealVectorValue> InertialForceBeam::_local_moment
private

Definition at line 160 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_accel_0

RealVectorValue InertialForceBeam::_local_rot_accel_0
private

Definition at line 154 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_accel_1

RealVectorValue InertialForceBeam::_local_rot_accel_1
private

Definition at line 154 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_vel_0

RealVectorValue InertialForceBeam::_local_rot_vel_0
private

Definition at line 148 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_vel_1

RealVectorValue InertialForceBeam::_local_rot_vel_1
private

Definition at line 148 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_vel_old_0

RealVectorValue InertialForceBeam::_local_rot_vel_old_0
private

Definition at line 142 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_rot_vel_old_1

RealVectorValue InertialForceBeam::_local_rot_vel_old_1
private

Definition at line 142 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_vel_0

RealVectorValue InertialForceBeam::_local_vel_0
private

Current translational and rotational velocities at the two nodes of the beam in the initial beam local coordinate system.

Definition at line 148 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_vel_1

RealVectorValue InertialForceBeam::_local_vel_1
private

Definition at line 148 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_vel_old_0

RealVectorValue InertialForceBeam::_local_vel_old_0
private

Old translational and rotational velocities at the two nodes of the beam in the initial beam local coordinate system.

Definition at line 142 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _local_vel_old_1

RealVectorValue InertialForceBeam::_local_vel_old_1
private

Definition at line 142 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _ndisp

unsigned int InertialForceBeam::_ndisp
private

Number of coupled displacement variables.

Definition at line 43 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), computeResidual(), and InertialForceBeam().

◆ _nrot

unsigned int InertialForceBeam::_nrot
private

Number of coupled rotational variables.

Definition at line 40 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), and InertialForceBeam().

◆ _original_length

const MaterialProperty<Real>& InertialForceBeam::_original_length
private

Initial length of beam.

Definition at line 115 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _original_local_config

const MaterialProperty<RankTwoTensor>& InertialForceBeam::_original_local_config
private

Rotational transformation from global to initial beam local coordinate system.

Definition at line 112 of file InertialForceBeam.h.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

◆ _rot_accel_0

RealVectorValue InertialForceBeam::_rot_accel_0
private

Definition at line 136 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _rot_accel_1

RealVectorValue InertialForceBeam::_rot_accel_1
private

Definition at line 136 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _rot_accel_num

std::vector<unsigned int> InertialForceBeam::_rot_accel_num
private

Variable numbers corresponding to rotational acceleration aux variables.

Definition at line 61 of file InertialForceBeam.h.

Referenced by computeResidual(), and InertialForceBeam().

◆ _rot_num

std::vector<unsigned int> InertialForceBeam::_rot_num
private

Variable numbers corresponding to rotational variables.

Definition at line 46 of file InertialForceBeam.h.

Referenced by computeOffDiagJacobian(), computeResidual(), and InertialForceBeam().

◆ _rot_vel_0

RealVectorValue InertialForceBeam::_rot_vel_0
private

Definition at line 130 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _rot_vel_1

RealVectorValue InertialForceBeam::_rot_vel_1
private

Definition at line 130 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _rot_vel_num

std::vector<unsigned int> InertialForceBeam::_rot_vel_num
private

Variable numbers corresponding to rotational velocity aux variables.

Definition at line 58 of file InertialForceBeam.h.

Referenced by computeResidual(), and InertialForceBeam().

◆ _rot_vel_old_0

RealVectorValue InertialForceBeam::_rot_vel_old_0
private

Definition at line 124 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _rot_vel_old_1

RealVectorValue InertialForceBeam::_rot_vel_old_1
private

Definition at line 124 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _vel_0

RealVectorValue InertialForceBeam::_vel_0
private

Current translational and rotational velocities at the two nodes of the beam in the global coordinate system.

Definition at line 130 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _vel_1

RealVectorValue InertialForceBeam::_vel_1
private

Definition at line 130 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _vel_num

std::vector<unsigned int> InertialForceBeam::_vel_num
private

Variable numbers corresponding to velocity aux variables.

Definition at line 52 of file InertialForceBeam.h.

Referenced by computeResidual(), and InertialForceBeam().

◆ _vel_old_0

RealVectorValue InertialForceBeam::_vel_old_0
private

Old translational and rotational velocities at the two nodes of the beam in the global coordinate system.

Definition at line 124 of file InertialForceBeam.h.

Referenced by computeResidual().

◆ _vel_old_1

RealVectorValue InertialForceBeam::_vel_old_1
private

Definition at line 124 of file InertialForceBeam.h.

Referenced by computeResidual().


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