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
00060 MeshEntity * getMeshEntity(MeshElement * me);
00061
00067 void destroyMeshElement(MeshElement* e);
00068
00073 enum ValueType {
00075 SCALAR,
00077 VECTOR,
00079 MATRIX,
00081 PACKED,
00083 VALUE_TYPES
00084 };
00085
00093 Field* createLagrangeField(Mesh* m, const char* name, int valueType, int order);
00094
00105 Field* createStepField(Mesh* m, const char* name, int valueType);
00106
00114 Field* createIPField(Mesh* m, const char* name, int valueType, int order);
00115
00118 Field* createField(Mesh* m, const char* name, int valueType, FieldShape* shape);
00119
00121 Field* createFieldOn(Mesh* m, const char* name, int valueType);
00122
00129 Field* createPackedField(Mesh* m, const char* name, int components,
00130 apf::FieldShape* shape = 0);
00131
00133 Field* createGeneralField(
00134 Mesh* m,
00135 const char* name,
00136 int valueType,
00137 int components,
00138 FieldShape* shape);
00139
00143 Field* cloneField(Field* f, Mesh* onto);
00144
00147 Mesh* getMesh(Field* f);
00148
00151 bool hasEntity(Field* f, MeshEntity* e);
00152
00158 const char* getName(Field* f);
00159
00162 int getValueType(Field* f);
00163
00169 void destroyField(Field* f);
00170
00177 void setScalar(Field* f, MeshEntity* e, int node, double value);
00178
00183 double getScalar(Field* f, MeshEntity* e, int node);
00184
00189 void setVector(Field* f, MeshEntity* e, int node, Vector3 const& value);
00190
00193 void getVector(Field* f, MeshEntity* e, int node, Vector3& value);
00194
00199 void setMatrix(Field* f, MeshEntity* e, int node, Matrix3x3 const& value);
00200
00203 void getMatrix(Field* f, MeshEntity* e, int node, Matrix3x3& value);
00204
00206 void setComponents(Field* f, MeshEntity* e, int node, double const* components);
00207
00212 void getComponents(Field* f, MeshEntity* e, int node, double* components);
00213
00224 Element* createElement(Field* f, MeshElement* e);
00225
00230 Element* createElement(Field* f, MeshEntity* e);
00231
00234 void destroyElement(Element* e);
00235
00243 MeshElement* getMeshElement(Element* e);
00244
00246 MeshEntity* getMeshEntity(Element* e);
00247
00253 double getScalar(Element* e, Vector3 const& param);
00254
00260 void getGrad(Element* e, Vector3 const& param, Vector3& grad);
00261
00267 void getVector(Element* e, Vector3 const& param, Vector3& value);
00268
00274 double getDiv(Element* e, Vector3 const& param);
00275
00281 void getCurl(Element* e, Vector3 const& param, Vector3& curl);
00282
00288 void getVectorGrad(Element* e, Vector3 const& param, Matrix3x3& deriv);
00289
00295 void getMatrix(Element* e, Vector3 const& param, Matrix3x3& value);
00296
00298 void getComponents(Element* e, Vector3 const& param, double* components);
00299
00305 int countIntPoints(MeshElement* e, int order);
00306
00313 void getIntPoint(MeshElement* e, int order, int point, Vector3& param);
00314
00321 double getIntWeight(MeshElement* e, int order, int point);
00322
00325 void mapLocalToGlobal(MeshElement* e, Vector3 const& local, Vector3& global);
00326
00336 double getDV(MeshElement* e, Vector3 const& param);
00337
00350 class Integrator
00351 {
00352 public:
00354 Integrator(int o);
00355 virtual ~Integrator();
00357 void process(Mesh* m);
00359 void process(MeshElement* e);
00366 virtual void inElement(MeshElement*);
00373 virtual void outElement();
00384 virtual void atPoint(Vector3 const& p, double w, double dV) = 0;
00391 virtual void parallelReduce();
00392 protected:
00393 int order;
00394 int ipnode;
00395 };
00396
00405 double measure(MeshElement* e);
00406
00411 double measure(Mesh* m, MeshEntity* e);
00412
00415 int getOrder(MeshElement* e);
00416
00419 void getJacobian(MeshElement* e, Vector3 const& local, Matrix3x3& j);
00420
00428 void getJacobianInv(MeshElement* e, Vector3 const& local, Matrix3x3& jinv);
00429
00435 int countNodes(Element* e);
00436
00439 void getScalarNodes(Element* e, NewArray<double>& values);
00440
00443 void getVectorNodes(Element* e, NewArray<Vector3>& values);
00444
00447 void getMatrixNodes(Element* e, NewArray<Matrix3x3>& values);
00448
00451 void getShapeValues(Element* e, Vector3 const& local,
00452 NewArray<double>& values);
00453
00458 void getShapeGrads(Element* e, Vector3 const& local,
00459 NewArray<Vector3>& grads);
00460
00461
00464 FieldShape* getShape(Field* f);
00465
00468 int countComponents(Field* f);
00469
00478 #define APF_ITERATE(t,w,i) \
00479 for (t::iterator (i) = (w).begin(); \
00480 (i) != (w).end(); ++(i))
00481
00483 #define APF_CONST_ITERATE(t,w,i) \
00484 for (t::const_iterator (i) = (w).begin(); \
00485 (i) != (w).end(); ++(i))
00486
00492 void writeVtkFiles(const char* prefix, Mesh* m);
00493
00500 void writeVtkFiles(const char* prefix, Mesh* m,
00501 std::vector<std::string> writeFields);
00502
00506 void writeOneVtkFile(const char* prefix, Mesh* m);
00507
00513 void writeASCIIVtkFiles(const char* prefix, Mesh* m);
00514
00521 void writeASCIIVtkFiles(const char* prefix, Mesh* m,
00522 std::vector<std::string> writeFields);
00523
00530 void getGaussPoint(int type, int order, int point, Vector3& param);
00531
00533 int countGaussPoints(int type, int order);
00534
00542 double getJacobianDeterminant(Matrix3x3 const& J, int dimension);
00543
00545 int getDimension(MeshElement* me);
00546
00552 void synchronize(Field* f, Sharing* shr = 0);
00553
00560 void accumulate(Field* f, Sharing* shr = 0);
00561
00568 void fail(const char* why) __attribute__((noreturn));
00569
00571 void freeze(Field* f);
00572
00574 void unfreeze(Field* f);
00575
00577 bool isFrozen(Field* f);
00578
00584 double* getArrayData(Field* f);
00585
00587 void zeroField(Field* f);
00588
00590 struct Function
00591 {
00593 virtual ~Function();
00600 virtual void eval(MeshEntity* e, double* result) = 0;
00601 };
00602
00608 Field* createUserField(Mesh* m, const char* name, int valueType, FieldShape* s,
00609 Function* f);
00610
00615 Field* recoverGradientByVolume(Field* f);
00616
00617 void copyData(Field* to, Field* from);
00618
00620 void projectField(Field* to, Field* from);
00621
00622 void axpy(double a, Field* x, Field* y);
00623
00624 void renameField(Field* f, const char* name);
00625
00633 void getBF(FieldShape* s, MeshElement* e, Vector3 const& p,
00634 NewArray<double>& BF);
00635
00640 void getGradBF(FieldShape* s, MeshElement* e, Vector3 const& p,
00641 NewArray<Vector3>& gradBF);
00642
00643 }
00644
00645 #endif