00001
00002
00003
00004
00005
00006
00007
00008 #ifndef APF_H
00009 #define APF_H
00010
00011 #include "apfMatrix.h"
00012 #include "apfNew.h"
00013 #include "apfDynamicArray.h"
00014
00015 #include <vector>
00029 namespace apf {
00030
00031 class Field;
00032 class Element;
00033 class Mesh;
00034 class MeshEntity;
00035 class VectorElement;
00037 typedef VectorElement MeshElement;
00038 class FieldShape;
00039 struct Sharing;
00040
00047 void destroyMesh(Mesh* m);
00048
00056 MeshElement* createMeshElement(Mesh* m, MeshEntity* e);
00057
00067 MeshElement* createMeshElement(Field* c, MeshEntity* e);
00068
00071 MeshEntity * getMeshEntity(MeshElement * me);
00072
00078 void destroyMeshElement(MeshElement* e);
00079
00084 enum ValueType {
00086 SCALAR,
00088 VECTOR,
00090 MATRIX,
00092 PACKED,
00094 VALUE_TYPES
00095 };
00096
00104 Field* createLagrangeField(Mesh* m, const char* name, int valueType, int order);
00105
00116 Field* createStepField(Mesh* m, const char* name, int valueType);
00117
00125 Field* createIPField(Mesh* m, const char* name, int valueType, int order);
00126
00129 Field* createField(Mesh* m, const char* name, int valueType, FieldShape* shape);
00130
00132 Field* createFieldOn(Mesh* m, const char* name, int valueType);
00133
00140 Field* createPackedField(Mesh* m, const char* name, int components,
00141 apf::FieldShape* shape = 0);
00142
00144 Field* createGeneralField(
00145 Mesh* m,
00146 const char* name,
00147 int valueType,
00148 int components,
00149 FieldShape* shape);
00150
00154 Field* cloneField(Field* f, Mesh* onto);
00155
00158 Mesh* getMesh(Field* f);
00159
00162 bool hasEntity(Field* f, MeshEntity* e);
00163
00169 const char* getName(Field* f);
00170
00173 int getValueType(Field* f);
00174
00180 void destroyField(Field* f);
00181
00188 void setScalar(Field* f, MeshEntity* e, int node, double value);
00189
00194 double getScalar(Field* f, MeshEntity* e, int node);
00195
00200 void setVector(Field* f, MeshEntity* e, int node, Vector3 const& value);
00201
00204 void getVector(Field* f, MeshEntity* e, int node, Vector3& value);
00205
00210 void setMatrix(Field* f, MeshEntity* e, int node, Matrix3x3 const& value);
00211
00214 void getMatrix(Field* f, MeshEntity* e, int node, Matrix3x3& value);
00215
00217 void setComponents(Field* f, MeshEntity* e, int node, double const* components);
00218
00223 void getComponents(Field* f, MeshEntity* e, int node, double* components);
00224
00235 Element* createElement(Field* f, MeshElement* e);
00236
00241 Element* createElement(Field* f, MeshEntity* e);
00242
00245 void destroyElement(Element* e);
00246
00254 MeshElement* getMeshElement(Element* e);
00255
00257 MeshEntity* getMeshEntity(Element* e);
00258
00264 double getScalar(Element* e, Vector3 const& param);
00265
00271 void getGrad(Element* e, Vector3 const& param, Vector3& grad);
00272
00278 void getVector(Element* e, Vector3 const& param, Vector3& value);
00279
00285 double getDiv(Element* e, Vector3 const& param);
00286
00292 void getCurl(Element* e, Vector3 const& param, Vector3& curl);
00293
00303 void getVectorGrad(Element* e, Vector3 const& param, Matrix3x3& deriv);
00304
00310 void getMatrix(Element* e, Vector3 const& param, Matrix3x3& value);
00311
00313 void getComponents(Element* e, Vector3 const& param, double* components);
00314
00320 int countIntPoints(MeshElement* e, int order);
00321
00328 void getIntPoint(MeshElement* e, int order, int point, Vector3& param);
00329
00336 double getIntWeight(MeshElement* e, int order, int point);
00337
00340 void mapLocalToGlobal(MeshElement* e, Vector3 const& local, Vector3& global);
00341
00351 double getDV(MeshElement* e, Vector3 const& param);
00352
00365 class Integrator
00366 {
00367 public:
00369 Integrator(int o);
00370 virtual ~Integrator();
00372 void process(Mesh* m);
00374 void process(MeshElement* e);
00381 virtual void inElement(MeshElement*);
00388 virtual void outElement();
00399 virtual void atPoint(Vector3 const& p, double w, double dV) = 0;
00406 virtual void parallelReduce();
00407 protected:
00408 int order;
00409 int ipnode;
00410 };
00411
00420 double measure(MeshElement* e);
00421
00426 double measure(Mesh* m, MeshEntity* e);
00427
00430 int getOrder(MeshElement* e);
00431
00434 void getJacobian(MeshElement* e, Vector3 const& local, Matrix3x3& j);
00435
00443 void getJacobianInv(MeshElement* e, Vector3 const& local, Matrix3x3& jinv);
00444
00460 double computeCosAngle(Mesh* m, MeshEntity* pe, MeshEntity* e1, MeshEntity* e2,
00461 const Matrix3x3& Q);
00462
00463 double computeCosAngleInTri(Mesh* m, MeshEntity* tri,
00464 MeshEntity* e1, MeshEntity* e2, const Matrix3x3& Q);
00465
00466 double computeCosAngleInTet(Mesh* m, MeshEntity* tet,
00467 MeshEntity* e1, MeshEntity* e2, const Matrix3x3& Q);
00468
00469 Vector3 computeFaceNormalAtVertex(Mesh* m, MeshEntity* face,
00470 MeshEntity* vert, const Matrix3x3& Q);
00471
00477 int countNodes(Element* e);
00478
00481 void getScalarNodes(Element* e, NewArray<double>& values);
00482
00485 void getVectorNodes(Element* e, NewArray<Vector3>& values);
00486
00489 void getMatrixNodes(Element* e, NewArray<Matrix3x3>& values);
00490
00493 void getShapeValues(Element* e, Vector3 const& local,
00494 NewArray<double>& values);
00495
00500 void getShapeGrads(Element* e, Vector3 const& local,
00501 NewArray<Vector3>& grads);
00502
00503
00506 FieldShape* getShape(Field* f);
00507
00510 int countComponents(Field* f);
00511
00520 #define APF_ITERATE(t,w,i) \
00521 for (t::iterator (i) = (w).begin(); \
00522 (i) != (w).end(); ++(i))
00523
00525 #define APF_CONST_ITERATE(t,w,i) \
00526 for (t::const_iterator (i) = (w).begin(); \
00527 (i) != (w).end(); ++(i))
00528
00534 void writeVtkFiles(const char* prefix, Mesh* m, int cellDim = -1);
00535
00542 void writeVtkFiles(const char* prefix, Mesh* m,
00543 std::vector<std::string> writeFields, int cellDim = -1);
00544
00548 void writeOneVtkFile(const char* prefix, Mesh* m);
00549
00555 void writeASCIIVtkFiles(const char* prefix, Mesh* m);
00556
00563 void writeASCIIVtkFiles(const char* prefix, Mesh* m,
00564 std::vector<std::string> writeFields);
00565
00572 void getGaussPoint(int type, int order, int point, Vector3& param);
00573
00575 int countGaussPoints(int type, int order);
00576
00584 double getJacobianDeterminant(Matrix3x3 const& J, int dimension);
00585
00587 int getDimension(MeshElement* me);
00588
00594 void synchronize(Field* f, Sharing* shr = 0);
00595
00602 void accumulate(Field* f, Sharing* shr = 0);
00603
00610 void fail(const char* why) __attribute__((noreturn));
00611
00613 void freeze(Field* f);
00614
00616 void unfreeze(Field* f);
00617
00619 bool isFrozen(Field* f);
00620
00626 double* getArrayData(Field* f);
00627
00629 void zeroField(Field* f);
00630
00632 struct Function
00633 {
00635 virtual ~Function();
00642 virtual void eval(MeshEntity* e, double* result) = 0;
00643 };
00644
00650 Field* createUserField(Mesh* m, const char* name, int valueType, FieldShape* s,
00651 Function* f);
00652
00657 Field* recoverGradientByVolume(Field* f);
00658
00659 void copyData(Field* to, Field* from);
00660
00662 void projectField(Field* to, Field* from);
00663
00664 void axpy(double a, Field* x, Field* y);
00665
00666 void renameField(Field* f, const char* name);
00667
00675 void getBF(FieldShape* s, MeshElement* e, Vector3 const& p,
00676 NewArray<double>& BF);
00677
00682 void getGradBF(FieldShape* s, MeshElement* e, Vector3 const& p,
00683 NewArray<Vector3>& gradBF);
00684
00685 }
00686
00687 #endif