diff --git a/AMDiS/src/ElInfo2d.cc b/AMDiS/src/ElInfo2d.cc index 01360e1105af9e7ad1c60cdc273d7b72d93d507c..33993e4c8326387b83fae88124621998bb1f0213 100644 --- a/AMDiS/src/ElInfo2d.cc +++ b/AMDiS/src/ElInfo2d.cc @@ -32,6 +32,33 @@ namespace AMDiS { mtl::dense2D<double> ElInfo2d::mat_d1_right(mat_d1_right_val); + + + double ElInfo2d::mat_d2_val[6][6] = {{1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, + {0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, + {0.0, 0.0, 1.0, 0.0, 0.0, 0.0}, + {0.0, 0.0, 0.0, 1.0, 0.0, 0.0}, + {0.0, 0.0, 0.0, 0.0, 1.0, 0.0}, + {0.0, 0.0, 0.0, 0.0, 0.0, 1.0}}; + mtl::dense2D<double> ElInfo2d::mat_d2(mat_d2_val); + + double ElInfo2d::mat_d2_left_val[6][6] = {{0.0, 1.0, 0.0, 0.375, -0.125, 0.0}, + {0.0, 0.0, 0.0, -0.125, -0.125, 0.0}, + {1.0, 0.0, 0.0, 0.0, 0.0, 0.0}, + {0.0, 0.0, 0.0, 0.0, 0.5, 0.0}, + {0.0, 0.0, 0.0, 0.0, 0.5, 1.0}, + {0.0, 0.0, 1.0, 0.75, 0.25, 0.0}}; + mtl::dense2D<double> ElInfo2d::mat_d2_left(mat_d2_left_val); + + double ElInfo2d::mat_d2_right_val[6][6] = {{0.0, 0.0, 0.0, -0.125, -0.125, 0.0}, + {1.0, 0.0, 0.0, -0.125, 0.375, 0.0}, + {0.0, 1.0, 0.0, 0.0, 0.0, 0.0}, + {0.0, 0.0, 0.0, 0.5, 0.0, 1.0}, + {0.0, 0.0, 0.0, 0.5, 0.0, 0.0}, + {0.0, 0.0, 1.0, 0.25, 0.75, 0.0}}; + mtl::dense2D<double> ElInfo2d::mat_d2_right(mat_d2_right_val); + + ElInfo2d::ElInfo2d(Mesh *aMesh) : ElInfo(aMesh) { @@ -714,6 +741,24 @@ namespace AMDiS { subElemMatrices[1][refinementPath] = mat; } break; + case 2: + { + dense2D<double> mat(mat_d2); + dense2D<double> tmpMat(num_rows(mat), num_rows(mat)); + + for (int i = 0; i < refinementPathLength; i++) { + if (refinementPath & (1 << i)) { + tmpMat = mat_d2_right * mat; + mat = tmpMat; + } else { + tmpMat = mat_d2_left * mat; + mat = tmpMat; + } + } + + subElemMatrices[2][refinementPath] = mat; + } + break; default: ERROR_EXIT("Not supported for basis function degree: %d\n", degree); } diff --git a/AMDiS/src/ElInfo2d.h b/AMDiS/src/ElInfo2d.h index 07b39770dee65d98239686ad9d5c6e10e2ef7d75..54f506510fab398d872c414d24423dee9ec7c6ae 100644 --- a/AMDiS/src/ElInfo2d.h +++ b/AMDiS/src/ElInfo2d.h @@ -67,16 +67,22 @@ namespace AMDiS { WorldVector<double> *e1, *e2, *normal; static double mat_d1_val[3][3]; - static mtl::dense2D<double> mat_d1; static double mat_d1_left_val[3][3]; - static mtl::dense2D<double> mat_d1_left; static double mat_d1_right_val[3][3]; - static mtl::dense2D<double> mat_d1_right; + + static double mat_d2_val[6][6]; + static mtl::dense2D<double> mat_d2; + + static double mat_d2_left_val[6][6]; + static mtl::dense2D<double> mat_d2_left; + + static double mat_d2_right_val[6][6]; + static mtl::dense2D<double> mat_d2_right; }; } diff --git a/AMDiS/src/ElInfo3d.cc b/AMDiS/src/ElInfo3d.cc index e16c9f49e802aa6755a967ade7ec92c9080492d4..07c2fbbcae7f4fcbdad42d83624cc667035de7c0 100644 --- a/AMDiS/src/ElInfo3d.cc +++ b/AMDiS/src/ElInfo3d.cc @@ -14,6 +14,25 @@ namespace AMDiS { + double ElInfo3d::mat_d1_left_val[4][4] = {{1.0, 0.0, 0.0, 0.5}, + {0.0, 0.0, 0.0, 0.5}, + {0.0, 1.0, 0.0, 0.0}, + {0.0, 0.0, 1.0, 0.0}}; + mtl::dense2D<double> ElInfo3d::mat_d1_left(mat_d1_left_val); + + + double ElInfo3d::mat_d1_l0_right_val[4][4] = {{0.0, 0.0, 0.0, 0.5}, + {1.0, 0.0, 0.0, 0.5}, + {0.0, 0.0, 1.0, 0.0}, + {0.0, 1.0, 0.0, 0.0}}; + double ElInfo3d::mat_d1_l12_right_val[4][4] = {{0.0, 0.0, 0.0, 0.5}, + {1.0, 0.0, 0.0, 0.5}, + {0.0, 1.0, 0.0, 0.0}, + {0.0, 0.0, 1.0, 0.0}}; + mtl::dense2D<double> ElInfo3d::mat_d1_l0_right(mat_d1_l0_right_val); + mtl::dense2D<double> ElInfo3d::mat_d1_l12_right(mat_d1_l12_right_val); + + void ElInfo3d::fillMacroInfo(const MacroElement * mel) { FUNCNAME("ElInfo3d::fillMacroInfo()"); @@ -104,6 +123,7 @@ namespace AMDiS { } } + double ElInfo3d::calcGrdLambda(DimVec<WorldVector<double> >& grd_lam) { FUNCNAME("ElInfo3d::calcGrdLambda()"); @@ -169,6 +189,7 @@ namespace AMDiS { return adet; } + const int ElInfo3d::worldToCoord(const WorldVector<double>& xy, DimVec<double>* lambda) const { @@ -343,8 +364,6 @@ namespace AMDiS { } - - void ElInfo3d::fillElInfo(int ichild, const ElInfo *elInfoOld) { FUNCNAME("ElInfo3d::fillElInfo()"); @@ -595,7 +614,35 @@ namespace AMDiS { { FUNCNAME("ElInfo3d::getSubElemCoordsMat()"); - ERROR_EXIT("Not yet implemented!\n"); + using namespace mtl; + + if (subElemMatrices[degree].count(refinementPath) == 0) { + switch (degree) { + case 1: + { + dense2D<double> mat(4, 4), tmpMat(4, 4); + mat = 1; + + for (int i = 0; i < refinementPathLength; i++) { + if (refinementPath & (1 << i)) { + if ((level + i) % 3 == 0) + tmpMat = mat_d1_l0_right * mat; + else + tmpMat = mat_d1_l12_right * mat; + + mat = tmpMat; + } else { + tmpMat = mat_d1_left * mat; + mat = tmpMat; + } + } + + subElemMatrices[1][refinementPath] = mat; + } + default: + ERROR_EXIT("Not supported for basis function degree: %d\n", degree); + } + } return subElemMatrices[degree][refinementPath]; } diff --git a/AMDiS/src/ElInfo3d.h b/AMDiS/src/ElInfo3d.h index ec2102fc35da91e32741e7a99c4a95f0650dedf3..e3ac0a2e8ccc3a648b380ab2ea51d6944838d69f 100644 --- a/AMDiS/src/ElInfo3d.h +++ b/AMDiS/src/ElInfo3d.h @@ -93,6 +93,14 @@ namespace AMDiS { /// Tmp vectors used for calculations in calcGrdLambda and getNormal(). std::vector< WorldVector<double> > tmpWorldVecs; + + static double mat_d1_left_val[4][4]; + static mtl::dense2D<double> mat_d1_left; + + static double mat_d1_l0_right_val[4][4]; + static double mat_d1_l12_right_val[4][4]; + static mtl::dense2D<double> mat_d1_l0_right; + static mtl::dense2D<double> mat_d1_l12_right; }; }