escript Revision_
ReferenceElementSets.h
Go to the documentation of this file.
1
2/*****************************************************************************
3*
4* Copyright (c) 2003-2020 by The University of Queensland
5* http://www.uq.edu.au
6*
7* Primary Business: Queensland, Australia
8* Licensed under the Apache License, version 2.0
9* http://www.apache.org/licenses/LICENSE-2.0
10*
11* Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12* Development 2012-2013 by School of Earth Sciences
13* Development from 2014-2017 by Centre for Geoscience Computing (GeoComp)
14* Development from 2019 by School of Earth and Environmental Sciences
15**
16*****************************************************************************/
17
18#ifndef __FINLEY_REFERENCEELEMENTSETS_H__
19#define __FINLEY_REFERENCEELEMENTSETS_H__
20
21#include "ReferenceElements.h"
22
23namespace finley {
24
28 ReferenceElementSet(ElementTypeId id, int order, int reduced_order)
29 {
32 id_info->BasisFunctions);
33 if (order<0)
34 order=std::max(2*bf_info->numOrder, 0);
35
36 referenceElement.reset(new ReferenceElement(id, order));
37 if (reduced_order<0)
38 reduced_order=std::max(2*(bf_info->numOrder-1), 0);
40 reduced_order));
41
42 if (referenceElement->getNumNodes() != referenceElementReducedQuadrature->getNumNodes()) {
43 throw escript::ValueError("ReferenceElementSet: numNodes in referenceElement and referenceElementReducedQuadrature don't match.");
44 }
45 }
46
48 bool reducedIntegrationOrder) const
49 {
50 if (reducedShapefunction) {
51 return (reducedIntegrationOrder ?
52 referenceElementReducedQuadrature->LinearBasisFunctions :
53 referenceElement->LinearBasisFunctions);
54 }
55 return (reducedIntegrationOrder ?
56 referenceElementReducedQuadrature->BasisFunctions :
57 referenceElement->BasisFunctions);
58 }
59
60 const_ShapeFunction_ptr borrowParametrization(bool reducedIntegrationOrder) const
61 {
62 return (reducedIntegrationOrder ?
63 referenceElementReducedQuadrature->Parametrization :
64 referenceElement->Parametrization);
65 }
66
68 {
69 return (reducedIntOrder ? referenceElementReducedQuadrature :
71 }
72
73 inline int getNumNodes() const { return referenceElement->getNumNodes(); }
74
77};
78
79
80typedef boost::shared_ptr<const ReferenceElementSet> const_ReferenceElementSet_ptr;
81
82
83} // namespace finley
84
85#endif // __FINLEY_REFERENCEELEMENTSETS_H__
86
An exception class that signals an invalid argument value.
Definition: EsysException.h:100
A suite of factory methods for creating various finley domains.
Definition: finley/src/Assemble.h:32
boost::shared_ptr< const ReferenceElement > const_ReferenceElement_ptr
Definition: ReferenceElements.h:216
boost::shared_ptr< const ShapeFunction > const_ShapeFunction_ptr
Definition: ShapeFunctions.h:102
boost::shared_ptr< ReferenceElement > ReferenceElement_ptr
Definition: ReferenceElements.h:215
ElementTypeId
Definition: ReferenceElements.h:41
boost::shared_ptr< const ReferenceElementSet > const_ReferenceElementSet_ptr
Definition: ReferenceElementSets.h:80
this struct holds the definition of the reference element
Definition: ReferenceElements.h:123
ShapeFunctionTypeId BasisFunctions
shape function for the basis functions
Definition: ReferenceElements.h:148
Definition: ReferenceElementSets.h:27
const_ShapeFunction_ptr borrowParametrization(bool reducedIntegrationOrder) const
Definition: ReferenceElementSets.h:60
int getNumNodes() const
Definition: ReferenceElementSets.h:73
ReferenceElement_ptr referenceElementReducedQuadrature
Definition: ReferenceElementSets.h:75
ReferenceElement_ptr referenceElement
Definition: ReferenceElementSets.h:76
const_ShapeFunction_ptr borrowBasisFunctions(bool reducedShapefunction, bool reducedIntegrationOrder) const
Definition: ReferenceElementSets.h:47
ReferenceElementSet(ElementTypeId id, int order, int reduced_order)
Definition: ReferenceElementSets.h:28
const_ReferenceElement_ptr borrowReferenceElement(bool reducedIntOrder) const
Definition: ReferenceElementSets.h:67
this struct holds the realization of a reference element
Definition: ReferenceElements.h:179
static const ReferenceElementInfo * getInfo(ElementTypeId id)
returns the element information structure for the given type id
Definition: ReferenceElements.cpp:664
this struct holds the definition of the shape functions on an element
Definition: ShapeFunctions.h:60
int numOrder
order of the shape functions
Definition: ShapeFunctions.h:70
static const ShapeFunctionInfo * getInfo(ShapeFunctionTypeId id)
Definition: ShapeFunctions.cpp:105