www.mooseframework.org
Assembly.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "MooseArray.h"
13 #include "MooseTypes.h"
14 
15 #include "libmesh/dense_matrix.h"
16 #include "libmesh/dense_vector.h"
17 #include "libmesh/enum_quadrature_type.h"
18 #include "libmesh/fe_type.h"
19 #include "libmesh/point.h"
20 
21 #include "metaphysicl/numberarray.h"
22 #include "metaphysicl/dualnumber.h"
23 
24 // libMesh forward declarations
25 namespace libMesh
26 {
27 class DofMap;
28 class CouplingMatrix;
29 class Elem;
30 template <typename>
33 template <typename T>
37 class Node;
38 template <typename T>
40 template <typename T>
42 }
43 
44 // MOOSE Forward Declares
45 class MooseMesh;
47 class SystemBase;
49 class MooseVariableBase;
50 template <typename>
51 class MooseVariableFE;
55 class XFEMInterface;
56 class SubProblem;
57 
62 class Assembly
63 {
64 public:
65  Assembly(SystemBase & sys, THREAD_ID tid);
66  virtual ~Assembly();
67 
74  const FEBase * const & getFE(FEType type, unsigned int dim) const
75  {
76  buildFE(type);
77  return _const_fe[dim][type];
78  }
79 
86  const FEBase * const & getFENeighbor(FEType type, unsigned int dim) const
87  {
89  return _const_fe_neighbor[dim][type];
90  }
91 
98  const FEBase * const & getFEFace(FEType type, unsigned int dim) const
99  {
100  buildFaceFE(type);
101  return _const_fe_face[dim][type];
102  }
103 
110  const FEBase * const & getFEFaceNeighbor(FEType type, unsigned int dim) const
111  {
113  return _const_fe_face_neighbor[dim][type];
114  }
115 
122  const FEVectorBase * const & getVectorFE(FEType type, unsigned int dim) const
123  {
125  return _const_vector_fe[dim][type];
126  }
127 
134  const FEVectorBase * const & getVectorFENeighbor(FEType type, unsigned int dim) const
135  {
137  return _const_vector_fe_neighbor[dim][type];
138  }
139 
146  const FEVectorBase * const & getVectorFEFace(FEType type, unsigned int dim) const
147  {
149  return _const_vector_fe_face[dim][type];
150  }
151 
158  const FEVectorBase * const & getVectorFEFaceNeighbor(FEType type, unsigned int dim) const
159  {
162  }
163 
168  const QBase * const & qRule() const { return _const_current_qrule; }
169 
174  QBase * const & writeableQRule() { return _current_qrule; }
175 
180  const MooseArray<Point> & qPoints() const { return _current_q_points; }
181 
187 
192  const MooseArray<Real> & JxW() const { return _current_JxW; }
193 
194  template <ComputeStage compute_stage>
195  const MooseArray<ADReal> & adJxW() const
196  {
197  return _ad_JxW;
198  }
199 
200  template <ComputeStage compute_stage>
202  {
203  return _current_JxW_face;
204  }
205 
206  template <ComputeStage compute_stage>
207  const MooseArray<ADReal> & adCurvatures() const;
208 
213  const MooseArray<Real> & coordTransformation() const { return _coord; }
214 
219  template <ComputeStage compute_stage>
221  {
222  return _coord;
223  }
224 
230 
235  const QBase * const & qRuleFace() const { return _const_current_qrule_face; }
236 
241  QBase * const & writeableQRuleFace() { return _current_qrule_face; }
242 
248 
253  const MooseArray<Real> & JxWFace() const { return _current_JxW_face; }
254 
259  const MooseArray<Point> & normals() const { return _current_normals; }
260 
261  template <ComputeStage compute_stage>
262  const ADPoint & adNormals() const
263  {
264  return _current_normals;
265  }
266 
267  template <ComputeStage compute_stage>
268  const ADPoint & adQPoints() const
269  {
270  return _current_q_points;
271  }
272 
273  template <ComputeStage compute_stage>
274  const ADPoint & adQPointsFace() const
275  {
276  return _current_q_points_face;
277  }
278 
283  const Elem * const & elem() const { return _current_elem; }
284 
289 
294 
299  const Real & elemVolume() { return _current_elem_volume; }
300 
305  const unsigned int & side() const { return _current_side; }
306 
311  const unsigned int & neighborSide() const { return _current_neighbor_side; }
312 
317  const Elem *& sideElem() { return _current_side_elem; }
318 
323  const Real & sideElemVolume() { return _current_side_volume; }
324 
329  const Elem * const & neighbor() const { return _current_neighbor_elem; }
330 
335  const Elem * const & lowerDElem() const { return _current_lower_d_elem; }
336 
341 
346 
351  const Real & neighborVolume() const;
352 
357  const QBase * const & qRuleNeighbor() const { return _const_current_qrule_neighbor; }
358 
363  QBase * const & writeableQRuleNeighbor() { return _current_qrule_neighbor; }
364 
370 
376 
381  const Node * const & node() const { return _current_node; }
382 
387  const Node * const & nodeNeighbor() const { return _current_neighbor_node; }
388 
392  void createQRules(QuadratureType type, Order order, Order volume_order, Order face_order);
393 
402  void setVolumeQRule(QBase * qrule, unsigned int dim);
403 
412  void setFaceQRule(QBase * qrule, unsigned int dim);
413 
422  void setNeighborQRule(QBase * qrule, unsigned int dim);
423 
429  void reinit(const Elem * elem);
430 
435  void reinitElemFaceRef(const Elem * elem,
436  unsigned int elem_side,
437  Real tolerance,
438  const std::vector<Point> * const pts = nullptr,
439  const std::vector<Real> * const weights = nullptr);
440 
445  void reinitNeighborFaceRef(const Elem * neighbor_elem,
446  unsigned int neighbor_side,
447  Real tolerance,
448  const std::vector<Point> * const pts,
449  const std::vector<Real> * const weights = nullptr);
450 
454  void reinitLowerDElemRef(const Elem * elem,
455  const std::vector<Point> * const pts,
456  const std::vector<Real> * const weights = nullptr);
457 
461  void reinitMortarElem(const Elem * elem);
462 
466  const std::vector<Real> & jxWMortar() const { return *_JxW_msm; }
467 
471  const QBase * const & qRuleMortar() const { return _const_qrule_msm; }
472 
473 private:
477  void computeADFace(const Elem * elem, unsigned int side);
478 
479 public:
483  void reinitAtPhysical(const Elem * elem, const std::vector<Point> & physical_points);
484 
488  void reinit(const Elem * elem, const std::vector<Point> & reference_points);
489 
493  void reinit(const Elem * elem, unsigned int side);
494 
498  void reinit(const Elem * elem, unsigned int side, const std::vector<Point> & reference_points);
499 
508  void reinitElemAndNeighbor(const Elem * elem,
509  unsigned int side,
510  const Elem * neighbor,
511  unsigned int neighbor_side);
512 
516  void reinitNeighborAtPhysical(const Elem * neighbor,
517  unsigned int neighbor_side,
518  const std::vector<Point> & physical_points);
519 
523  void reinitNeighborAtPhysical(const Elem * neighbor, const std::vector<Point> & physical_points);
524 
525  void reinitNeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
526 
530  void reinit(const Node * node);
531 
535  void init(const CouplingMatrix * cm);
536 
538  void init();
539 
541  void initNonlocalCoupling();
542 
544  void prepareJacobianBlock();
545 
547  void prepareResidual();
548 
549  void prepare();
550  void prepareNonlocal();
551 
559  void prepareNeighbor();
560 
565  void prepareLowerD();
566 
567  void prepareBlock(unsigned int ivar, unsigned jvar, const std::vector<dof_id_type> & dof_indices);
568  void prepareBlockNonlocal(unsigned int ivar,
569  unsigned jvar,
570  const std::vector<dof_id_type> & idof_indices,
571  const std::vector<dof_id_type> & jdof_indices);
572  void prepareScalar();
573  void prepareOffDiagScalar();
574 
575  template <typename T>
576  void copyShapes(MooseVariableFE<T> & v);
577  void copyShapes(unsigned int var);
578 
579  template <typename T>
581  void copyFaceShapes(unsigned int var);
582 
583  template <typename T>
585  void copyNeighborShapes(unsigned int var);
586 
587  void addResidual(NumericVector<Number> & residual, TagID tag_id = 0);
588  void addResidual(const std::map<TagName, TagID> & tags);
589  void addResidualNeighbor(NumericVector<Number> & residual, TagID tag_id = 0);
590  void addResidualNeighbor(const std::map<TagName, TagID> & tags);
591  void addResidualScalar(TagID tag_id);
592  void addResidualScalar(const std::map<TagName, TagID> & tags);
593 
597  void cacheResidual();
598 
607  void cacheResidualContribution(dof_id_type dof, Real value, TagID tag_id);
608 
617  void cacheResidualContribution(dof_id_type dof, Real value, const std::set<TagID> & tags);
618 
622  void cacheResidualNodes(const DenseVector<Number> & res,
623  const std::vector<dof_id_type> & dof_index,
624  TagID tag = 0);
625 
629  void cacheResidualNeighbor();
630 
634  void cacheResidualLower();
635 
636  void addCachedResiduals();
637 
644  void addCachedResidual(NumericVector<Number> & residual, TagID tag_id);
645 
646  void setResidual(NumericVector<Number> & residual, TagID tag_id = 0);
647  void setResidualNeighbor(NumericVector<Number> & residual, TagID tag_id = 0);
648 
649  void addJacobian();
654  void addJacobianNonlocal();
655  void addJacobianBlock(SparseMatrix<Number> & jacobian,
656  unsigned int ivar,
657  unsigned int jvar,
658  const DofMap & dof_map,
659  std::vector<dof_id_type> & dof_indices);
660  void addJacobianBlockNonlocal(SparseMatrix<Number> & jacobian,
661  unsigned int ivar,
662  unsigned int jvar,
663  const DofMap & dof_map,
664  const std::vector<dof_id_type> & idof_indices,
665  const std::vector<dof_id_type> & jdof_indices);
666 
671  void addJacobianNeighbor();
672 
678  void addJacobianLower();
679 
680  void addJacobianNeighbor(SparseMatrix<Number> & jacobian,
681  unsigned int ivar,
682  unsigned int jvar,
683  const DofMap & dof_map,
684  std::vector<dof_id_type> & dof_indices,
685  std::vector<dof_id_type> & neighbor_dof_indices);
686  void addJacobianScalar();
687  void addJacobianOffDiagScalar(unsigned int ivar);
688 
692  void cacheJacobian();
693 
698 
702  void cacheJacobianNonlocal();
703 
708  void cacheJacobianNeighbor();
709 
716  void addCachedJacobian(SparseMatrix<Number> & jacobian);
717 
718  void addCachedJacobian();
719 
720  DenseVector<Number> & residualBlock(unsigned int var_num, TagID tag_id = 0)
721  {
722  return _sub_Re[tag_id][var_num];
723  }
724 
725  DenseVector<Number> & residualBlockNeighbor(unsigned int var_num, TagID tag_id = 0)
726  {
727  return _sub_Rn[tag_id][var_num];
728  }
729 
730  DenseVector<Number> & residualBlockLower(unsigned int var_num, TagID tag_id = 0)
731  {
732  return _sub_Rl[tag_id][var_num];
733  }
734 
735  DenseMatrix<Number> & jacobianBlock(unsigned int ivar, unsigned int jvar, TagID tag = 0);
736 
737  DenseMatrix<Number> & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar, TagID tag = 0);
738  DenseMatrix<Number> & jacobianBlockNeighbor(Moose::DGJacobianType type,
739  unsigned int ivar,
740  unsigned int jvar,
741  TagID tag = 0);
742 
747  unsigned int ivar,
748  unsigned int jvar,
749  TagID tag = 0);
750 
751  void cacheJacobianBlock(DenseMatrix<Number> & jac_block,
752  const std::vector<dof_id_type> & idof_indices,
753  const std::vector<dof_id_type> & jdof_indices,
754  Real scaling_factor,
755  TagID tag = 0);
756  void cacheJacobianBlockNonlocal(DenseMatrix<Number> & jac_block,
757  const std::vector<dof_id_type> & idof_indices,
758  const std::vector<dof_id_type> & jdof_indices,
759  Real scaling_factor,
760  TagID tag = 0);
761 
762  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & couplingEntries()
763  {
764  return _cm_ff_entry;
765  }
766  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & nonlocalCouplingEntries()
767  {
768  return _cm_nonlocal_entry;
769  }
770 
771  // Read-only references
772  const VariablePhiValue & phi() const { return _phi; }
773  template <typename T, ComputeStage compute_stage>
775  adGradPhi(const MooseVariableFE<T> & v) const
776  {
777  return gradPhi(v);
778  }
779  const VariablePhiValue & phi(const MooseVariable &) const { return _phi; }
780  const VariablePhiGradient & gradPhi() const { return _grad_phi; }
781  const VariablePhiGradient & gradPhi(const MooseVariable &) const { return _grad_phi; }
782  const VariablePhiSecond & secondPhi() const { return _second_phi; }
783  const VariablePhiSecond & secondPhi(const MooseVariable &) const { return _second_phi; }
784 
785  const VariablePhiValue & phiFace() const { return _phi_face; }
786  const VariablePhiValue & phiFace(const MooseVariable &) const { return _phi_face; }
787  const VariablePhiGradient & gradPhiFace() const { return _grad_phi_face; }
788  const VariablePhiGradient & gradPhiFace(const MooseVariable &) const { return _grad_phi_face; }
790 
791  const VariablePhiValue & phiNeighbor(const MooseVariable &) const { return _phi_neighbor; }
793  {
794  return _grad_phi_neighbor;
795  }
797  {
798  return _second_phi_neighbor;
799  }
800 
802  {
803  return _phi_face_neighbor;
804  }
806  {
808  }
810  {
812  }
813 
814  const VectorVariablePhiValue & phi(const VectorMooseVariable &) const { return _vector_phi; }
816  {
817  return _vector_grad_phi;
818  }
820  {
821  return _vector_second_phi;
822  }
824  {
825  return _vector_curl_phi;
826  }
827 
829  {
830  return _vector_phi_face;
831  }
833  {
834  return _vector_grad_phi_face;
835  }
837  {
839  }
841  {
842  return _vector_curl_phi_face;
843  }
844 
846  {
847  return _vector_phi_neighbor;
848  }
850  {
852  }
854  {
856  }
858  {
860  }
861 
863  {
865  }
867  {
869  }
871  {
873  }
875  {
877  }
878 
879  // Writeable references
880  VariablePhiValue & phi(const MooseVariable &) { return _phi; }
883 
887 
891 
894  {
896  }
898  {
900  }
901 
906 
909  {
910  return _vector_grad_phi_face;
911  }
913  {
915  }
917 
920  {
922  }
924  {
926  }
928  {
930  }
931 
933  {
935  }
937  {
939  }
941  {
943  }
945  {
947  }
948 
949  template <typename OutputType>
950  const typename OutputTools<OutputType>::VariablePhiValue & fePhi(FEType type) const
951  {
952  buildFE(type);
953  return _fe_shape_data[type]->_phi;
954  }
955 
956  template <typename OutputType>
958  {
959  buildFE(type);
960  return _fe_shape_data[type]->_grad_phi;
961  }
962 
963  template <typename OutputType>
965  feADGradPhi(FEType type) const
966  {
967  return _ad_grad_phi_data[type];
968  }
969 
970  template <typename OutputType>
972  {
974  buildFE(type);
975  return _fe_shape_data[type]->_second_phi;
976  }
977 
978  template <typename OutputType>
979  const typename OutputTools<OutputType>::VariablePhiValue & fePhiLower(FEType type) const;
980 
981  template <typename OutputType>
983 
984  template <typename OutputType>
986  {
987  buildFaceFE(type);
988  return _fe_shape_data_face[type]->_phi;
989  }
990 
991  template <typename OutputType>
993  {
994  buildFaceFE(type);
995  return _fe_shape_data_face[type]->_grad_phi;
996  }
997 
998  template <typename OutputType>
1000  feADGradPhiFace(FEType type) const
1001  {
1002  return _ad_grad_phi_data_face[type];
1003  }
1004 
1005  template <typename OutputType>
1007  {
1008  _need_second_derivative[type] = true;
1009  buildFaceFE(type);
1010  return _fe_shape_data_face[type]->_second_phi;
1011  }
1012 
1013  template <typename OutputType>
1015  {
1017  return _fe_shape_data_neighbor[type]->_phi;
1018  }
1019 
1020  template <typename OutputType>
1022  {
1024  return _fe_shape_data_neighbor[type]->_grad_phi;
1025  }
1026 
1027  template <typename OutputType>
1029  {
1032  return _fe_shape_data_neighbor[type]->_second_phi;
1033  }
1034 
1035  template <typename OutputType>
1037  {
1039  return _fe_shape_data_face_neighbor[type]->_phi;
1040  }
1041 
1042  template <typename OutputType>
1045  {
1047  return _fe_shape_data_face_neighbor[type]->_grad_phi;
1048  }
1049 
1050  template <typename OutputType>
1053  {
1056  return _fe_shape_data_face_neighbor[type]->_second_phi;
1057  }
1058 
1059  template <typename OutputType>
1061  {
1062  _need_curl[type] = true;
1063  buildFE(type);
1064  return _fe_shape_data[type]->_curl_phi;
1065  }
1066 
1067  template <typename OutputType>
1069  {
1070  _need_curl[type] = true;
1071  buildFaceFE(type);
1072  return _fe_shape_data_face[type]->_curl_phi;
1073  }
1074 
1075  template <typename OutputType>
1077  {
1078  _need_curl[type] = true;
1080  return _fe_shape_data_neighbor[type]->_curl_phi;
1081  }
1082 
1083  template <typename OutputType>
1085  {
1086  _need_curl[type] = true;
1088  return _fe_shape_data_face_neighbor[type]->_curl_phi;
1089  }
1090 
1099  void
1100  cacheJacobianContribution(numeric_index_type i, numeric_index_type j, Real value, TagID tag = 0);
1101 
1102  void cacheJacobianContribution(numeric_index_type i,
1103  numeric_index_type j,
1104  Real value,
1105  const std::set<TagID> & tags);
1106 
1111 
1116 
1121 
1125  void setXFEM(std::shared_ptr<XFEMInterface> xfem) { _xfem = xfem; }
1126 
1127  void assignDisplacements(std::vector<unsigned> && disp_numbers) { _displacements = disp_numbers; }
1128 
1129 protected:
1135  void reinitFE(const Elem * elem);
1136 
1143  void reinitFEFace(const Elem * elem, unsigned int side);
1144 
1145  void computeFaceMap(unsigned dim, const std::vector<Real> & qw, const Elem * side);
1146 
1147  void reinitFEFaceNeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
1148 
1149  void reinitFENeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
1150 
1151  template <ComputeStage compute_stage>
1152  void setCoordinateTransformation(const QBase * qrule,
1153  const ADPoint & q_points,
1154  MooseArray<ADReal> & coord);
1155 
1156  void computeCurrentElemVolume();
1157 
1158  void computeCurrentFaceVolume();
1159 
1161 
1162  void addResidualBlock(NumericVector<Number> & residual,
1163  DenseVector<Number> & res_block,
1164  const std::vector<dof_id_type> & dof_indices,
1165  Real scaling_factor);
1166  void cacheResidualBlock(std::vector<Real> & cached_residual_values,
1167  std::vector<dof_id_type> & cached_residual_rows,
1168  DenseVector<Number> & res_block,
1169  const std::vector<dof_id_type> & dof_indices,
1170  Real scaling_factor);
1171 
1172  void setResidualBlock(NumericVector<Number> & residual,
1173  DenseVector<Number> & res_block,
1174  const std::vector<dof_id_type> & dof_indices,
1175  Real scaling_factor);
1176 
1177  void addJacobianBlock(SparseMatrix<Number> & jacobian,
1178  DenseMatrix<Number> & jac_block,
1179  const std::vector<dof_id_type> & idof_indices,
1180  const std::vector<dof_id_type> & jdof_indices,
1181  Real scaling_factor);
1182 
1190 
1196  void modifyWeightsDueToXFEM(const Elem * elem);
1197 
1204  void modifyFaceWeightsDueToXFEM(const Elem * elem, unsigned int side = 0);
1205 
1206  template <typename OutputType>
1207  void computeGradPhiAD(
1208  const Elem * elem,
1209  unsigned int n_qp,
1211  FEGenericBase<OutputType> * fe);
1212  void resizeMappingObjects(unsigned int n_qp, unsigned int dim);
1213  void computeAffineMapAD(const Elem * elem,
1214  const std::vector<Real> & qw,
1215  unsigned int n_qp,
1216  FEBase * fe);
1217  void
1218  computeSinglePointMapAD(const Elem * elem, const std::vector<Real> & qw, unsigned p, FEBase * fe);
1219 
1220 private:
1225  void buildFE(FEType type) const;
1226 
1231  void buildFaceFE(FEType type) const;
1232 
1237  void buildNeighborFE(FEType type) const;
1238 
1243  void buildFaceNeighborFE(FEType type) const;
1244 
1249  void buildLowerDFE(FEType type) const;
1250 
1255  void buildVectorFE(FEType type) const;
1256 
1261  void buildVectorFaceFE(FEType type) const;
1262 
1267  void buildVectorNeighborFE(FEType type) const;
1268 
1273  void buildVectorFaceNeighborFE(FEType type) const;
1274 
1279  void buildVectorLowerDFE(FEType type) const;
1280 
1281 private:
1284 
1285  const bool _displaced;
1286 
1288  const CouplingMatrix * _cm;
1289  const CouplingMatrix & _nonlocal_cm;
1290 
1291  const bool & _computing_jacobian;
1292 
1294  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> _cm_ff_entry;
1296  std::vector<std::pair<MooseVariableFEBase *, MooseVariableScalar *>> _cm_fs_entry;
1298  std::vector<std::pair<MooseVariableScalar *, MooseVariableFEBase *>> _cm_sf_entry;
1300  std::vector<std::pair<MooseVariableScalar *, MooseVariableScalar *>> _cm_ss_entry;
1302  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> _cm_nonlocal_entry;
1304  std::vector<std::vector<std::vector<unsigned char>>> _jacobian_block_used;
1305  std::vector<std::vector<std::vector<unsigned char>>> _jacobian_block_nonlocal_used;
1307  std::vector<std::vector<std::vector<unsigned char>>> _jacobian_block_neighbor_used;
1309  std::vector<std::vector<std::vector<unsigned char>>> _jacobian_block_lower_used;
1311  const DofMap & _dof_map;
1314 
1316 
1317  unsigned int _mesh_dimension;
1318 
1320  std::shared_ptr<XFEMInterface> _xfem;
1321 
1323  std::map<FEType, FEBase *> _current_fe;
1325  std::map<FEType, FEBase *> _current_fe_face;
1327  std::map<FEType, FEBase *> _current_fe_neighbor;
1329  std::map<FEType, FEBase *> _current_fe_face_neighbor;
1330 
1332  std::map<FEType, FEVectorBase *> _current_vector_fe;
1334  std::map<FEType, FEVectorBase *> _current_vector_fe_face;
1336  std::map<FEType, FEVectorBase *> _current_vector_fe_neighbor;
1338  std::map<FEType, FEVectorBase *> _current_vector_fe_face_neighbor;
1339 
1340  /**** Volume Stuff ****/
1341 
1343  mutable std::map<unsigned int, std::map<FEType, FEBase *>> _fe;
1345  mutable std::map<unsigned int, std::map<FEType, const FEBase *>> _const_fe;
1347  mutable std::map<unsigned int, std::map<FEType, FEVectorBase *>> _vector_fe;
1349  mutable std::map<unsigned int, std::map<FEType, const FEVectorBase *>> _const_vector_fe;
1351  std::map<unsigned int, FEBase **> _holder_fe_helper;
1355  const QBase * _const_current_qrule;
1374 
1376  std::map<unsigned int, QBase *> _holder_qrule_volume;
1378  std::map<unsigned int, ArbitraryQuadrature *> _holder_qrule_arbitrary;
1380  std::map<unsigned int, ArbitraryQuadrature *> _holder_qrule_arbitrary_face;
1382  std::map<unsigned int, const std::vector<Point> *> _holder_q_points;
1384  std::map<unsigned int, const std::vector<Real> *> _holder_JxW;
1385 
1386  /**** Face Stuff ****/
1387 
1389  mutable std::map<unsigned int, std::map<FEType, FEBase *>> _fe_face;
1391  mutable std::map<unsigned int, std::map<FEType, const FEBase *>> _const_fe_face;
1393  mutable std::map<unsigned int, std::map<FEType, FEVectorBase *>> _vector_fe_face;
1395  mutable std::map<unsigned int, std::map<FEType, const FEVectorBase *>> _const_vector_fe_face;
1397  std::map<unsigned int, FEBase **> _holder_fe_face_helper;
1413  std::map<unsigned int, QBase *> _holder_qrule_face;
1415  std::map<unsigned int, const std::vector<Point> *> _holder_q_points_face;
1417  std::map<unsigned int, const std::vector<Real> *> _holder_JxW_face;
1419  std::map<unsigned int, const std::vector<Point> *> _holder_normals;
1420 
1421  /**** Neighbor Stuff ****/
1422 
1424  mutable std::map<unsigned int, std::map<FEType, FEBase *>> _fe_neighbor;
1425  mutable std::map<unsigned int, std::map<FEType, FEBase *>> _fe_face_neighbor;
1426  mutable std::map<unsigned int, std::map<FEType, FEVectorBase *>> _vector_fe_neighbor;
1427  mutable std::map<unsigned int, std::map<FEType, FEVectorBase *>> _vector_fe_face_neighbor;
1428  mutable std::map<unsigned int, std::map<FEType, const FEBase *>> _const_fe_neighbor;
1429  mutable std::map<unsigned int, std::map<FEType, const FEBase *>> _const_fe_face_neighbor;
1430  mutable std::map<unsigned int, std::map<FEType, const FEVectorBase *>> _const_vector_fe_neighbor;
1431  mutable std::map<unsigned int, std::map<FEType, const FEVectorBase *>>
1433 
1435  std::map<unsigned int, FEBase **> _holder_fe_neighbor_helper;
1436  std::map<unsigned int, FEBase **> _holder_fe_face_neighbor_helper;
1437 
1439  mutable std::map<unsigned int, std::map<FEType, FEBase *>> _fe_lower;
1441  mutable std::map<unsigned int, std::map<FEType, const FEBase *>> _const_fe_lower;
1443  mutable std::map<unsigned int, std::map<FEType, FEVectorBase *>> _vector_fe_lower;
1445  mutable std::map<unsigned int, std::map<FEType, const FEVectorBase *>> _const_vector_fe_lower;
1446 
1454  std::map<unsigned int, ArbitraryQuadrature *> _holder_qrule_neighbor;
1459 
1460  /********** mortar stuff *************/
1461 
1463  const std::vector<Real> * _JxW_msm;
1465  std::unique_ptr<FEBase> _fe_msm;
1470  QBase * _qrule_msm;
1472  const QBase * _const_qrule_msm;
1473 
1475  const Elem * _current_elem;
1481  unsigned int _current_side;
1483  const Elem * _current_side_elem;
1499  const Node * _current_node;
1506 
1509 
1512 
1514  std::vector<std::vector<DenseVector<Number>>> _sub_Re;
1516  std::vector<std::vector<DenseVector<Number>>> _sub_Rn;
1518  std::vector<std::vector<DenseVector<Number>>> _sub_Rl;
1520  DenseVector<Number> _tmp_Re;
1521 
1523  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kee;
1524  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Keg;
1525 
1527  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Ken;
1529  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kne;
1531  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Knn;
1533  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kll;
1535  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kle;
1537  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kln;
1539  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Kel;
1541  std::vector<std::vector<std::vector<DenseMatrix<Number>>>> _sub_Knl;
1542 
1544  DenseMatrix<Number> _tmp_Ke;
1545 
1546  // Shape function values, gradients. second derivatives
1550 
1554 
1558 
1562 
1563  // Shape function values, gradients, second derivatives
1568 
1573 
1578 
1583 
1585  {
1586  public:
1591  };
1592 
1594  {
1595  public:
1600  };
1601 
1603  mutable std::map<FEType, FEShapeData *> _fe_shape_data;
1604  mutable std::map<FEType, FEShapeData *> _fe_shape_data_face;
1605  mutable std::map<FEType, FEShapeData *> _fe_shape_data_neighbor;
1606  mutable std::map<FEType, FEShapeData *> _fe_shape_data_face_neighbor;
1607  mutable std::map<FEType, FEShapeData *> _fe_shape_data_lower;
1608 
1610  mutable std::map<FEType, VectorFEShapeData *> _vector_fe_shape_data;
1611  mutable std::map<FEType, VectorFEShapeData *> _vector_fe_shape_data_face;
1612  mutable std::map<FEType, VectorFEShapeData *> _vector_fe_shape_data_neighbor;
1613  mutable std::map<FEType, VectorFEShapeData *> _vector_fe_shape_data_face_neighbor;
1614  mutable std::map<FEType, VectorFEShapeData *> _vector_fe_shape_data_lower;
1615 
1618  mutable std::map<FEType,
1623  mutable std::map<FEType,
1626 
1628  std::vector<std::vector<Real>> _cached_residual_values;
1629 
1631  std::vector<std::vector<dof_id_type>> _cached_residual_rows;
1632 
1634 
1636  std::vector<std::vector<Real>> _cached_jacobian_values;
1638  std::vector<std::vector<dof_id_type>> _cached_jacobian_rows;
1640  std::vector<std::vector<dof_id_type>> _cached_jacobian_cols;
1641 
1643 
1646 
1648  std::vector<dof_id_type> _temp_dof_indices;
1649 
1651  std::vector<Point> _temp_reference_points;
1652 
1656  std::vector<std::vector<Real>> _cached_jacobian_contribution_vals;
1657  std::vector<std::vector<numeric_index_type>> _cached_jacobian_contribution_rows;
1658  std::vector<std::vector<numeric_index_type>> _cached_jacobian_contribution_cols;
1659 
1661  std::vector<VectorValue<DualReal>> _ad_dxyzdxi_map;
1662  std::vector<VectorValue<DualReal>> _ad_dxyzdeta_map;
1663  std::vector<VectorValue<DualReal>> _ad_dxyzdzeta_map;
1664  std::vector<VectorValue<DualReal>> _ad_d2xyzdxi2_map;
1665  std::vector<VectorValue<DualReal>> _ad_d2xyzdxideta_map;
1666  std::vector<VectorValue<DualReal>> _ad_d2xyzdeta2_map;
1667  std::vector<DualReal> _ad_jac;
1670  std::vector<DualReal> _ad_dxidx_map;
1671  std::vector<DualReal> _ad_dxidy_map;
1672  std::vector<DualReal> _ad_dxidz_map;
1673  std::vector<DualReal> _ad_detadx_map;
1674  std::vector<DualReal> _ad_detady_map;
1675  std::vector<DualReal> _ad_detadz_map;
1676  std::vector<DualReal> _ad_dzetadx_map;
1677  std::vector<DualReal> _ad_dzetady_map;
1678  std::vector<DualReal> _ad_dzetadz_map;
1679 
1685 
1686  std::vector<unsigned> _displacements;
1687 
1688  mutable bool _calculate_xyz;
1689  mutable bool _calculate_face_xyz;
1691 
1692  mutable std::map<FEType, bool> _need_second_derivative;
1693  mutable std::map<FEType, bool> _need_second_derivative_neighbor;
1694  mutable std::map<FEType, bool> _need_curl;
1695 };
1696 
1697 template <typename OutputType>
1700 {
1702  return _fe_shape_data_lower[type]->_phi;
1703 }
1704 
1705 template <typename OutputType>
1708 {
1710  return _fe_shape_data_lower[type]->_grad_phi;
1711 }
1712 
1713 template <>
1715 Assembly::feADGradPhi<RealVectorValue>(FEType type) const
1716 {
1717  return _ad_vector_grad_phi_data[type];
1718 }
1719 
1720 template <>
1722 Assembly::feADGradPhiFace<RealVectorValue>(FEType type) const
1723 {
1724  return _ad_vector_grad_phi_data_face[type];
1725 }
1726 
1727 template <>
1728 inline void
1729 Assembly::computeGradPhiAD<RealVectorValue>(
1730  const Elem *,
1731  unsigned int,
1733  FEGenericBase<RealVectorValue> *)
1734 {
1735  mooseError("Not implemented");
1736 }
1737 
1738 template <>
1740 Assembly::adCurvatures<ComputeStage::JACOBIAN>() const
1741 {
1742  _calculate_curvatures = true;
1743  return _ad_curvatures;
1744 }
1745 
1746 template <>
1747 inline const MooseArray<VectorValue<DualReal>> &
1748 Assembly::adNormals<ComputeStage::JACOBIAN>() const
1749 {
1750  return _ad_normals;
1751 }
1752 
1753 template <>
1754 inline const typename PointType<ComputeStage::JACOBIAN>::type &
1755 Assembly::adQPoints<ComputeStage::JACOBIAN>() const
1756 {
1757  _calculate_xyz = true;
1758  return _ad_q_points;
1759 }
1760 
1761 template <>
1762 inline const typename PointType<ComputeStage::JACOBIAN>::type &
1763 Assembly::adQPointsFace<ComputeStage::JACOBIAN>() const
1764 {
1765  _calculate_face_xyz = true;
1766  return _ad_q_points_face;
1767 }
1768 
1769 template <>
1770 inline const MooseArray<DualReal> &
1771 Assembly::adJxWFace<ComputeStage::JACOBIAN>() const
1772 {
1773  return _ad_JxW_face;
1774 }
1775 
1776 template <>
1777 inline const MooseArray<DualReal> &
1778 Assembly::adCoordTransformation<ComputeStage::JACOBIAN>() const
1779 {
1780  _calculate_xyz = _calculate_face_xyz = true;
1781  return _ad_coord;
1782 }
1783 
1784 template <>
1786 Assembly::fePhi<VectorValue<Real>>(FEType type) const;
1787 
1788 template <>
1790 Assembly::feGradPhi<VectorValue<Real>>(FEType type) const;
1791 
1792 template <>
1794 Assembly::feSecondPhi<VectorValue<Real>>(FEType type) const;
1795 
1796 template <>
1798 Assembly::fePhiLower<VectorValue<Real>>(FEType type) const;
1799 
1800 template <>
1802 Assembly::feGradPhiLower<VectorValue<Real>>(FEType type) const;
1803 
1804 template <>
1806 Assembly::fePhiFace<VectorValue<Real>>(FEType type) const;
1807 
1808 template <>
1810 Assembly::feGradPhiFace<VectorValue<Real>>(FEType type) const;
1811 
1812 template <>
1814 Assembly::feSecondPhiFace<VectorValue<Real>>(FEType type) const;
1815 
1816 template <>
1818 Assembly::fePhiNeighbor<VectorValue<Real>>(FEType type) const;
1819 
1820 template <>
1822 Assembly::feGradPhiNeighbor<VectorValue<Real>>(FEType type) const;
1823 
1824 template <>
1826 Assembly::feSecondPhiNeighbor<VectorValue<Real>>(FEType type) const;
1827 
1828 template <>
1830 Assembly::fePhiFaceNeighbor<VectorValue<Real>>(FEType type) const;
1831 
1832 template <>
1834 Assembly::feGradPhiFaceNeighbor<VectorValue<Real>>(FEType type) const;
1835 
1836 template <>
1838 Assembly::feSecondPhiFaceNeighbor<VectorValue<Real>>(FEType type) const;
1839 
1840 template <>
1842 Assembly::feCurlPhi<VectorValue<Real>>(FEType type) const;
1843 
1844 template <>
1846 Assembly::feCurlPhiFace<VectorValue<Real>>(FEType type) const;
1847 
1848 template <>
1850 Assembly::feCurlPhiNeighbor<VectorValue<Real>>(FEType type) const;
1851 
1852 template <>
1854 Assembly::feCurlPhiFaceNeighbor<VectorValue<Real>>(FEType type) const;
1855 
1856 template <>
1858 Assembly::adGradPhi<Real, ComputeStage::JACOBIAN>(const MooseVariableFE<Real> & v) const;
1859 
1860 template <>
1862 Assembly::adGradPhi<RealVectorValue, ComputeStage::JACOBIAN>(
1863  const MooseVariableFE<RealVectorValue> & v) const;
1864 
1865 template <>
1866 const MooseArray<typename Moose::RealType<RESIDUAL>::type> & Assembly::adJxW<RESIDUAL>() const;
1867 
std::vector< std::pair< MooseVariableFEBase *, MooseVariableScalar * > > _cm_fs_entry
Entries in the coupling matrix for field variables vs scalar variables.
Definition: Assembly.h:1296
VariablePhiSecond _second_phi_face
Definition: Assembly.h:1553
const Elem *const & elem() const
Return the current element.
Definition: Assembly.h:283
VariablePhiSecond _second_phi_face_neighbor
Definition: Assembly.h:1561
const VectorVariablePhiGradient & gradPhiNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:849
const FEVectorBase *const & getVectorFEFace(FEType type, unsigned int dim) const
GetVector a reference to a pointer that will contain the current "face" FE.
Definition: Assembly.h:146
const ADPoint & adQPoints() const
Definition: Assembly.h:268
MooseArray< DualReal > _ad_JxW
Definition: Assembly.h:1668
const FEVectorBase *const & getVectorFE(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current volume FEVector.
Definition: Assembly.h:122
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_face
types of finite elements
Definition: Assembly.h:1389
std::map< unsigned int, std::map< FEType, const FEBase * > > _const_fe_lower
FE objects for lower dimensional elements.
Definition: Assembly.h:1441
bool _need_neighbor_elem_volume
true is apps need to compute neighbor element volume
Definition: Assembly.h:1495
const unsigned int & neighborSide() const
Returns the current neighboring side.
Definition: Assembly.h:311
VariablePhiSecond & secondPhiNeighbor(const MooseVariable &)
Definition: Assembly.h:890
VectorVariablePhiSecond & secondPhi(const VectorMooseVariable &)
Definition: Assembly.h:904
ArbitraryQuadrature * _current_qrule_arbitrary
The current arbitrary quadrature rule used within the element interior.
Definition: Assembly.h:1361
std::map< FEType, FEBase * > _current_fe
The "volume" fe object that matches the current elem.
Definition: Assembly.h:1323
MooseArray< Real > _curvatures
Definition: Assembly.h:1683
VectorVariablePhiValue _phi
Definition: Assembly.h:1596
std::map< unsigned int, FEBase ** > _holder_fe_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:1351
SystemBase & _sys
Definition: Assembly.h:1282
const MooseArray< ADReal > & adJxW() const
Definition: Assembly.h:195
std::vector< std::vector< dof_id_type > > _cached_residual_rows
Where the cached values should go (the first vector is for TIME vs NONTIME)
Definition: Assembly.h:1631
std::vector< std::vector< dof_id_type > > _cached_jacobian_rows
Row where the corresponding cached value should go.
Definition: Assembly.h:1638
VariablePhiSecond & secondPhiFace(const MooseVariable &)
Definition: Assembly.h:886
unsigned int _max_cached_residuals
Definition: Assembly.h:1633
MooseArray< VectorValue< DualReal > > _ad_q_points_face
Definition: Assembly.h:1682
const OutputTools< OutputType >::VariablePhiCurl & feCurlPhi(FEType type) const
Definition: Assembly.h:1060
void addJacobianBlock(SparseMatrix< Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, std::vector< dof_id_type > &dof_indices)
Definition: Assembly.C:3321
std::vector< DualReal > _ad_detadz_map
Definition: Assembly.h:1675
MooseArray< Point > _current_physical_points
This will be filled up with the physical points passed into reinitAtPhysical() if it is called...
Definition: Assembly.h:1511
const VectorVariablePhiCurl & curlPhiNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:857
std::vector< DualReal > _ad_dxidy_map
Definition: Assembly.h:1671
const VariablePhiSecond & secondPhi(const MooseVariable &) const
Definition: Assembly.h:783
std::map< unsigned int, ArbitraryQuadrature * > _holder_qrule_arbitrary_face
Holds arbitrary qrules for each dimension for faces.
Definition: Assembly.h:1380
MooseArray< VectorValue< DualReal > > _ad_normals
Definition: Assembly.h:1681
const QBase * _const_current_qrule
The current current quadrature rule being used (could be either volumetric or arbitrary - for dirac k...
Definition: Assembly.h:1355
std::map< unsigned int, std::map< FEType, const FEVectorBase * > > _const_vector_fe_face_neighbor
Definition: Assembly.h:1432
void buildNeighborFE(FEType type) const
Build FEs for a neighbor with a type.
Definition: Assembly.C:277
void cacheResidualNodes(const DenseVector< Number > &res, const std::vector< dof_id_type > &dof_index, TagID tag=0)
Lets an external class cache residual at a set of nodes.
Definition: Assembly.C:2875
const VariablePhiSecond & secondPhi() const
Definition: Assembly.h:782
std::map< unsigned int, std::map< FEType, FEVectorBase * > > _vector_fe
Each dimension&#39;s actual vector fe objects indexed on type.
Definition: Assembly.h:1347
void buildVectorLowerDFE(FEType type) const
Build Vector FEs for a lower dimensional element with a type.
Definition: Assembly.C:339
std::vector< VectorValue< DualReal > > _ad_d2xyzdxi2_map
Definition: Assembly.h:1664
MooseArray< Real > _coord_neighbor
The current coordinate transformation coefficients.
Definition: Assembly.h:1458
const FEBase *const & getFENeighbor(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current &#39;neighbor&#39; FE.
Definition: Assembly.h:86
void buildFE(FEType type) const
Build FEs with a type.
Definition: Assembly.C:234
QBase * _qrule_msm
A qrule object for working on mortar segement elements.
Definition: Assembly.h:1470
QBase * _current_qrule
The current current quadrature rule being used (could be either volumetric or arbitrary - for dirac k...
Definition: Assembly.h:1357
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kle
dlower/dslave (or dlower/delement)
Definition: Assembly.h:1535
std::map< unsigned int, FEBase ** > _holder_fe_face_neighbor_helper
Definition: Assembly.h:1436
const VectorVariablePhiCurl & curlPhiFace(const VectorMooseVariable &) const
Definition: Assembly.h:840
const OutputTools< OutputType >::VariablePhiValue & fePhi(FEType type) const
Definition: Assembly.h:950
MooseArray< DualReal > _ad_JxW_face
Definition: Assembly.h:1680
void computeGradPhiAD(const Elem *elem, unsigned int n_qp, typename VariableTestGradientType< OutputType, ComputeStage::JACOBIAN >::type &grad_phi, FEGenericBase< OutputType > *fe)
Definition: Assembly.C:675
VectorVariablePhiSecond & secondPhiNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:923
Keeps track of stuff related to assembling.
Definition: Assembly.h:62
Assembly(SystemBase &sys, THREAD_ID tid)
Definition: Assembly.C:36
Class for stuff related to variables.
Definition: Adaptivity.h:29
VectorVariablePhiGradient & gradPhiFace(const VectorMooseVariable &)
Definition: Assembly.h:908
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:725
unsigned int TagID
Definition: MooseTypes.h:162
void init()
Deprecated init method.
VariablePhiValue & phiFace(const MooseVariable &)
Definition: Assembly.h:884
VectorVariablePhiValue _vector_phi
Definition: Assembly.h:1564
const VariablePhiSecond & secondPhiFaceNeighbor(const MooseVariable &) const
Definition: Assembly.h:809
MooseArray< Real > _current_JxW_neighbor
The current transformed jacobian weights on a neighbor&#39;s face.
Definition: Assembly.h:1456
std::shared_ptr< XFEMInterface > _xfem
The XFEM controller.
Definition: Assembly.h:1320
virtual ~Assembly()
Definition: Assembly.C:120
void prepareJacobianBlock()
Sizes and zeroes the Jacobian blocks used for the current element.
Definition: Assembly.C:2302
void setResidualNeighbor(NumericVector< Number > &residual, TagID tag_id=0)
Definition: Assembly.C:2970
VariablePhiValue & phiNeighbor(const MooseVariable &)
Definition: Assembly.h:888
VectorValue< Real > RealVectorValue
Definition: Assembly.h:31
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:207
std::vector< DualReal > _ad_dzetady_map
Definition: Assembly.h:1677
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kel
dslave/dlower (or delement/dlower)
Definition: Assembly.h:1539
MooseArray< Real > _coord
The current coordinate transformation coefficients.
Definition: Assembly.h:1371
const OutputTools< OutputType >::VariablePhiValue & fePhiLower(FEType type) const
Definition: Assembly.h:1699
MooseArray< VectorValue< DualReal > > _ad_q_points
Definition: Assembly.h:1669
void reinitFE(const Elem *elem)
Just an internal helper function to reinit the volume FE objects.
Definition: Assembly.C:550
std::map< FEType, VectorFEShapeData * > _vector_fe_shape_data_neighbor
Definition: Assembly.h:1612
const MooseArray< Point > & physicalPoints() const
The current points in physical space where we have reinited through reinitAtPhysical() ...
Definition: Assembly.h:186
void prepareNonlocal()
Definition: Assembly.C:2341
VariablePhiGradient _grad_phi
Definition: Assembly.h:1548
std::unique_ptr< FEBase > _fe_msm
A FE object for working on mortar segement elements.
Definition: Assembly.h:1465
VariablePhiSecond _second_phi_neighbor
Definition: Assembly.h:1557
void setXFEM(std::shared_ptr< XFEMInterface > xfem)
Set the pointer to the XFEM controller object.
Definition: Assembly.h:1125
const VariablePhiGradient & gradPhiFaceNeighbor(const MooseVariable &) const
Definition: Assembly.h:805
std::vector< DualReal > _ad_dzetadz_map
Definition: Assembly.h:1678
void reinitAtPhysical(const Elem *elem, const std::vector< Point > &physical_points)
Reinitialize the assembly data at specific physical point in the given element.
Definition: Assembly.C:1546
std::map< unsigned int, ArbitraryQuadrature * > _holder_qrule_neighbor
Holds arbitrary qrules for each dimension.
Definition: Assembly.h:1454
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries()
Definition: Assembly.h:762
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2019
bool _current_elem_volume_computed
Boolean to indicate whether current element volumes has been computed.
Definition: Assembly.h:1503
MooseVariableFE< Real > MooseVariable
Definition: Assembly.h:52
std::map< unsigned int, std::map< FEType, const FEVectorBase * > > _const_vector_fe_lower
Vector FE objects for lower dimensional elements.
Definition: Assembly.h:1445
VectorVariablePhiGradient _vector_grad_phi_neighbor
Definition: Assembly.h:1575
const MooseArray< Point > & qPointsFaceNeighbor() const
Returns the reference to the current quadrature points being used on the neighbor face...
Definition: Assembly.h:375
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kll
dlower/dlower
Definition: Assembly.h:1533
VariablePhiValue _phi
Definition: Assembly.h:1587
Real _current_neighbor_volume
Volume of the current neighbor.
Definition: Assembly.h:1497
const VariableTestGradientType< T, compute_stage >::type & adGradPhi(const MooseVariableFE< T > &v) const
Definition: Assembly.h:775
std::map< FEType, FEBase * > _current_fe_face
The "face" fe object that matches the current elem.
Definition: Assembly.h:1325
void addResidualNeighbor(NumericVector< Number > &residual, TagID tag_id=0)
Definition: Assembly.C:2772
VariablePhiSecond & secondPhiFaceNeighbor(const MooseVariable &)
Definition: Assembly.h:897
void computeAffineMapAD(const Elem *elem, const std::vector< Real > &qw, unsigned int n_qp, FEBase *fe)
Definition: Assembly.C:772
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kln
dlower/dmaster (or dlower/dneighbor)
Definition: Assembly.h:1537
MooseVariableFE< RealVectorValue > VectorMooseVariable
Definition: Assembly.h:54
const Elem * _current_neighbor_elem
The current neighbor "element".
Definition: Assembly.h:1487
std::vector< std::vector< DenseVector< Number > > > _sub_Rn
residual contributions for each variable from the neighbor
Definition: Assembly.h:1516
VariablePhiSecond & secondPhi(const MooseVariable &)
Definition: Assembly.h:882
MooseMesh & _mesh
Definition: Assembly.h:1315
const VariablePhiValue & phi() const
Definition: Assembly.h:772
void cacheJacobianCoupledVarPair(MooseVariableBase *ivar, MooseVariableBase *jvar)
Caches element matrix for ivar rows and jvar columns.
Definition: Assembly.C:3239
const VariablePhiGradient & gradPhi() const
Definition: Assembly.h:780
MooseArray< Real > _current_JxW_face
The current transformed jacobian weights on a face.
Definition: Assembly.h:1409
void modifyWeightsDueToXFEM(const Elem *elem)
Update the integration weights for XFEM partial elements.
Definition: Assembly.C:3534
DenseMatrix< Number > _tmp_Ke
auxiliary matrix for scaling jacobians (optimization to avoid expensive construction/destruction) ...
Definition: Assembly.h:1544
const VariablePhiSecond & secondPhiFace(const MooseVariable &) const
Definition: Assembly.h:789
const Elem * _current_elem
The current "element" we are currently on.
Definition: Assembly.h:1475
VectorVariablePhiCurl _vector_curl_phi_neighbor
Definition: Assembly.h:1577
const FEBase *const & getFEFace(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current "face" FE.
Definition: Assembly.h:98
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & nonlocalCouplingEntries()
Definition: Assembly.h:766
THREAD_ID _tid
Thread number (id)
Definition: Assembly.h:1313
void reinitFEFaceNeighbor(const Elem *neighbor, const std::vector< Point > &reference_points)
Definition: Assembly.C:1318
VectorVariablePhiGradient _vector_grad_phi
Definition: Assembly.h:1565
void reinitNeighborAtPhysical(const Elem *neighbor, unsigned int neighbor_side, const std::vector< Point > &physical_points)
Reinitializes the neighbor at the physical coordinates on neighbor side given.
Definition: Assembly.C:1976
const MooseArray< ADReal > & adCurvatures() const
Definition: Assembly.C:3750
std::vector< std::vector< std::vector< unsigned char > > > _jacobian_block_nonlocal_used
Definition: Assembly.h:1305
VectorVariablePhiValue _vector_phi_neighbor
Definition: Assembly.h:1574
void prepareNeighbor()
Definition: Assembly.C:2419
OutputTools< Real >::VariablePhiValue VariablePhiValue
Definition: MooseTypes.h:201
std::map< FEType, typename VariableTestGradientType< Real, ComputeStage::JACOBIAN >::type > _ad_grad_phi_data
Definition: Assembly.h:1617
std::map< unsigned int, std::map< FEType, FEVectorBase * > > _vector_fe_face_neighbor
Definition: Assembly.h:1427
std::map< FEType, bool > _need_second_derivative
Definition: Assembly.h:1692
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > _cm_ff_entry
Entries in the coupling matrix for field variables.
Definition: Assembly.h:1294
Real _current_elem_volume
Volume of the current element.
Definition: Assembly.h:1479
OutputTools< RealVectorValue >::VariablePhiGradient VectorVariablePhiGradient
Definition: MooseTypes.h:216
VariablePhiGradient & gradPhiFaceNeighbor(const MooseVariable &)
Definition: Assembly.h:893
std::map< FEType, FEBase * > _current_fe_face_neighbor
The "neighbor face" fe object that matches the current elem.
Definition: Assembly.h:1329
const MooseArray< Point > & qPoints() const
Returns the reference to the quadrature points.
Definition: Assembly.h:180
std::map< FEType, FEVectorBase * > _current_vector_fe_face
The "face" vector fe object that matches the current elem.
Definition: Assembly.h:1334
QBase *const & writeableQRuleNeighbor()
Returns the reference to the current quadrature being used on a current neighbor. ...
Definition: Assembly.h:363
VectorVariablePhiValue & phiFaceNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:932
std::vector< DualReal > _ad_dzetadx_map
Definition: Assembly.h:1676
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
std::map< FEType, FEShapeData * > _fe_shape_data_face_neighbor
Definition: Assembly.h:1606
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Knn
jacobian contributions from the neighbor <Tag, ivar, jvar>
Definition: Assembly.h:1531
Base class for a system (of equations)
Definition: SystemBase.h:92
void setVolumeQRule(QBase *qrule, unsigned int dim)
Set the qrule to be used for volume integration.
Definition: Assembly.C:514
void zeroCachedJacobianContributions()
Zero out previously-cached Jacobian rows.
Definition: Assembly.C:3476
unsigned int _block_diagonal_matrix
Will be true if our preconditioning matrix is a block-diagonal matrix. Which means that we can take s...
Definition: Assembly.h:1645
const VectorVariablePhiValue & phiFaceNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:862
const VariableTestGradientType< OutputType, ComputeStage::JACOBIAN >::type & feADGradPhiFace(FEType type) const
Definition: Assembly.h:1000
void addJacobianNonlocal()
Definition: Assembly.C:3137
std::vector< std::vector< Real > > _cached_jacobian_values
Values cached by calling cacheJacobian()
Definition: Assembly.h:1636
VariablePhiValue _phi
Definition: Assembly.h:1547
const OutputTools< OutputType >::VariablePhiSecond & feSecondPhiFaceNeighbor(FEType type) const
Definition: Assembly.h:1052
std::map< FEType, FEVectorBase * > _current_vector_fe
The "volume" vector fe object that matches the current elem.
Definition: Assembly.h:1332
VariablePhiGradient _grad_phi
Definition: Assembly.h:1588
void prepareScalar()
Definition: Assembly.C:2550
void copyNeighborShapes(MooseVariableFE< T > &v)
Definition: Assembly.C:2660
VectorVariablePhiValue & phiNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:918
void modifyFaceWeightsDueToXFEM(const Elem *elem, unsigned int side=0)
Update the face integration weights for XFEM partial elements.
Definition: Assembly.C:3554
const CouplingMatrix * _cm
Coupling matrices.
Definition: Assembly.h:1288
std::map< FEType, bool > _need_second_derivative_neighbor
Definition: Assembly.h:1693
const VectorVariablePhiSecond & secondPhiFaceNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:870
const OutputTools< OutputType >::VariablePhiGradient & feGradPhiFace(FEType type) const
Definition: Assembly.h:992
void createQRules(QuadratureType type, Order order, Order volume_order, Order face_order)
Creates the volume, face and arbitrary qrules based on the orders passed in.
Definition: Assembly.C:487
OutputTools< Real >::VariablePhiSecond VariablePhiSecond
Definition: MooseTypes.h:203
DenseVector< Number > _tmp_Re
auxiliary vector for scaling residuals (optimization to avoid expensive construction/destruction) ...
Definition: Assembly.h:1520
unsigned int _mesh_dimension
Definition: Assembly.h:1317
void reinit(const Elem *elem)
Reinitialize objects (JxW, q_points, ...) for an elements.
Definition: Assembly.C:1567
QBase * _current_qrule_volume
The current volumetric quadrature for the element.
Definition: Assembly.h:1359
VectorVariablePhiSecond & secondPhiFace(const VectorMooseVariable &)
Definition: Assembly.h:912
const DofMap & _dof_map
DOF map.
Definition: Assembly.h:1311
std::vector< std::vector< Real > > _cached_jacobian_contribution_vals
Storage for cached Jacobian entries.
Definition: Assembly.h:1656
std::vector< DualReal > _ad_dxidx_map
Definition: Assembly.h:1670
unsigned int _max_cached_jacobians
Definition: Assembly.h:1642
const VariablePhiValue & phiFaceNeighbor(const MooseVariable &) const
Definition: Assembly.h:801
VectorVariablePhiCurl & curlPhiNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:927
void cacheResidual()
Takes the values that are currently in _sub_Re and appends them to the cached values.
Definition: Assembly.C:2809
DenseVector< Number > & residualBlockLower(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:730
std::map< FEType, VectorFEShapeData * > _vector_fe_shape_data_face_neighbor
Definition: Assembly.h:1613
const MooseArray< ADReal > & adJxWFace() const
Definition: Assembly.h:201
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Ken
jacobian contributions from the element and neighbor <Tag, ivar, jvar>
Definition: Assembly.h:1527
const QBase * _const_qrule_msm
A pointer to const qrule_msm.
Definition: Assembly.h:1472
const CouplingMatrix & _nonlocal_cm
Definition: Assembly.h:1289
VectorVariablePhiValue _vector_phi_face
Definition: Assembly.h:1569
std::vector< Point > _temp_reference_points
Temporary work data for reinitAtPhysical()
Definition: Assembly.h:1651
std::vector< std::pair< MooseVariableScalar *, MooseVariableScalar * > > _cm_ss_entry
Entries in the coupling matrix for scalar variables.
Definition: Assembly.h:1300
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kne
jacobian contributions from the neighbor and element <Tag, ivar, jvar>
Definition: Assembly.h:1529
const VectorVariablePhiGradient & gradPhi(const VectorMooseVariable &) const
Definition: Assembly.h:815
unsigned int _current_neighbor_side
The current side of the selected neighboring element (valid only when working with sides) ...
Definition: Assembly.h:1491
void prepareLowerD()
Prepare the Jacobians and residuals for a lower dimensional element.
Definition: Assembly.C:2459
QBase *const & writeableQRuleFace()
Returns the reference to the current quadrature being used on a current face.
Definition: Assembly.h:241
const VectorVariablePhiValue & phiNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:845
void buildFaceNeighborFE(FEType type) const
Build FEs for a neighbor face with a type.
Definition: Assembly.C:297
void computeCurrentFaceVolume()
Definition: Assembly.C:1527
std::vector< unsigned > _displacements
Definition: Assembly.h:1686
std::vector< std::vector< DenseVector< Number > > > _sub_Rl
residual contributions for each variable from the lower dimensional element
Definition: Assembly.h:1518
VariablePhiCurl _curl_phi
Definition: Assembly.h:1590
void computeADFace(const Elem *elem, unsigned int side)
compute AD things on an element face
Definition: Assembly.C:1782
VectorVariablePhiGradient _vector_grad_phi_face_neighbor
Definition: Assembly.h:1580
DenseMatrix< Number > & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2026
std::map< unsigned int, QBase * > _holder_qrule_face
Holds face qrules for each dimension.
Definition: Assembly.h:1413
void cacheJacobianBlockNonlocal(DenseMatrix< Number > &jac_block, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, Real scaling_factor, TagID tag=0)
Definition: Assembly.C:3032
void prepareOffDiagScalar()
Definition: Assembly.C:2578
std::map< unsigned int, std::map< FEType, const FEBase * > > _const_fe_neighbor
Definition: Assembly.h:1428
VectorVariablePhiSecond _second_phi
Definition: Assembly.h:1598
std::map< unsigned int, std::map< FEType, const FEBase * > > _const_fe_face_neighbor
Definition: Assembly.h:1429
std::map< unsigned int, std::map< FEType, const FEBase * > > _const_fe_face
types of finite elements
Definition: Assembly.h:1391
Implements a fake quadrature rule where you can specify the locations (in the reference domain) of th...
void addCachedResidual(NumericVector< Number > &residual, TagID tag_id)
Adds the values that have been cached by calling cacheResidual(), cacheResidualNeighbor(), and/or cacheResidualLower() to the residual.
Definition: Assembly.C:2905
void copyShapes(MooseVariableFE< T > &v)
Definition: Assembly.C:2606
void addJacobian()
Definition: Assembly.C:3124
VectorVariablePhiGradient & gradPhiFaceNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:936
void setCoordinateTransformation(const QBase *qrule, const ADPoint &q_points, MooseArray< ADReal > &coord)
Definition: Assembly.C:1473
void setResidualBlock(NumericVector< Number > &residual, DenseVector< Number > &res_block, const std::vector< dof_id_type > &dof_indices, Real scaling_factor)
Definition: Assembly.C:2939
std::map< FEType, FEShapeData * > _fe_shape_data_neighbor
Definition: Assembly.h:1605
QBase * _current_qrule_neighbor
quadrature rule used on neighbors
Definition: Assembly.h:1450
std::vector< DualReal > _ad_dxidz_map
Definition: Assembly.h:1672
VectorVariablePhiCurl & curlPhiFaceNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:944
std::map< unsigned int, std::map< FEType, const FEBase * > > _const_fe
Each dimension&#39;s actual fe objects indexed on type.
Definition: Assembly.h:1345
std::vector< std::vector< numeric_index_type > > _cached_jacobian_contribution_cols
Definition: Assembly.h:1658
VectorVariablePhiGradient _grad_phi
Definition: Assembly.h:1597
ArbitraryQuadrature * _current_qrule_arbitrary_face
The current arbitrary quadrature rule used on the element face.
Definition: Assembly.h:1363
const OutputTools< OutputType >::VariablePhiCurl & feCurlPhiFaceNeighbor(FEType type) const
Definition: Assembly.h:1084
FEBase * _current_fe_face_helper
helper object for transforming coordinates
Definition: Assembly.h:1399
const std::vector< Real > * _JxW_msm
A JxW for working on mortar segement elements.
Definition: Assembly.h:1463
const MooseArray< ADReal > & adCoordTransformation() const
Returns the reference to the AD version of the coordinate transformation coefficients.
Definition: Assembly.h:220
const MooseArray< Point > & normals() const
Returns the array of normals for quadrature points on a current side.
Definition: Assembly.h:259
const QBase *const & qRuleNeighbor() const
Returns the reference to the current quadrature being used on a current neighbor. ...
Definition: Assembly.h:357
void prepareResidual()
Sizes and zeroes the residual for the current element.
Definition: Assembly.C:2322
void prepareVariableNonlocal(MooseVariableFEBase *var)
Definition: Assembly.C:2393
void addJacobianOffDiagScalar(unsigned int ivar)
Definition: Assembly.C:3421
void assignDisplacements(std::vector< unsigned > &&disp_numbers)
Definition: Assembly.h:1127
VariablePhiValue & phiFaceNeighbor(const MooseVariable &)
Definition: Assembly.h:892
VectorVariablePhiCurl _curl_phi
Definition: Assembly.h:1599
std::map< FEType, FEBase * > _current_fe_neighbor
The "neighbor" fe object that matches the current elem.
Definition: Assembly.h:1327
void cacheJacobianContribution(numeric_index_type i, numeric_index_type j, Real value, TagID tag=0)
Caches the Jacobian entry &#39;value&#39;, to eventually be added/set in the (i,j) location of the matrix...
Definition: Assembly.C:3430
SubdomainID _current_neighbor_subdomain_id
The current neighbor subdomain ID.
Definition: Assembly.h:1489
std::map< unsigned int, FEBase ** > _holder_fe_neighbor_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:1435
QBase *const & writeableQRule()
Returns the reference to the current quadrature being used.
Definition: Assembly.h:174
VectorVariablePhiValue _vector_phi_face_neighbor
Definition: Assembly.h:1579
const ADPoint & adQPointsFace() const
Definition: Assembly.h:274
std::vector< std::vector< DenseVector< Number > > > _sub_Re
residual contributions for each variable from the element
Definition: Assembly.h:1514
VectorVariablePhiValue & phi(const VectorMooseVariable &)
Definition: Assembly.h:902
OutputTools< RealVectorValue >::VariablePhiSecond VectorVariablePhiSecond
Definition: MooseTypes.h:217
void reinitNeighbor(const Elem *neighbor, const std::vector< Point > &reference_points)
Definition: Assembly.C:1413
DenseMatrix< Number > & jacobianBlockLower(Moose::ConstraintJacobianType type, unsigned int ivar, unsigned int jvar, TagID tag=0)
Returns the jacobian block for the given mortar Jacobian type.
Definition: Assembly.C:2072
const ADPoint & adNormals() const
Definition: Assembly.h:262
const MooseArray< Real > & coordTransformation() const
Returns the reference to the coordinate transformation coefficients.
Definition: Assembly.h:213
const OutputTools< OutputType >::VariablePhiCurl & feCurlPhiNeighbor(FEType type) const
Definition: Assembly.h:1076
VariablePhiGradient & gradPhi(const MooseVariable &)
Definition: Assembly.h:881
std::vector< std::vector< std::vector< unsigned char > > > _jacobian_block_used
Flag that indicates if the jacobian block was used.
Definition: Assembly.h:1304
void prepare()
Definition: Assembly.C:2334
DofMap & dof_map
void setFaceQRule(QBase *qrule, unsigned int dim)
Set the qrule to be used for face integration.
Definition: Assembly.C:528
VectorVariablePhiCurl & curlPhiFace(const VectorMooseVariable &)
Definition: Assembly.h:916
void reinitElemAndNeighbor(const Elem *elem, unsigned int side, const Elem *neighbor, unsigned int neighbor_side)
Reinitialize an element and its neighbor along a particular side.
Definition: Assembly.C:1679
std::vector< DualReal > _ad_jac
Definition: Assembly.h:1667
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
void addJacobianLower()
Add LowerLower, LowerSlave (LowerElement), LowerMaster (LowerNeighbor), SlaveLower (ElementLower)...
Definition: Assembly.C:3193
SubProblem & _subproblem
Definition: Assembly.h:1283
void buildFaceFE(FEType type) const
Build FEs for a face with a type.
Definition: Assembly.C:257
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_neighbor
types of finite elements
Definition: Assembly.h:1424
bool _calculate_xyz
Definition: Assembly.h:1688
bool _calculate_curvatures
Definition: Assembly.h:1690
VariablePhiGradient _grad_phi_face_neighbor
Definition: Assembly.h:1560
std::map< FEType, FEVectorBase * > _current_vector_fe_face_neighbor
The "neighbor face" vector fe object that matches the current elem.
Definition: Assembly.h:1338
const OutputTools< OutputType >::VariablePhiCurl & feCurlPhiFace(FEType type) const
Definition: Assembly.h:1068
MooseArray< Real > _current_JxW
The current list of transformed jacobian weights.
Definition: Assembly.h:1367
void prepareVariable(MooseVariableFEBase *var)
Used for preparing the dense residual and jacobian blocks for one particular variable.
Definition: Assembly.C:2364
std::vector< dof_id_type > _temp_dof_indices
Temporary work vector to keep from reallocating it.
Definition: Assembly.h:1648
OutputTools< Real >::VariablePhiCurl VariablePhiCurl
Definition: MooseTypes.h:204
const VectorVariablePhiGradient & gradPhiFaceNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:866
const QBase *const & qRule() const
Returns the reference to the current quadrature being used.
Definition: Assembly.h:168
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Kee
jacobian contributions <Tag, ivar, jvar>
Definition: Assembly.h:1523
void cacheResidualNeighbor()
Takes the values that are currently in _sub_Rn and appends them to the cached values.
Definition: Assembly.C:2839
const VectorVariablePhiValue & phi(const VectorMooseVariable &) const
Definition: Assembly.h:814
std::map< unsigned int, ArbitraryQuadrature * > _holder_qrule_arbitrary
Holds arbitrary qrules for each dimension.
Definition: Assembly.h:1378
std::map< unsigned int, const std::vector< Point > * > _holder_q_points_face
Holds pointers to the dimension&#39;s q_points on a face.
Definition: Assembly.h:1415
VariablePhiGradient _grad_phi_neighbor
Definition: Assembly.h:1556
const Real & elemVolume()
Returns the reference to the current element volume.
Definition: Assembly.h:299
ArbitraryQuadrature * _current_qface_arbitrary
The current arbitrary quadrature rule used on element faces.
Definition: Assembly.h:1405
std::map< unsigned int, std::map< FEType, FEVectorBase * > > _vector_fe_face
types of vector finite elements
Definition: Assembly.h:1393
VariablePhiValue _phi_face
Definition: Assembly.h:1551
VectorVariablePhiSecond _vector_second_phi
Definition: Assembly.h:1566
subdomain_id_type SubdomainID
std::map< FEType, typename VariableTestGradientType< RealVectorValue, ComputeStage::JACOBIAN >::type > _ad_vector_grad_phi_data_face
Definition: Assembly.h:1625
VectorVariablePhiValue & phiFace(const VectorMooseVariable &)
Definition: Assembly.h:907
const OutputTools< OutputType >::VariablePhiGradient & feGradPhi(FEType type) const
Definition: Assembly.h:957
const VariablePhiSecond & secondPhiNeighbor(const MooseVariable &) const
Definition: Assembly.h:796
VariablePhiGradient & gradPhiNeighbor(const MooseVariable &)
Definition: Assembly.h:889
const Elem * _current_lower_d_elem
The current lower dimensional element.
Definition: Assembly.h:1508
std::map< FEType, typename VariableTestGradientType< RealVectorValue, ComputeStage::JACOBIAN >::type > _ad_vector_grad_phi_data
Definition: Assembly.h:1620
const VariablePhiValue & phiFace(const MooseVariable &) const
Definition: Assembly.h:786
std::vector< VectorValue< DualReal > > _ad_dxyzdeta_map
Definition: Assembly.h:1662
const VariablePhiValue & phiFace() const
Definition: Assembly.h:785
VectorVariablePhiGradient & gradPhiNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:919
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_lower
FE objects for lower dimensional elements.
Definition: Assembly.h:1439
VectorVariablePhiSecond & secondPhiFaceNeighbor(const VectorMooseVariable &)
Definition: Assembly.h:940
const VectorVariablePhiCurl & curlPhi(const VectorMooseVariable &) const
Definition: Assembly.h:823
FEGenericBase< Real > FEBase
Definition: Assembly.h:34
std::map< unsigned int, FEBase ** > _holder_fe_face_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:1397
const OutputTools< OutputType >::VariablePhiGradient & feGradPhiLower(FEType type) const
Definition: Assembly.h:1707
VectorVariablePhiGradient & gradPhi(const VectorMooseVariable &)
Definition: Assembly.h:903
VariablePhiValue _phi_face_neighbor
Definition: Assembly.h:1559
const OutputTools< OutputType >::VariablePhiValue & fePhiNeighbor(FEType type) const
Definition: Assembly.h:1014
std::map< unsigned int, std::map< FEType, FEBase * > > _fe
Each dimension&#39;s actual fe objects indexed on type.
Definition: Assembly.h:1343
VariablePhiValue _phi_neighbor
Definition: Assembly.h:1555
const Node *const & nodeNeighbor() const
Returns the reference to the neighboring node.
Definition: Assembly.h:387
const MooseArray< Real > & JxW() const
Returns the reference to the transformed jacobian weights.
Definition: Assembly.h:192
void reinitMortarElem(const Elem *elem)
reinitialize a mortar segment mesh element in order to get a proper JxW
Definition: Assembly.C:1967
std::map< unsigned int, const std::vector< Real > * > _holder_JxW_face
Holds pointers to the dimension&#39;s transformed jacobian weights on a face.
Definition: Assembly.h:1417
const OutputTools< OutputType >::VariablePhiGradient & feGradPhiFaceNeighbor(FEType type) const
Definition: Assembly.h:1044
Real _current_side_volume
Volume of the current side element.
Definition: Assembly.h:1485
CoordinateSystemType
Definition: MooseTypes.h:556
VectorVariablePhiCurl _vector_curl_phi
Definition: Assembly.h:1567
std::vector< std::vector< Real > > _cached_residual_values
Values cached by calling cacheResidual() (the first vector is for TIME vs NONTIME) ...
Definition: Assembly.h:1628
void addJacobianCoupledVarPair(MooseVariableBase *ivar, MooseVariableBase *jvar)
Adds element matrix for ivar rows and jvar columns.
Definition: Assembly.C:3110
const bool _displaced
Definition: Assembly.h:1285
void addJacobianNeighbor()
Add ElementNeighbor, NeighborElement, and NeighborNeighbor portions of the Jacobian for compute objec...
Definition: Assembly.C:3158
void cacheJacobianBlock(DenseMatrix< Number > &jac_block, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, Real scaling_factor, TagID tag=0)
Definition: Assembly.C:3003
std::map< FEType, VectorFEShapeData * > _vector_fe_shape_data_lower
Definition: Assembly.h:1614
void buildVectorFE(FEType type) const
Build Vector FEs with a type.
Definition: Assembly.C:361
VectorVariablePhiCurl _vector_curl_phi_face_neighbor
Definition: Assembly.h:1582
MooseArray< DualReal > _ad_curvatures
Definition: Assembly.h:1684
std::map< unsigned int, std::map< FEType, FEVectorBase * > > _vector_fe_lower
Vector FE objects for lower dimensional elements.
Definition: Assembly.h:1443
const MooseArray< Real > & JxWFace() const
Returns the reference to the transformed jacobian weights on a current face.
Definition: Assembly.h:253
std::map< unsigned int, std::map< FEType, const FEVectorBase * > > _const_vector_fe
Each dimension&#39;s actual vector fe objects indexed on type.
Definition: Assembly.h:1349
std::vector< VectorValue< DualReal > > _ad_d2xyzdxideta_map
Definition: Assembly.h:1665
const Node * _current_node
The current node we are working with.
Definition: Assembly.h:1499
bool _calculate_face_xyz
Definition: Assembly.h:1689
This is the XFEMInterface class.
Definition: XFEMInterface.h:37
std::map< FEType, typename VariableTestGradientType< Real, ComputeStage::JACOBIAN >::type > _ad_grad_phi_data_face
Definition: Assembly.h:1622
std::vector< std::pair< MooseVariableScalar *, MooseVariableFEBase * > > _cm_sf_entry
Entries in the coupling matrix for scalar variables vs field variables.
Definition: Assembly.h:1298
DGJacobianType
Definition: MooseTypes.h:515
void addCachedJacobian()
Definition: Assembly.C:3070
const VectorVariablePhiValue & phiFace(const VectorMooseVariable &) const
Definition: Assembly.h:828
VectorVariablePhiCurl & curlPhi(const VectorMooseVariable &)
Definition: Assembly.h:905
VariablePhiValue & phi(const MooseVariable &)
Definition: Assembly.h:880
OutputTools< Real >::VariablePhiGradient VariablePhiGradient
Definition: MooseTypes.h:202
MooseArray< Point > _current_normals
The current Normal vectors at the quadrature points.
Definition: Assembly.h:1411
MatType type
void reinitLowerDElemRef(const Elem *elem, const std::vector< Point > *const pts, const std::vector< Real > *const weights=nullptr)
Reinitialize FE data for a lower dimenesional element with a given set of reference points...
Definition: Assembly.C:1935
void buildVectorFaceFE(FEType type) const
Build Vector FEs for a face with a type.
Definition: Assembly.C:393
const QBase *const & qRuleFace() const
Returns the reference to the current quadrature being used on a current face.
Definition: Assembly.h:235
void cacheJacobian()
Takes the values that are currently in _sub_Kee and appends them to the cached values.
Definition: Assembly.C:3253
void setCurrentSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:293
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > _cm_nonlocal_entry
Entries in the coupling matrix for field variables for nonlocal calculations.
Definition: Assembly.h:1302
VectorVariablePhiGradient _vector_grad_phi_face
Definition: Assembly.h:1570
const MooseArray< Real > & JxWNeighbor() const
Returns the reference to the transformed jacobian weights on a current face.
Definition: Assembly.h:369
const VariablePhiGradient & gradPhi(const MooseVariable &) const
Definition: Assembly.h:781
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:59
const QBase * _const_current_qrule_face
quadrature rule used on faces
Definition: Assembly.h:1401
OutputTools< RealVectorValue >::VariablePhiValue VectorVariablePhiValue
Definition: MooseTypes.h:215
void setNeighborQRule(QBase *qrule, unsigned int dim)
Set the qrule to be used for neighbor integration.
Definition: Assembly.C:539
std::map< FEType, VectorFEShapeData * > _vector_fe_shape_data
Shape function values, gradients, second derivatives for each vector FE type.
Definition: Assembly.h:1610
std::vector< VectorValue< DualReal > > _ad_dxyzdzeta_map
Definition: Assembly.h:1663
void reinitNeighborFaceRef(const Elem *neighbor_elem, unsigned int neighbor_side, Real tolerance, const std::vector< Point > *const pts, const std::vector< Real > *const weights=nullptr)
Reinitialize FE data for the given neighbor_element on the given side with a given set of reference p...
Definition: Assembly.C:1866
void setResidual(NumericVector< Number > &residual, TagID tag_id=0)
Definition: Assembly.C:2961
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Keg
Definition: Assembly.h:1524
void resizeMappingObjects(unsigned int n_qp, unsigned int dim)
Definition: Assembly.C:745
const MooseArray< Point > & qPointsFace() const
Returns the reference to the current quadrature being used.
Definition: Assembly.h:247
std::map< unsigned int, std::map< FEType, const FEVectorBase * > > _const_vector_fe_neighbor
Definition: Assembly.h:1430
void buildVectorNeighborFE(FEType type) const
Build Vector FEs for a neighbor with a type.
Definition: Assembly.C:422
void reinitFENeighbor(const Elem *neighbor, const std::vector< Point > &reference_points)
Definition: Assembly.C:1366
void initNonlocalCoupling()
Create pair of variables requiring nonlocal jacobian contributions.
Definition: Assembly.C:2283
void cacheJacobianNeighbor()
Takes the values that are currently in the neighbor Dense Matrices and appends them to the cached val...
Definition: Assembly.C:3287
ConstraintJacobianType
Definition: MooseTypes.h:543
const QBase *const & qRuleMortar() const
Returns a reference to the quadrature rule for the mortar segments.
Definition: Assembly.h:471
void addResidual(NumericVector< Number > &residual, TagID tag_id=0)
Definition: Assembly.C:2755
const bool & _computing_jacobian
Definition: Assembly.h:1291
void computeCurrentNeighborVolume()
const Elem *const & neighbor() const
Return the neighbor element.
Definition: Assembly.h:329
std::map< unsigned int, std::map< FEType, FEVectorBase * > > _vector_fe_neighbor
Definition: Assembly.h:1426
void addResidualScalar(TagID tag_id)
Definition: Assembly.C:2789
const OutputTools< OutputType >::VariablePhiValue & fePhiFaceNeighbor(FEType type) const
Definition: Assembly.h:1036
Class for scalar variables (they are different).
void computeFaceMap(unsigned dim, const std::vector< Real > &qw, const Elem *side)
Definition: Assembly.C:1107
const VectorVariablePhiSecond & secondPhi(const VectorMooseVariable &) const
Definition: Assembly.h:819
void computeSinglePointMapAD(const Elem *elem, const std::vector< Real > &qw, unsigned p, FEBase *fe)
Definition: Assembly.C:829
void addJacobianScalar()
Definition: Assembly.C:3414
const QBase * _const_current_qrule_neighbor
quadrature rule used on neighbors
Definition: Assembly.h:1448
VectorVariablePhiSecond _vector_second_phi_face_neighbor
Definition: Assembly.h:1581
std::vector< DualReal > _ad_detadx_map
Definition: Assembly.h:1673
void buildLowerDFE(FEType type) const
Build FEs for a lower dimensional element with a type.
Definition: Assembly.C:317
unsigned int _current_side
The current side of the selected element (valid only when working with sides)
Definition: Assembly.h:1481
VariablePhiSecond _second_phi
Definition: Assembly.h:1589
const OutputTools< OutputType >::VariablePhiSecond & feSecondPhiNeighbor(FEType type) const
Definition: Assembly.h:1028
const FEVectorBase *const & getVectorFEFaceNeighbor(FEType type, unsigned int dim) const
GetVector a reference to a pointer that will contain the current "neighbor" FE.
Definition: Assembly.h:158
const OutputTools< OutputType >::VariablePhiSecond & feSecondPhiFace(FEType type) const
Definition: Assembly.h:1006
const Elem * _current_neighbor_side_elem
The current side element of the ncurrent neighbor element.
Definition: Assembly.h:1493
std::map< unsigned int, QBase * > _holder_qrule_volume
Holds volume qrules for each dimension.
Definition: Assembly.h:1376
void addCachedJacobianContributions()
Adds previously-cached Jacobian values via SparseMatrix::add() calls.
Definition: Assembly.C:3488
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2033
void reinitElemFaceRef(const Elem *elem, unsigned int elem_side, Real tolerance, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)
Reinitialize FE data for the given element on the given side, optionally with a given set of referenc...
Definition: Assembly.C:1699
const Real & neighborVolume() const
Returns the reference to the current neighbor volume.
Definition: Assembly.C:480
bool _current_side_volume_computed
Boolean to indicate whether current element side volumes has been computed.
Definition: Assembly.h:1505
void cacheResidualContribution(dof_id_type dof, Real value, TagID tag_id)
Cache individual residual contributions.
Definition: Assembly.C:2825
std::vector< VectorValue< DualReal > > _ad_d2xyzdeta2_map
Definition: Assembly.h:1666
const Real & sideElemVolume()
Returns the reference to the volume of current side element.
Definition: Assembly.h:323
const VariablePhiValue & phiNeighbor(const MooseVariable &) const
Definition: Assembly.h:791
const SubdomainID & currentSubdomainID() const
Return the current subdomain ID.
Definition: Assembly.h:288
void prepareBlock(unsigned int ivar, unsigned jvar, const std::vector< dof_id_type > &dof_indices)
Definition: Assembly.C:2515
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_face_neighbor
Definition: Assembly.h:1425
std::vector< std::vector< std::vector< unsigned char > > > _jacobian_block_lower_used
Flag that indicates if the jacobian block for the lower dimensional element was used.
Definition: Assembly.h:1309
OutputTools< RealVectorValue >::VariablePhiCurl VectorVariablePhiCurl
Definition: MooseTypes.h:218
const VectorVariablePhiGradient & gradPhiFace(const VectorMooseVariable &) const
Definition: Assembly.h:832
void cacheResidualBlock(std::vector< Real > &cached_residual_values, std::vector< dof_id_type > &cached_residual_rows, DenseVector< Number > &res_block, const std::vector< dof_id_type > &dof_indices, Real scaling_factor)
Definition: Assembly.C:2719
const Moose::CoordinateSystemType & coordSystem()
Get the coordinate system type.
Definition: Assembly.h:229
const unsigned int & side() const
Returns the current side.
Definition: Assembly.h:305
const FEBase *const & getFE(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current volume FE.
Definition: Assembly.h:74
void computeCurrentElemVolume()
Definition: Assembly.C:1510
std::map< FEType, FEShapeData * > _fe_shape_data
Shape function values, gradients, second derivatives for each FE type.
Definition: Assembly.h:1603
const VariableTestGradientType< OutputType, ComputeStage::JACOBIAN >::type & feADGradPhi(FEType type) const
Definition: Assembly.h:965
std::map< FEType, FEShapeData * > _fe_shape_data_lower
Definition: Assembly.h:1607
DenseVector< Number > & residualBlock(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:720
MooseArray< Point > _current_q_points
The current list of quadrature points.
Definition: Assembly.h:1365
const VectorVariablePhiCurl & curlPhiFaceNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:874
Moose::CoordinateSystemType _coord_type
The coordinate system.
Definition: Assembly.h:1369
void cacheResidualLower()
Takes the values that are currently in _sub_Rl and appends them to the cached values.
Definition: Assembly.C:2857
void buildVectorFaceNeighborFE(FEType type) const
Build Vector FEs for a neighbor face with a type.
Definition: Assembly.C:451
const SubdomainID & currentNeighborSubdomainID() const
Return the current subdomain ID.
Definition: Assembly.h:340
VariablePhiGradient & gradPhiFace(const MooseVariable &)
Definition: Assembly.h:885
VectorVariablePhiSecond _vector_second_phi_neighbor
Definition: Assembly.h:1576
const Node *const & node() const
Returns the reference to the node.
Definition: Assembly.h:381
std::vector< std::vector< numeric_index_type > > _cached_jacobian_contribution_rows
Definition: Assembly.h:1657
VariablePhiSecond _second_phi
Definition: Assembly.h:1549
std::vector< std::vector< std::vector< unsigned char > > > _jacobian_block_neighbor_used
Flag that indicates if the jacobian block for neighbor was used.
Definition: Assembly.h:1307
const OutputTools< OutputType >::VariablePhiGradient & feGradPhiNeighbor(FEType type) const
Definition: Assembly.h:1021
std::map< unsigned int, const std::vector< Point > * > _holder_q_points
Holds pointers to the dimension&#39;s q_points.
Definition: Assembly.h:1382
void setCurrentNeighborSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:345
std::map< FEType, FEShapeData * > _fe_shape_data_face
Definition: Assembly.h:1604
std::vector< DualReal > _ad_detady_map
Definition: Assembly.h:1674
std::vector< std::vector< dof_id_type > > _cached_jacobian_cols
Column where the corresponding cached value should go.
Definition: Assembly.h:1640
void clearCachedJacobianContributions()
Clear any currently cached jacobian contributions.
Definition: Assembly.C:3508
const VariablePhiGradient & gradPhiNeighbor(const MooseVariable &) const
Definition: Assembly.h:792
MooseArray< Point > _current_q_points_face
The current quadrature points on a face.
Definition: Assembly.h:1407
std::map< FEType, VectorFEShapeData * > _vector_fe_shape_data_face
Definition: Assembly.h:1611
void prepareBlockNonlocal(unsigned int ivar, unsigned jvar, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices)
Definition: Assembly.C:2534
const Elem *const & lowerDElem() const
Return the lower dimensional element.
Definition: Assembly.h:335
const Elem * _current_side_elem
The current "element" making up the side we are currently on.
Definition: Assembly.h:1483
MooseArray< DualReal > _ad_coord
The AD version of the current coordinate transformation coefficients.
Definition: Assembly.h:1373
const OutputTools< OutputType >::VariablePhiSecond & feSecondPhi(FEType type) const
Definition: Assembly.h:971
FEBase * _current_fe_helper
The current helper object for transforming coordinates.
Definition: Assembly.h:1353
const VariablePhiGradient & gradPhiFace(const MooseVariable &) const
Definition: Assembly.h:788
std::map< FEType, bool > _need_curl
Definition: Assembly.h:1694
const VectorVariablePhiSecond & secondPhiNeighbor(const VectorMooseVariable &) const
Definition: Assembly.h:853
VectorVariablePhiSecond _vector_second_phi_face
Definition: Assembly.h:1571
void cacheJacobianNonlocal()
Takes the values that are currently in _sub_Keg and appends them to the cached values.
Definition: Assembly.C:3266
VariablePhiGradient _grad_phi_face
Definition: Assembly.h:1552
MooseArray< Point > _current_q_points_face_neighbor
The current quadrature points on the neighbor face.
Definition: Assembly.h:1452
FEGenericBase< VectorValue< Real > > FEVectorBase
Definition: Assembly.h:36
void copyFaceShapes(MooseVariableFE< T > &v)
Definition: Assembly.C:2633
void addCachedResiduals()
Definition: Assembly.C:2890
VectorVariablePhiCurl _vector_curl_phi_face
Definition: Assembly.h:1572
std::map< unsigned int, std::map< FEType, const FEVectorBase * > > _const_vector_fe_face
types of vector finite elements
Definition: Assembly.h:1395
const std::vector< Real > & jxWMortar() const
Returns a reference to JxW for mortar segment elements.
Definition: Assembly.h:466
const FEVectorBase *const & getVectorFENeighbor(FEType type, unsigned int dim) const
GetVector a reference to a pointer that will contain the current &#39;neighbor&#39; FE.
Definition: Assembly.h:134
SubdomainID _current_subdomain_id
The current subdomain ID.
Definition: Assembly.h:1477
void addResidualBlock(NumericVector< Number > &residual, DenseVector< Number > &res_block, const std::vector< dof_id_type > &dof_indices, Real scaling_factor)
Definition: Assembly.C:2695
const Elem *& sideElem()
Returns the side element.
Definition: Assembly.h:317
const FEBase *const & getFEFaceNeighbor(FEType type, unsigned int dim) const
Get a reference to a pointer that will contain the current "neighbor" FE.
Definition: Assembly.h:110
std::map< unsigned int, const std::vector< Point > * > _holder_normals
Holds pointers to the dimension&#39;s normal vectors.
Definition: Assembly.h:1419
QBase * _current_qrule_face
quadrature rule used on faces
Definition: Assembly.h:1403
const VariablePhiValue & phi(const MooseVariable &) const
Definition: Assembly.h:779
std::map< FEType, FEVectorBase * > _current_vector_fe_neighbor
The "neighbor" vector fe object that matches the current elem.
Definition: Assembly.h:1336
unsigned int THREAD_ID
Definition: MooseTypes.h:161
const Node * _current_neighbor_node
The current neighboring node we are working with.
Definition: Assembly.h:1501
const VectorVariablePhiSecond & secondPhiFace(const VectorMooseVariable &) const
Definition: Assembly.h:836
std::map< unsigned int, const std::vector< Real > * > _holder_JxW
Holds pointers to the dimension&#39;s transformed jacobian weights.
Definition: Assembly.h:1384
const VariablePhiGradient & gradPhiFace() const
Definition: Assembly.h:787
std::vector< std::vector< std::vector< DenseMatrix< Number > > > > _sub_Knl
dmaster/dlower (or dneighbor/dlower)
Definition: Assembly.h:1541
const OutputTools< OutputType >::VariablePhiValue & fePhiFace(FEType type) const
Definition: Assembly.h:985
void addJacobianBlockNonlocal(SparseMatrix< Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices)
Definition: Assembly.C:3345
void setCachedJacobianContributions()
Sets previously-cached Jacobian values via SparseMatrix::set() calls.
Definition: Assembly.C:3452
std::vector< VectorValue< DualReal > > _ad_dxyzdxi_map
AD quantities.
Definition: Assembly.h:1661
void reinitFEFace(const Elem *elem, unsigned int side)
Just an internal helper function to reinit the face FE objects.
Definition: Assembly.C:1046