33 #ifndef __H__UG__LIB_DISC__REFERENCE_ELEMENT__REFERENCE_ELEMENT_MAPPING__
34 #define __H__UG__LIB_DISC__REFERENCE_ELEMENT__REFERENCE_ELEMENT_MAPPING__
63 template <
typename TRefElem,
int TWorldDim>
110 const size_t maxIter = 1000,
111 const number tol = 1e-10)
const;
116 const size_t maxIter = 1000,
117 const number tol = 1e-10)
const;
122 const size_t maxIter = 1000,
123 const number tol = 1e-10)
const;
168 std::vector<number>& vDet,
190 template <
int dim,
int worldDim,
bool isLinear,
typename TImpl>
201 for(
size_t ip = 0; ip < n; ++ip)
202 getImpl().local_to_global(vGlobPos[ip], vLocPos[ip]);
209 const size_t n = vLocPos.size();
219 const size_t maxIter = 1000,
235 for (
size_t i = 0; i < maxIter; ++i) {
238 getImpl().local_to_global(compGlobPos, locPos);
257 "reference_mapping.h: global_to_local() Newton iteration: Iter# "
258 << i <<
"; fabs(VecTwoNorm(dist)) = " << fabs(
VecTwoNorm(dist)) <<
259 "; dist = " << dist <<
"; locPos: " << locPos << std::endl);
287 UG_COND_THROW(!maxIter,
"Without a single iteration, local-to-global "
288 "mapping can never converge.");
291 for(
int i = 0; i < 3; ++i)
294 JInv(i, 0) <<
"; " << JInv(i, 1) <<
"; " << JInv(i, 2) << std::endl);
297 UG_THROW(
"ReferenceMapping::global_to_local: Newton method did not"
298 " reach a tolerance "<<tol<<
" after "<<maxIter<<
" steps. "
299 "Global Pos: "<<globPos<<
", dim: "<<
dim<<
", worldDim: "<<
300 worldDim<<
", last newton defect: "<<fabs(
VecTwoNorm(dist)));
306 const size_t maxIter = 1000,
307 const number tol = 1e-10)
const
317 getImpl().jacobian(J, vLocPos[0]);
322 for(
size_t ip = 0; ip < n; ++ip)
325 getImpl().local_to_global(compGlobPos, vLocPos[ip]);
332 for(
size_t ip = 0; ip < n; ++ip)
333 getImpl().global_to_local(vLocPos[ip], vGlobPos[ip], maxIter, tol);
340 const size_t maxIter = 1000,
341 const number tol = 1e-10)
const
343 const size_t n = vGlobPos.size();
353 getImpl().jacobian_transposed(JT, locPos);
363 getImpl().jacobian(vJ[0], vLocPos[0]);
364 for(
size_t ip = 1; ip < n; ++ip) vJ[ip] = vJ[0];
367 for(
size_t ip = 0; ip < n; ++ip)
368 getImpl().jacobian(vJ[ip], vLocPos[ip]);
376 const size_t n = vLocPos.size();
388 getImpl().jacobian_transposed(vJT[0], vLocPos[0]);
389 for(
size_t ip = 1; ip < n; ++ip) vJT[ip] = vJT[0];
392 for(
size_t ip = 0; ip < n; ++ip)
393 getImpl().jacobian_transposed(vJT[ip], vLocPos[ip]);
401 const size_t n = vLocPos.size();
411 getImpl().jacobian_transposed(JT, locPos);
412 UG_DLOG(
DID_REFERENCE_MAPPING, 2,
">>OCT_DISC_DEBUG:: " <<
"reference_mapping.h: " <<
"jacobian_transposed_inverse(): JT " << std::endl);
413 for(
int i = 0; i < 3; ++i)
425 vDet[0] =
getImpl().jacobian_transposed_inverse(vJTInv[0], vLocPos[0]);
426 for(
size_t ip = 1; ip < n; ++ip) vJTInv[ip] = vJTInv[0];
427 for(
size_t ip = 1; ip < n; ++ip) vDet[ip] = vDet[0];
430 for(
size_t ip = 0; ip < n; ++ip)
431 vDet[ip] =
getImpl().jacobian_transposed_inverse(vJTInv[ip], vLocPos[ip]);
441 getImpl().jacobian_transposed_inverse(vJTInv[0], vLocPos[0]);
442 for(
size_t ip = 1; ip < n; ++ip) vJTInv[ip] = vJTInv[0];
445 for(
size_t ip = 0; ip < n; ++ip)
446 getImpl().jacobian_transposed_inverse(vJTInv[ip], vLocPos[ip]);
452 std::vector<number>& vDet,
455 const size_t n = vLocPos.size();
456 vJTInv.resize(n); vDet.resize(n);
464 const size_t n = vLocPos.size();
473 getImpl().jacobian_transposed(JT, locPos);
483 for(
size_t ip = 1; ip < n; ++ip) vDet[ip] = vDet[0];
486 for(
size_t ip = 0; ip < n; ++ip)
495 const size_t n = vLocPos.size();
502 TImpl&
getImpl() {
return static_cast<TImpl&
>(*this);}
505 const TImpl&
getImpl()
const {
return static_cast<const TImpl&
>(*this);}
512 template <
int TWorldDim>
515 ReferenceMapping<ReferenceVertex, TWorldDim> >
530 using base_type::local_to_global;
531 using base_type::jacobian;
532 using base_type::jacobian_transposed;
533 using base_type::jacobian_transposed_inverse;
534 using base_type::sqrt_gram_det;
550 "ReferenceMapping: to few Corner Coordinates.");
557 co0 = vCornerCoord[0];
581 template <
int TWorldDim>
584 ReferenceMapping<ReferenceEdge, TWorldDim> >
599 using base_type::local_to_global;
600 using base_type::jacobian;
601 using base_type::jacobian_transposed;
602 using base_type::jacobian_transposed_inverse;
603 using base_type::sqrt_gram_det;
619 "ReferenceMapping: to few Corner Coordinates.");
626 co0 = vCornerCoord[0];
627 VecSubtract(a10, vCornerCoord[1], vCornerCoord[0]);
641 for(
int i = 0; i <
worldDim; ++i) JT(0,i) = a10[i];
651 template <
int TWorldDim>
654 ReferenceMapping<ReferenceTriangle, TWorldDim> >
669 using base_type::local_to_global;
670 using base_type::jacobian;
671 using base_type::jacobian_transposed;
672 using base_type::jacobian_transposed_inverse;
673 using base_type::sqrt_gram_det;
689 "ReferenceMapping: to few Corner Coordinates.");
696 co0 = vCornerCoord[0];
697 VecSubtract(a10, vCornerCoord[1], vCornerCoord[0]);
698 VecSubtract(a20, vCornerCoord[2], vCornerCoord[0]);
705 VecScaleAdd(globPos, 1.0, co0, locPos[0], a10, locPos[1], a20);
712 for(
int i = 0; i <
worldDim; ++i) JT(0, i) = a10[i];
713 for(
int i = 0; i <
worldDim; ++i) JT(1, i) = a20[i];
725 template <
int TWorldDim>
728 ReferenceMapping<ReferenceQuadrilateral, TWorldDim> >
743 using base_type::local_to_global;
744 using base_type::jacobian;
745 using base_type::jacobian_transposed;
746 using base_type::jacobian_transposed_inverse;
747 using base_type::sqrt_gram_det;
763 "ReferenceMapping: to few Corner Coordinates.");
771 x[co] = vCornerCoord[co];
778 const number a = (1.-locPos[0]);
779 const number b = (1.-locPos[1]);
783 locPos[0]*locPos[1], x[2],
791 const number a = 1. - locPos[1];
792 const number b = 1. - locPos[0];
795 JT(0, i) = a*(x[1][i] - x[0][i]) + locPos[1]*(x[2][i] - x[3][i]);
798 JT(1, i) = b*(x[3][i] - x[0][i]) + locPos[0]*(x[2][i] - x[1][i]);
809 template <
int TWorldDim>
812 ReferenceMapping<ReferenceTetrahedron, TWorldDim> >
827 using base_type::local_to_global;
828 using base_type::jacobian;
829 using base_type::jacobian_transposed;
830 using base_type::jacobian_transposed_inverse;
831 using base_type::sqrt_gram_det;
847 "ReferenceMapping: to few Corner Coordinates.");
854 co0 = vCornerCoord[0];
855 VecSubtract(a10, vCornerCoord[1], vCornerCoord[0]);
856 VecSubtract(a20, vCornerCoord[2], vCornerCoord[0]);
857 VecSubtract(a30, vCornerCoord[3], vCornerCoord[0]);
874 for(
int i = 0; i <
worldDim; ++i) JT[0][i] = a10[i];
875 for(
int i = 0; i <
worldDim; ++i) JT[1][i] = a20[i];
876 for(
int i = 0; i <
worldDim; ++i) JT[2][i] = a30[i];
886 template <
int TWorldDim>
889 ReferenceMapping<ReferencePyramid, TWorldDim> >
904 using base_type::local_to_global;
905 using base_type::jacobian;
906 using base_type::jacobian_transposed;
907 using base_type::jacobian_transposed_inverse;
908 using base_type::sqrt_gram_det;
924 "ReferenceMapping: to few Corner Coordinates.");
932 x[co] = vCornerCoord[co];
939 const number a = 1.0 - locPos[0];
940 const number b = 1.0 - locPos[1];
941 if (locPos[0] > locPos[1])
943 const number a0 = a * b - locPos[2] * b;
944 const number a1 = locPos[0] * b - locPos[2]*locPos[1];
945 const number a2 = locPos[0] * locPos[1] + locPos[2]*locPos[1];
946 const number a3 = a * locPos[1] - locPos[2] * locPos[1];
949 globPos[d] = a0*x[0][d]+a1*x[1][d]+a2*x[2][d]
950 +a3*x[3][d]+locPos[2]*x[4][d];
954 const number a0 = a * b - locPos[2] * a;
955 const number a1 = locPos[0] * b - locPos[2]*locPos[0];
956 const number a2 = locPos[0] * locPos[1] + locPos[2]*locPos[0];
957 const number a3 = a * locPos[1] - locPos[2] * locPos[0];
959 globPos[d] = a0*x[0][d]+a1*x[1][d]+
960 a2*x[2][d]+a3*x[3][d]+locPos[2]*x[4][d];
970 a[d] = x[0][d]-x[1][d]+x[2][d]-x[3][d];
972 if (locPos[0] > locPos[1])
975 JT(0,d) = x[1][d]-x[0][d]+locPos[1]*a[d];
978 JT(1,d) = x[3][d]-x[0][d]+(locPos[0]+locPos[2])*a[d];
981 JT(2,d) = x[4][d]-x[0][d]+locPos[1]*a[d];
986 JT(0,d) = x[1][d]-x[0][d]+(locPos[1]+locPos[2])*a[d];
989 JT(1,d) = x[3][d]-x[0][d]+locPos[0]*a[d];
992 JT(2,d) = x[4][d]-x[0][d]+locPos[0]*a[d];
1004 template <
int TWorldDim>
1007 ReferenceMapping<ReferencePrism, TWorldDim> >
1022 using base_type::local_to_global;
1023 using base_type::jacobian;
1024 using base_type::jacobian_transposed;
1025 using base_type::jacobian_transposed_inverse;
1026 using base_type::sqrt_gram_det;
1042 "ReferenceMapping: to few Corner Coordinates.");
1043 update(&vCornerCoord[0]);
1050 x[co] = vCornerCoord[co];
1057 const number a = 1.0 - locPos[0] - locPos[1];
1058 const number b = 1.0 - locPos[2];
1060 const number a1 = locPos[0] * b;
1061 const number a2 = locPos[1] * b;
1062 const number a3 = a * locPos[2];
1063 const number a4 = locPos[0] * locPos[2];
1064 const number a5 = locPos[1] * locPos[2];
1067 globPos[d] = a0*x[0][d]+a1*x[1][d]+a2*x[2][d]+a3*x[3][d]+a4*x[4][d]+a5*x[5][d];
1077 a[d] = x[0][d]-x[1][d]-x[3][d]+x[4][d];
1080 b[d] = x[0][d]-x[2][d]-x[3][d]+x[5][d];
1083 JT(0,d) = x[1][d]-x[0][d]+locPos[2]*a[d];
1084 JT(1,d) = x[2][d]-x[0][d]+locPos[2]*b[d];
1085 JT(2,d) = x[3][d]-x[0][d]+locPos[0]*a[d]+locPos[1]*b[d];
1098 template <
int TWorldDim>
1101 ReferenceMapping<ReferenceHexahedron, TWorldDim> >
1116 using base_type::local_to_global;
1117 using base_type::jacobian;
1118 using base_type::jacobian_transposed;
1119 using base_type::jacobian_transposed_inverse;
1120 using base_type::sqrt_gram_det;
1136 "ReferenceMapping: to few Corner Coordinates.");
1137 update(&vCornerCoord[0]);
1144 x[co] = vCornerCoord[co];
1151 const number a = 1.0 - locPos[0];
1152 const number b = 1.0 - locPos[1];
1153 const number c = 1.0 - locPos[2];
1154 const number a0 = a * b * c;
1155 const number a1 = locPos[0] * b * c;
1156 const number a2 = locPos[0] * locPos[1] * c;
1157 const number a3 = a * locPos[1] * c;
1158 const number a4 = a * b * locPos[2];
1159 const number a5 = locPos[0] * b * locPos[2];
1160 const number a6 = locPos[0] * locPos[1] * locPos[2];
1161 const number a7 = a * locPos[1] * locPos[2];
1164 globPos[d] = a0*x[0][d]+a1*x[1][d]+a2*x[2][d]+a3*x[3][d]+
1165 a4*x[4][d]+a5*x[5][d]+a6*x[6][d]+a7*x[7][d];
1174 const number a = 1.0 - locPos[0];
1175 const number b = 1.0 - locPos[1];
1176 const number c = 1.0 - locPos[2];
1179 a2 = locPos[1] * locPos[2];
1182 JT(0,d) = a0*(x[1][d]-x[0][d])+a1*(x[2][d]-x[3][d])
1183 + a2*(x[6][d]-x[7][d])+a3*(x[5][d]-x[4][d]);
1187 a2 = locPos[0] * locPos[2];
1190 JT(1,d) = a0*(x[3][d]-x[0][d])+a1*(x[2][d]-x[1][d])
1191 + a2*(x[6][d]-x[5][d])+a3*(x[7][d]-x[4][d]);
1195 a2 = locPos[0] * locPos[1];
1198 JT(2,d) = a0*(x[4][d]-x[0][d])+a1*(x[5][d]-x[1][d])
1199 + a2*(x[6][d]-x[2][d])+a3*(x[7][d]-x[3][d]);
1209 template <
int TWorldDim>
1212 ReferenceMapping<ReferenceOctahedron, TWorldDim> >
1227 using base_type::local_to_global;
1228 using base_type::jacobian;
1229 using base_type::jacobian_transposed;
1230 using base_type::jacobian_transposed_inverse;
1231 using base_type::sqrt_gram_det;
1247 "ReferenceMapping: to few Corner Coordinates.");
1248 update(&vCornerCoord[0]);
1255 x[co] = vCornerCoord[co];
1265 const number z_sgn = (locPos[2] < 0) ? -1.0 : 1.0;
1266 const number a0 = (locPos[2] < 0) ? -locPos[2] : 0.0;
1267 const number a5 = (locPos[2] < 0) ? 0.0 : locPos[2];
1269 if (locPos[0] > locPos[1])
1271 const number a1 = 1.0 - locPos[0] - z_sgn * locPos[2];
1272 const number a2 = locPos[0] - locPos[1];
1273 const number a3 = locPos[1];
1277 globPos[d] = a0*x[0][d]+a1*x[1][d]+a2*x[2][d]
1278 +a3*x[3][d]+a4*x[4][d]+a5*x[5][d];
1282 const number a1 = 1.0 - locPos[1] - z_sgn * locPos[2];
1284 const number a3 = locPos[0];
1285 const number a4 = locPos[1] - locPos[0];
1288 globPos[d] = a0*x[0][d]+a1*x[1][d]+a2*x[2][d]
1289 +a3*x[3][d]+a4*x[4][d]+a5*x[5][d];
1300 const number z_sgn = (locPos[2] < 0) ? -1.0 : 1.0;
1301 const number Da0 = (locPos[2] < 0) ? -1.0 : 0.0;
1302 const number Da5 = (locPos[2] < 0) ? 0.0 : 1.0;
1304 if (locPos[0] > locPos[1])
1307 JT(0,d) = -x[1][d]+x[2][d];
1310 JT(1,d) = -x[2][d]+x[3][d];
1313 JT(2,d) = Da0*x[0][d] - z_sgn*x[1][d] + Da5*x[5][d];
1318 JT(0,d) = x[3][d]-x[4][d];
1321 JT(1,d) = -x[1][d]+x[4][d];
1324 JT(2,d) = Da0*x[0][d] - z_sgn*x[1][d] + Da5*x[5][d];
Base class for Reference mappings helping to implement interface.
Definition: reference_mapping.h:192
void local_to_global(std::vector< MathVector< worldDim > > &vGlobPos, const std::vector< MathVector< dim > > &vLocPos) const
map local coordinate to global coordinate for a vector of local positions
Definition: reference_mapping.h:206
void global_to_local(std::vector< MathVector< dim > > &vLocPos, const std::vector< MathVector< worldDim > > &vGlobPos, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate for a vector of local positions
Definition: reference_mapping.h:338
TImpl & getImpl()
access to implementation
Definition: reference_mapping.h:502
void jacobian_transposed_inverse(std::vector< MathMatrix< worldDim, dim > > &vJTInv, std::vector< number > &vDet, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of the inverse of the jacobian for a vector of positions
Definition: reference_mapping.h:451
number sqrt_gram_det(const MathVector< dim > &locPos) const
returns the determinate of the jacobian
Definition: reference_mapping.h:470
void global_to_local(MathVector< dim > &locPos, const MathVector< worldDim > &globPos, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate
Definition: reference_mapping.h:216
void sqrt_gram_det(std::vector< number > &vDet, const std::vector< MathVector< dim > > &vLocPos) const
returns the determinate of the jacobian for a vector of local positions
Definition: reference_mapping.h:492
void jacobian_transposed(MathMatrix< dim, worldDim > *vJT, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of jacobian for n local positions
Definition: reference_mapping.h:383
void jacobian(MathMatrix< worldDim, dim > *vJ, const MathVector< dim > *vLocPos, size_t n) const
returns jacobian for n local positions
Definition: reference_mapping.h:358
void jacobian_transposed_inverse(MathMatrix< worldDim, dim > *vJTInv, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of the inverse of the jacobian for n local positions
Definition: reference_mapping.h:436
void local_to_global(MathVector< worldDim > *vGlobPos, const MathVector< dim > *vLocPos, size_t n) const
map local coordinate to global coordinate for n local positions
Definition: reference_mapping.h:198
void jacobian(MathMatrix< worldDim, dim > &J, const MathVector< dim > &locPos) const
returns jacobian
Definition: reference_mapping.h:349
void jacobian_transposed_inverse(std::vector< MathMatrix< worldDim, dim > > &vJTInv, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of the inverse of the jacobian for a vector of positions
Definition: reference_mapping.h:461
bool is_linear() const
returns if mapping is affine
Definition: reference_mapping.h:195
void jacobian_transposed(std::vector< MathMatrix< dim, worldDim > > &vJT, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of jacobian for a vector of positions
Definition: reference_mapping.h:398
void jacobian_transposed_inverse(MathMatrix< worldDim, dim > *vJTInv, number *vDet, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of the inverse of the jacobian for n local positions
Definition: reference_mapping.h:419
const TImpl & getImpl() const
const access to implementation
Definition: reference_mapping.h:505
number jacobian_transposed_inverse(MathMatrix< worldDim, dim > &JTInv, const MathVector< dim > &locPos) const
returns transposed of the inverse of the jacobian
Definition: reference_mapping.h:407
void global_to_local(MathVector< dim > *vLocPos, const MathVector< worldDim > *vGlobPos, size_t n, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate for n local positions
Definition: reference_mapping.h:304
void jacobian(std::vector< MathMatrix< worldDim, dim > > &vJ, const std::vector< MathVector< dim > > &vLocPos) const
returns jacobian for a vector of local positions
Definition: reference_mapping.h:373
void sqrt_gram_det(number *vDet, const MathVector< dim > *vLocPos, size_t n) const
returns the determinate of the jacobian for n local positions
Definition: reference_mapping.h:478
Definition: reference_element.h:345
static const int dim
dimension of reference element
Definition: reference_element.h:351
static const int numCorners
number of corners
Definition: reference_element.h:354
reference element for a hexahedron
Definition: reference_element.h:648
static const int dim
dimension of reference element
Definition: reference_element.h:654
static const int numCorners
number of corners
Definition: reference_element.h:657
Definition: reference_mapping.h:585
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:638
BaseReferenceMapping< ReferenceEdge::dim, TWorldDim, true, ReferenceMapping< ReferenceEdge, TWorldDim > > base_type
Definition: reference_mapping.h:598
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:624
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:607
MathVector< worldDim > a10
Definition: reference_mapping.h:645
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:616
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:612
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:631
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:611
Definition: reference_mapping.h:1102
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:1141
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:1128
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:1124
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:1170
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:1148
BaseReferenceMapping< ReferenceHexahedron::dim, TWorldDim, false, ReferenceMapping< ReferenceHexahedron, TWorldDim > > base_type
Definition: reference_mapping.h:1115
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:1133
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:1129
Definition: reference_mapping.h:1213
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:1244
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:1239
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:1252
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:1294
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:1259
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:1240
BaseReferenceMapping< ReferenceOctahedron::dim, TWorldDim, false, ReferenceMapping< ReferenceOctahedron, TWorldDim > > base_type
Definition: reference_mapping.h:1226
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:1235
Definition: reference_mapping.h:1008
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:1035
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:1034
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:1071
BaseReferenceMapping< ReferencePrism::dim, TWorldDim, false, ReferenceMapping< ReferencePrism, TWorldDim > > base_type
Definition: reference_mapping.h:1021
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:1054
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:1047
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:1030
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:1039
Definition: reference_mapping.h:890
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:917
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:921
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:929
BaseReferenceMapping< ReferencePyramid::dim, TWorldDim, false, ReferenceMapping< ReferencePyramid, TWorldDim > > base_type
Definition: reference_mapping.h:903
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:936
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:912
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:916
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:965
Definition: reference_mapping.h:729
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:760
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:768
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:775
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:751
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:788
BaseReferenceMapping< ReferenceQuadrilateral::dim, TWorldDim, false, ReferenceMapping< ReferenceQuadrilateral, TWorldDim > > base_type
Definition: reference_mapping.h:742
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:756
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:755
Definition: reference_mapping.h:813
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:852
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:835
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:871
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:840
BaseReferenceMapping< ReferenceTetrahedron::dim, TWorldDim, true, ReferenceMapping< ReferenceTetrahedron, TWorldDim > > base_type
Definition: reference_mapping.h:826
MathVector< worldDim > a10
Definition: reference_mapping.h:880
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:861
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:839
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:844
Definition: reference_mapping.h:655
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:694
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:677
MathVector< worldDim > a10
Definition: reference_mapping.h:717
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:702
BaseReferenceMapping< ReferenceTriangle::dim, TWorldDim, true, ReferenceMapping< ReferenceTriangle, TWorldDim > > base_type
Definition: reference_mapping.h:668
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:686
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:681
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:709
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:682
Definition: reference_mapping.h:516
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
Definition: reference_mapping.h:561
MathVector< worldDim > co0
Definition: reference_mapping.h:575
BaseReferenceMapping< ReferenceVertex::dim, TWorldDim, true, ReferenceMapping< ReferenceVertex, TWorldDim > > base_type
Definition: reference_mapping.h:529
void update(const MathVector< worldDim > *vCornerCoord)
update the mapping for a new set of corners
Definition: reference_mapping.h:555
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_mapping.h:568
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Definition: reference_mapping.h:542
ReferenceMapping()
Default Constructor.
Definition: reference_mapping.h:538
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
Definition: reference_mapping.h:547
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Definition: reference_mapping.h:543
Definition: reference_mapping.h:65
number sqrt_gram_det(const MathVector< dim > &locPos) const
returns the determinate of the jacobian
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
void global_to_local(MathVector< dim > &locPos, const MathVector< worldDim > &globPos, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate
ReferenceMapping(const MathVector< worldDim > *vCornerCoord)
Constructor setting the corners of the element.
static const int worldDim
world dimension (range space dimension)
Definition: reference_mapping.h:68
void jacobian_transposed(MathMatrix< dim, worldDim > *vJT, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of jacobian for n local positions
bool is_linear() const
returns if mapping is affine
void local_to_global(std::vector< MathVector< worldDim > > &vGlobPos, const std::vector< MathVector< dim > > &vLocPos) const
map local coordinate to global coordinate for a vector of local positions
void jacobian(std::vector< MathMatrix< worldDim, dim > > &J, const std::vector< MathVector< dim > > &vLocPos) const
returns jacobian for a vector of local positions
void jacobian_transposed_inverse(MathMatrix< worldDim, dim > *vJTInv, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of the inverse of the jacobian for n local positions
void jacobian_transposed(std::vector< MathMatrix< dim, worldDim > > &vJT, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of jacobian for a vector of positions
void local_to_global(MathVector< worldDim > *vGlobPos, const MathVector< dim > *vLocPos, size_t n) const
map local coordinate to global coordinate for n local positions
void sqrt_gram_det(std::vector< number > vDet, const std::vector< MathVector< dim > > &vLocPos) const
returns the determinate of the jacobian for a vector of local positions
void jacobian_transposed_inverse(std::vector< MathMatrix< worldDim, dim > > &vJTInv, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of the inverse of the jacobian for a vector of positions
void jacobian(MathMatrix< worldDim, dim > *vJ, const MathVector< dim > *vLocPos, size_t n) const
returns jacobian for n local positions
ReferenceMapping()
Default Constructor.
void jacobian(MathMatrix< worldDim, dim > &J, const MathVector< dim > &locPos) const
returns jacobian
number jacobian_transposed_inverse(MathMatrix< worldDim, dim > &JTInv, const MathVector< dim > &locPos) const
returns transposed of the inverse of the jacobian and sqrt of gram determinante
void update(const MathVector< worldDim > *vCornerCoord)
refresh mapping for new set of corners
void global_to_local(MathVector< dim > *vLocPos, const MathVector< worldDim > *vGlobPos, size_t n, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate for n local positions
void jacobian_transposed_inverse(std::vector< MathMatrix< worldDim, dim > > &vJTInv, std::vector< number > &vDet, const std::vector< MathVector< dim > > &vLocPos) const
returns transposed of the inverse of the jacobian for a vector of positions
void sqrt_gram_det(number *vDet, const MathVector< dim > *vLocPos, size_t n) const
returns the determinate of the jacobian for n local positions
static const int dim
reference dimension (domain space dimension)
Definition: reference_mapping.h:71
void update(const std::vector< MathVector< worldDim > > &vCornerCoord)
refresh mapping for new set of corners
ReferenceMapping(const std::vector< MathVector< worldDim > > &vCornerCoord)
Constructor setting the corners of the element.
void global_to_local(std::vector< MathVector< dim > > &vLocPos, const std::vector< MathVector< worldDim > > &vGlobPos, const size_t maxIter=1000, const number tol=1e-10) const
map global coordinate to local coordinate for a vector of local positions
static const bool isLinear
flag if mapping is linear (i.e. Jacobian does not depend on x)
Definition: reference_mapping.h:74
void jacobian_transposed_inverse(MathMatrix< worldDim, dim > *vJTInv, number *vDet, const MathVector< dim > *vLocPos, size_t n) const
returns transposed of the inverse of the jacobian for n local positions
void jacobian_transposed(MathMatrix< dim, worldDim > &JT, const MathVector< dim > &locPos) const
returns transposed of jacobian
Definition: reference_element.h:699
static const int dim
dimension of reference element
Definition: reference_element.h:705
static const int numCorners
number of corners
Definition: reference_element.h:708
Definition: reference_element.h:595
static const int numCorners
number of corners
Definition: reference_element.h:604
static const int dim
dimension of reference element
Definition: reference_element.h:601
Definition: reference_element.h:546
static const int dim
dimension of reference element
Definition: reference_element.h:552
static const int numCorners
number of corners
Definition: reference_element.h:555
Definition: reference_element.h:445
static const int numCorners
number of corners
Definition: reference_element.h:454
static const int dim
dimension of reference element
Definition: reference_element.h:451
Definition: reference_element.h:494
static const int dim
dimension of reference element
Definition: reference_element.h:500
static const int numCorners
number of corners
Definition: reference_element.h:503
Definition: reference_element.h:394
static const int dim
dimension of reference element
Definition: reference_element.h:400
static const int numCorners
number of corners
Definition: reference_element.h:403
Definition: reference_element.h:303
static const int dim
dimension of reference element
Definition: reference_element.h:309
MathMatrix< N, M, T >::value_type SqrtGramDeterminant(const MathMatrix< N, M, T > &m)
Square root of Gram Determinant of a matrix.
Definition: math_matrix_functions_common_impl.hpp:508
MathMatrix< N, M, T >::value_type RightInverse(MathMatrix< N, M, T > &mOut, const MathMatrix< M, N, T > &m)
Right-Inverse of a Matrix.
Definition: math_matrix_functions_common_impl.hpp:680
void Transpose(MathMatrix< N, M, T > &mOut, const MathMatrix< M, N, T > &m)
transpose a matrix
Definition: math_matrix_functions_common_impl.hpp:345
MathMatrix< N, M, T >::value_type LeftInverse(MathMatrix< N, M, T > &mOut, const MathMatrix< M, N, T > &m)
Left-Inverse of a Matrix.
Definition: math_matrix_functions_common_impl.hpp:716
#define UG_ASSERT(expr, msg)
Definition: assert.h:70
#define UG_THROW(msg)
Definition: error.h:57
#define UG_DLOG(__debugID__, level, msg)
Definition: log.h:298
#define UG_COND_THROW(cond, msg)
UG_COND_THROW(cond, msg) : performs a UG_THROW(msg) if cond == true.
Definition: error.h:61
double number
Definition: types.h:124
void MatVecScaleMultAppend(vector_t_out &vOut, typename vector_t_out::value_type s, const matrix_t &m, const vector_t_in &v)
Matrix - Vector Multiplication added scaled to a second vector.
Definition: math_matrix_vector_functions_common_impl.hpp:90
void MatVecMult(vector_t_out &vOut, const matrix_t &m, const vector_t_in &v)
Matrix - Vector Multiplication.
Definition: math_matrix_vector_functions_common_impl.hpp:49
void VecScaleAdd(vector_t &vOut, typename vector_t::value_type s1, const vector_t &v1, typename vector_t::value_type s2, const vector_t &v2)
Scales two Vectors, adds them and returns the sum in a third vector.
Definition: math_vector_functions_common_impl.hpp:265
void VecSubtract(vector_t &vOut, const vector_t &v1, const vector_t &v2)
subtracts v2 from v1 and stores the result in a vOut
Definition: math_vector_functions_common_impl.hpp:226
void VecScale(vector_t &vOut, const vector_t &v, typename vector_t::value_type s)
scales a MathVector<N>
Definition: math_vector_functions_common_impl.hpp:252
vector_t::value_type VecTwoNorm(const vector_t &v)
Definition: math_vector_functions_common_impl.hpp:577
DebugID DID_REFERENCE_MAPPING
DebugID DID_REFERENCE_MAPPING_GLOB_TO_LOC
bool VecAbsIsLess(const vector_t &v1, const vector_t &v2)
component-wise comparison of two vectors (in the absolute values)
Definition: math_vector_functions_common_impl.hpp:677