www.mooseframework.org
Assembly.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #ifndef ASSEMBLY_H
16 #define ASSEMBLY_H
17 
18 #include "MooseArray.h"
19 #include "MooseTypes.h"
20 
21 // libMesh
22 #include "libmesh/dense_matrix.h"
23 #include "libmesh/dense_vector.h"
24 #include "libmesh/enum_quadrature_type.h"
25 #include "libmesh/fe_type.h"
26 
27 // libMesh forward declarations
28 namespace libMesh
29 {
30 class DofMap;
31 class CouplingMatrix;
32 class Elem;
33 template <typename T>
36 class Node;
37 template <typename T>
38 class NumericVector;
39 template <typename T>
41 
42 typedef VectorValue<Real> RealVectorValue;
43 typedef RealVectorValue RealGradient;
44 
45 typedef TensorValue<Real> RealTensorValue;
46 typedef RealTensorValue RealTensor;
47 }
48 
49 // MOOSE Forward Declares
50 class MooseMesh;
52 class SystemBase;
53 class MooseVariable;
58 
63 class Assembly
64 {
65 public:
66  Assembly(SystemBase & sys, THREAD_ID tid);
67  virtual ~Assembly();
68 
73  void buildFE(FEType type);
74 
79  void buildFaceFE(FEType type);
80 
85  void buildNeighborFE(FEType type);
86 
91  void buildFaceNeighborFE(FEType type);
92 
99  FEBase *& getFE(FEType type, unsigned int dim);
100 
107  FEBase *& getFEFace(FEType type, unsigned int dim);
108 
115  FEBase *& getFENeighbor(FEType type, unsigned int dim);
116 
123  FEBase *& getFEFaceNeighbor(FEType type, unsigned int dim);
124 
129  QBase *& qRule() { return _current_qrule; }
130 
135  const MooseArray<Point> & qPoints() { return _current_q_points; }
136 
141  const MooseArray<Point> & physicalPoints() { return _current_physical_points; }
142 
147  const MooseArray<Real> & JxW() { return _current_JxW; }
148 
153  const MooseArray<Real> & coordTransformation() { return _coord; }
154 
159  const Moose::CoordinateSystemType & coordSystem() { return _coord_type; }
160 
165  QBase *& qRuleFace() { return _current_qrule_face; }
166 
171  const MooseArray<Point> & qPointsFace() { return _current_q_points_face; }
172 
177  const MooseArray<Real> & JxWFace() { return _current_JxW_face; }
178 
183  const MooseArray<Point> & normals() { return _current_normals; }
184 
189  const Elem *& elem() { return _current_elem; }
190 
194  const SubdomainID & currentSubdomainID() const { return _current_subdomain_id; }
195 
199  void setCurrentSubdomainID(SubdomainID i) { _current_subdomain_id = i; }
200 
205  const Real & elemVolume() { return _current_elem_volume; }
206 
211  unsigned int & side() { return _current_side; }
212 
217  unsigned int & neighborSide() { return _current_neighbor_side; }
218 
223  const Elem *& sideElem() { return _current_side_elem; }
224 
229  const Real & sideElemVolume() { return _current_side_volume; }
230 
235  const Elem *& neighbor() { return _current_neighbor_elem; }
236 
240  const SubdomainID & currentNeighborSubdomainID() const { return _current_neighbor_subdomain_id; }
241 
245  void setCurrentNeighborSubdomainID(SubdomainID i) { _current_neighbor_subdomain_id = i; }
246 
251  const Real & neighborVolume();
252 
257  QBase *& qRuleNeighbor() { return _current_qrule_neighbor; }
258 
263  const MooseArray<Real> & JxWNeighbor() { return _current_JxW_neighbor; }
264 
269  const Node *& node() { return _current_node; }
270 
275  const Node *& nodeNeighbor() { return _current_neighbor_node; }
276 
280  void createQRules(QuadratureType type, Order order, Order volume_order, Order face_order);
281 
290  void setVolumeQRule(QBase * qrule, unsigned int dim);
291 
300  void setFaceQRule(QBase * qrule, unsigned int dim);
301 
310  void setNeighborQRule(QBase * qrule, unsigned int dim);
311 
317  void reinit(const Elem * elem);
318 
322  void reinitAtPhysical(const Elem * elem, const std::vector<Point> & physical_points);
323 
327  void reinit(const Elem * elem, const std::vector<Point> & reference_points);
328 
332  void reinit(const Elem * elem, unsigned int side);
333 
342  void reinitElemAndNeighbor(const Elem * elem,
343  unsigned int side,
344  const Elem * neighbor,
345  unsigned int neighbor_side);
346 
350  void reinitNeighborAtPhysical(const Elem * neighbor,
351  unsigned int neighbor_side,
352  const std::vector<Point> & physical_points);
353 
357  void reinitNeighborAtPhysical(const Elem * neighbor, const std::vector<Point> & physical_points);
358 
359  void reinitNeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
360 
364  void reinit(const Node * node);
365 
369  void reinitNodeNeighbor(const Node * node);
370 
374  void init(const CouplingMatrix * cm);
375 
377  void init();
378 
380  void initNonlocalCoupling();
381 
387  void useFECache(bool fe_cache) { _should_use_fe_cache = fe_cache; }
388 
389  void prepare();
390  void prepareNonlocal();
391 
397  void prepareVariable(MooseVariable * var);
398  void prepareVariableNonlocal(MooseVariable * var);
399  void prepareNeighbor();
400  void prepareBlock(unsigned int ivar, unsigned jvar, const std::vector<dof_id_type> & dof_indices);
401  void prepareBlockNonlocal(unsigned int ivar,
402  unsigned jvar,
403  const std::vector<dof_id_type> & idof_indices,
404  const std::vector<dof_id_type> & jdof_indices);
405  void prepareScalar();
406  void prepareOffDiagScalar();
407 
408  void copyShapes(unsigned int var);
409  void copyFaceShapes(unsigned int var);
410  void copyNeighborShapes(unsigned int var);
411 
412  void addResidual(NumericVector<Number> & residual, Moose::KernelType type = Moose::KT_NONTIME);
413  void addResidualNeighbor(NumericVector<Number> & residual,
415  void addResidualScalar(NumericVector<Number> & residual,
417 
421  void cacheResidual();
422 
431  void cacheResidualContribution(dof_id_type dof, Real value, Moose::KernelType type);
432 
436  void cacheResidualNodes(DenseVector<Number> & res, std::vector<dof_id_type> & dof_index);
437 
441  void cacheResidualNeighbor();
442 
443  void addCachedResiduals();
444 
451  void addCachedResidual(NumericVector<Number> & residual, Moose::KernelType type);
452 
453  void setResidual(NumericVector<Number> & residual, Moose::KernelType type = Moose::KT_NONTIME);
454  void setResidualNeighbor(NumericVector<Number> & residual,
456 
457  void addJacobian(SparseMatrix<Number> & jacobian);
458  void addJacobianNonlocal(SparseMatrix<Number> & jacobian);
459  void addJacobianBlock(SparseMatrix<Number> & jacobian,
460  unsigned int ivar,
461  unsigned int jvar,
462  const DofMap & dof_map,
463  std::vector<dof_id_type> & dof_indices);
464  void addJacobianBlockNonlocal(SparseMatrix<Number> & jacobian,
465  unsigned int ivar,
466  unsigned int jvar,
467  const DofMap & dof_map,
468  const std::vector<dof_id_type> & idof_indices,
469  const std::vector<dof_id_type> & jdof_indices);
470  void addJacobianNeighbor(SparseMatrix<Number> & jacobian);
471  void addJacobianNeighbor(SparseMatrix<Number> & jacobian,
472  unsigned int ivar,
473  unsigned int jvar,
474  const DofMap & dof_map,
475  std::vector<dof_id_type> & dof_indices,
476  std::vector<dof_id_type> & neighbor_dof_indices);
477  void addJacobianScalar(SparseMatrix<Number> & jacobian);
478  void addJacobianOffDiagScalar(SparseMatrix<Number> & jacobian, unsigned int ivar);
479 
483  void cacheJacobian();
484 
488  void cacheJacobianNonlocal();
489 
494  void cacheJacobianNeighbor();
495 
502  void addCachedJacobian(SparseMatrix<Number> & jacobian);
503 
504  DenseVector<Number> & residualBlock(unsigned int var_num,
506  {
507  return _sub_Re[static_cast<unsigned int>(type)][var_num];
508  }
509  DenseVector<Number> & residualBlockNeighbor(unsigned int var_num,
511  {
512  return _sub_Rn[static_cast<unsigned int>(type)][var_num];
513  }
514 
515  DenseMatrix<Number> & jacobianBlock(unsigned int ivar, unsigned int jvar);
516  DenseMatrix<Number> & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar);
517  DenseMatrix<Number> &
518  jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar);
519  void cacheJacobianBlock(DenseMatrix<Number> & jac_block,
520  std::vector<dof_id_type> & idof_indices,
521  std::vector<dof_id_type> & jdof_indices,
522  Real scaling_factor);
523  void cacheJacobianBlockNonlocal(DenseMatrix<Number> & jac_block,
524  const std::vector<dof_id_type> & idof_indices,
525  const std::vector<dof_id_type> & jdof_indices,
526  Real scaling_factor);
527 
528  std::vector<std::pair<MooseVariable *, MooseVariable *>> & couplingEntries() { return _cm_entry; }
529  std::vector<std::pair<MooseVariable *, MooseVariable *>> & nonlocalCouplingEntries()
530  {
531  return _cm_nonlocal_entry;
532  }
533 
534  const VariablePhiValue & phi() { return _phi; }
535  const VariablePhiGradient & gradPhi() { return _grad_phi; }
536  const VariablePhiSecond & secondPhi() { return _second_phi; }
537 
538  const VariablePhiValue & phiFace() { return _phi_face; }
539  const VariablePhiGradient & gradPhiFace() { return _grad_phi_face; }
540  const VariablePhiSecond & secondPhiFace() { return _second_phi_face; }
541 
542  const VariablePhiValue & phiNeighbor() { return _phi_neighbor; }
543  const VariablePhiGradient & gradPhiNeighbor() { return _grad_phi_neighbor; }
544  const VariablePhiSecond & secondPhiNeighbor() { return _second_phi_neighbor; }
545 
546  const VariablePhiValue & phiFaceNeighbor() { return _phi_face_neighbor; }
547  const VariablePhiGradient & gradPhiFaceNeighbor() { return _grad_phi_face_neighbor; }
548  const VariablePhiSecond & secondPhiFaceNeighbor() { return _second_phi_face_neighbor; }
549 
550  const VariablePhiValue & fePhi(FEType type);
551  const VariablePhiGradient & feGradPhi(FEType type);
552  const VariablePhiSecond & feSecondPhi(FEType type);
553 
554  const VariablePhiValue & fePhiFace(FEType type);
555  const VariablePhiGradient & feGradPhiFace(FEType type);
556  const VariablePhiSecond & feSecondPhiFace(FEType type);
557 
558  const VariablePhiValue & fePhiNeighbor(FEType type);
559  const VariablePhiGradient & feGradPhiNeighbor(FEType type);
560  const VariablePhiSecond & feSecondPhiNeighbor(FEType type);
561 
562  const VariablePhiValue & fePhiFaceNeighbor(FEType type);
563  const VariablePhiGradient & feGradPhiFaceNeighbor(FEType type);
564  const VariablePhiSecond & feSecondPhiFaceNeighbor(FEType type);
565 
569  void invalidateCache();
570 
571  std::map<FEType, bool> _need_second_derivative;
572 
581  void cacheJacobianContribution(numeric_index_type i, numeric_index_type j, Real value);
582 
586  void setCachedJacobianContributions(SparseMatrix<Number> & jacobian);
587 
591  void zeroCachedJacobianContributions(SparseMatrix<Number> & jacobian);
592 
596  void addCachedJacobianContributions(SparseMatrix<Number> & jacobian);
597 
601  void setXFEM(std::shared_ptr<XFEMInterface> xfem) { _xfem = xfem; }
602 
603 protected:
609  void reinitFE(const Elem * elem);
610 
617  void reinitFEFace(const Elem * elem, unsigned int side);
618 
619  void reinitFEFaceNeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
620 
621  void reinitFENeighbor(const Elem * neighbor, const std::vector<Point> & reference_points);
622 
623  void setCoordinateTransformation(const QBase * qrule, const MooseArray<Point> & q_points);
624 
625  void computeCurrentElemVolume();
626 
627  void computeCurrentFaceVolume();
628 
629  void computeCurrentNeighborVolume();
630 
631  void addResidualBlock(NumericVector<Number> & residual,
632  DenseVector<Number> & res_block,
633  const std::vector<dof_id_type> & dof_indices,
634  Real scaling_factor);
635  void cacheResidualBlock(std::vector<Real> & cached_residual_values,
636  std::vector<dof_id_type> & cached_residual_rows,
637  DenseVector<Number> & res_block,
638  std::vector<dof_id_type> & dof_indices,
639  Real scaling_factor);
640 
641  void setResidualBlock(NumericVector<Number> & residual,
642  DenseVector<Number> & res_block,
643  std::vector<dof_id_type> & dof_indices,
644  Real scaling_factor);
645 
646  void addJacobianBlock(SparseMatrix<Number> & jacobian,
647  DenseMatrix<Number> & jac_block,
648  const std::vector<dof_id_type> & idof_indices,
649  const std::vector<dof_id_type> & jdof_indices,
650  Real scaling_factor);
651 
658  void clearCachedJacobianContributions();
659 
665  void modifyWeightsDueToXFEM(const Elem * elem);
666 
668 
670  const CouplingMatrix * _cm;
671  const CouplingMatrix & _nonlocal_cm;
673  std::vector<std::pair<MooseVariable *, MooseVariable *>> _cm_entry;
674  std::vector<std::pair<MooseVariable *, MooseVariable *>> _cm_nonlocal_entry;
676  std::vector<std::vector<unsigned char>> _jacobian_block_used;
677  std::vector<std::vector<unsigned char>> _jacobian_block_nonlocal_used;
679  std::vector<std::vector<unsigned char>> _jacobian_block_neighbor_used;
681  const DofMap & _dof_map;
684 
686 
687  unsigned int _mesh_dimension;
688 
690  std::shared_ptr<XFEMInterface> _xfem;
691 
693  std::map<FEType, FEBase *> _current_fe;
695  std::map<FEType, FEBase *> _current_fe_face;
697  std::map<FEType, FEBase *> _current_fe_neighbor;
699  std::map<FEType, FEBase *> _current_fe_face_neighbor;
700 
701  /**** Volume Stuff ****/
702 
704  std::map<unsigned int, std::map<FEType, FEBase *>> _fe;
706  std::map<unsigned int, FEBase **> _holder_fe_helper;
710  QBase * _current_qrule;
724  std::map<unsigned int, QBase *> _holder_qrule_volume;
726  std::map<unsigned int, ArbitraryQuadrature *> _holder_qrule_arbitrary;
728  std::map<unsigned int, const std::vector<Point> *> _holder_q_points;
730  std::map<unsigned int, const std::vector<Real> *> _holder_JxW;
731 
732  /**** Face Stuff ****/
733 
735  std::map<unsigned int, std::map<FEType, FEBase *>> _fe_face;
737  std::map<unsigned int, FEBase **> _holder_fe_face_helper;
751  std::map<unsigned int, QBase *> _holder_qrule_face;
753  std::map<unsigned int, ArbitraryQuadrature *> _holder_qface_arbitrary;
755  std::map<unsigned int, const std::vector<Point> *> _holder_q_points_face;
757  std::map<unsigned int, const std::vector<Real> *> _holder_JxW_face;
759  std::map<unsigned int, const std::vector<Point> *> _holder_normals;
760 
761  /**** Neighbor Stuff ****/
762 
764  std::map<unsigned int, std::map<FEType, FEBase *>> _fe_neighbor;
765  std::map<unsigned int, std::map<FEType, FEBase *>> _fe_face_neighbor;
767  std::map<unsigned int, FEBase **> _holder_fe_neighbor_helper;
768  std::map<unsigned int, FEBase **> _holder_fe_face_neighbor_helper;
769 
773  std::map<unsigned int, ArbitraryQuadrature *> _holder_qrule_neighbor;
778 
780  const Elem * _current_elem;
786  unsigned int _current_side;
788  const Elem * _current_side_elem;
804  const Node * _current_node;
811 
814 
816  std::vector<std::vector<DenseVector<Number>>> _sub_Re;
818  std::vector<std::vector<DenseVector<Number>>> _sub_Rn;
820  DenseVector<Number> _tmp_Re;
821 
823  std::vector<std::vector<DenseMatrix<Number>>> _sub_Kee;
824  std::vector<std::vector<DenseMatrix<Number>>> _sub_Keg;
825 
827  std::vector<std::vector<DenseMatrix<Number>>> _sub_Ken;
829  std::vector<std::vector<DenseMatrix<Number>>> _sub_Kne;
831  std::vector<std::vector<DenseMatrix<Number>>> _sub_Knn;
832 
834  DenseMatrix<Number> _tmp_Ke;
835 
836  // Shape function values, gradients. second derivatives
840 
844 
848 
852 
854  {
855  public:
859  };
860 
871  {
872  public:
874  std::map<FEType, FEShapeData *> _shape_data;
875 
878 
881 
884  };
885 
887  std::map<dof_id_type, ElementFEShapeData *> _element_fe_shape_data_cache;
888 
891 
894 
895  // Shape function values, gradients. second derivatives for each FE type
896  std::map<FEType, FEShapeData *> _fe_shape_data;
897  std::map<FEType, FEShapeData *> _fe_shape_data_face;
898  std::map<FEType, FEShapeData *> _fe_shape_data_neighbor;
899  std::map<FEType, FEShapeData *> _fe_shape_data_face_neighbor;
900 
902  std::vector<std::vector<Real>> _cached_residual_values;
903 
905  std::vector<std::vector<dof_id_type>> _cached_residual_rows;
906 
907  unsigned int _max_cached_residuals;
908 
910  std::vector<Real> _cached_jacobian_values;
912  std::vector<dof_id_type> _cached_jacobian_rows;
914  std::vector<dof_id_type> _cached_jacobian_cols;
915 
916  unsigned int _max_cached_jacobians;
917 
920 
922  std::vector<dof_id_type> _temp_dof_indices;
923 
925  std::vector<Point> _temp_reference_points;
926 
931  std::vector<numeric_index_type> _cached_jacobian_contribution_rows;
932  std::vector<numeric_index_type> _cached_jacobian_contribution_cols;
933 };
934 
935 #endif /* ASSEMBLY_H */
VariablePhiSecond _second_phi_face
Definition: Assembly.h:843
VariablePhiSecond _second_phi_face_neighbor
Definition: Assembly.h:851
Ok - here&#39;s the design.
Definition: Assembly.h:870
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_face
types of finite elements
Definition: Assembly.h:735
bool _need_neighbor_elem_volume
true is apps need to compute neighbor element volume
Definition: Assembly.h:800
QBase *& qRuleFace()
Returns the reference to the current quadrature being used on a current face.
Definition: Assembly.h:165
ArbitraryQuadrature * _current_qrule_arbitrary
The current arbitrary quadrature rule used within the element interior.
Definition: Assembly.h:714
std::map< FEType, FEBase * > _current_fe
The "volume" fe object that matches the current elem.
Definition: Assembly.h:693
RealVectorValue RealGradient
Definition: Assembly.h:43
std::map< unsigned int, FEBase ** > _holder_fe_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:706
SystemBase & _sys
Definition: Assembly.h:667
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:905
unsigned int _max_cached_residuals
Definition: Assembly.h:907
MooseArray< Point > _current_physical_points
This will be filled up with the physical points passed into reinitAtPhysical() if it is called...
Definition: Assembly.h:813
std::vector< std::vector< unsigned char > > _jacobian_block_nonlocal_used
Definition: Assembly.h:677
const VariablePhiValue & phi()
Definition: Assembly.h:534
MooseArray< Real > _coord_neighbor
The current coordinate transformation coefficients.
Definition: Assembly.h:777
QBase * _current_qrule
The current current quadrature rule being used (could be either volumetric or arbitrary - for dirac k...
Definition: Assembly.h:710
std::map< unsigned int, FEBase ** > _holder_fe_face_neighbor_helper
Definition: Assembly.h:768
const MooseArray< Point > & qPointsFace()
Returns the reference to the current quadrature being used.
Definition: Assembly.h:171
VectorValue< Real > RealVectorValue
Definition: Assembly.h:40
Keeps track of stuff related to assembling.
Definition: Assembly.h:63
Class for stuff related to variables.
Definition: MooseVariable.h:43
FEGenericBase< Real > FEBase
Definition: Assembly.h:34
MooseArray< Real > _current_JxW_neighbor
The current transformed jacobian weights on a neighbor&#39;s face.
Definition: Assembly.h:775
std::shared_ptr< XFEMInterface > _xfem
The XFEM controller.
Definition: Assembly.h:690
MooseArray< std::vector< Real > > _phi
Definition: Assembly.h:856
subdomain_id_type SubdomainID
Definition: MooseTypes.h:77
MooseArray< Real > _coord
The current coordinate transformation coefficients.
Definition: Assembly.h:722
VariablePhiGradient _grad_phi
Definition: Assembly.h:838
VariablePhiSecond _second_phi_neighbor
Definition: Assembly.h:847
void setXFEM(std::shared_ptr< XFEMInterface > xfem)
Set the pointer to the XFEM controller object.
Definition: Assembly.h:601
std::map< unsigned int, ArbitraryQuadrature * > _holder_qrule_neighbor
Holds arbitrary qrules for each dimension.
Definition: Assembly.h:773
bool _current_elem_volume_computed
Boolean to indicate whether current element volumes has been computed.
Definition: Assembly.h:808
bool _currently_fe_caching
Whether or not fe should currently be cached - This will be false if something funky is going on with...
Definition: Assembly.h:893
Real _current_neighbor_volume
Volume of the current neighbor.
Definition: Assembly.h:802
std::map< FEType, FEBase * > _current_fe_face
The "face" fe object that matches the current elem.
Definition: Assembly.h:695
std::vector< std::pair< MooseVariable *, MooseVariable * > > _cm_entry
Entries in the coupling matrix (only for field variables)
Definition: Assembly.h:673
const Elem * _current_neighbor_elem
The current neighbor "element".
Definition: Assembly.h:792
std::vector< std::vector< DenseVector< Number > > > _sub_Rn
residual contributions for each variable from the neighbor
Definition: Assembly.h:818
MooseMesh & _mesh
Definition: Assembly.h:685
MooseArray< Real > _current_JxW_face
The current transformed jacobian weights on a face.
Definition: Assembly.h:747
bool _invalidated
Whether or not this data is invalid (needs to be recached) note that there is no constructor so the v...
Definition: Assembly.h:877
DenseMatrix< Number > _tmp_Ke
auxiliary matrix for scaling jacobians (optimization to avoid expensive construction/destruction) ...
Definition: Assembly.h:834
const Elem * _current_elem
The current "element" we are currently on.
Definition: Assembly.h:780
THREAD_ID _tid
Thread number (id)
Definition: Assembly.h:683
const Elem *& neighbor()
Return the neighbor element.
Definition: Assembly.h:235
const MooseArray< Real > & coordTransformation()
Returns the reference to the coordinate transformation coefficients.
Definition: Assembly.h:153
const MooseArray< Real > & JxW()
Returns the reference to the transformed jacobian weights.
Definition: Assembly.h:147
RealTensorValue RealTensor
Definition: Assembly.h:46
const MooseArray< Point > & normals()
Returns the array of normals for quadrature points on a current side.
Definition: Assembly.h:183
std::map< FEType, bool > _need_second_derivative
Definition: Assembly.h:571
std::vector< numeric_index_type > _cached_jacobian_contribution_cols
Definition: Assembly.h:932
Real _current_elem_volume
Volume of the current element.
Definition: Assembly.h:784
std::map< FEType, FEBase * > _current_fe_face_neighbor
The "neighbor face" fe object that matches the current elem.
Definition: Assembly.h:699
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:899
Base class for a system (of equations)
Definition: SystemBase.h:91
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:919
const MooseArray< Point > & qPoints()
Returns the reference to the quadrature points.
Definition: Assembly.h:135
VariablePhiValue _phi
Definition: Assembly.h:837
MooseArray< Point > _q_points
Cached xyz positions of quadrature points.
Definition: Assembly.h:883
const VariablePhiSecond & secondPhi()
Definition: Assembly.h:536
QBase *& qRule()
Returns the reference to the current quadrature being used.
Definition: Assembly.h:129
const CouplingMatrix * _cm
Coupling matrices.
Definition: Assembly.h:670
std::vector< std::vector< DenseMatrix< Number > > > _sub_Kee
jacobian contributions
Definition: Assembly.h:823
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
DenseVector< Number > _tmp_Re
auxiliary vector for scaling residuals (optimization to avoid expensive construction/destruction) ...
Definition: Assembly.h:820
unsigned int _mesh_dimension
Definition: Assembly.h:687
QBase * _current_qrule_volume
The current volumetric quadrature for the element.
Definition: Assembly.h:712
const DofMap & _dof_map
DOF map.
Definition: Assembly.h:681
MooseArray< std::vector< RealGradient > > _grad_phi
Definition: Assembly.h:857
unsigned int _max_cached_jacobians
Definition: Assembly.h:916
std::vector< std::vector< DenseMatrix< Number > > > _sub_Kne
jacobian contributions from the neighbor and element
Definition: Assembly.h:829
const VariablePhiGradient & gradPhi()
Definition: Assembly.h:535
std::map< dof_id_type, ElementFEShapeData * > _element_fe_shape_data_cache
Cached shape function values stored by element.
Definition: Assembly.h:887
const CouplingMatrix & _nonlocal_cm
Definition: Assembly.h:671
std::vector< Point > _temp_reference_points
Temporary work data for reinitAtPhysical()
Definition: Assembly.h:925
unsigned int _current_neighbor_side
The current side of the selected neighboring element (valid only when working with sides) ...
Definition: Assembly.h:796
unsigned int & side()
Returns the current side.
Definition: Assembly.h:211
const Node *& node()
Returns the reference to the node.
Definition: Assembly.h:269
const SubdomainID & currentNeighborSubdomainID() const
Return the current subdomain ID.
Definition: Assembly.h:240
std::map< unsigned int, QBase * > _holder_qrule_face
Holds face qrules for each dimension.
Definition: Assembly.h:751
const VariablePhiSecond & secondPhiNeighbor()
Definition: Assembly.h:544
Implements a fake quadrature rule where you can specify the locations (in the reference domain) of th...
std::map< FEType, FEShapeData * > _fe_shape_data_neighbor
Definition: Assembly.h:898
QBase * _current_qrule_neighbor
quadrature rule used on neighbors
Definition: Assembly.h:771
const SubdomainID & currentSubdomainID() const
Return the current subdomain ID.
Definition: Assembly.h:194
FEBase * _current_fe_face_helper
helper object for transforming coordinates
Definition: Assembly.h:739
MooseArray< std::vector< Real > > VariablePhiValue
Definition: Assembly.h:54
std::vector< dof_id_type > _cached_jacobian_cols
Column where the corresponding cached value should go.
Definition: Assembly.h:914
const VariablePhiGradient & gradPhiFace()
Definition: Assembly.h:539
std::map< FEType, FEBase * > _current_fe_neighbor
The "neighbor" fe object that matches the current elem.
Definition: Assembly.h:697
SubdomainID _current_neighbor_subdomain_id
The current neighbor subdomain ID.
Definition: Assembly.h:794
std::map< unsigned int, FEBase ** > _holder_fe_neighbor_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:767
std::vector< std::vector< DenseVector< Number > > > _sub_Re
residual contributions for each variable from the element
Definition: Assembly.h:816
DofMap & dof_map
std::vector< std::vector< DenseMatrix< Number > > > _sub_Knn
jacobian contributions from the neighbor
Definition: Assembly.h:831
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
std::vector< std::vector< DenseMatrix< Number > > > _sub_Ken
jacobian contributions from the element and neighbor
Definition: Assembly.h:827
std::vector< numeric_index_type > _cached_jacobian_contribution_rows
Definition: Assembly.h:931
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_neighbor
types of finite elements
Definition: Assembly.h:764
VariablePhiGradient _grad_phi_face_neighbor
Definition: Assembly.h:850
MooseArray< Real > _current_JxW
The current list of transformed jacobian weights.
Definition: Assembly.h:718
std::vector< dof_id_type > _temp_dof_indices
Temporary work vector to keep from reallocating it.
Definition: Assembly.h:922
std::map< unsigned int, ArbitraryQuadrature * > _holder_qrule_arbitrary
Holds arbitrary qrules for each dimension.
Definition: Assembly.h:726
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:755
const MooseArray< Real > & JxWNeighbor()
Returns the reference to the transformed jacobian weights on a current face.
Definition: Assembly.h:263
VariablePhiGradient _grad_phi_neighbor
Definition: Assembly.h:846
const Real & elemVolume()
Returns the reference to the current element volume.
Definition: Assembly.h:205
ArbitraryQuadrature * _current_qface_arbitrary
The current arbitrary quadrature rule used on element faces.
Definition: Assembly.h:743
VariablePhiValue _phi_face
Definition: Assembly.h:841
const VariablePhiGradient & gradPhiNeighbor()
Definition: Assembly.h:543
std::vector< std::pair< MooseVariable *, MooseVariable * > > & nonlocalCouplingEntries()
Definition: Assembly.h:529
std::vector< Real > _cached_jacobian_values
Values cached by calling cacheJacobian()
Definition: Assembly.h:910
KernelType
Definition: MooseTypes.h:162
std::map< unsigned int, FEBase ** > _holder_fe_face_helper
Each dimension&#39;s helper objects.
Definition: Assembly.h:737
VariablePhiValue _phi_face_neighbor
Definition: Assembly.h:849
std::map< unsigned int, std::map< FEType, FEBase * > > _fe
Each dimension&#39;s actual fe objects indexed on type.
Definition: Assembly.h:704
std::vector< dof_id_type > _cached_jacobian_rows
Row where the corresponding cached value should go.
Definition: Assembly.h:912
VariablePhiValue _phi_neighbor
Definition: Assembly.h:845
std::map< unsigned int, ArbitraryQuadrature * > _holder_qface_arbitrary
Holds arbitrary face qrules for each dimension.
Definition: Assembly.h:753
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:757
Real _current_side_volume
Volume of the current side element.
Definition: Assembly.h:790
CoordinateSystemType
Definition: MooseTypes.h:212
std::vector< std::vector< Real > > _cached_residual_values
Values cached by calling cacheResidual() (the first vector is for TIME vs NONTIME) ...
Definition: Assembly.h:902
const MooseArray< Real > & JxWFace()
Returns the reference to the transformed jacobian weights on a current face.
Definition: Assembly.h:177
const VariablePhiGradient & gradPhiFaceNeighbor()
Definition: Assembly.h:547
const Node * _current_node
The current node we are working with.
Definition: Assembly.h:804
const Node *& nodeNeighbor()
Returns the reference to the neighboring node.
Definition: Assembly.h:275
This is the XFEMInterface class.
Definition: XFEMInterface.h:43
DGJacobianType
Definition: MooseTypes.h:190
std::vector< Real > _cached_jacobian_contribution_vals
Storage for cached Jacobian entries.
Definition: Assembly.h:930
MooseArray< Point > _current_normals
The current Normal vectors at the quadrature points.
Definition: Assembly.h:749
MatType type
void setCurrentSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:199
std::vector< std::pair< MooseVariable *, MooseVariable * > > _cm_nonlocal_entry
Definition: Assembly.h:674
const Elem *& elem()
Return the current element.
Definition: Assembly.h:189
MooseArray< std::vector< RealTensor > > _second_phi
Definition: Assembly.h:858
MooseArray< Real > _JxW
Cached JxW.
Definition: Assembly.h:880
MooseArray< std::vector< RealGradient > > VariablePhiGradient
Definition: Assembly.h:56
QBase *& qRuleNeighbor()
Returns the reference to the current quadrature being used on a current neighbor. ...
Definition: Assembly.h:257
const VariablePhiValue & phiFace()
Definition: Assembly.h:538
std::map< FEType, FEShapeData * > _shape_data
This is where the cached shape functions will be held.
Definition: Assembly.h:874
const MooseArray< Point > & physicalPoints()
The current points in physical space where we have reinited through reinitAtPhysical() ...
Definition: Assembly.h:141
const VariablePhiSecond & secondPhiFace()
Definition: Assembly.h:540
bool _should_use_fe_cache
Whether or not fe cache should be built at all.
Definition: Assembly.h:890
unsigned int & neighborSide()
Returns the current neighboring side.
Definition: Assembly.h:217
unsigned int _current_side
The current side of the selected element (valid only when working with sides)
Definition: Assembly.h:786
const Elem * _current_neighbor_side_elem
The current side element of the ncurrent neighbor element.
Definition: Assembly.h:798
std::map< unsigned int, QBase * > _holder_qrule_volume
Holds volume qrules for each dimension.
Definition: Assembly.h:724
std::vector< std::vector< DenseMatrix< Number > > > _sub_Keg
Definition: Assembly.h:824
bool _current_side_volume_computed
Boolean to indicate whether current element side volumes has been computed.
Definition: Assembly.h:810
const Real & sideElemVolume()
Returns the reference to the volume of current side element.
Definition: Assembly.h:229
std::vector< std::vector< unsigned char > > _jacobian_block_used
Flag that indicates if the jacobian block was used.
Definition: Assembly.h:676
std::map< unsigned int, std::map< FEType, FEBase * > > _fe_face_neighbor
Definition: Assembly.h:765
MooseArray< std::vector< RealTensor > > VariablePhiSecond
Definition: Assembly.h:57
const Moose::CoordinateSystemType & coordSystem()
Get the coordinate system type.
Definition: Assembly.h:159
std::map< FEType, FEShapeData * > _fe_shape_data
Definition: Assembly.h:896
MooseArray< Point > _current_q_points
The current list of quadrature points.
Definition: Assembly.h:716
Moose::CoordinateSystemType _coord_type
The coordinate system.
Definition: Assembly.h:720
VariablePhiSecond _second_phi
Definition: Assembly.h:839
const VariablePhiValue & phiNeighbor()
Definition: Assembly.h:542
std::map< unsigned int, const std::vector< Point > * > _holder_q_points
Holds pointers to the dimension&#39;s q_points.
Definition: Assembly.h:728
void setCurrentNeighborSubdomainID(SubdomainID i)
set the current subdomain ID
Definition: Assembly.h:245
const VariablePhiValue & phiFaceNeighbor()
Definition: Assembly.h:546
TensorValue< Real > RealTensorValue
Definition: Assembly.h:45
std::map< FEType, FEShapeData * > _fe_shape_data_face
Definition: Assembly.h:897
MooseArray< Point > _current_q_points_face
The current quadrature points on a face.
Definition: Assembly.h:745
std::vector< std::vector< unsigned char > > _jacobian_block_neighbor_used
Flag that indicates if the jacobian block for neighbor was used.
Definition: Assembly.h:679
const Elem * _current_side_elem
The current "element" making up the side we are currently on.
Definition: Assembly.h:788
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:509
FEBase * _current_fe_helper
The current helper object for transforming coordinates.
Definition: Assembly.h:708
VariablePhiGradient _grad_phi_face
Definition: Assembly.h:842
SubdomainID _current_subdomain_id
The current subdomain ID.
Definition: Assembly.h:782
const Elem *& sideElem()
Returns the side element.
Definition: Assembly.h:223
std::map< unsigned int, const std::vector< Point > * > _holder_normals
Holds pointers to the dimension&#39;s normal vectors.
Definition: Assembly.h:759
QBase * _current_qrule_face
quadrature rule used on faces
Definition: Assembly.h:741
unsigned int THREAD_ID
Definition: MooseTypes.h:79
const Node * _current_neighbor_node
The current neighboring node we are working with.
Definition: Assembly.h:806
const VariablePhiSecond & secondPhiFaceNeighbor()
Definition: Assembly.h:548
std::map< unsigned int, const std::vector< Real > * > _holder_JxW
Holds pointers to the dimension&#39;s transformed jacobian weights.
Definition: Assembly.h:730
void useFECache(bool fe_cache)
Whether or not this assembly should utilize FE shape function caching.
Definition: Assembly.h:387
std::vector< std::pair< MooseVariable *, MooseVariable * > > & couplingEntries()
Definition: Assembly.h:528