From 9c12ccd2ccdeaa8674d4a5a103659c6e473ca443 Mon Sep 17 00:00:00 2001 From: Jonathan Lambrechts <jonathan.lambrechts@uclouvain.be> Date: Thu, 17 Mar 2011 14:35:13 +0000 Subject: [PATCH] spetialize xyz2uvw for p1Triangle --- Geo/MTriangle.cpp | 55 ++++++++++++++++------------------------------- Geo/MTriangle.h | 2 +- 2 files changed, 19 insertions(+), 38 deletions(-) diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp index 71c4e44fdc..7555c2b259 100644 --- a/Geo/MTriangle.cpp +++ b/Geo/MTriangle.cpp @@ -91,43 +91,24 @@ double MTriangle::gammaShapeMeasure() #endif } -// void MTriangle::xyz2uvw(double xyz[3], double uvw[3]){ - -// double X[3] = {_v[0]->x(), _v[1]->x(), _v[2]->x()}; -// double Y[3] = {_v[0]->y(), _v[1]->y(), _v[2]->y()}; -// double Z[3] = {_v[0]->z(), _v[1]->z(), _v[2]->z()}; - - // uvw[0] = uvw[1] = uvw[2] = 0.; - - // double jac[3][3] = {{X[1]-X[0], X[2]-X[Ø], 0},{}, {}}; - - - // if(!getJacobian(uvw[0], uvw[1], uvw[2], jac)) break; - - // double xn = 0., yn = 0., zn = 0.; - // double sf[1256]; - // getShapeFunctions(uvw[0], uvw[1], uvw[2], sf); - // for (int i = 0; i < getNumShapeFunctions(); i++) { - // MVertex *v = getShapeFunctionNode(i); - // xn += v->x() * sf[i]; - // yn += v->y() * sf[i]; - // zn += v->z() * sf[i]; - // } - // double inv[3][3]; - // inv3x3(jac, inv); - - // double un = uvw[0] + inv[0][0] * (xyz[0] - xn) + - // inv[1][0] * (xyz[1] - yn) + inv[2][0] * (xyz[2] - zn); - // double vn = uvw[1] + inv[0][1] * (xyz[0] - xn) + - // inv[1][1] * (xyz[1] - yn) + inv[2][1] * (xyz[2] - zn); - // double wn = uvw[2] + inv[0][2] * (xyz[0] - xn) + - // inv[1][2] * (xyz[1] - yn) + inv[2][2] * (xyz[2] - zn); - - // uvw[0] = un; - // uvw[1] = vn; - // uvw[2] = wn; - - //} + void MTriangle::xyz2uvw(double xyz[3], double uvw[3]) + { + const double O[3] = {_v[0]->x(), _v[0]->y(), _v[0]->z()}; + const double d[3] = {xyz[0] - O[0], xyz[1] - O[1], xyz[2] - O[2]}; + const double d1[3] = {_v[1]->x() - O[0], _v[1]->y() - O[1], _v[1]->z() - O[2]}; + const double d2[3] = {_v[2]->x() - O[0], _v[2]->y() - O[1], _v[2]->z() - O[2]}; + const double Jxy = d1[0] * d2[1] - d1[1] * d2[0]; + const double Jxz = d1[0] * d2[2] - d1[2] * d2[0]; + if (fabs(Jxy) > fabs(Jxz)) { + uvw[0] = (d[0] * d2[1] - d[1] * d2[0]) / Jxy; + uvw[1] = (d[1] * d1[0] - d[0] * d1[1]) / Jxy; + } + else { + uvw[0] = (d[0] * d2[2] - d[2] * d2[0]) / Jxz; + uvw[1] = (d[2] * d1[0] - d[0] * d1[2]) / Jxz; + } + uvw[2] = 0.; +} const polynomialBasis* MTriangle::getFunctionSpace(int o) const { diff --git a/Geo/MTriangle.h b/Geo/MTriangle.h index c2f5b9088d..f749d2933e 100644 --- a/Geo/MTriangle.h +++ b/Geo/MTriangle.h @@ -63,7 +63,7 @@ class MTriangle : public MElement { static const int map[3] = {0, 2, 1}; return getVertex(map[num]); } - //virtual void xyz2uvw(double xyz[3], double uvw[3]); + virtual void xyz2uvw(double xyz[3], double uvw[3]); virtual MVertex *getOtherVertex(MVertex *v1, MVertex *v2) { if(_v[0] != v1 && _v[0] != v2) return _v[0]; -- GitLab