00001
00002
00003
00004
00005
00006
00007
00008 #ifndef CRVBEZIER_H
00009 #define CRVBEZIER_H
00010
00011 #include "crv.h"
00012 #include "mth.h"
00013
00017 namespace crv {
00018
00020 double Bij(const int i, const int j,const double u, const double v);
00021 double Bijk(const int i, const int j, const int k, const double u,
00022 const double v, const double w);
00023 double Bijkl(const int i, const int j, const int k, const int l,
00024 const double u, const double v, const double w, const double t);
00025
00027 double Bij(const int ij[], const double xi[]);
00028 double Bijk(const int ijk[], const double xi[]);
00029 double Bijkl(const int ijkl[], const double xi[]);
00030
00032 int computeTriNodeIndex(int P, int i, int j);
00034 int computeTetNodeIndex(int P, int i, int j, int k);
00035
00040 int getNumControlPoints(int type, int order);
00045 int getNumInternalControlPoints(int type, int order);
00046
00049 void getTriNodesFromTetNodes(int f, int P,
00050 apf::NewArray<apf::Vector3>& tetNodes,
00051 apf::NewArray<apf::Vector3>& triNodes);
00053 void getTriDetJacNodesFromTetDetJacNodes(int f, int P,
00054 apf::NewArray<double>& tetNodes,
00055 apf::NewArray<double>& triNodes);
00059 void getFullRepFromBlended(int type,
00060 apf::NewArray<double>& transformCoefficients,
00061 apf::NewArray<apf::Vector3>& elemNodes);
00066 double computeTriJacobianDetFromBezierFormulation(apf::Mesh* m,
00067 apf::MeshEntity* e, apf::Vector3& xi);
00072 double computeTetJacobianDetFromBezierFormulation(apf::Mesh* m,
00073 apf::MeshEntity* e, apf::Vector3& xi);
00074
00076 void getBezierNodeXi(int type, int P, int node, apf::Vector3& xi);
00078 void collectNodeXi(int parentType, int childType, int P,
00079 const apf::Vector3* range, apf::NewArray<apf::Vector3>& xi);
00081 void elevateBezierEdge(int P, int r, apf::NewArray<apf::Vector3>& nodes,
00082 apf::NewArray<apf::Vector3>& elevatedNodes);
00083 void elevateBezierTriangle(int P, int r, apf::NewArray<apf::Vector3>& nodes,
00084 apf::NewArray<apf::Vector3>& elevatedNodes);
00085 void elevateBezierTet(int P, int r, apf::NewArray<apf::Vector3>& nodes,
00086 apf::NewArray<apf::Vector3>& elevatedNodes);
00087 void elevateBezier(int type, int P, int r,
00088 apf::NewArray<apf::Vector3>& nodes,
00089 apf::NewArray<apf::Vector3>& elevatedNodes);
00090
00095 void elevateBezierCurve(apf::Mesh2* m, apf::MeshEntity* edge, int n, int r);
00096
00098 void subdivideBezierEdge(int P, double t, apf::NewArray<apf::Vector3>& nodes,
00099 apf::NewArray<apf::Vector3> (&subNodes)[2]);
00100 void subdivideBezierTriangle(int P, apf::Vector3& p,
00101 apf::NewArray<apf::Vector3>& nodes,
00102 apf::NewArray<apf::Vector3> (&subNodes)[3]);
00103 void subdivideBezierTriangle(int P, apf::NewArray<apf::Vector3>& nodes,
00104 apf::NewArray<apf::Vector3> (&subNodes)[4]);
00105 void subdivideBezierTet(int P, apf::Vector3& p,
00106 apf::NewArray<apf::Vector3>& nodes,
00107 apf::NewArray<apf::Vector3> (&subNodes)[4]);
00108
00114 void convertInterpolationPoints(int n, int ne,
00115 apf::NewArray<apf::Vector3>& nodes,
00116 apf::NewArray<double>& c,
00117 apf::NewArray<apf::Vector3>& newNodes);
00118
00119 void convertInterpolationPoints(apf::Mesh2* m, apf::MeshEntity* e,
00120 int n, int ne, apf::NewArray<double>& c);
00121
00125 void getBezierTransformationCoefficients(int P, int type,
00126 apf::NewArray<double>& c);
00127 void getInternalBezierTransformationCoefficients(apf::Mesh* m, int P, int blend,
00128 int type, apf::NewArray<double>& c);
00129 void getGregoryTransformationCoefficients(int type, apf::NewArray<double>& c);
00130 void getGregoryBlendedTransformationCoefficients(int blend, int type,
00131 apf::NewArray<double>& c);
00132
00134 void snapToInterpolate(apf::Mesh2* m, apf::MeshEntity* e);
00135
00143 void getTransformationMatrix(apf::Mesh* m, apf::MeshEntity* e,
00144 mth::Matrix<double>& A, const apf::Vector3 *range);
00145
00146 }
00147
00148 #endif