//
#include "MEDCouplingCMesh.hxx"
#include "MEDCouplingMemArray.hxx"
+#include "MEDCouplingFieldDouble.hxx"
using namespace ParaMEDMEM;
return 0;
}
+MEDCouplingFieldDouble *MEDCouplingCMesh::getMeasureFieldOnNode(bool isAbs) const
+{
+ //not implemented yet !
+ return 0;
+}
+
+MEDCouplingFieldDouble *MEDCouplingCMesh::buildOrthogonalField() const
+{
+ if(getMeshDimension()!=2)
+ throw INTERP_KERNEL::Exception("Expected a cmesh with meshDim == 2 !");
+ MEDCouplingFieldDouble *ret=MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+ DataArrayDouble *array=DataArrayDouble::New();
+ int nbOfCells=getNumberOfCells();
+ array->alloc(nbOfCells,3);
+ double *vals=array->getPointer();
+ for(int i=0;i<nbOfCells;i++)
+ { vals[3*i]=1.; vals[3*i+1]=1.; vals[3*i+2]=1.; }
+ ret->setArray(array);
+ array->decrRef();
+ ret->setMesh(this);
+ return ret;
+}
+
void MEDCouplingCMesh::rotate(const double *center, const double *vector, double angle)
{
throw INTERP_KERNEL::Exception("No rotation available on CMesh : Traduce it to StructuredMesh to apply it !");
// tools
void getBoundingBox(double *bbox) const;
MEDCouplingFieldDouble *getMeasureField(bool isAbs) const;
+ MEDCouplingFieldDouble *getMeasureFieldOnNode(bool isAbs) const;
+ MEDCouplingFieldDouble *buildOrthogonalField() const;
void rotate(const double *center, const double *vector, double angle);
void translate(const double *vector);
MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const;
return 0;
}
+MEDCouplingFieldDouble *MEDCouplingExtrudedMesh::getMeasureFieldOnNode(bool) const
+{
+ //not implemented yet
+ return 0;
+}
+
+MEDCouplingFieldDouble *MEDCouplingExtrudedMesh::buildOrthogonalField() const
+{
+ //not implemented yet
+ throw INTERP_KERNEL::Exception("MEDCouplingExtrudedMesh::buildOrthogonalField not implemented yet !");
+}
+
MEDCouplingExtrudedMesh::~MEDCouplingExtrudedMesh()
{
if(_mesh2D)
MEDCouplingUMesh *getMesh1D() const { return _mesh1D; }
DataArrayInt *getMesh3DIds() const { return _mesh3D_ids; }
MEDCouplingFieldDouble *getMeasureField(bool) const;
+ MEDCouplingFieldDouble *getMeasureFieldOnNode(bool) const;
+ MEDCouplingFieldDouble *buildOrthogonalField() const;
static int findCorrespCellByNodalConn(const std::vector<int>& nodalConnec,
const int *revNodalPtr, const int *revNodalIndxPtr) throw(INTERP_KERNEL::Exception);
void rotate(const double *center, const double *vector, double angle);
MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationP1::getWeightingField(const MEDCouplingMesh *mesh, bool isAbs) const
{
- //not implemented yet.
- //Dual mesh to build
- return 0;
+ return mesh->getMeasureFieldOnNode(isAbs);
}
void MEDCouplingFieldDiscretizationP1::renumberValuesOnNodes(const DataArrayInt *old2New, DataArrayDouble *arr) const
return ret;
}
+void MEDCouplingFieldDouble::accumulate(double *res) const
+{
+ const double *ptr=getArray()->getConstPointer();
+ int nbTuple=getArray()->getNumberOfTuples();
+ int nbComps=getArray()->getNumberOfComponents();
+ std::fill(res,res+nbComps,0.);
+ for(int i=0;i<nbTuple;i++)
+ std::transform(ptr+i*nbComps,ptr+(i+1)*nbComps,res,res,std::plus<double>());
+}
+
double MEDCouplingFieldDouble::measureAccumulate(int compId, bool isWAbs) const
{
if(!_mesh)
return ret;
}
+void MEDCouplingFieldDouble::measureAccumulate(bool isWAbs, double *res) const
+{
+ if(!_mesh)
+ throw INTERP_KERNEL::Exception("No mesh underlying this field to perform measureAccumulate2");
+ MEDCouplingFieldDouble *weight=_type->getWeightingField(_mesh,isWAbs);
+ const double *ptr=weight->getArray()->getConstPointer();
+ int nbOfValues=weight->getArray()->getNbOfElems();
+ int nbComps=getArray()->getNumberOfComponents();
+ const double *vals=getArray()->getConstPointer();
+ std::fill(res,res+nbComps,0.);
+ double *tmp=new double[nbComps];
+ for (int i=0; i<nbOfValues; i++)
+ {
+ std::transform(vals+i*nbComps,vals+(i+1)*nbComps,tmp,std::bind2nd(std::multiplies<double>(),ptr[i]));
+ std::transform(tmp,tmp+nbComps,res,res,std::plus<double>());
+ }
+ weight->decrRef();
+ delete [] tmp;
+}
+
void MEDCouplingFieldDouble::getValueOn(const double *spaceLoc, double *res) const throw(INTERP_KERNEL::Exception)
{
_time_discr->checkNoTimePresence();
void setArray(DataArrayDouble *array);
DataArrayDouble *getArray() const { return _time_discr->getArray(); }
double accumulate(int compId) const;
+ void accumulate(double *res) const;
double measureAccumulate(int compId, bool isWAbs) const;
+ void measureAccumulate(bool isWAbs, double *res) const;
void getValueOn(const double *spaceLoc, double *res) const throw(INTERP_KERNEL::Exception);
void getValueOn(const double *spaceLoc, double time, double *res) const throw(INTERP_KERNEL::Exception);
//! \b temporary
// tools
virtual void getBoundingBox(double *bbox) const = 0;
virtual MEDCouplingFieldDouble *getMeasureField(bool isAbs) const = 0;
+ virtual MEDCouplingFieldDouble *getMeasureFieldOnNode(bool isAbs) const = 0;
virtual MEDCouplingFieldDouble *fillFromAnalytic(TypeOfField t, int nbOfComp, FunctionToEvaluate func) const;
+ virtual MEDCouplingFieldDouble *buildOrthogonalField() const = 0;
virtual void rotate(const double *center, const double *vector, double angle) = 0;
virtual void translate(const double *vector) = 0;
virtual MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const = 0;
#include "MEDCouplingFieldDouble.hxx"
#include "CellModel.hxx"
#include "VolSurfUser.txx"
+#include "InterpolationUtils.hxx"
#include <sstream>
#include <limits>
return field;
}
+MEDCouplingFieldDouble *MEDCouplingUMesh::getMeasureFieldOnNode(bool isAbs) const
+{
+ MEDCouplingFieldDouble *tmp=getMeasureField(abs);
+ std::string name="MeasureOnNodeOfMesh_";
+ name+=getName();
+ int nbNodes=getNumberOfNodes();
+ MEDCouplingFieldDouble *ret=MEDCouplingFieldDouble::New(ON_NODES);
+ double cst=1./((double)getMeshDimension()+1.);
+ DataArrayDouble* array=DataArrayDouble::New();
+ array->alloc(nbNodes,1);
+ double *valsToFill=array->getPointer();
+ std::fill(valsToFill,valsToFill+nbNodes,0.);
+ const double *values=tmp->getArray()->getConstPointer();
+ DataArrayInt *da=DataArrayInt::New();
+ DataArrayInt *daInd=DataArrayInt::New();
+ getReverseNodalConnectivity(da,daInd);
+ const int *daPtr=da->getConstPointer();
+ const int *daIPtr=daInd->getConstPointer();
+ for(int i=0;i<nbNodes;i++)
+ for(const int *cell=daPtr+daIPtr[i];cell!=daPtr+daIPtr[i+1];cell++)
+ valsToFill[i]+=cst*values[*cell];
+ ret->setMesh(this);
+ da->decrRef();
+ daInd->decrRef();
+ ret->setArray(array);
+ array->decrRef();
+ tmp->decrRef();
+ return ret;
+}
+
+MEDCouplingFieldDouble *MEDCouplingUMesh::buildOrthogonalField() const
+{
+ if(getMeshDimension()!=2)
+ throw INTERP_KERNEL::Exception("Expected a umesh with meshDim == 2 !");
+ MEDCouplingFieldDouble *ret=MEDCouplingFieldDouble::New(ON_CELLS,NO_TIME);
+ DataArrayDouble *array=DataArrayDouble::New();
+ int nbOfCells=getNumberOfCells();
+ array->alloc(nbOfCells,3);
+ double *vals=array->getPointer();
+ const int *connI=_nodal_connec_index->getConstPointer();
+ const int *conn=_nodal_connec->getConstPointer();
+ const double *coords=_coords->getConstPointer();
+ if(getSpaceDimension()==3)
+ {
+ for(int i=0;i<nbOfCells;i++,vals+=3)
+ {
+ int offset=connI[i];
+ INTERP_KERNEL::crossprod<3>(coords+3*conn[offset+1],coords+3*conn[offset+2],coords+3*conn[offset+3],vals);
+ double n=INTERP_KERNEL::norm<3>(vals);
+ std::transform(vals,vals+3,vals,std::bind2nd(std::multiplies<double>(),1./n));
+ }
+ }
+ else
+ {
+ for(int i=0;i<nbOfCells;i++)
+ { vals[3*i]=1.; vals[3*i+1]=1.; vals[3*i+2]=1.; }
+ }
+ ret->setArray(array);
+ array->decrRef();
+ ret->setMesh(this);
+ return ret;
+}
+
/*!
* This method is only available for a mesh with meshDim==2 and spaceDim==2||spaceDim==3.
* This method returns a vector 'cells' where all detected butterfly cells have been added to cells.
void renumberNodes(const int *newNodeNumbers, int newNbOfNodes);
void giveElemsInBoundingBox(const double *bbox, double eps, std::vector<int>& elems);
MEDCouplingFieldDouble *getMeasureField(bool isAbs) const;
+ MEDCouplingFieldDouble *getMeasureFieldOnNode(bool isAbs) const;
+ MEDCouplingFieldDouble *buildOrthogonalField() const;
void checkButterflyCells(std::vector<int>& cells) const;
bool checkConsecutiveCellTypes() const;
MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const;
return 0;
}
+MEDCouplingFieldDouble *MEDCouplingUMeshDesc::getMeasureFieldOnNode(bool isAbs) const
+{
+ //not implemented yet.
+ return 0;
+}
+
+MEDCouplingFieldDouble *MEDCouplingUMeshDesc::buildOrthogonalField() const
+{
+ if(getMeshDimension()!=2)
+ throw INTERP_KERNEL::Exception("Expected a cmesh with meshDim == 2 !");
+ //not implemented yet !
+ return 0;
+}
+
DataArrayInt *MEDCouplingUMeshDesc::zipCoordsTraducer()
{
//not implemented yet.
MEDCouplingPointSet *buildBoundaryMesh(bool keepCoords) const;
void renumberNodes(const int *newNodeNumbers, int newNbOfNodes);
MEDCouplingFieldDouble *getMeasureField(bool isAbs) const;
+ MEDCouplingFieldDouble *getMeasureFieldOnNode(bool isAbs) const;
+ MEDCouplingFieldDouble *buildOrthogonalField() const;
DataArrayInt *zipCoordsTraducer();
MEDCouplingMesh *mergeMyselfWith(const MEDCouplingMesh *other) const;
DataArrayDouble *getBarycenterAndOwner() const;
std::transform(values3,values3+18,values3,std::ptr_fun<double,double>(fabs));
max=*std::max_element(values3,values3+18);
CPPUNIT_ASSERT_DOUBLES_EQUAL(0.,max,1.e-12);
+ double values4[2];
+ f1->accumulate(values4);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(3.6,values4[0],1.e-12);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(7.2,values4[1],1.e-12);
+ f1->measureAccumulate(true,values4);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(0.5,values4[0],1.e-12);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(1.,values4[1],1.e-12);
f1->decrRef();
//
CPPUNIT_ASSERT_THROW(f1=m->fillFromAnalytic(ON_NODES,1,func3),INTERP_KERNEL::Exception);
sourceMesh->decrRef();
}
+void MEDCouplingBasicsTest::testBuildOrthogonalField()
+{
+ MEDCouplingUMesh *targetMesh=build3DSurfTargetMesh_1();
+ MEDCouplingFieldDouble *field=targetMesh->buildOrthogonalField();
+ double expected[3]={0.70710678118654746,0.,-0.70710678118654746};
+ CPPUNIT_ASSERT_EQUAL(5,field->getNumberOfTuples());
+ CPPUNIT_ASSERT_EQUAL(3,field->getNumberOfComponents());
+ const double *vals=field->getArray()->getConstPointer();
+ for(int i=0;i<15;i++)
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(expected[i%3],vals[i],1e-12);
+ field->decrRef();
+ targetMesh->decrRef();
+}
+
void MEDCouplingBasicsTest::test2DInterpP0P0_1()
{
MEDCouplingUMesh *sourceMesh=build2DSourceMesh_1();
CPPUNIT_TEST( testOperationsOnFields );
CPPUNIT_TEST( testMergeNodesOnField );
CPPUNIT_TEST( testCheckConsecutiveCellTypes );
+ CPPUNIT_TEST( testBuildOrthogonalField );
CPPUNIT_TEST( test2DInterpP0P0_1 );
CPPUNIT_TEST( test2DInterpP0P0PL_1 );
CPPUNIT_TEST( test2DInterpP0P0PL_2 );
void testOperationsOnFields();
void testMergeNodesOnField();
void testCheckConsecutiveCellTypes();
+ void testBuildOrthogonalField();
void test2DInterpP0P0_1();
void test2DInterpP0P0PL_1();
void test2DInterpP0P0PL_2();
