connINew->alloc(nbOfCells+1,1);
int *connINewPtr=connINew->getPointer(); *connINewPtr++=0;
MCAuto<DataArrayInt> connNew=DataArrayInt::New(); connNew->alloc(0,1);
+ MCAuto<DataArrayInt> E_Fi(DataArrayInt::New()), E_F(DataArrayInt::New()), F_Ei(DataArrayInt::New()), F_E(DataArrayInt::New());
+ MCAuto<MEDCouplingUMesh> m_faces(buildDescendingConnectivity(E_F, E_Fi, F_E, F_Ei));
bool changed=false;
for(int i=0;i<nbOfCells;i++,connINewPtr++)
{
if(conn[index[i]]==(int)INTERP_KERNEL::NORM_POLYHED)
{
- SimplifyPolyhedronCell(eps,coords,conn+index[i],conn+index[i+1],connNew);
+ SimplifyPolyhedronCell(eps,coords, i,connNew, m_faces, E_Fi, E_F, F_Ei, F_E);
changed=true;
}
else
* \param [in] end end of nodal connectivity of a single polyhedron cell (excluded)
* \param [out] res the result is put at the end of the vector without any alteration of the data.
*/
-void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res)
+void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, int index, DataArrayInt *res, MEDCouplingUMesh *faces,
+ DataArrayInt *E_Fi, DataArrayInt *E_F, DataArrayInt *F_Ei, DataArrayInt *F_E)
{
- int nbFaces=std::count(begin+1,end,-1)+1;
+ int nbFaces = E_Fi->getIJ(index + 1, 0) - E_Fi->getIJ(index, 0);
MCAuto<DataArrayDouble> v=DataArrayDouble::New(); v->alloc(nbFaces,3);
double *vPtr=v->getPointer();
- MCAuto<DataArrayDouble> p=DataArrayDouble::New(); p->alloc(nbFaces,1);
+ MCAuto<DataArrayDouble> p=DataArrayDouble::New(); p->alloc(nbFaces,2);
double *pPtr=p->getPointer();
- const int *stFaceConn=begin+1;
+ int *e_fi = E_Fi->getPointer(), *e_f = E_F->getPointer(), *f_ei = F_Ei->getPointer(), *f_e = F_E->getPointer();
+ const int *f_idx = faces->getNodalConnectivityIndex()->getPointer(), *f_cnn = faces->getNodalConnectivity()->getPointer();
for(int i=0;i<nbFaces;i++,vPtr+=3,pPtr++)
{
- const int *endFaceConn=std::find(stFaceConn,end,-1);
- ComputeVecAndPtOfFace(eps,coords->begin(),stFaceConn,endFaceConn,vPtr,pPtr);
- stFaceConn=endFaceConn+1;
+ int face = e_f[e_fi[index] + i];
+ ComputeVecAndPtOfFace(eps, coords->begin(), f_cnn + f_idx[face] + 1, f_cnn + f_idx[face + 1], vPtr, pPtr);
+ // to differentiate faces going to different cells:
+ pPtr++, *pPtr = 0;
+ for (int j = f_ei[face]; j < f_ei[face + 1]; j++)
+ *pPtr += f_e[j];
}
pPtr=p->getPointer(); vPtr=v->getPointer();
DataArrayInt *comm1=0,*commI1=0;
v->findCommonTuples(eps,-1,comm1,commI1);
+ for (int i = 0; i < nbFaces; i++)
+ if (comm1->findIdFirstEqual(i) < 0)
+ {
+ comm1->pushBackSilent(i);
+ commI1->pushBackSilent(comm1->getNumberOfTuples());
+ }
MCAuto<DataArrayInt> comm1Auto(comm1),commI1Auto(commI1);
const int *comm1Ptr=comm1->begin();
const int *commI1Ptr=commI1->begin();
int nbOfGrps1=commI1Auto->getNumberOfTuples()-1;
res->pushBackSilent((int)INTERP_KERNEL::NORM_POLYHED);
//
- MCAuto<MEDCouplingUMesh> mm=MEDCouplingUMesh::New("",3);
- mm->setCoords(const_cast<DataArrayDouble *>(coords)); mm->allocateCells(1); mm->insertNextCell(INTERP_KERNEL::NORM_POLYHED,(int)std::distance(begin+1,end),begin+1);
- mm->finishInsertingCells();
- //
for(int i=0;i<nbOfGrps1;i++)
{
int vecId=comm1Ptr[commI1Ptr[i]];
MCAuto<DataArrayDouble> tmpgrp2=p->selectByTupleId(comm1Ptr+commI1Ptr[i],comm1Ptr+commI1Ptr[i+1]);
DataArrayInt *comm2=0,*commI2=0;
tmpgrp2->findCommonTuples(eps,-1,comm2,commI2);
+ for (int j = 0; j < commI1Ptr[i+1] - commI1Ptr[i]; j++)
+ if (comm2->findIdFirstEqual(j) < 0)
+ {
+ comm2->pushBackSilent(j);
+ commI2->pushBackSilent(comm2->getNumberOfTuples());
+ }
MCAuto<DataArrayInt> comm2Auto(comm2),commI2Auto(commI2);
const int *comm2Ptr=comm2->begin();
const int *commI2Ptr=commI2->begin();
int nbOfGrps2=commI2Auto->getNumberOfTuples()-1;
for(int j=0;j<nbOfGrps2;j++)
{
- if(commI2Ptr[j+1]-commI2Ptr[j]<=1)
+ if(commI2Ptr[j+1] == commI2Ptr[j] + 1)
{
- res->insertAtTheEnd(begin,end);
+ int face = e_f[e_fi[index] + comm1Ptr[commI1Ptr[i] + comm2Ptr[commI2Ptr[j]]]]; //hmmm
+ res->insertAtTheEnd(f_cnn + f_idx[face] + 1, f_cnn + f_idx[face + 1]);
res->pushBackSilent(-1);
}
else
int pointId=comm1Ptr[commI1Ptr[i]+comm2Ptr[commI2Ptr[j]]];
MCAuto<DataArrayInt> ids2=comm2->selectByTupleIdSafeSlice(commI2Ptr[j],commI2Ptr[j+1],1);
ids2->transformWithIndArr(comm1Ptr+commI1Ptr[i],comm1Ptr+commI1Ptr[i+1]);
- DataArrayInt *tmp0=DataArrayInt::New(),*tmp1=DataArrayInt::New(),*tmp2=DataArrayInt::New(),*tmp3=DataArrayInt::New();
- MCAuto<MEDCouplingUMesh> mm2=mm->buildDescendingConnectivity(tmp0,tmp1,tmp2,tmp3); tmp0->decrRef(); tmp1->decrRef(); tmp2->decrRef(); tmp3->decrRef();
- MCAuto<MEDCouplingUMesh> mm3=static_cast<MEDCouplingUMesh *>(mm2->buildPartOfMySelf(ids2->begin(),ids2->end(),true));
+ ids2->transformWithIndArr(e_f + e_fi[index], e_f + e_fi[index + 1]);
+ MCAuto<MEDCouplingUMesh> mm3=static_cast<MEDCouplingUMesh *>(faces->buildPartOfMySelf(ids2->begin(),ids2->end(),true));
MCAuto<DataArrayInt> idsNodeTmp=mm3->zipCoordsTraducer();
MCAuto<DataArrayInt> idsNode=idsNodeTmp->invertArrayO2N2N2O(mm3->getNumberOfNodes());
const int *idsNodePtr=idsNode->begin();
- double center[3]; center[0]=pPtr[pointId]*vPtr[3*vecId]; center[1]=pPtr[pointId]*vPtr[3*vecId+1]; center[2]=pPtr[pointId]*vPtr[3*vecId+2];
+ double center[3]; center[0]=pPtr[2*pointId]*vPtr[3*vecId]; center[1]=pPtr[2*pointId]*vPtr[3*vecId+1]; center[2]=pPtr[2*pointId]*vPtr[3*vecId+2];
double vec[3]; vec[0]=vPtr[3*vecId+1]; vec[1]=-vPtr[3*vecId]; vec[2]=0.;
double norm=vec[0]*vec[0]+vec[1]*vec[1]+vec[2]*vec[2];
if(std::abs(norm)>eps)
#
# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
#
-
from MEDCoupling import *
import unittest
from math import pi,e,sqrt,cos,sin
mea=fieldOnCell.getMesh().getMeasureField(True).getArray(); mea.rearrange(48); mea2 = mea.sumPerTuple()
self.assertTrue(mea2.isEqual(meaRef2,1e-9))
pass
-
+
def testVoronoi3DSurf_1(self):
tmp=MEDCouplingCMesh("mesh")
arr=DataArrayDouble(5) ; arr.iota()
self.assertTrue(f3.getMesh().getMeasureField(False).getArray().isEqual(ref,1e-12))
self.assertTrue(f3.getArray().isEqual(DataArrayDouble([0,1,2,3]),1e-12))
pass
-
+
def testVoronoi3D_4(self):
"""Idem testVoronoi3D_3 except that here quadratic cells are considered"""
coo=DataArrayDouble([0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,1.0,0.0,0.0,0.0,0.5,0.0,0.0,0.0,0.5,0.0,0.5,0.5,0.5,0.5,0.0,0.5,0.0,0.0,0.5,0.0,0.5],10,3)
m.insertNextCell(NORM_POLYGON,[0,1,2,3,4,5])
self.assertTrue(m.computePlaneEquationOf3DFaces().isEqual(DataArrayDouble([0,1,0,-1],1,4),1e-12))
pass
-
+
+ def testSimplifyPolyhedra(self):
+ mesh = MEDCouplingUMesh('mesh', 3)
+ coo = DataArrayDouble([(-0.01225,-0.0212176,0.02),(-0.00634107,-0.0236652,0.02),(1.50019e-18,-0.0245,0.02),(0.00634107,-0.0236652,0.02),(0.01225,-0.0212176,0.02),(-0.0153864,-0.02665,0),(-0.00714085,-0.02665,0),(1.63184e-18,-0.02665,0),(0.00714085,-0.02665,0),(0.0153864,-0.02665,0),(-0.00714085,-0.02665,0.0101475),(1.63184e-18,-0.02665,0.013145),(0.00714085,-0.02665,0.0101475),(-0.013,-0.0225167,0.02),(-0.0067293,-0.0251141,0.02),(1.59204e-18,-0.026,0.02),(0.0067293,-0.0251141,0.02),(0.013,-0.0225167,0.02),(-0.0161658,-0.028,0),(-0.00750258,-0.028,0),(1.71451e-18,-0.028,0),(0.00750258,-0.028,0),(0.0161658,-0.028,0),(-0.00750258,-0.028,0.0105625),(1.71451e-18,-0.028,0.0136825),(0.00750258,-0.028,0.0105625)])
+ mesh.setCoords(coo)
+ c = DataArrayInt([31, 13, 14, 15, 16, 17, 4, 3, 2, 1, 0, -1, 18, 5, 6, 7, 8, 9, 22, 21, 20, 19, -1, 19, 23, 18, -1, 23, 14, 13, 18, -1, 20, 24, 23, 19, -1, 24, 15, 14, 23, -1, 21, 25, 24, 20, -1, 25, 16, 15, 24, -1, 22, 25, 21, -1, 22, 17, 16, 25, -1, 9, 4, 17, 22, -1, 8, 12, 9, -1, 12, 3, 4, 9, -1, 7, 11, 12, 8, -1, 11, 2, 3, 12, -1, 6, 10, 11, 7, -1, 10, 1, 2, 11, -1, 5, 10, 6, -1, 5, 0, 1, 10, -1, 18, 13, 0, 5])
+ cI = DataArrayInt([0, 108])
+ mesh.setConnectivity(c, cI)
+ mesh.simplifyPolyhedra(1.0e-8)
+ c, cI = mesh.getNodalConnectivity(), mesh.getNodalConnectivityIndex()
+ tgt_c = DataArrayInt([31, 23, 18, 19, 20, 21, 22, 25, 24, -1, 12, 9, 8, 7, 6, 5, 10, 11, -1, 13, 14, 15, 16, 17, 4, 3, 2, 1, 0, -1, 18, 5, 6, 7, 8, 9, 22, 21, 20, 19, -1, 23, 14, 13, 18, -1, 24, 15, 14, 23, -1, 25, 16, 15, 24, -1, 22, 17, 16, 25, -1, 9, 4, 17, 22, -1, 12, 3, 4, 9, -1, 11, 2, 3, 12, -1, 10, 1, 2, 11, -1, 5, 0, 1, 10, -1, 18, 13, 0, 5])
+ tgt_cI = DataArrayInt([0, 90])
+ self.assertEqual(c.getValues(), tgt_c.getValues())
+ self.assertEqual(cI.getValues(), tgt_cI.getValues())
+ pass
+
pass
if __name__ == '__main__':
