const CellModel& cellModel(CellModel::GetCellModel(_my_geometry));
switch( _my_geometry )
{
+ case NORM_POINT1:
+ _my_local_ref_dim = 0;
+ _my_local_nb_ref = 1;
+ point1Init();
+ aSatify = isSatisfy();
+ CHECK_MACRO;
+ break;
+
case NORM_SEG2:
_my_local_ref_dim = 1;
_my_local_nb_ref = 2;
}
break;
+ case NORM_TRI7:
+ _my_local_ref_dim = 2;
+ _my_local_nb_ref = 7;
+ tria7aInit();
+ aSatify = isSatisfy();
+ CHECK_MACRO;
+ break;
+
case NORM_QUAD4:
{
_my_local_ref_dim = 2;
break;
case NORM_PENTA15:
- _my_local_ref_dim = 3;
- _my_local_nb_ref = 15;
- penta15aInit();
- aSatify = isSatisfy();
-
- if(!aSatify)
- {
- penta15bInit();
- aSatify = isSatisfy();
- CHECK_MACRO;
- }
- break;
+ {
+ _my_local_ref_dim = 3;
+ _my_local_nb_ref = 15;
+ MapToShapeFunction PENTA15PTR[]={Penta15aInit,Penta15bInit};
+ std::size_t NB_OF_PENTA15PTR(sizeof(PENTA15PTR)/sizeof(MapToShapeFunction));
+ for(std::size_t i=0;i<NB_OF_PENTA15PTR && !aSatify;i++)
+ {
+ (PENTA15PTR[i])(*this);
+ aSatify = isSatisfy();
+ }
+ CHECK_MACRO;
+ break;
+ }
case NORM_HEXA8:
{
return &_my_function_value[ _my_nb_ref*theGaussId ];
}
+void GaussInfo::point1Init()
+{
+ double *funValue(&_my_function_value[0]);
+ funValue[0] = 1. ;
+}
+
/*!
* Init Segment 2 Reference coordinates ans Shape function.
*/
SHAPE_FUN_MACRO_END;
}
+void GaussInfo::tria7aInit()
+{
+ LOCAL_COORD_MACRO_BEGIN;
+ case 0:
+ coords[0] = 0.0;
+ coords[1] = 0.0;
+ break;
+ case 1:
+ coords[0] = 1.0;
+ coords[1] = 0.0;
+ break;
+ case 2:
+ coords[0] = 0.0;
+ coords[1] = 1.0;
+ break;
+ case 3:
+ coords[0] = 0.5;
+ coords[1] = 0.0;
+ break;
+ case 4:
+ coords[0] = 0.5;
+ coords[1] = 0.5;
+ break;
+ case 5:
+ coords[0] = 0.0;
+ coords[1] = 0.5;
+ break;
+ case 6:
+ coords[0] = 0.3333333333333333;
+ coords[1] = 0.3333333333333333;
+ break;
+
+ LOCAL_COORD_MACRO_END;
+
+ SHAPE_FUN_MACRO_BEGIN;
+ funValue[0]=1-3*(gc[0]+gc[1])+2*(gc[0]*gc[0]+gc[1]*gc[1])+7*gc[0]*gc[1]-3*gc[0]*gc[1]*(gc[0]+gc[1]);
+ funValue[1]=gc[0]*(-1+2*gc[0]+3*gc[1]-3*gc[1]*(gc[0]+gc[1]));
+ funValue[2]=gc[1]*(-1.+3.*gc[0]+2.*gc[1]-3.*gc[0]*(gc[0]+gc[1]));
+ funValue[3]=4*gc[0]*(1-gc[0]-4*gc[1]+3*gc[1]*(gc[0]+gc[1]));
+ funValue[4]=4*gc[0]*gc[1]*(-2+3*(gc[0]+gc[1]));
+ funValue[5]=4*gc[1]*(1-4*gc[0]-gc[1]+3*gc[0]*(gc[0]+gc[1]));
+ funValue[6]=27*gc[0]*gc[1]*(1-gc[0]-gc[1]);
+ SHAPE_FUN_MACRO_END;
+}
+
/*!
* Init Quadrangle Reference coordinates ans Shape function.
* Case A.
const TypeOfField MEDCouplingFieldDiscretizationKriging::TYPE=ON_NODES_KR;
// doc is here http://www.code-aster.org/V2/doc/default/fr/man_r/r3/r3.01.01.pdf
+const double MEDCouplingFieldDiscretizationGaussNE::FGP_POINT1[1]={0.};
const double MEDCouplingFieldDiscretizationGaussNE::FGP_SEG2[2]={1.,1.};
const double MEDCouplingFieldDiscretizationGaussNE::FGP_SEG3[3]={0.5555555555555556,0.8888888888888888,0.5555555555555556};
const double MEDCouplingFieldDiscretizationGaussNE::FGP_SEG4[4]={0.347854845137454,0.347854845137454,0.652145154862546,0.652145154862546};
{
switch(geoType)
{
+ case INTERP_KERNEL::NORM_POINT1:
+ lgth=(int)sizeof(FGP_POINT1)/sizeof(double);
+ return FGP_POINT1;
case INTERP_KERNEL::NORM_SEG2:
lgth=(int)sizeof(FGP_SEG2)/sizeof(double);
return FGP_SEG2;
{
switch(geoType)
{
+ case INTERP_KERNEL::NORM_POINT1:
+ lgth=0;
+ return 0;
case INTERP_KERNEL::NORM_SEG2:
lgth=(int)sizeof(REF_SEG2)/sizeof(double);
return REF_SEG2;
{
switch(geoType)
{
+ case INTERP_KERNEL::NORM_POINT1:
+ {
+ lgth=0;
+ return 0;
+ }
case INTERP_KERNEL::NORM_SEG2:
{
lgth=(int)sizeof(LOC_SEG2)/sizeof(double);
self.assertTrue(arrOfDisc2.isEqual(coo,1e-10)) # be less exigent 1e-10 instead of 1e-12 due to shape function sensitivity arount 0.,0.,1. !
pass
+ def testSwig2Tri7GP1(self):
+ coo=DataArrayDouble([[0,0],[0,2],[2,0],[0,1],[1,1],[1,0],[0.6666666666666667,0.6666666666666667]])
+ m=MEDCouplingUMesh("mesh",2) ; m.setCoords(coo)
+ m.allocateCells()
+ # the cell description is exactly those described in the description of TRI7 in MED file 3.0.7 documentation
+ m.insertNextCell(NORM_TRI7,range(7))
+ refCoords=[0.,0.,1.,0.,0.,1.,0.5,0.,0.5,0.5,0.,0.5,0.3333333333333333,0.3333333333333333]
+ gaussCoords=[0.3333333333333333,0.3333333333333333,0.470142064105115,0.470142064105115,0.05971587178977,0.470142064105115,0.470142064105115,0.05971587178977,0.101286507323456,0.101286507323456,0.797426985353088,0.101286507323456,0.101286507323456,0.797426985353088]
+ weights=[0.062969590272413,0.062969590272413,0.062969590272413,0.066197076394253,0.066197076394253,0.066197076394253,0.1125]
+ fGauss=MEDCouplingFieldDouble(ON_GAUSS_PT) ; fGauss.setName("fGauss")
+ fGauss.setMesh(m)
+ fGauss.setGaussLocalizationOnType(NORM_TRI7,refCoords,gaussCoords,weights)
+ arr=DataArrayDouble(fGauss.getNumberOfTuplesExpected()) ; arr.iota()
+ fGauss.setArray(arr)
+ arrOfDisc=fGauss.getLocalizationOfDiscr()
+ self.assertTrue(arrOfDisc.isEqual(DataArrayDouble([0.666666666666667,0.666666666666667,0.9402841282102293,0.9402841282102293,0.9402841282102299,0.11943174357954002,0.11943174357953992,0.9402841282102299,0.20257301464691194,0.20257301464691196,0.20257301464691205,1.5948539707061757,1.5948539707061757,0.20257301464691202],7,2),1e-12))
+ #
+ weights=7*[1]
+ gaussCoords=refCoords
+ fGauss.setGaussLocalizationOnType(NORM_TRI7,refCoords,gaussCoords,weights)
+ arrOfDisc2=fGauss.getLocalizationOfDiscr()
+ self.assertTrue(arrOfDisc2.isEqual(coo,1e-12))
+ pass
+
def setUp(self):
pass
pass