00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #ifndef SUNDANCE_ELEMENTINTEGRAL_H
00032 #define SUNDANCE_ELEMENTINTEGRAL_H
00033
00034 #include "SundanceDefs.hpp"
00035 #include "SundanceCellJacobianBatch.hpp"
00036 #include "SundanceQuadratureFamily.hpp"
00037 #include "SundanceBasisFamily.hpp"
00038 #include "SundanceParametrizedCurve.hpp"
00039 #include "SundanceMesh.hpp"
00040 #include "Teuchos_Array.hpp"
00041
00042
00043 namespace Sundance
00044 {
00045 using namespace Sundance;
00046 using namespace Teuchos;
00047
00048
00049
00050
00051
00052 class ElementIntegral
00053 {
00054 public:
00055
00056 ElementIntegral(int spatialDim,
00057 const CellType& maxCellType,
00058 int dim,
00059 const CellType& cellType,
00060 bool isInternalBdry,
00061 const ParametrizedCurve& globalCurve,
00062 const Mesh& mesh,
00063 int verb);
00064
00065
00066 ElementIntegral(int spatialDim,
00067 const CellType& maxCellType,
00068 int dim,
00069 const CellType& cellType,
00070 const BasisFamily& testBasis,
00071 int alpha,
00072 int testDerivOrder,
00073 bool isInternalBdry,
00074 const ParametrizedCurve& globalCurve,
00075 const Mesh& mesh,
00076 int verb);
00077
00078
00079 ElementIntegral(int spatialDim,
00080 const CellType& maxCellType,
00081 int dim,
00082 const CellType& cellType,
00083 const BasisFamily& testBasis,
00084 int alpha,
00085 int testDerivOrder,
00086 const BasisFamily& unkBasis,
00087 int beta,
00088 int unkDerivOrder,
00089 bool isInternalBdry,
00090 const ParametrizedCurve& globalCurve,
00091 const Mesh& mesh,
00092 int verb);
00093
00094
00095 virtual ~ElementIntegral(){;}
00096
00097
00098
00099 int order() const {return order_;}
00100
00101
00102 int nNodesTest() const {return nNodesTest_;}
00103
00104
00105 int nNodesUnk() const {return nNodesUnk_;}
00106
00107
00108
00109 int nNodes() const {return nNodes_;}
00110
00111
00112
00113
00114 int nFacetCases() const {return nFacetCases_;}
00115
00116
00117 bool isInternalBdry() const {return isInternalBdry_;}
00118
00119
00120 static bool& alwaysUseCofacets();
00121
00122
00123 void setVerb(
00124 int integrationVerb,
00125 int transformVerb);
00126
00127
00128 int setupVerb() const {return setupVerb_;}
00129
00130
00131 int integrationVerb() const {return integrationVerb_;}
00132
00133
00134 int transformVerb() const {return transformVerb_;}
00135
00136
00137 void describe(std::ostream& os) const ;
00138
00139
00140 static int& transformationMatrixIsValid(int alpha);
00141
00142
00143
00144 static int& transformationMatrixIsValid(int alpha, int beta);
00145
00146
00147 static void invalidateTransformationMatrices();
00148
00149
00150 static double& totalFlops() {static double rtn = 0; return rtn;}
00151
00152 BasisFamily &getTestBasis() { return testBasis_; }
00153 BasisFamily &getUnknownBasis() { return unkBasis_; }
00154
00155 protected:
00156
00157
00158 void assertBilinearForm() const ;
00159
00160
00161 void assertLinearForm() const ;
00162
00163
00164 static void addFlops(const double& flops) {totalFlops() += flops;}
00165
00166
00167 int dim() const {return dim_;}
00168
00169
00170 int spatialDim() const {return spatialDim_;}
00171
00172
00173
00174
00175 int nRefDerivTest() const {return nRefDerivTest_;}
00176
00177
00178
00179
00180 int nRefDerivUnk() const {return nRefDerivUnk_;}
00181
00182
00183
00184 int testDerivOrder() const {return testDerivOrder_;}
00185
00186
00187
00188 int unkDerivOrder() const {return unkDerivOrder_;}
00189
00190
00191 int alpha() const {return alpha_;}
00192
00193
00194 int beta() const {return beta_;}
00195
00196
00197 const CellType& cellType() const {return cellType_;}
00198
00199
00200 const CellType& maxCellType() const {return maxCellType_;}
00201
00202
00203 const CellType& evalCellType() const {return evalCellType_;}
00204
00205
00206 const BasisFamily& testBasis() const {return testBasis_;}
00207
00208
00209 const BasisFamily& unkBasis() const {return unkBasis_;}
00210
00211
00212 static Array<double>& G(int gamma) ;
00213
00214
00215 static Array<double>& G(int gamma, int delta) ;
00216
00217
00218 static int ipow(int base, int power);
00219
00220
00221 static double chopVal() {static double rtn=1.0e-14; return rtn;}
00222
00223
00224
00225 static double chop(const double& x)
00226 {
00227 if (::fabs(x) > chopVal()) return x;
00228 else return 0.0;
00229 }
00230
00231
00232 void getQuad(const QuadratureFamily& quad, int evalCase,
00233 Array<Point>& quadPts, Array<double>& quadWeights) const ;
00234
00235
00236 void createTwoFormTransformationMatrix(const CellJacobianBatch& JTrans,
00237 const CellJacobianBatch& JVol) const;
00238
00239 void createOneFormTransformationMatrix(const CellJacobianBatch& JTrans,
00240 const CellJacobianBatch& JVol) const;
00241
00242
00243 const Mesh& mesh() const {return mesh_;}
00244
00245
00246 const ParametrizedCurve& globalCurve() const {return globalCurve_;}
00247
00248 private:
00249
00250 int setupVerb_;
00251 int integrationVerb_;
00252 int transformVerb_;
00253
00254 int spatialDim_;
00255
00256 int dim_;
00257
00258 bool isInternalBdry_;
00259
00260 int nFacetCases_;
00261
00262 int testDerivOrder_;
00263
00264 int nRefDerivTest_;
00265
00266 int nNodesTest_;
00267
00268 int unkDerivOrder_;
00269
00270 int nRefDerivUnk_;
00271
00272 int nNodesUnk_;
00273
00274 int nNodes_;
00275
00276 int order_;
00277
00278 int alpha_;
00279
00280 int beta_;
00281
00282 CellType cellType_;
00283
00284 CellType maxCellType_;
00285
00286 CellType evalCellType_;
00287
00288 BasisFamily testBasis_;
00289
00290 BasisFamily unkBasis_;
00291
00292
00293 const ParametrizedCurve globalCurve_;
00294
00295
00296 const Mesh mesh_;
00297
00298 };
00299 }
00300
00301
00302 #endif