1 // Copyright (C) 2007-2014 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay (CEA/DEN)
21 #include "MEDCouplingIMesh.hxx"
22 #include "MEDCouplingCMesh.hxx"
23 #include "MEDCouplingMemArray.hxx"
24 #include "MEDCouplingFieldDouble.hxx"
31 using namespace ParaMEDMEM;
33 MEDCouplingIMesh::MEDCouplingIMesh():_space_dim(-1)
35 _origin[0]=0.; _origin[1]=0.; _origin[2]=0.;
36 _dxyz[0]=0.; _dxyz[1]=0.; _dxyz[2]=0.;
37 _structure[0]=0; _structure[1]=0; _structure[2]=0;
40 MEDCouplingIMesh::MEDCouplingIMesh(const MEDCouplingIMesh& other, bool deepCopy):MEDCouplingStructuredMesh(other,deepCopy),_space_dim(other._space_dim),_axis_unit(other._axis_unit)
42 _origin[0]=other._origin[0]; _origin[1]=other._origin[1]; _origin[2]=other._origin[2];
43 _dxyz[0]=other._dxyz[0]; _dxyz[1]=other._dxyz[1]; _dxyz[2]=other._dxyz[2];
44 _structure[0]=other._structure[0]; _structure[1]=other._structure[1]; _structure[2]=other._structure[2];
47 MEDCouplingIMesh::~MEDCouplingIMesh()
51 MEDCouplingIMesh *MEDCouplingIMesh::New()
53 return new MEDCouplingIMesh;
56 MEDCouplingIMesh *MEDCouplingIMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
57 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
59 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(new MEDCouplingIMesh);
60 ret->setName(meshName);
61 ret->setSpaceDimension(spaceDim);
62 ret->setNodeStruct(nodeStrctStart,nodeStrctStop);
63 ret->setOrigin(originStart,originStop);
64 ret->setDXYZ(dxyzStart,dxyzStop);
68 MEDCouplingMesh *MEDCouplingIMesh::deepCpy() const
73 MEDCouplingIMesh *MEDCouplingIMesh::clone(bool recDeepCpy) const
75 return new MEDCouplingIMesh(*this,recDeepCpy);
78 void MEDCouplingIMesh::setNodeStruct(const int *nodeStrctStart, const int *nodeStrctStop)
80 checkSpaceDimension();
81 int sz((int)std::distance(nodeStrctStart,nodeStrctStop));
83 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setNodeStruct : input vector of node structure has not the right size ! Or change space dimension before calling it !");
84 std::copy(nodeStrctStart,nodeStrctStop,_structure);
88 std::vector<int> MEDCouplingIMesh::getNodeStruct() const
90 checkSpaceDimension();
91 return std::vector<int>(_structure,_structure+_space_dim);
94 void MEDCouplingIMesh::setOrigin(const double *originStart, const double *originStop)
96 checkSpaceDimension();
97 int sz((int)std::distance(originStart,originStop));
99 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setOrigin : input vector of origin vector has not the right size ! Or change space dimension before calling it !");
100 std::copy(originStart,originStop,_origin);
104 std::vector<double> MEDCouplingIMesh::getOrigin() const
106 checkSpaceDimension();
107 return std::vector<double>(_origin,_origin+_space_dim);
110 void MEDCouplingIMesh::setDXYZ(const double *dxyzStart, const double *dxyzStop)
112 checkSpaceDimension();
113 int sz((int)std::distance(dxyzStart,dxyzStop));
115 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setDXYZ : input vector of dxyz vector has not the right size ! Or change space dimension before calling it !");
116 std::copy(dxyzStart,dxyzStop,_dxyz);
120 std::vector<double> MEDCouplingIMesh::getDXYZ() const
122 checkSpaceDimension();
123 return std::vector<double>(_dxyz,_dxyz+_space_dim);
126 void MEDCouplingIMesh::setAxisUnit(const std::string& unitName)
132 std::string MEDCouplingIMesh::getAxisUnit() const
138 * This method returns the measure of any cell in \a this.
139 * This specific method of image grid mesh utilizes the fact that any cell in \a this have the same measure.
140 * The value returned by this method is those used to feed the returned field in the MEDCouplingIMesh::getMeasureField.
142 * \sa getMeasureField
144 double MEDCouplingIMesh::getMeasureOfAnyCell() const
147 int dim(getSpaceDimension());
149 for(int i=0;i<dim;i++)
155 * This method is allows to convert \a this into MEDCouplingCMesh instance.
156 * This method is the middle level between MEDCouplingIMesh and the most general MEDCouplingUMesh.
157 * This method is useful for MED writers that do not have still the image grid support.
159 * \sa MEDCouplingMesh::buildUnstructured
161 MEDCouplingCMesh *MEDCouplingIMesh::convertToCartesian() const
164 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCMesh> ret(MEDCouplingCMesh::New());
166 { ret->copyTinyInfoFrom(this); }
167 catch(INTERP_KERNEL::Exception& ) { }
168 int spaceDim(getSpaceDimension());
169 std::vector<std::string> infos(buildInfoOnComponents());
170 for(int i=0;i<spaceDim;i++)
172 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(_structure[i],1); arr->setInfoOnComponent(0,infos[i]);
173 arr->iota(); arr->applyLin(_dxyz[i],_origin[i]);
174 ret->setCoordsAt(i,arr);
180 * This method refines \a this uniformaly along all of its dimensions. In case of success the space covered by \a this will remain
181 * the same before the invocation except that the number of cells will be multiplied by \a factor ^ this->getMeshDimension().
182 * The origin of \a this will be not touched only spacing and node structure will be changed.
183 * This method can be useful for AMR users.
185 void MEDCouplingIMesh::refineWithFactor(int factor)
188 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::refineWithFactor : refinement factor must be != 0 !");
190 int factAbs(std::abs(factor));
191 double fact2(1./(double)factor);
192 std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::plus<int>(),-1));
193 std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::multiplies<int>(),factAbs));
194 std::transform(_structure,_structure+_space_dim,_structure,std::bind2nd(std::plus<int>(),1));
195 std::transform(_dxyz,_dxyz+_space_dim,_dxyz,std::bind2nd(std::multiplies<double>(),fact2));
200 * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
201 * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
202 * to a coarse image mesh. Only tuples ( deduced from \a fineLocInCoarse ) of \a coarseDA will be modified. Other tuples of \a coarseDA will be let unchanged.
204 * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
205 * \param [in] coarseSt The cell structure of coarse mesh.
206 * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
207 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
208 * \param [in] facts The refinement coefficient per axis.
209 * \sa SpreadCoarseToFine
211 void MEDCouplingIMesh::CondenseFineToCoarse(DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts)
213 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
214 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the parameters 1 or 3 are NULL or not allocated !");
215 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
216 int nbCompo(fineDA->getNumberOfComponents());
217 if(coarseDA->getNumberOfComponents()!=nbCompo)
218 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the number of components of fine DA and coarse one mismatches !");
219 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
220 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
221 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
223 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
224 throw INTERP_KERNEL::Exception(oss.str().c_str());
226 int nbTuplesFine(fineDA->getNumberOfTuples());
227 if(nbTuplesFine%nbOfTuplesInFineExp!=0)
228 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Invalid nb of tuples in fine DataArray regarding its structure !");
229 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
230 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
232 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
233 throw INTERP_KERNEL::Exception(oss.str().c_str());
235 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
236 double *outPtr(coarseDA->getPointer());
237 const double *inPtr(fineDA->begin());
239 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
244 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
245 for(int i=0;i<dims[0];i++)
247 double *loc(outPtr+(offset+i)*nbCompo);
248 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
251 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
253 std::copy(inPtr,inPtr+nbCompo,loc);
260 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact1(facts[1]),fact0(facts[0]);
261 for(int j=0;j<dims[1];j++)
263 for(int jfact=0;jfact<fact1;jfact++)
265 for(int i=0;i<dims[0];i++)
267 double *loc(outPtr+(kk+i)*nbCompo);
268 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
270 if(jfact!=0 || ifact!=0)
271 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
273 std::copy(inPtr,inPtr+nbCompo,loc);
283 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first+coarseSt[0]*coarseSt[1]*fineLocInCoarse[2].first),fact2(facts[2]),fact1(facts[1]),fact0(facts[0]);
284 for(int k=0;k<dims[2];k++)
286 for(int kfact=0;kfact<fact2;kfact++)
288 for(int j=0;j<dims[1];j++)
290 for(int jfact=0;jfact<fact1;jfact++)
292 for(int i=0;i<dims[0];i++)
294 double *loc(outPtr+(kk+i+j*coarseSt[0])*nbCompo);
295 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
297 if(kfact!=0 || jfact!=0 || ifact!=0)
298 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
300 std::copy(inPtr,inPtr+nbCompo,loc);
306 kk+=coarseSt[0]*coarseSt[1];
311 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : only dimensions 1, 2 and 3 supported !");
316 * This method spreads the values of coarse data \a coarseDA into \a fineDA.
318 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
319 * \param [in] coarseSt The cell structure of coarse mesh.
320 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
321 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
322 * \param [in] facts The refinement coefficient per axis.
323 * \sa CondenseFineToCoarse
325 void MEDCouplingIMesh::SpreadCoarseToFine(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts)
327 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
328 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the parameters 1 or 3 are NULL or not allocated !");
329 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
330 int nbCompo(fineDA->getNumberOfComponents());
331 if(coarseDA->getNumberOfComponents()!=nbCompo)
332 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the number of components of fine DA and coarse one mismatches !");
333 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
334 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
335 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
337 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
338 throw INTERP_KERNEL::Exception(oss.str().c_str());
340 int nbTuplesFine(fineDA->getNumberOfTuples());
341 if(nbTuplesFine%nbOfTuplesInFineExp!=0)
342 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : Invalid nb of tuples in fine DataArray regarding its structure !");
343 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
344 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
346 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
347 throw INTERP_KERNEL::Exception(oss.str().c_str());
349 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
350 double *outPtr(fineDA->getPointer());
351 const double *inPtr(coarseDA->begin());
353 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
358 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
359 for(int i=0;i<dims[0];i++)
361 const double *loc(inPtr+(offset+i)*nbCompo);
362 for(int ifact=0;ifact<fact0;ifact++)
363 outPtr=std::copy(loc,loc+nbCompo,outPtr);
369 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact0(facts[0]),fact1(facts[1]);
370 for(int j=0;j<dims[1];j++)
372 for(int jfact=0;jfact<fact1;jfact++)
374 for(int i=0;i<dims[0];i++)
376 const double *loc(inPtr+(kk+i)*nbCompo);
377 for(int ifact=0;ifact<fact0;ifact++)
378 outPtr=std::copy(loc,loc+nbCompo,outPtr);
387 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first+coarseSt[0]*coarseSt[1]*fineLocInCoarse[2].first),fact0(facts[0]),fact1(facts[2]),fact2(facts[2]);
388 for(int k=0;k<dims[2];k++)
390 for(int kfact=0;kfact<fact2;kfact++)
392 for(int j=0;j<dims[1];j++)
394 for(int jfact=0;jfact<fact1;jfact++)
396 for(int i=0;i<dims[0];i++)
398 const double *loc(inPtr+(kk+i+j*coarseSt[0])*nbCompo);
399 for(int ifact=0;ifact<fact0;ifact++)
400 outPtr=std::copy(loc,loc+nbCompo,outPtr);
405 kk+=coarseSt[0]*coarseSt[1];
410 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : only dimensions 1, 2 and 3 supported !");
414 void MEDCouplingIMesh::setSpaceDimension(int spaceDim)
416 if(spaceDim==_space_dim)
418 CheckSpaceDimension(spaceDim);
423 void MEDCouplingIMesh::updateTime() const
427 std::size_t MEDCouplingIMesh::getHeapMemorySizeWithoutChildren() const
429 return MEDCouplingStructuredMesh::getHeapMemorySizeWithoutChildren();
432 std::vector<const BigMemoryObject *> MEDCouplingIMesh::getDirectChildren() const
434 return std::vector<const BigMemoryObject *>();
438 * This method copyies all tiny strings from other (name and components name).
439 * @throw if other and this have not same mesh type.
441 void MEDCouplingIMesh::copyTinyStringsFrom(const MEDCouplingMesh *other)
443 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
445 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::copyTinyStringsFrom : meshes have not same type !");
446 MEDCouplingStructuredMesh::copyTinyStringsFrom(other);
450 bool MEDCouplingIMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
453 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::isEqualIfNotWhy : input other pointer is null !");
454 const MEDCouplingIMesh *otherC(dynamic_cast<const MEDCouplingIMesh *>(other));
457 reason="mesh given in input is not castable in MEDCouplingIMesh !";
460 if(!MEDCouplingStructuredMesh::isEqualIfNotWhy(other,prec,reason))
462 if(!isEqualWithoutConsideringStrInternal(otherC,prec,reason))
464 if(_axis_unit!=otherC->_axis_unit)
466 reason="The units of axis are not the same !";
472 bool MEDCouplingIMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
474 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
478 return isEqualWithoutConsideringStrInternal(other,prec,tmp);
481 bool MEDCouplingIMesh::isEqualWithoutConsideringStrInternal(const MEDCouplingMesh *other, double prec, std::string& reason) const
483 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
486 if(_space_dim!=otherC->_space_dim)
488 std::ostringstream oss;
489 oss << "The spaceDimension of this (" << _space_dim << ") is not equal to those of other (" << otherC->_space_dim << ") !";
492 checkSpaceDimension();
493 for(int i=0;i<_space_dim;i++)
495 if(fabs(_origin[i]-otherC->_origin[i])>prec)
497 std::ostringstream oss;
498 oss << "The origin of this and other differs at " << i << " !";
503 for(int i=0;i<_space_dim;i++)
505 if(fabs(_dxyz[i]-otherC->_dxyz[i])>prec)
507 std::ostringstream oss;
508 oss << "The delta of this and other differs at " << i << " !";
513 for(int i=0;i<_space_dim;i++)
515 if(_structure[i]!=otherC->_structure[i])
517 std::ostringstream oss;
518 oss << "The structure of this and other differs at " << i << " !";
526 void MEDCouplingIMesh::checkDeepEquivalWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
527 DataArrayInt *&cellCor, DataArrayInt *&nodeCor) const
529 if(!isEqualWithoutConsideringStr(other,prec))
530 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalWith : Meshes are not the same !");
534 * Nothing is done here (except to check that the other is a ParaMEDMEM::MEDCouplingIMesh instance too).
535 * The user intend that the nodes are the same, so by construction of ParaMEDMEM::MEDCouplingIMesh, \a this and \a other are the same !
537 void MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
538 DataArrayInt *&cellCor) const
540 if(!isEqualWithoutConsideringStr(other,prec))
541 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith : Meshes are not the same !");
544 void MEDCouplingIMesh::checkCoherency() const
546 checkSpaceDimension();
547 for(int i=0;i<_space_dim;i++)
550 std::ostringstream oss; oss << "MEDCouplingIMesh::checkCoherency : On axis " << i << "/" << _space_dim << ", number of nodes is equal to " << _structure[i] << " ! must be >=1 !";
551 throw INTERP_KERNEL::Exception(oss.str().c_str());
555 void MEDCouplingIMesh::checkCoherency1(double eps) const
560 void MEDCouplingIMesh::checkCoherency2(double eps) const
562 checkCoherency1(eps);
565 void MEDCouplingIMesh::getNodeGridStructure(int *res) const
567 checkSpaceDimension();
568 std::copy(_structure,_structure+_space_dim,res);
571 std::vector<int> MEDCouplingIMesh::getNodeGridStructure() const
573 checkSpaceDimension();
574 std::vector<int> ret(_structure,_structure+_space_dim);
578 MEDCouplingStructuredMesh *MEDCouplingIMesh::buildStructuredSubPart(const std::vector< std::pair<int,int> >& cellPart) const
581 int dim(getSpaceDimension());
582 if(dim!=(int)cellPart.size())
584 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : the space dimension is " << dim << " and cell part size is " << cellPart.size() << " !";
585 throw INTERP_KERNEL::Exception(oss.str().c_str());
587 double retOrigin[3]={0.,0.,0.};
588 int retStruct[3]={0,0,0};
589 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(dynamic_cast<MEDCouplingIMesh *>(deepCpy()));
590 for(int i=0;i<dim;i++)
592 int startNode(cellPart[i].first),endNode(cellPart[i].second+1);
593 int myDelta(endNode-startNode);
594 if(startNode<0 || startNode>=_structure[i])
596 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : At dimension #" << i << " the start node id is " << startNode << " it should be in [0," << _structure[i] << ") !";
597 throw INTERP_KERNEL::Exception(oss.str().c_str());
599 if(myDelta<0 || myDelta>_structure[i])
601 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : Along dimension #" << i << " the number of nodes is " << _structure[i] << ", and you are requesting for " << myDelta << " nodes wide range !" << std::endl;
602 throw INTERP_KERNEL::Exception(oss.str().c_str());
604 retOrigin[i]=_origin[i]+startNode*_dxyz[i];
605 retStruct[i]=myDelta;
607 ret->setNodeStruct(retStruct,retStruct+dim);
608 ret->setOrigin(retOrigin,retOrigin+dim);
609 ret->checkCoherency();
614 * Return the space dimension of \a this.
616 int MEDCouplingIMesh::getSpaceDimension() const
621 void MEDCouplingIMesh::getCoordinatesOfNode(int nodeId, std::vector<double>& coo) const
624 int spaceDim(getSpaceDimension());
625 getSplitNodeValues(tmp);
627 GetPosFromId(nodeId,spaceDim,tmp,tmp2);
628 for(int j=0;j<spaceDim;j++)
629 coo.push_back(_origin[j]+_dxyz[j]*tmp2[j]);
632 std::string MEDCouplingIMesh::simpleRepr() const
634 std::ostringstream ret;
635 ret << "Image grid with name : \"" << getName() << "\"\n";
636 ret << "Description of mesh : \"" << getDescription() << "\"\n";
638 double tt(getTime(tmpp1,tmpp2));
639 int spaceDim(_space_dim);
640 ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
641 ret << "Iteration : " << tmpp1 << " Order : " << tmpp2 << "\n";
642 ret << "Space dimension : " << spaceDim << "\n";
643 if(spaceDim<0 || spaceDim>3)
645 ret << "The nodal structure is : "; std::copy(_structure,_structure+spaceDim,std::ostream_iterator<int>(ret," ")); ret << "\n";
646 ret << "The origin position is [" << _axis_unit << "]: ";
647 std::copy(_origin,_origin+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
648 ret << "The intervals along axis are : ";
649 std::copy(_dxyz,_dxyz+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
653 std::string MEDCouplingIMesh::advancedRepr() const
658 void MEDCouplingIMesh::getBoundingBox(double *bbox) const
661 int dim(getSpaceDimension());
662 for(int idim=0; idim<dim; idim++)
664 bbox[2*idim]=_origin[idim];
665 bbox[2*idim+1]=_origin[idim]+_dxyz[idim]*_structure[idim];
670 * Returns a new MEDCouplingFieldDouble containing volumes of cells constituting \a this
672 * For 1D cells, the returned field contains lengths.<br>
673 * For 2D cells, the returned field contains areas.<br>
674 * For 3D cells, the returned field contains volumes.
675 * \param [in] isAbs - a not used parameter.
676 * \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble on cells
677 * and one time . The caller is to delete this field using decrRef() as it is no
680 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureField(bool isAbs) const
683 std::string name="MeasureOfMesh_";
685 int nbelem(getNumberOfCells());
686 MEDCouplingFieldDouble *field(MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME));
687 field->setName(name);
688 DataArrayDouble* array(DataArrayDouble::New());
689 array->alloc(nbelem,1);
690 array->fillWithValue(getMeasureOfAnyCell());
691 field->setArray(array) ;
693 field->setMesh(const_cast<MEDCouplingIMesh *>(this));
694 field->synchronizeTimeWithMesh();
699 * not implemented yet !
701 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureFieldOnNode(bool isAbs) const
703 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::getMeasureFieldOnNode : not implemented yet !");
707 int MEDCouplingIMesh::getCellContainingPoint(const double *pos, double eps) const
709 int dim(getSpaceDimension()),ret(0),coeff(1);
710 for(int i=0;i<dim;i++)
712 int nbOfCells(_structure[i]-1);
714 int tmp((int)((ref-_origin[i])/_dxyz[i]));
715 if(tmp>=0 && tmp<nbOfCells)
726 void MEDCouplingIMesh::rotate(const double *center, const double *vector, double angle)
728 throw INTERP_KERNEL::Exception("No rotation available on IMesh : Traduce it to unstructured mesh to apply it !");
732 * Translates all nodes of \a this mesh by a given vector. Actually, it adds each
733 * component of the \a vector to all node coordinates of a corresponding axis.
734 * \param [in] vector - the translation vector whose size must be not less than \a
735 * this->getSpaceDimension().
737 void MEDCouplingIMesh::translate(const double *vector)
739 checkSpaceDimension();
740 int dim(getSpaceDimension());
741 std::transform(_origin,_origin+dim,vector,_origin,std::plus<double>());
746 * Applies scaling transformation to all nodes of \a this mesh.
747 * \param [in] point - coordinates of a scaling center. This array is to be of
748 * size \a this->getSpaceDimension() at least.
749 * \param [in] factor - a scale factor.
751 void MEDCouplingIMesh::scale(const double *point, double factor)
753 checkSpaceDimension();
754 int dim(getSpaceDimension());
755 std::transform(_origin,_origin+dim,point,_origin,std::minus<double>());
756 std::transform(_origin,_origin+dim,_origin,std::bind2nd(std::multiplies<double>(),factor));
757 std::transform(_dxyz,_dxyz+dim,_dxyz,std::bind2nd(std::multiplies<double>(),factor));
758 std::transform(_origin,_origin+dim,point,_origin,std::plus<double>());
762 MEDCouplingMesh *MEDCouplingIMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
764 //not implemented yet !
769 * Returns a new DataArrayDouble holding coordinates of all nodes of \a this mesh.
770 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
771 * this->getNumberOfNodes() tuples per \a this->getSpaceDimension()
772 * components. The caller is to delete this array using decrRef() as it is
775 DataArrayDouble *MEDCouplingIMesh::getCoordinatesAndOwner() const
778 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
779 int spaceDim(getSpaceDimension()),nbNodes(getNumberOfNodes());
780 ret->alloc(nbNodes,spaceDim);
781 double *pt(ret->getPointer());
782 ret->setInfoOnComponents(buildInfoOnComponents());
784 getSplitNodeValues(tmp);
785 for(int i=0;i<nbNodes;i++)
787 GetPosFromId(i,spaceDim,tmp,tmp2);
788 for(int j=0;j<spaceDim;j++)
789 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+_origin[j];
795 * Returns a new DataArrayDouble holding barycenters of all cells. The barycenter is
796 * computed by averaging coordinates of cell nodes.
797 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
798 * this->getNumberOfCells() tuples per \a this->getSpaceDimension()
799 * components. The caller is to delete this array using decrRef() as it is
802 DataArrayDouble *MEDCouplingIMesh::getBarycenterAndOwner() const
805 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
806 int spaceDim(getSpaceDimension()),nbCells(getNumberOfCells()),tmp[3],tmp2[3];
807 ret->alloc(nbCells,spaceDim);
808 double *pt(ret->getPointer()),shiftOrigin[3];
809 std::transform(_dxyz,_dxyz+spaceDim,shiftOrigin,std::bind2nd(std::multiplies<double>(),0.5));
810 std::transform(_origin,_origin+spaceDim,shiftOrigin,shiftOrigin,std::plus<double>());
811 getSplitCellValues(tmp);
812 ret->setInfoOnComponents(buildInfoOnComponents());
813 for(int i=0;i<nbCells;i++)
815 GetPosFromId(i,spaceDim,tmp,tmp2);
816 for(int j=0;j<spaceDim;j++)
817 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+shiftOrigin[j];
822 DataArrayDouble *MEDCouplingIMesh::computeIsoBarycenterOfNodesPerCell() const
824 return MEDCouplingIMesh::getBarycenterAndOwner();
827 void MEDCouplingIMesh::renumberCells(const int *old2NewBg, bool check)
829 throw INTERP_KERNEL::Exception("Functionnality of renumbering cell not available for IMesh !");
832 void MEDCouplingIMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
835 double time(getTime(it,order));
838 littleStrings.clear();
839 littleStrings.push_back(getName());
840 littleStrings.push_back(getDescription());
841 littleStrings.push_back(getTimeUnit());
842 littleStrings.push_back(getAxisUnit());
843 tinyInfo.push_back(it);
844 tinyInfo.push_back(order);
845 tinyInfo.push_back(_space_dim);
846 tinyInfo.insert(tinyInfo.end(),_structure,_structure+3);
847 tinyInfoD.push_back(time);
848 tinyInfoD.insert(tinyInfoD.end(),_dxyz,_dxyz+3);
849 tinyInfoD.insert(tinyInfoD.end(),_origin,_origin+3);
852 void MEDCouplingIMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
858 void MEDCouplingIMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
860 a1=DataArrayInt::New();
862 a2=DataArrayDouble::New();
866 void MEDCouplingIMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
867 const std::vector<std::string>& littleStrings)
869 setName(littleStrings[0]);
870 setDescription(littleStrings[1]);
871 setTimeUnit(littleStrings[2]);
872 setAxisUnit(littleStrings[3]);
873 setTime(tinyInfoD[0],tinyInfo[0],tinyInfo[1]);
874 _space_dim=tinyInfo[2];
875 _structure[0]=tinyInfo[3]; _structure[1]=tinyInfo[4]; _structure[2]=tinyInfo[5];
876 _dxyz[0]=tinyInfoD[1]; _dxyz[1]=tinyInfoD[2]; _dxyz[2]=tinyInfoD[3];
877 _origin[0]=tinyInfoD[4]; _origin[1]=tinyInfoD[5]; _origin[2]=tinyInfoD[6];
881 void MEDCouplingIMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
884 std::ostringstream extent,origin,spacing;
888 { extent << "0 " << _structure[i]-1 << " "; origin << _origin[i] << " "; spacing << _dxyz[i] << " "; }
890 { extent << "0 0 "; origin << "0 "; spacing << "0 "; }
892 ofs << " <" << getVTKDataSetType() << " WholeExtent=\"" << extent.str() << "\" Origin=\"" << origin.str() << "\" Spacing=\"" << spacing.str() << "\">\n";
893 ofs << " <Piece Extent=\"" << extent.str() << "\">\n";
894 ofs << " <PointData>\n" << pointData << std::endl;
895 ofs << " </PointData>\n";
896 ofs << " <CellData>\n" << cellData << std::endl;
897 ofs << " </CellData>\n";
898 ofs << " <Coordinates>\n";
899 ofs << " </Coordinates>\n";
900 ofs << " </Piece>\n";
901 ofs << " </" << getVTKDataSetType() << ">\n";
904 void MEDCouplingIMesh::reprQuickOverview(std::ostream& stream) const
906 stream << "MEDCouplingIMesh C++ instance at " << this << ". Name : \"" << getName() << "\". Space dimension : " << _space_dim << ".";
907 if(_space_dim<0 || _space_dim>3)
910 std::ostringstream stream0,stream1;
911 int nbNodes(1),nbCells(0);
913 for(int i=0;i<_space_dim;i++)
916 int tmpNodes(_structure[i]);
917 stream1 << "- Axis " << tmp << " : " << tmpNodes << " nodes (orig=" << _origin[i] << ", inter=" << _dxyz[i] << ").";
919 stream1 << std::endl;
925 nbCells=nbCells==0?tmpNodes-1:nbCells*(tmpNodes-1);
929 stream0 << "Number of cells : " << nbCells << ", Number of nodes : " << nbNodes;
930 stream << stream0.str();
934 stream << stream1.str();
937 std::string MEDCouplingIMesh::getVTKDataSetType() const
939 return std::string("ImageData");
942 std::vector<std::string> MEDCouplingIMesh::buildInfoOnComponents() const
944 checkSpaceDimension();
945 int dim(getSpaceDimension());
946 std::vector<std::string> ret(dim);
947 for(int i=0;i<dim;i++)
949 std::ostringstream oss;
950 char tmp('X'+i); oss << tmp;
951 ret[i]=DataArray::BuildInfoFromVarAndUnit(oss.str(),_axis_unit);
956 void MEDCouplingIMesh::checkSpaceDimension() const
958 CheckSpaceDimension(_space_dim);
961 void MEDCouplingIMesh::CheckSpaceDimension(int spaceDim)
963 if(spaceDim<0 || spaceDim>3)
964 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CheckSpaceDimension : input spaceDim must be in [0,1,2,3] !");
967 int MEDCouplingIMesh::FindIntRoot(int val, int order)
972 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : input val is < 0 ! Not possible to compute a root !");
975 if(order!=2 && order!=3)
976 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the order available are 0,1,2 or 3 !");
977 double valf((double)val);
980 double retf(sqrt(valf));
983 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect square root !");
988 double retf(std::pow(val,0.3333333333333333));
989 int ret((int)retf),ret2(ret+1);
990 if(ret*ret*ret!=val && ret2*ret2*ret2!=val)
991 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect cublic root !");