-// Copyright (C) 2007-2014 CEA/DEN, EDF R&D
+// Copyright (C) 2007-2016 CEA/DEN, EDF R&D
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
//
// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
-// Author : Anthony Geay (CEA/DEN)
+// Author : Anthony Geay (EDF R&D)
#include "MEDCouplingFieldDiscretization.hxx"
#include "MEDCouplingCMesh.hxx"
#include "MEDCouplingUMesh.hxx"
#include "MEDCouplingFieldDouble.hxx"
-#include "MEDCouplingAutoRefCountObjectPtr.hxx"
+#include "MCAuto.hxx"
#include "CellModel.hxx"
#include "InterpolationUtils.hxx"
#include <algorithm>
#include <functional>
-using namespace ParaMEDMEM;
+using namespace MEDCoupling;
const double MEDCouplingFieldDiscretization::DFLT_PRECISION=1.e-12;
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretization::New(TypeOfField type)
{
switch(type)
- {
+ {
case MEDCouplingFieldDiscretizationP0::TYPE:
return new MEDCouplingFieldDiscretizationP0;
case MEDCouplingFieldDiscretizationP1::TYPE:
case MEDCouplingFieldDiscretizationKriging::TYPE:
return new MEDCouplingFieldDiscretizationKriging;
default:
- throw INTERP_KERNEL::Exception("Choosen discretization is not implemented yet.");
- }
+ throw INTERP_KERNEL::Exception("Chosen discretization is not implemented yet.");
+ }
}
TypeOfField MEDCouplingFieldDiscretization::GetTypeOfFieldFromStringRepr(const std::string& repr)
throw INTERP_KERNEL::Exception("Representation does not match with any field discretization !");
}
+std::string MEDCouplingFieldDiscretization::GetTypeOfFieldRepr(TypeOfField type)
+{
+ if(type==MEDCouplingFieldDiscretizationP0::TYPE)
+ return MEDCouplingFieldDiscretizationP0::REPR;
+ if(type==MEDCouplingFieldDiscretizationP1::TYPE)
+ return MEDCouplingFieldDiscretizationP1::REPR;
+ if(type==MEDCouplingFieldDiscretizationGauss::TYPE)
+ return MEDCouplingFieldDiscretizationGauss::REPR;
+ if(type==MEDCouplingFieldDiscretizationGaussNE::TYPE)
+ return MEDCouplingFieldDiscretizationGaussNE::REPR;
+ if(type==MEDCouplingFieldDiscretizationKriging::TYPE)
+ return MEDCouplingFieldDiscretizationKriging::REPR;
+ throw INTERP_KERNEL::Exception("GetTypeOfFieldRepr : Representation does not match with any field discretization !");
+}
+
bool MEDCouplingFieldDiscretization::isEqual(const MEDCouplingFieldDiscretization *other, double eps) const
{
std::string reason;
* This method is an alias of MEDCouplingFieldDiscretization::clone. It is only here for coherency with all the remaining of MEDCoupling.
* \sa MEDCouplingFieldDiscretization::clone.
*/
-MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretization::deepCpy() const
+MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretization::deepCopy() const
{
return clone();
}
return 0;
}
-std::vector<const BigMemoryObject *> MEDCouplingFieldDiscretization::getDirectChildren() const
+std::vector<const BigMemoryObject *> MEDCouplingFieldDiscretization::getDirectChildrenWithNull() const
{
return std::vector<const BigMemoryObject *>();
}
*/
void MEDCouplingFieldDiscretization::normL1(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const
{
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=getMeasureField(mesh,true);
+ MCAuto<MEDCouplingFieldDouble> vol=getMeasureField(mesh,true);
int nbOfCompo=arr->getNumberOfComponents();
int nbOfElems=getNumberOfTuples(mesh);
std::fill(res,res+nbOfCompo,0.);
*/
void MEDCouplingFieldDiscretization::normL2(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, double *res) const
{
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=getMeasureField(mesh,true);
+ MCAuto<MEDCouplingFieldDouble> vol=getMeasureField(mesh,true);
int nbOfCompo=arr->getNumberOfComponents();
int nbOfElems=getNumberOfTuples(mesh);
std::fill(res,res+nbOfCompo,0.);
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::integral : mesh is NULL !");
if(!arr)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::integral : input array is NULL !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=getMeasureField(mesh,isWAbs);
- int nbOfCompo=arr->getNumberOfComponents();
- int nbOfElems=getNumberOfTuples(mesh);
+ MCAuto<MEDCouplingFieldDouble> vol=getMeasureField(mesh,isWAbs);
+ std::size_t nbOfCompo(arr->getNumberOfComponents()),nbOfElems(getNumberOfTuples(mesh));
if(nbOfElems!=arr->getNumberOfTuples())
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretization::integral : field is not correct ! number of tuples in array is " << arr->getNumberOfTuples();
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
std::fill(res,res+nbOfCompo,0.);
- const double *arrPtr=arr->getConstPointer();
- const double *volPtr=vol->getArray()->getConstPointer();
+ const double *arrPtr(arr->begin()),*volPtr(vol->getArray()->begin());
INTERP_KERNEL::AutoPtr<double> tmp=new double[nbOfCompo];
- for (int i=0;i<nbOfElems;i++)
+ for(std::size_t i=0;i<nbOfElems;i++)
{
std::transform(arrPtr+i*nbOfCompo,arrPtr+(i+1)*nbOfCompo,(double *)tmp,std::bind2nd(std::multiplies<double>(),volPtr[i]));
std::transform((double *)tmp,(double *)tmp+nbOfCompo,res,res,std::plus<double>());
*/
MEDCouplingMesh *MEDCouplingFieldDiscretization::buildSubMeshDataRange(const MEDCouplingMesh *mesh, int beginCellIds, int endCellIds, int stepCellIds, int& beginOut, int& endOut, int& stepOut, DataArrayInt *&di) const
{
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> da=DataArrayInt::Range(beginCellIds,endCellIds,stepCellIds);
+ MCAuto<DataArrayInt> da=DataArrayInt::Range(beginCellIds,endCellIds,stepCellIds);
return buildSubMeshData(mesh,da->begin(),da->end(),di);
}
arr=0;
}
+/*!
+ * Empty : Not a bug
+ */
+void MEDCouplingFieldDiscretization::checkForUnserialization(const std::vector<int>& tinyInfo, const DataArrayInt *arr)
+{
+}
+
/*!
* Empty : Not a bug
*/
/*!
* This method is typically the first step of renumbering. The implementation is empty it is not a bug only gauss is impacted
- * virtualy by this method.
+ * virtually by this method.
*/
void MEDCouplingFieldDiscretization::renumberCells(const int *old2NewBg, bool check)
{
}
void MEDCouplingFieldDiscretization::setGaussLocalizationOnType(const MEDCouplingMesh *m, INTERP_KERNEL::NormalizedCellType type, const std::vector<double>& refCoo,
- const std::vector<double>& gsCoo, const std::vector<double>& wg) throw(INTERP_KERNEL::Exception)
+ const std::vector<double>& gsCoo, const std::vector<double>& wg)
{
throw INTERP_KERNEL::Exception("Invalid method for the corresponding field discretization : available only for GaussPoint discretization !");
}
void MEDCouplingFieldDiscretization::setGaussLocalizationOnCells(const MEDCouplingMesh *m, const int *begin, const int *end, const std::vector<double>& refCoo,
- const std::vector<double>& gsCoo, const std::vector<double>& wg) throw(INTERP_KERNEL::Exception)
+ const std::vector<double>& gsCoo, const std::vector<double>& wg)
{
throw INTERP_KERNEL::Exception("Invalid method for the corresponding field discretization : available only for GaussPoint discretization !");
}
int oldNbOfElems=arr->getNumberOfTuples();
int nbOfComp=arr->getNumberOfComponents();
int newNbOfTuples=newNbOfEntity;
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arrCpy=arr->deepCpy();
+ MCAuto<DataArrayDouble> arrCpy=arr->deepCopy();
const double *ptSrc=arrCpy->getConstPointer();
arr->reAlloc(newNbOfTuples);
double *ptToFill=arr->getPointer();
if(newNb>=0)//if newNb<0 the node is considered as out.
{
if(std::find_if(ptToFill+newNb*nbOfComp,ptToFill+(newNb+1)*nbOfComp,std::bind2nd(std::not_equal_to<double>(),std::numeric_limits<double>::max()))
- ==ptToFill+(newNb+1)*nbOfComp)
+ ==ptToFill+(newNb+1)*nbOfComp)
std::copy(ptSrc+i*nbOfComp,ptSrc+(i+1)*nbOfComp,ptToFill+newNb*nbOfComp);
else
{
void MEDCouplingFieldDiscretization::RenumberEntitiesFromN2OArr(const int *new2OldPtr, int new2OldSz, DataArrayDouble *arr, const std::string& msg)
{
int nbOfComp=arr->getNumberOfComponents();
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arrCpy=arr->deepCpy();
+ MCAuto<DataArrayDouble> arrCpy=arr->deepCopy();
const double *ptSrc=arrCpy->getConstPointer();
arr->reAlloc(new2OldSz);
double *ptToFill=arr->getPointer();
}
/*!
- * This method is simply called by MEDCouplingFieldDiscretization::deepCpy. It performs the deep copy of \a this.
+ * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this.
*
- * \sa MEDCouplingFieldDiscretization::deepCpy.
+ * \sa MEDCouplingFieldDiscretization::deepCopy.
*/
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationP0::clone() const
{
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
const DataArrayInt *ids(idsPerType[pos]);
- if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
+ if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || (int)ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretizationP0::getNumberOfTuplesExpectedRegardingCode : input pfl chunck at pos " << pos << " should have " << i << " tuples and one component and with ids all >=0 !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
void MEDCouplingFieldDiscretizationP0::renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector<DataArray *>& arrays,
- const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
+ const int *old2NewBg, bool check)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::renumberArraysForCell : NULL input mesh !");
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::getLocalizationOfDiscValues : NULL input mesh !");
- return mesh->getBarycenterAndOwner();
+ return mesh->computeCellCenterOfMass();
}
void MEDCouplingFieldDiscretizationP0::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *tupleIdsBg, const int *tupleIdsEnd,
- DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const throw(INTERP_KERNEL::Exception)
+ DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::computeMeshRestrictionFromTupleIds : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp=DataArrayInt::New();
+ MCAuto<DataArrayInt> tmp=DataArrayInt::New();
tmp->alloc((int)std::distance(tupleIdsBg,tupleIdsEnd),1);
std::copy(tupleIdsBg,tupleIdsEnd,tmp->getPointer());
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp2(tmp->deepCpy());
+ MCAuto<DataArrayInt> tmp2(tmp->deepCopy());
cellRestriction=tmp.retn();
trueTupleRestriction=tmp2.retn();
}
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::getValueOnMulti : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsArr,eltsIndexArr;
+ MCAuto<DataArrayInt> eltsArr,eltsIndexArr;
mesh->getCellsContainingPoints(loc,nbOfPoints,_precision,eltsArr,eltsIndexArr);
const int *elts(eltsArr->begin()),*eltsIndex(eltsIndexArr->begin());
int spaceDim=mesh->getSpaceDimension();
int nbOfComponents=arr->getNumberOfComponents();
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> ret=DataArrayDouble::New();
ret->alloc(nbOfPoints,nbOfComponents);
double *ptToFill=ret->getPointer();
for(int i=0;i<nbOfPoints;i++,ptToFill+=nbOfComponents)
*/
DataArrayInt *MEDCouplingFieldDiscretizationP0::computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const int *startCellIds, const int *endCellIds) const
{
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
+ MCAuto<DataArrayInt> ret=DataArrayInt::New();
ret->alloc((int)std::distance(startCellIds,endCellIds),1);
std::copy(startCellIds,endCellIds,ret->getPointer());
return ret.retn();
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::buildSubMeshData : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPart(start,end);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> diSafe=DataArrayInt::New();
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPart(start,end);
+ MCAuto<DataArrayInt> diSafe=DataArrayInt::New();
diSafe->alloc((int)std::distance(start,end),1);
std::copy(start,end,diSafe->getPointer());
di=diSafe.retn();
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP0::buildSubMeshDataRange : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
di=0; beginOut=beginCellIds; endOut=endCellIds; stepOut=stepCellIds;
return ret.retn();
}
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
const DataArrayInt *ids(idsPerType[pos]);
- if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
+ if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || (int)ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretizationOnNodes::getNumberOfTuplesExpectedRegardingCode : input pfl chunck at pos " << pos << " should have " << i << " tuples and one component and with ids all >=0 !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
* Nothing to do here.
*/
void MEDCouplingFieldDiscretizationOnNodes::renumberArraysForCell(const MEDCouplingMesh *, const std::vector<DataArray *>& arrays,
- const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
+ const int *old2NewBg, bool check)
{
}
}
void MEDCouplingFieldDiscretizationOnNodes::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *tupleIdsBg, const int *tupleIdsEnd,
- DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const throw(INTERP_KERNEL::Exception)
+ DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationOnNodes::computeMeshRestrictionFromTupleIds : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret1=mesh->getCellIdsFullyIncludedInNodeIds(tupleIdsBg,tupleIdsEnd);
+ MCAuto<DataArrayInt> ret1=mesh->getCellIdsFullyIncludedInNodeIds(tupleIdsBg,tupleIdsEnd);
const MEDCouplingUMesh *meshc=dynamic_cast<const MEDCouplingUMesh *>(mesh);
if(!meshc)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationOnNodes::computeMeshRestrictionFromTupleIds : trying to subpart field on nodes by node ids ! Your mesh has to be unstructured !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> meshPart=static_cast<MEDCouplingUMesh *>(meshc->buildPartOfMySelf(ret1->begin(),ret1->end(),true));
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret2=meshPart->computeFetchedNodeIds();
+ MCAuto<MEDCouplingUMesh> meshPart=static_cast<MEDCouplingUMesh *>(meshc->buildPartOfMySelf(ret1->begin(),ret1->end(),true));
+ MCAuto<DataArrayInt> ret2=meshPart->computeFetchedNodeIds();
cellRestriction=ret1.retn();
trueTupleRestriction=ret2.retn();
}
{
if(!mesh || !da)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationNodes::checkCoherencyBetween : NULL input mesh or DataArray !");
- if(mesh->getNumberOfNodes()!=da->getNumberOfTuples())
+ if(mesh->getNumberOfNodes()!=(int)da->getNumberOfTuples())
{
std::ostringstream message;
message << "Field on nodes invalid because there are " << mesh->getNumberOfNodes();
/*!
* This method returns a submesh of 'mesh' instance constituting cell ids contained in array defined as an interval [start;end).
-* @param di is an array returned that specifies entity ids (here nodes ids) in mesh 'mesh' of entity in returned submesh.
+ * @param di is an array returned that specifies entity ids (here nodes ids) in mesh 'mesh' of entity in returned submesh.
* Example : The first node id of returned mesh has the (*di)[0] id in 'mesh'
*/
MEDCouplingMesh *MEDCouplingFieldDiscretizationOnNodes::buildSubMeshData(const MEDCouplingMesh *mesh, const int *start, const int *end, DataArrayInt *&di) const
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationNodes::buildSubMeshData : NULL input mesh !");
DataArrayInt *diTmp=0;
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPartAndReduceNodes(start,end,diTmp);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> diTmpSafe(diTmp);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> di2=diTmpSafe->invertArrayO2N2N2O(ret->getNumberOfNodes());
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPartAndReduceNodes(start,end,diTmp);
+ MCAuto<DataArrayInt> diTmpSafe(diTmp);
+ MCAuto<DataArrayInt> di2=diTmpSafe->invertArrayO2N2N2O(ret->getNumberOfNodes());
di=di2.retn();
return ret.retn();
}
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationOnNodes::buildSubMeshDataRange : NULL input mesh !");
DataArrayInt *diTmp=0;
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPartRangeAndReduceNodes(beginCellIds,endCellIds,stepCellIds,beginOut,endOut,stepOut,diTmp);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPartRangeAndReduceNodes(beginCellIds,endCellIds,stepCellIds,beginOut,endOut,stepOut,diTmp);
if(diTmp)
{
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> diTmpSafe(diTmp);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> di2=diTmpSafe->invertArrayO2N2N2O(ret->getNumberOfNodes());
+ MCAuto<DataArrayInt> diTmpSafe(diTmp);
+ MCAuto<DataArrayInt> di2=diTmpSafe->invertArrayO2N2N2O(ret->getNumberOfNodes());
di=di2.retn();
}
return ret.retn();
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP1::computeTupleIdsToSelectFromCellIds : NULL input mesh !");
- const MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh=mesh->buildUnstructured();
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh2=static_cast<MEDCouplingUMesh *>(umesh->buildPartOfMySelf(startCellIds,endCellIds,true));
+ const MCAuto<MEDCouplingUMesh> umesh=mesh->buildUnstructured();
+ MCAuto<MEDCouplingUMesh> umesh2=static_cast<MEDCouplingUMesh *>(umesh->buildPartOfMySelf(startCellIds,endCellIds,true));
return umesh2->computeFetchedNodeIds();
}
}
/*!
- * This method is simply called by MEDCouplingFieldDiscretization::deepCpy. It performs the deep copy of \a this.
+ * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this.
*
- * \sa MEDCouplingFieldDiscretization::deepCpy.
+ * \sa MEDCouplingFieldDiscretization::deepCopy.
*/
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationP1::clone() const
{
void MEDCouplingFieldDiscretizationP1::checkCompatibilityWithNature(NatureOfField nat) const
{
- if(nat!=ConservativeVolumic)
- throw INTERP_KERNEL::Exception("Invalid nature for P1 field : expected ConservativeVolumic !");
+ if(nat!=IntensiveMaximum)
+ throw INTERP_KERNEL::Exception("Invalid nature for P1 field : expected IntensiveMaximum !");
}
MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationP1::getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationP1::getValueOnMulti : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsArr,eltsIndexArr;
+ MCAuto<DataArrayInt> eltsArr,eltsIndexArr;
mesh->getCellsContainingPoints(loc,nbOfPoints,_precision,eltsArr,eltsIndexArr);
const int *elts(eltsArr->begin()),*eltsIndex(eltsIndexArr->begin());
int spaceDim=mesh->getSpaceDimension();
int nbOfComponents=arr->getNumberOfComponents();
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> ret=DataArrayDouble::New();
ret->alloc(nbOfPoints,nbOfComponents);
double *ptToFill=ret->getPointer();
for(int i=0;i<nbOfPoints;i++)
}
/*!
- * This constructor deep copies ParaMEDMEM::DataArrayInt instance from other (if any).
+ * This constructor deep copies MEDCoupling::DataArrayInt instance from other (if any).
*/
MEDCouplingFieldDiscretizationPerCell::MEDCouplingFieldDiscretizationPerCell(const MEDCouplingFieldDiscretizationPerCell& other, const int *startCellIds, const int *endCellIds):_discr_per_cell(0)
{
if(arr)
{
if(startCellIds==0 && endCellIds==0)
- _discr_per_cell=arr->deepCpy();
+ _discr_per_cell=arr->deepCopy();
else
_discr_per_cell=arr->selectByTupleIdSafe(startCellIds,endCellIds);
}
DataArrayInt *arr=other._discr_per_cell;
if(arr)
{
- _discr_per_cell=arr->selectByTupleId2(beginCellIds,endCellIds,stepCellIds);
+ _discr_per_cell=arr->selectByTupleIdSafeSlice(beginCellIds,endCellIds,stepCellIds);
}
}
return ret;
}
-std::vector<const BigMemoryObject *> MEDCouplingFieldDiscretizationPerCell::getDirectChildren() const
+std::vector<const BigMemoryObject *> MEDCouplingFieldDiscretizationPerCell::getDirectChildrenWithNull() const
{
- std::vector<const BigMemoryObject *> ret(MEDCouplingFieldDiscretization::getDirectChildren());
- if(_discr_per_cell)
- ret.push_back(_discr_per_cell);
+ std::vector<const BigMemoryObject *> ret(MEDCouplingFieldDiscretization::getDirectChildrenWithNull());
+ ret.push_back(_discr_per_cell);
return ret;
}
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationPerCell has no discretization per cell !");
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationPerCell::checkCoherencyBetween : NULL input mesh or DataArray !");
- int nbOfTuples=_discr_per_cell->getNumberOfTuples();
+ std::size_t nbOfTuples(_discr_per_cell->getNumberOfTuples());
if(nbOfTuples!=mesh->getNumberOfCells())
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationPerCell has a discretization per cell but it's not matching the underlying mesh !");
}
/*!
* This method is typically the first step of renumbering. The impact on _discr_per_cell is necessary here.
- * virtualy by this method.
+ * virtually by this method.
*/
void MEDCouplingFieldDiscretizationPerCell::renumberCells(const int *old2NewBg, bool check)
{
{
if(!_discr_per_cell)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationPerCell::checkNoOrphanCells : no discretization defined !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> test=_discr_per_cell->getIdsEqual(DFT_INVALID_LOCID_VALUE);
+ MCAuto<DataArrayInt> test=_discr_per_cell->findIdsEqual(DFT_INVALID_LOCID_VALUE);
if(test->getNumberOfTuples()!=0)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationPerCell::checkNoOrphanCells : presence of orphan cells !");
}
}
/*!
- * This method is simply called by MEDCouplingFieldDiscretization::deepCpy. It performs the deep copy of \a this.
+ * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this.
*
- * \sa MEDCouplingFieldDiscretization::deepCpy.
+ * \sa MEDCouplingFieldDiscretization::deepCopy.
*/
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationGauss::clone() const
{
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getNumberOfTuplesExpectedRegardingCode : invalid input code !");
int nbOfSplit=(int)idsPerType.size();
int nbOfTypes=(int)code.size()/3;
- int ret=0;
+ std::size_t ret(0);
for(int i=0;i<nbOfTypes;i++)
{
int nbOfEltInChunk=code[3*i+1];
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
const DataArrayInt *ids(idsPerType[pos]);
- if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
+ if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || (int)ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretizationGauss::getNumberOfTuplesExpectedRegardingCode : input pfl chunck at pos " << pos << " should have " << i << " tuples and one component and with ids all >=0 !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
/*!
* This method is redevelopped for performance reasons, but it is equivalent to a call to MEDCouplingFieldDiscretizationGauss::buildNbOfGaussPointPerCellField
- * and a call to DataArrayDouble::computeOffsets2 on the returned array.
+ * and a call to DataArrayDouble::computeOffsetsFull on the returned array.
*/
DataArrayInt *MEDCouplingFieldDiscretizationGauss::getOffsetArr(const MEDCouplingMesh *mesh) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getOffsetArr : NULL input mesh !");
- int nbOfTuples=mesh->getNumberOfCells();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
+ std::size_t nbOfTuples(mesh->getNumberOfCells());
+ MCAuto<DataArrayInt> ret=DataArrayInt::New();
ret->alloc(nbOfTuples+1,1);
- int *retPtr=ret->getPointer();
- const int *start=_discr_per_cell->getConstPointer();
+ int *retPtr(ret->getPointer());
+ const int *start(_discr_per_cell->begin());
if(_discr_per_cell->getNumberOfTuples()!=nbOfTuples)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getOffsetArr : mismatch between the mesh and the discretization ids array length !");
int maxPossible=(int)_loc.size();
retPtr[0]=0;
- for(int i=0;i<nbOfTuples;i++,start++)
+ for(std::size_t i=0;i<nbOfTuples;i++,start++)
{
if(*start>=0 && *start<maxPossible)
retPtr[i+1]=retPtr[i]+_loc[*start].getNumberOfGaussPt();
}
void MEDCouplingFieldDiscretizationGauss::renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector<DataArray *>& arrays,
- const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
+ const int *old2NewBg, bool check)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::renumberArraysForCell : NULL input mesh !");
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getLocalizationOfDiscValues : NULL input mesh !");
checkNoOrphanCells();
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh=mesh->buildUnstructured();//in general do nothing
+ MCAuto<MEDCouplingUMesh> umesh=mesh->buildUnstructured();//in general do nothing
int nbOfTuples=getNumberOfTuples(mesh);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> ret=DataArrayDouble::New();
int spaceDim=mesh->getSpaceDimension();
ret->alloc(nbOfTuples,spaceDim);
std::vector< int > locIds;
std::vector<DataArrayInt *> parts=splitIntoSingleGaussDicrPerCellType(locIds);
- std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > parts2(parts.size());
+ std::vector< MCAuto<DataArrayInt> > parts2(parts.size());
std::copy(parts.begin(),parts.end(),parts2.begin());
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> offsets=buildNbOfGaussPointPerCellField();
+ MCAuto<DataArrayInt> offsets=buildNbOfGaussPointPerCellField();
offsets->computeOffsets();
const int *ptrOffsets=offsets->getConstPointer();
const double *coords=umesh->getCoords()->getConstPointer();
{
INTERP_KERNEL::GaussCoords calculator;
//
- const MEDCouplingGaussLocalization& cli=_loc[locIds[i]];//curLocInfo
- INTERP_KERNEL::NormalizedCellType typ=cli.getType();
- const std::vector<double>& wg=cli.getWeights();
+ const MEDCouplingGaussLocalization& cli(_loc[locIds[i]]);//curLocInfo
+ INTERP_KERNEL::NormalizedCellType typ(cli.getType());
+ const std::vector<double>& wg(cli.getWeights());
calculator.addGaussInfo(typ,INTERP_KERNEL::CellModel::GetCellModel(typ).getDimension(),
- &cli.getGaussCoords()[0],(int)wg.size(),&cli.getRefCoords()[0],
- INTERP_KERNEL::CellModel::GetCellModel(typ).getNumberOfNodes());
+ &cli.getGaussCoords()[0],(int)wg.size(),&cli.getRefCoords()[0],
+ INTERP_KERNEL::CellModel::GetCellModel(typ).getNumberOfNodes());
//
- int nbt=parts2[i]->getNumberOfTuples();
- for(const int *w=parts2[i]->getConstPointer();w!=parts2[i]->getConstPointer()+nbt;w++)
+ for(const int *w=parts2[i]->begin();w!=parts2[i]->end();w++)
calculator.calculateCoords(cli.getType(),coords,spaceDim,conn+connI[*w]+1,valsToFill+spaceDim*(ptrOffsets[*w]));
}
ret->copyStringInfoFrom(*umesh->getCoords());
}
void MEDCouplingFieldDiscretizationGauss::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *tupleIdsBg, const int *tupleIdsEnd,
- DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const throw(INTERP_KERNEL::Exception)
+ DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeMeshRestrictionFromTupleIds : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp=DataArrayInt::New(); tmp->alloc((int)std::distance(tupleIdsBg,tupleIdsEnd),1);
+ MCAuto<DataArrayInt> tmp=DataArrayInt::New(); tmp->alloc((int)std::distance(tupleIdsBg,tupleIdsEnd),1);
std::copy(tupleIdsBg,tupleIdsEnd,tmp->getPointer());
tmp->sort(true);
tmp=tmp->buildUnique();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=buildNbOfGaussPointPerCellField();
- nbOfNodesPerCell->computeOffsets2();
- nbOfNodesPerCell->searchRangesInListOfIds(tmp,cellRestriction,trueTupleRestriction);
+ MCAuto<DataArrayInt> nbOfNodesPerCell=buildNbOfGaussPointPerCellField();
+ nbOfNodesPerCell->computeOffsetsFull();
+ nbOfNodesPerCell->findIdsRangesInListOfIds(tmp,cellRestriction,trueTupleRestriction);
}
/*!
else
_discr_per_cell=0;
arr=_discr_per_cell;
- int nbOfLoc=tinyInfo[1];
- _loc.clear();
- int dim=tinyInfo[2];
- int delta=-1;
- if(nbOfLoc>0)
- delta=((int)tinyInfo.size()-3)/nbOfLoc;
- for(int i=0;i<nbOfLoc;i++)
+ commonUnserialization(tinyInfo);
+}
+
+void MEDCouplingFieldDiscretizationGauss::checkForUnserialization(const std::vector<int>& tinyInfo, const DataArrayInt *arr)
+{
+ static const char MSG[]="MEDCouplingFieldDiscretizationGauss::checkForUnserialization : expect to have one not null DataArrayInt !";
+ int val=tinyInfo[0];
+ if(val>=0)
{
- std::vector<int> tmp(tinyInfo.begin()+3+i*delta,tinyInfo.begin()+3+(i+1)*delta);
- MEDCouplingGaussLocalization elt=MEDCouplingGaussLocalization::BuildNewInstanceFromTinyInfo(dim,tmp);
- _loc.push_back(elt);
+ if(!arr)
+ throw INTERP_KERNEL::Exception(MSG);
+ arr->checkNbOfTuplesAndComp(val,1,MSG);
+ _discr_per_cell=const_cast<DataArrayInt *>(arr);
+ _discr_per_cell->incrRef();
}
+ else
+ _discr_per_cell=0;
+ commonUnserialization(tinyInfo);
}
void MEDCouplingFieldDiscretizationGauss::finishUnserialization(const std::vector<double>& tinyInfo)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::checkCoherencyBetween : NULL input mesh or DataArray !");
MEDCouplingFieldDiscretizationPerCell::checkCoherencyBetween(mesh,da);
for(std::vector<MEDCouplingGaussLocalization>::const_iterator iter=_loc.begin();iter!=_loc.end();iter++)
- (*iter).checkCoherency();
+ (*iter).checkConsistencyLight();
int nbOfDesc=(int)_loc.size();
int nbOfCells=mesh->getNumberOfCells();
const int *dc=_discr_per_cell->getConstPointer();
{
if(dc[i]>=nbOfDesc)
{
- std::ostringstream oss; oss << "Cell # " << i << " of mesh \"" << mesh->getName() << "\" has an undefined gauss location ! Should never happend !";
+ std::ostringstream oss; oss << "Cell # " << i << " of mesh \"" << mesh->getName() << "\" has an undefined gauss location ! Should never happened !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
if(dc[i]<0)
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
}
- int nbOfTuples=getNumberOfTuples(mesh);
+ std::size_t nbOfTuples(getNumberOfTuples(mesh));
if(nbOfTuples!=da->getNumberOfTuples())
{
- std::ostringstream oss; oss << "Invalid number of tuples in the array : expecting " << nbOfTuples << " !";
+ std::ostringstream oss; oss << "Invalid number of tuples in the array : expecting " << nbOfTuples << " having " << da->getNumberOfTuples() << " !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
}
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : mesh instance specified is NULL !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isAbs);
+ MCAuto<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isAbs);
const double *volPtr=vol->getArray()->begin();
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=MEDCouplingFieldDouble::New(ON_GAUSS_PT);
+ MCAuto<MEDCouplingFieldDouble> ret=MEDCouplingFieldDouble::New(ON_GAUSS_PT);
ret->setMesh(mesh);
ret->setDiscretization(const_cast<MEDCouplingFieldDiscretizationGauss *>(this));
if(!_discr_per_cell)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : no discr per cell array defined but with nb of components different from 1 !");
if(_discr_per_cell->getNumberOfTuples()!=vol->getNumberOfTuples())
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::getMeasureField : no discr per cell array defined but mismatch between nb of cells of mesh and size of spatial disr array !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> offset=getOffsetArr(mesh);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=DataArrayDouble::New(); arr->alloc(getNumberOfTuples(mesh),1);
+ MCAuto<DataArrayInt> offset=getOffsetArr(mesh);
+ MCAuto<DataArrayDouble> arr=DataArrayDouble::New(); arr->alloc(getNumberOfTuples(mesh),1);
ret->setArray(arr);
double *arrPtr=arr->getPointer();
const int *offsetPtr=offset->getConstPointer();
int maxGaussLoc=(int)_loc.size();
std::vector<int> locIds;
std::vector<DataArrayInt *> ids=splitIntoSingleGaussDicrPerCellType(locIds);
- std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > ids2(ids.size()); std::copy(ids.begin(),ids.end(),ids2.begin());
+ std::vector< MCAuto<DataArrayInt> > ids2(ids.size()); std::copy(ids.begin(),ids.end(),ids2.begin());
for(std::size_t i=0;i<locIds.size();i++)
{
const DataArrayInt *curIds=ids[i];
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::buildSubMeshData : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> diSafe=computeTupleIdsToSelectFromCellIds(mesh,start,end);
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPart(start,end);
+ MCAuto<DataArrayInt> diSafe=computeTupleIdsToSelectFromCellIds(mesh,start,end);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPart(start,end);
di=diSafe.retn();
return ret.retn();
}
else
{ std::ostringstream oss; oss << msg << i << " is detected as orphan !"; throw INTERP_KERNEL::Exception(oss.str().c_str()); }
}
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
return ret.retn();
}
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds : null mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=buildNbOfGaussPointPerCellField();//check of _discr_per_cell not NULL pointer
- int nbOfCells=mesh->getNumberOfCells();
+ MCAuto<DataArrayInt> nbOfNodesPerCell=buildNbOfGaussPointPerCellField();//check of _discr_per_cell not NULL pointer
+ std::size_t nbOfCells(mesh->getNumberOfCells());
if(_discr_per_cell->getNumberOfTuples()!=nbOfCells)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::computeTupleIdsToSelectFromCellIds : mismatch of nb of tuples of cell ids array and number of cells !");
- nbOfNodesPerCell->computeOffsets2();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1);
+ nbOfNodesPerCell->computeOffsetsFull();
+ MCAuto<DataArrayInt> sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1);
return sel->buildExplicitArrByRanges(nbOfNodesPerCell);
}
}
void MEDCouplingFieldDiscretizationGauss::setGaussLocalizationOnType(const MEDCouplingMesh *mesh, INTERP_KERNEL::NormalizedCellType type, const std::vector<double>& refCoo,
- const std::vector<double>& gsCoo, const std::vector<double>& wg) throw(INTERP_KERNEL::Exception)
+ const std::vector<double>& gsCoo, const std::vector<double>& wg)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::setGaussLocalizationOnType : NULL input mesh !");
}
void MEDCouplingFieldDiscretizationGauss::setGaussLocalizationOnCells(const MEDCouplingMesh *mesh, const int *begin, const int *end, const std::vector<double>& refCoo,
- const std::vector<double>& gsCoo, const std::vector<double>& wg) throw(INTERP_KERNEL::Exception)
+ const std::vector<double>& gsCoo, const std::vector<double>& wg)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::setGaussLocalizationOnCells : NULL input mesh !");
if(!_discr_per_cell)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGauss::buildNbOfGaussPointPerCellField : no discretization array set !");
int nbOfTuples=_discr_per_cell->getNumberOfTuples();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
+ MCAuto<DataArrayInt> ret=DataArrayInt::New();
const int *w=_discr_per_cell->begin();
ret->alloc(nbOfTuples,1);
int *valsToFill=ret->getPointer();
_loc=tmpLoc;
}
+void MEDCouplingFieldDiscretizationGauss::commonUnserialization(const std::vector<int>& tinyInfo)
+{
+ int nbOfLoc=tinyInfo[1];
+ _loc.clear();
+ int dim=tinyInfo[2];
+ int delta=-1;
+ if(nbOfLoc>0)
+ delta=((int)tinyInfo.size()-3)/nbOfLoc;
+ for(int i=0;i<nbOfLoc;i++)
+ {
+ std::vector<int> tmp(tinyInfo.begin()+3+i*delta,tinyInfo.begin()+3+(i+1)*delta);
+ MEDCouplingGaussLocalization elt=MEDCouplingGaussLocalization::BuildNewInstanceFromTinyInfo(dim,tmp);
+ _loc.push_back(elt);
+ }
+}
+
MEDCouplingFieldDiscretizationGaussNE::MEDCouplingFieldDiscretizationGaussNE()
{
}
}
/*!
- * This method is simply called by MEDCouplingFieldDiscretization::deepCpy. It performs the deep copy of \a this.
+ * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this.
*
- * \sa MEDCouplingFieldDiscretization::deepCpy.
+ * \sa MEDCouplingFieldDiscretization::deepCopy.
*/
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationGaussNE::clone() const
{
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
const DataArrayInt *ids(idsPerType[pos]);
- if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
+ if(!ids || !ids->isAllocated() || ids->getNumberOfComponents()!=1 || (int)ids->getNumberOfTuples()!=nbOfEltInChunk || ids->getMinValueInArray()<0)
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretizationGaussNE::getNumberOfTuplesExpectedRegardingCode : input pfl chunck at pos " << pos << " should have " << i << " tuples and one component and with ids all >=0 !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
void MEDCouplingFieldDiscretizationGaussNE::renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector<DataArray *>& arrays,
- const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
+ const int *old2NewBg, bool check)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::renumberArraysForCell : NULL input mesh !");
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::getLocalizationOfDiscValues : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> umesh=mesh->buildUnstructured();//in general do nothing
+ MCAuto<DataArrayDouble> ret=DataArrayDouble::New();
+ MCAuto<MEDCouplingUMesh> umesh=mesh->buildUnstructured();//in general do nothing
int nbOfTuples=getNumberOfTuples(umesh);
int spaceDim=mesh->getSpaceDimension();
ret->alloc(nbOfTuples,spaceDim);
int nbOfCompo=arr->getNumberOfComponents();
std::fill(res,res+nbOfCompo,0.);
//
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isWAbs);
+ MCAuto<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isWAbs);
std::set<INTERP_KERNEL::NormalizedCellType> types=mesh->getAllGeoTypes();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
- nbOfNodesPerCell->computeOffsets2();
+ MCAuto<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
+ nbOfNodesPerCell->computeOffsetsFull();
const double *arrPtr=arr->begin(),*volPtr=vol->getArray()->begin();
for(std::set<INTERP_KERNEL::NormalizedCellType>::const_iterator it=types.begin();it!=types.end();it++)
{
INTERP_KERNEL::AutoPtr<double> wArr2=new double[wArrSz];
double sum=std::accumulate(wArr,wArr+wArrSz,0.);
std::transform(wArr,wArr+wArrSz,(double *)wArr2,std::bind2nd(std::multiplies<double>(),1./sum));
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids=mesh->giveCellsWithType(*it);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids2=ids->buildExplicitArrByRanges(nbOfNodesPerCell);
+ MCAuto<DataArrayInt> ids=mesh->giveCellsWithType(*it);
+ MCAuto<DataArrayInt> ids2=ids->buildExplicitArrByRanges(nbOfNodesPerCell);
const int *ptIds2=ids2->begin(),*ptIds=ids->begin();
int nbOfCellsWithCurGeoType=ids->getNumberOfTuples();
for(int i=0;i<nbOfCellsWithCurGeoType;i++,ptIds++,ptIds2+=wArrSz)
const double *MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(INTERP_KERNEL::NormalizedCellType geoType, std::size_t& lgth)
{
switch(geoType)
- {
+ {
case INTERP_KERNEL::NORM_POINT1:
lgth=(int)sizeof(FGP_POINT1)/sizeof(double);
return FGP_POINT1;
return FGP_PYRA13;
default:
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType : only SEG[2,3,4], TRI[3,6,7], QUAD[4,9], TETRA[4,10], PENTA[6,15], HEXA[8,20,27], PYRA[5,13] supported !");
- }
+ }
}
const double *MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(INTERP_KERNEL::NormalizedCellType geoType, std::size_t& lgth)
{
switch(geoType)
- {
+ {
case INTERP_KERNEL::NORM_POINT1:
lgth=0;
return 0;
return REF_PYRA13;
default:
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType : only SEG[2,3,4], TRI[3,6,7], QUAD[4,8,9], TETRA[4,10], PENTA[6,15], HEXA[8,20,27], PYRA[5,13] supported !");
- }
+ }
}
const double *MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(INTERP_KERNEL::NormalizedCellType geoType, std::size_t& lgth)
{
switch(geoType)
- {
+ {
case INTERP_KERNEL::NORM_POINT1:
{
lgth=0;
}
default:
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType : only SEG[2,3,4], TRI[3,6,7], QUAD[4,8,9], TETRA[4,10], PENTA[6,15], HEXA[8,20,27], PYRA[5,13] supported !");
- }
+ }
}
void MEDCouplingFieldDiscretizationGaussNE::computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const int *tupleIdsBg, const int *tupleIdsEnd,
- DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const throw(INTERP_KERNEL::Exception)
+ DataArrayInt *&cellRestriction, DataArrayInt *&trueTupleRestriction) const
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::computeMeshRestrictionFromTupleIds : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp=DataArrayInt::New(); tmp->alloc((int)std::distance(tupleIdsBg,tupleIdsEnd),1);
+ MCAuto<DataArrayInt> tmp=DataArrayInt::New(); tmp->alloc((int)std::distance(tupleIdsBg,tupleIdsEnd),1);
std::copy(tupleIdsBg,tupleIdsEnd,tmp->getPointer());
tmp->sort(true);
tmp=tmp->buildUnique();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
- nbOfNodesPerCell->computeOffsets2();
- nbOfNodesPerCell->searchRangesInListOfIds(tmp,cellRestriction,trueTupleRestriction);
+ MCAuto<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
+ nbOfNodesPerCell->computeOffsetsFull();
+ nbOfNodesPerCell->findIdsRangesInListOfIds(tmp,cellRestriction,trueTupleRestriction);
}
void MEDCouplingFieldDiscretizationGaussNE::checkCompatibilityWithNature(NatureOfField nat) const
void MEDCouplingFieldDiscretizationGaussNE::checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArray *da) const
{
- int nbOfTuples=getNumberOfTuples(mesh);
+ std::size_t nbOfTuples(getNumberOfTuples(mesh));
if(nbOfTuples!=da->getNumberOfTuples())
{
std::ostringstream oss; oss << "Invalid number of tuples in the array : expecting " << nbOfTuples << " !";
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::getMeasureField : mesh instance specified is NULL !");
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isAbs);
+ MCAuto<MEDCouplingFieldDouble> vol=mesh->getMeasureField(isAbs);
const double *volPtr=vol->getArray()->begin();
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=MEDCouplingFieldDouble::New(ON_GAUSS_NE);
+ MCAuto<MEDCouplingFieldDouble> ret=MEDCouplingFieldDouble::New(ON_GAUSS_NE);
ret->setMesh(mesh);
//
std::set<INTERP_KERNEL::NormalizedCellType> types=mesh->getAllGeoTypes();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
+ MCAuto<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
int nbTuples=nbOfNodesPerCell->accumulate(0);
- nbOfNodesPerCell->computeOffsets2();
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr=DataArrayDouble::New(); arr->alloc(nbTuples,1);
+ nbOfNodesPerCell->computeOffsetsFull();
+ MCAuto<DataArrayDouble> arr=DataArrayDouble::New(); arr->alloc(nbTuples,1);
ret->setArray(arr);
double *arrPtr=arr->getPointer();
for(std::set<INTERP_KERNEL::NormalizedCellType>::const_iterator it=types.begin();it!=types.end();it++)
INTERP_KERNEL::AutoPtr<double> wArr2=new double[wArrSz];
double sum=std::accumulate(wArr,wArr+wArrSz,0.);
std::transform(wArr,wArr+wArrSz,(double *)wArr2,std::bind2nd(std::multiplies<double>(),1./sum));
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids=mesh->giveCellsWithType(*it);
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids2=ids->buildExplicitArrByRanges(nbOfNodesPerCell);
+ MCAuto<DataArrayInt> ids=mesh->giveCellsWithType(*it);
+ MCAuto<DataArrayInt> ids2=ids->buildExplicitArrByRanges(nbOfNodesPerCell);
const int *ptIds2=ids2->begin(),*ptIds=ids->begin();
int nbOfCellsWithCurGeoType=ids->getNumberOfTuples();
for(int i=0;i<nbOfCellsWithCurGeoType;i++,ptIds++)
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::buildSubMeshData : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> diSafe=computeTupleIdsToSelectFromCellIds(mesh,start,end);
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPart(start,end);
+ MCAuto<DataArrayInt> diSafe=computeTupleIdsToSelectFromCellIds(mesh,start,end);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPart(start,end);
di=diSafe.retn();
return ret.retn();
}
if(i>=endCellIds)
break;
}
- MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
+ MCAuto<MEDCouplingMesh> ret=mesh->buildPartRange(beginCellIds,endCellIds,stepCellIds);
return ret.retn();
}
{
if(!mesh)
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationGaussNE::computeTupleIdsToSelectFromCellIds : null mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
- nbOfNodesPerCell->computeOffsets2();
- MEDCouplingAutoRefCountObjectPtr<DataArrayInt> sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1);
+ MCAuto<DataArrayInt> nbOfNodesPerCell=mesh->computeNbOfNodesPerCell();
+ nbOfNodesPerCell->computeOffsetsFull();
+ MCAuto<DataArrayInt> sel=DataArrayInt::New(); sel->useArray(startCellIds,false,CPP_DEALLOC,(int)std::distance(startCellIds,endCellIds),1);
return sel->buildExplicitArrByRanges(nbOfNodesPerCell);
}
}
/*!
- * This method is simply called by MEDCouplingFieldDiscretization::deepCpy. It performs the deep copy of \a this.
+ * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this.
*
- * \sa MEDCouplingFieldDiscretization::deepCpy.
+ * \sa MEDCouplingFieldDiscretization::deepCopy.
*/
MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationKriging::clone() const
{
void MEDCouplingFieldDiscretizationKriging::checkCompatibilityWithNature(NatureOfField nat) const
{
- if(nat!=ConservativeVolumic)
- throw INTERP_KERNEL::Exception("Invalid nature for Kriging field : expected ConservativeVolumic !");
+ if(nat!=IntensiveMaximum)
+ throw INTERP_KERNEL::Exception("Invalid nature for Kriging field : expected IntensiveMaximum !");
}
bool MEDCouplingFieldDiscretizationKriging::isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const
void MEDCouplingFieldDiscretizationKriging::getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const
{
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> res2=MEDCouplingFieldDiscretizationKriging::getValueOnMulti(arr,mesh,loc,1);
+ MCAuto<DataArrayDouble> res2=MEDCouplingFieldDiscretizationKriging::getValueOnMulti(arr,mesh,loc,1);
std::copy(res2->begin(),res2->end(),res);
}
{
if(!arr || !arr->isAllocated())
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::getValueOnMulti : input array is null or not allocated !");
- int nbOfRows(getNumberOfMeshPlaces(mesh));
+ std::size_t nbOfRows(getNumberOfMeshPlaces(mesh));
if(arr->getNumberOfTuples()!=nbOfRows)
{
std::ostringstream oss; oss << "MEDCouplingFieldDiscretizationKriging::getValueOnMulti : input array does not have correct number of tuples ! Excepted " << nbOfRows << " having " << arr->getNumberOfTuples() << " !";
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
int nbCols(-1),nbCompo(arr->getNumberOfComponents());
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> m(computeEvaluationMatrixOnGivenPts(mesh,loc,nbOfTargetPoints,nbCols));
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
+ MCAuto<DataArrayDouble> m(computeEvaluationMatrixOnGivenPts(mesh,loc,nbOfTargetPoints,nbCols));
+ MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
ret->alloc(nbOfTargetPoints,nbCompo);
INTERP_KERNEL::matrixProduct(m->begin(),nbOfTargetPoints,nbCols,arr->begin(),nbOfRows,nbCompo,ret->getPointer());
return ret.retn();
DataArrayDouble *MEDCouplingFieldDiscretizationKriging::computeEvaluationMatrixOnGivenPts(const MEDCouplingMesh *mesh, const double *loc, int nbOfTargetPoints, int& nbCols) const
{
int isDrift(-1),nbRows(-1);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrixInv(computeInverseMatrix(mesh,isDrift,nbRows));
+ MCAuto<DataArrayDouble> matrixInv(computeInverseMatrix(mesh,isDrift,nbRows));
//
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=getLocalizationOfDiscValues(mesh);
+ MCAuto<DataArrayDouble> coords=getLocalizationOfDiscValues(mesh);
int nbOfPts(coords->getNumberOfTuples()),dimension(coords->getNumberOfComponents());
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> locArr=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> locArr=DataArrayDouble::New();
locArr->useArray(loc,false,CPP_DEALLOC,nbOfTargetPoints,dimension);
nbCols=nbOfPts;
//
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrix2=coords->buildEuclidianDistanceDenseMatrixWith(locArr);
+ MCAuto<DataArrayDouble> matrix2=coords->buildEuclidianDistanceDenseMatrixWith(locArr);
operateOnDenseMatrix(mesh->getSpaceDimension(),nbOfTargetPoints*nbOfPts,matrix2->getPointer());
//
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrix3=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> matrix3=DataArrayDouble::New();
matrix3->alloc(nbOfTargetPoints*nbRows,1);
double *work=matrix3->getPointer();
const double *workCst(matrix2->begin()),*workCst2(loc);
for(int j=0;j<isDrift-1;j++)
work[i*nbRows+(nbOfPts+1+j)]=workCst2[j];
}
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
+ MCAuto<DataArrayDouble> ret(DataArrayDouble::New());
ret->alloc(nbOfTargetPoints,nbRows);
INTERP_KERNEL::matrixProduct(matrix3->begin(),nbOfTargetPoints,nbRows,matrixInv->begin(),nbRows,nbRows,ret->getPointer());
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret2(DataArrayDouble::New());
+ MCAuto<DataArrayDouble> ret2(DataArrayDouble::New());
ret2->alloc(nbOfTargetPoints*nbOfPts,1);
workCst=ret->begin(); work=ret2->getPointer();
for(int i=0;i<nbOfTargetPoints;i++,workCst+=nbRows)
* \param [out] isDrift return if drift coefficients are present in the returned vector of coefficients. If different from 0 there is presence of drift coefficients.
* \param [out] matSz the size of returned square matrix
* \return the new result matrix to be deallocated by the caller.
+ * \sa computeMatrix
*/
DataArrayDouble *MEDCouplingFieldDiscretizationKriging::computeInverseMatrix(const MEDCouplingMesh *mesh, int& isDrift, int& matSz) const
{
- if(!mesh)
- throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::computeVectorOfCoefficients : NULL input mesh !");
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=getLocalizationOfDiscValues(mesh);
- int nbOfPts=coords->getNumberOfTuples();
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrix=coords->buildEuclidianDistanceDenseMatrix();
- operateOnDenseMatrix(mesh->getSpaceDimension(),nbOfPts*nbOfPts,matrix->getPointer());
- // Drift
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrixWithDrift=performDrift(matrix,coords,isDrift);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrixInv=DataArrayDouble::New();
- matSz=nbOfPts+isDrift;
+ MCAuto<DataArrayDouble> matrixWithDrift(computeMatrix(mesh,isDrift,matSz));
+ MCAuto<DataArrayDouble> matrixInv(DataArrayDouble::New());
matrixInv->alloc(matSz*matSz,1);
INTERP_KERNEL::inverseMatrix(matrixWithDrift->getConstPointer(),matSz,matrixInv->getPointer());
return matrixInv.retn();
}
+/*!
+ * This method computes the kriging matrix.
+ * \return the new result matrix to be deallocated by the caller.
+ * \sa computeInverseMatrix
+ */
+DataArrayDouble *MEDCouplingFieldDiscretizationKriging::computeMatrix(const MEDCouplingMesh *mesh, int& isDrift, int& matSz) const
+{
+ if(!mesh)
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::computeMatrix : NULL input mesh !");
+ MCAuto<DataArrayDouble> coords(getLocalizationOfDiscValues(mesh));
+ int nbOfPts(coords->getNumberOfTuples());
+ MCAuto<DataArrayDouble> matrix(coords->buildEuclidianDistanceDenseMatrix());
+ operateOnDenseMatrix(mesh->getSpaceDimension(),nbOfPts*nbOfPts,matrix->getPointer());
+ // Drift
+ MCAuto<DataArrayDouble> matrixWithDrift(performDrift(matrix,coords,isDrift));
+ matSz=nbOfPts+isDrift;
+ return matrixWithDrift.retn();
+}
+
/*!
* This method computes coefficients to apply to each representing points of \a mesh, that is to say the nodes of \a mesh given a field array \a arr whose
* number of tuples should be equal to the number of representing points in \a mesh.
DataArrayDouble *MEDCouplingFieldDiscretizationKriging::computeVectorOfCoefficients(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, int& isDrift) const
{
int nbRows(-1);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> matrixInv(computeInverseMatrix(mesh,isDrift,nbRows));
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> KnewiK=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> matrixInv(computeInverseMatrix(mesh,isDrift,nbRows));
+ MCAuto<DataArrayDouble> KnewiK(DataArrayDouble::New());
KnewiK->alloc(nbRows*1,1);
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr2=DataArrayDouble::New();
- arr2->alloc(nbRows*1,1);
- double *work=std::copy(arr->begin(),arr->end(),arr2->getPointer());
- std::fill(work,work+isDrift,0.);
- INTERP_KERNEL::matrixProduct(matrixInv->getConstPointer(),nbRows,nbRows,arr2->getConstPointer(),nbRows,1,KnewiK->getPointer());
+ MCAuto<DataArrayDouble> arr2(PerformDriftOfVec(arr,isDrift));
+ INTERP_KERNEL::matrixProduct(matrixInv->getConstPointer(),nbRows,nbRows,arr2->getConstPointer(),arr2->getNumberOfTuples(),1,KnewiK->getPointer());
return KnewiK.retn();
}
void MEDCouplingFieldDiscretizationKriging::operateOnDenseMatrix(int spaceDimension, int nbOfElems, double *matrixPtr) const
{
switch(spaceDimension)
- {
+ {
case 1:
{
- for(int i=0;i<nbOfElems;i++)
- {
- double val=matrixPtr[i];
- matrixPtr[i]=val*val*val;
- }
+ OperateOnDenseMatrixH3(nbOfElems,matrixPtr);
break;
}
case 2:
{
- for(int i=0;i<nbOfElems;i++)
- {
- double val=matrixPtr[i];
- if(val!=0.)
- matrixPtr[i]=val*val*log(val);
- }
+ OperateOnDenseMatrixH2Ln(nbOfElems,matrixPtr);
break;
}
case 3:
}
default:
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::operateOnDenseMatrix : only dimension 1, 2 and 3 implemented !");
+ }
+}
+
+void MEDCouplingFieldDiscretizationKriging::OperateOnDenseMatrixH3(int nbOfElems, double *matrixPtr)
+{
+ for(int i=0;i<nbOfElems;i++)
+ {
+ double val=matrixPtr[i];
+ matrixPtr[i]=val*val*val;
}
}
+void MEDCouplingFieldDiscretizationKriging::OperateOnDenseMatrixH2Ln(int nbOfElems, double *matrixPtr)
+{
+ for(int i=0;i<nbOfElems;i++)
+ {
+ double val=matrixPtr[i];
+ if(val!=0.)
+ matrixPtr[i]=val*val*log(val);
+ }
+}
+
+/*!
+ * Performs a drift to the rectangular input matrix \a matr.
+ * This method generate a dense matrix starting from an input dense matrix \a matr and input array \a arr.
+ * \param [in] matr The rectangular dense matrix (with only one component). The number of rows of \a matr must be equal to the number of tuples of \a arr
+ * \param [in] arr The array of coords to be appended in the input dense matrix \a matr. Typically arr is an array of coordinates.
+ * \param [out] delta the delta of number of columns between returned dense matrix and input dense matrix \a matr. \a delta is equal to number of components of \a arr + 1.
+ * \sa performDrift
+ */
+DataArrayDouble *MEDCouplingFieldDiscretizationKriging::PerformDriftRect(const DataArrayDouble *matr, const DataArrayDouble *arr, int& delta)
+{
+ if(!matr || !matr->isAllocated() || matr->getNumberOfComponents()!=1)
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::PerformDriftRect : invalid input dense matrix ! Must be allocated not NULL and with exactly one component !");
+ if(!arr || !arr->isAllocated())
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::PerformDriftRect : invalid input array of coordiantes ! Must be allocated and not NULL !");
+ int spaceDimension(arr->getNumberOfComponents()),nbOfPts(arr->getNumberOfTuples()),nbOfEltInMatrx(matr->getNumberOfTuples());
+ delta=spaceDimension+1;
+ int nbOfCols(nbOfEltInMatrx/nbOfPts);
+ if(nbOfEltInMatrx%nbOfPts!=0)
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::PerformDriftRect : size of input dense matrix and input arrays mismatch ! NbOfElems in matrix % nb of tuples in array must be equal to 0 !");
+ MCAuto<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfPts*(nbOfCols+delta));
+ double *retPtr(ret->getPointer());
+ const double *mPtr(matr->begin()),*aPtr(arr->begin());
+ for(int i=0;i<nbOfPts;i++,aPtr+=spaceDimension,mPtr+=nbOfCols)
+ {
+ retPtr=std::copy(mPtr,mPtr+nbOfCols,retPtr);
+ *retPtr++=1.;
+ retPtr=std::copy(aPtr,aPtr+spaceDimension,retPtr);
+ }
+ return ret.retn();
+}
+
+/*!
+ * \return a newly allocated array having \a isDrift more tuples than \a arr.
+ * \sa computeVectorOfCoefficients
+ */
+DataArrayDouble *MEDCouplingFieldDiscretizationKriging::PerformDriftOfVec(const DataArrayDouble *arr, int isDrift)
+{
+ if(!arr || !arr->isAllocated() || arr->getNumberOfComponents()!=1)
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::PerformDriftOfVec : input array must be not NULL allocated and with one component !");
+ if(isDrift<0)
+ throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::PerformDriftOfVec : isDrift parameter must be >=0 !");
+ MCAuto<DataArrayDouble> arr2(DataArrayDouble::New());
+ arr2->alloc((arr->getNumberOfTuples()+isDrift)*1,1);
+ double *work(std::copy(arr->begin(),arr->end(),arr2->getPointer()));
+ std::fill(work,work+isDrift,0.);
+ return arr2.retn();
+}
+
/*!
* Starting from a square matrix \a matr, this method returns a newly allocated dense square matrix whose \a matr is included in returned matrix
* in the top left corner, and in the remaining returned matrix the parameters to take into account about the kriging drift.
* \param [in] matr input matrix whose drift part will be added
* \param [out] delta the difference between the size of the output matrix and the input matrix \a matr.
* \return a newly allocated matrix bigger than input matrix \a matr.
+ * \sa MEDCouplingFieldDiscretizationKriging::PerformDriftRect
*/
DataArrayDouble *MEDCouplingFieldDiscretizationKriging::performDrift(const DataArrayDouble *matr, const DataArrayDouble *arr, int& delta) const
{
- int spaceDimension=arr->getNumberOfComponents();
+ std::size_t spaceDimension(arr->getNumberOfComponents());
delta=spaceDimension+1;
- int szOfMatrix=arr->getNumberOfTuples();
+ std::size_t szOfMatrix(arr->getNumberOfTuples());
if(szOfMatrix*szOfMatrix!=matr->getNumberOfTuples())
throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationKriging::performDrift : invalid size");
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
+ MCAuto<DataArrayDouble> ret=DataArrayDouble::New();
ret->alloc((szOfMatrix+delta)*(szOfMatrix+delta),1);
const double *srcWork=matr->getConstPointer();
const double *srcWork2=arr->getConstPointer();
double *destWork=ret->getPointer();
- for(int i=0;i<szOfMatrix;i++)
+ for(std::size_t i=0;i<szOfMatrix;i++)
{
destWork=std::copy(srcWork,srcWork+szOfMatrix,destWork);
srcWork+=szOfMatrix;
}
std::fill(destWork,destWork+szOfMatrix,1.); destWork+=szOfMatrix;
std::fill(destWork,destWork+spaceDimension+1,0.); destWork+=spaceDimension+1;
- MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arrNoI=arr->toNoInterlace();
+ MCAuto<DataArrayDouble> arrNoI=arr->toNoInterlace();
srcWork2=arrNoI->getConstPointer();
- for(int i=0;i<spaceDimension;i++)
+ for(std::size_t i=0;i<spaceDimension;i++)
{
destWork=std::copy(srcWork2,srcWork2+szOfMatrix,destWork);
srcWork2+=szOfMatrix;