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);
79 * This method creates a copy of \a this enlarged by \a ghostLev cells on each axis.
80 * If \a ghostLev equal to 0 this method behaves as MEDCouplingIMesh::clone.
82 * \param [in] ghostLev - the ghost level expected
83 * \return MEDCouplingIMesh * - a newly alloacted object to be managed by the caller.
84 * \throw if \a ghostLev < 0.
86 MEDCouplingIMesh *MEDCouplingIMesh::buildWithGhost(int ghostLev) const
89 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::buildWithGhost : the ghostLev must be >= 0 !");
91 int spaceDim(getSpaceDimension());
92 double origin[3],dxyz[3];
94 for(int i=0;i<spaceDim;i++)
96 origin[i]=_origin[i]-ghostLev*_dxyz[i];
98 structure[i]=_structure[i]+2*ghostLev;
100 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(MEDCouplingIMesh::New(getName(),spaceDim,structure,structure+spaceDim,origin,origin+spaceDim,dxyz,dxyz+spaceDim));
101 ret->copyTinyInfoFrom(this);
105 void MEDCouplingIMesh::setNodeStruct(const int *nodeStrctStart, const int *nodeStrctStop)
107 checkSpaceDimension();
108 int sz((int)std::distance(nodeStrctStart,nodeStrctStop));
110 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setNodeStruct : input vector of node structure has not the right size ! Or change space dimension before calling it !");
111 std::copy(nodeStrctStart,nodeStrctStop,_structure);
115 std::vector<int> MEDCouplingIMesh::getNodeStruct() const
117 checkSpaceDimension();
118 return std::vector<int>(_structure,_structure+_space_dim);
121 void MEDCouplingIMesh::setOrigin(const double *originStart, const double *originStop)
123 checkSpaceDimension();
124 int sz((int)std::distance(originStart,originStop));
126 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setOrigin : input vector of origin vector has not the right size ! Or change space dimension before calling it !");
127 std::copy(originStart,originStop,_origin);
131 std::vector<double> MEDCouplingIMesh::getOrigin() const
133 checkSpaceDimension();
134 return std::vector<double>(_origin,_origin+_space_dim);
137 void MEDCouplingIMesh::setDXYZ(const double *dxyzStart, const double *dxyzStop)
139 checkSpaceDimension();
140 int sz((int)std::distance(dxyzStart,dxyzStop));
142 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::setDXYZ : input vector of dxyz vector has not the right size ! Or change space dimension before calling it !");
143 std::copy(dxyzStart,dxyzStop,_dxyz);
147 std::vector<double> MEDCouplingIMesh::getDXYZ() const
149 checkSpaceDimension();
150 return std::vector<double>(_dxyz,_dxyz+_space_dim);
153 void MEDCouplingIMesh::setAxisUnit(const std::string& unitName)
159 std::string MEDCouplingIMesh::getAxisUnit() const
165 * This method returns the measure of any cell in \a this.
166 * This specific method of image grid mesh utilizes the fact that any cell in \a this have the same measure.
167 * The value returned by this method is those used to feed the returned field in the MEDCouplingIMesh::getMeasureField.
169 * \sa getMeasureField
171 double MEDCouplingIMesh::getMeasureOfAnyCell() const
174 int dim(getSpaceDimension());
176 for(int i=0;i<dim;i++)
182 * This method is allows to convert \a this into MEDCouplingCMesh instance.
183 * This method is the middle level between MEDCouplingIMesh and the most general MEDCouplingUMesh.
184 * This method is useful for MED writers that do not have still the image grid support.
186 * \sa MEDCouplingMesh::buildUnstructured
188 MEDCouplingCMesh *MEDCouplingIMesh::convertToCartesian() const
191 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCMesh> ret(MEDCouplingCMesh::New());
193 { ret->copyTinyInfoFrom(this); }
194 catch(INTERP_KERNEL::Exception& ) { }
195 int spaceDim(getSpaceDimension());
196 std::vector<std::string> infos(buildInfoOnComponents());
197 for(int i=0;i<spaceDim;i++)
199 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr(DataArrayDouble::New()); arr->alloc(_structure[i],1); arr->setInfoOnComponent(0,infos[i]);
200 arr->iota(); arr->applyLin(_dxyz[i],_origin[i]);
201 ret->setCoordsAt(i,arr);
207 * This method refines \a this uniformaly along all of its dimensions. In case of success the space covered by \a this will remain
208 * the same before the invocation except that the number of cells will be multiplied by \a factor ^ this->getMeshDimension().
209 * The origin of \a this will be not touched only spacing and node structure will be changed.
210 * This method can be useful for AMR users.
212 void MEDCouplingIMesh::refineWithFactor(const std::vector<int>& factors)
214 if((int)factors.size()!=_space_dim)
215 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::refineWithFactor : refinement factors must have size equal to spaceDim !");
217 std::vector<int> structure(_structure,_structure+3);
218 std::vector<double> dxyz(_dxyz,_dxyz+3);
219 for(int i=0;i<_space_dim;i++)
223 std::ostringstream oss; oss << "MEDCouplingIMesh::refineWithFactor : factor for axis #" << i << " (" << factors[i] << ")is invalid ! Must be > 0 !";
224 throw INTERP_KERNEL::Exception(oss.str().c_str());
226 int factAbs(std::abs(factors[i]));
227 double fact2(1./(double)factors[i]);
228 structure[i]=(_structure[i]-1)*factAbs+1;
229 dxyz[i]=fact2*_dxyz[i];
231 std::copy(structure.begin(),structure.end(),_structure);
232 std::copy(dxyz.begin(),dxyz.end(),_dxyz);
237 * This method returns a newly created mesh containing a single cell in it. This returned cell covers exactly the space covered by \a this.
239 * \return MEDCouplingIMesh * - A newly created object (to be managed by the caller with decrRef) containing simply one cell.
241 * \throw if \a this does not pass the \c checkCoherency test.
243 MEDCouplingIMesh *MEDCouplingIMesh::asSingleCell() const
246 int spaceDim(getSpaceDimension()),nodeSt[3];
248 for(int i=0;i<spaceDim;i++)
253 dxyz[i]=(_structure[i]-1)*_dxyz[i];
257 nodeSt[i]=_structure[i];
261 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(MEDCouplingIMesh::New(getName(),getSpaceDimension(),nodeSt,nodeSt+spaceDim,_origin,_origin+spaceDim,dxyz,dxyz+spaceDim));
262 ret->copyTinyInfoFrom(this);
267 * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
268 * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
269 * 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.
271 * \param [in] coarseSt The cell structure of coarse mesh.
272 * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
273 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
274 * \param [in] facts The refinement coefficient per axis.
275 * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
277 * \sa CondenseFineToCoarseGhost,SpreadCoarseToFine
279 void MEDCouplingIMesh::CondenseFineToCoarse(const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, DataArrayDouble *coarseDA)
281 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
282 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : All input vectors (dimension) must have the same size !");
283 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
284 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the parameters 1 or 3 are NULL or not allocated !");
285 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
286 int nbCompo(fineDA->getNumberOfComponents());
287 if(coarseDA->getNumberOfComponents()!=nbCompo)
288 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : the number of components of fine DA and coarse one mismatches !");
289 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
290 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) !");
291 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
293 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
294 throw INTERP_KERNEL::Exception(oss.str().c_str());
296 int nbTuplesFine(fineDA->getNumberOfTuples());
297 if(nbOfTuplesInFineExp==0)
299 if(nbTuplesFine%nbOfTuplesInFineExp!=0)
300 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Invalid nb of tuples in fine DataArray regarding its structure !");
301 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
302 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
304 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
305 throw INTERP_KERNEL::Exception(oss.str().c_str());
307 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
308 double *outPtr(coarseDA->getPointer());
309 const double *inPtr(fineDA->begin());
311 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
316 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
317 for(int i=0;i<dims[0];i++)
319 double *loc(outPtr+(offset+i)*nbCompo);
320 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
323 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
325 std::copy(inPtr,inPtr+nbCompo,loc);
332 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact1(facts[1]),fact0(facts[0]);
333 for(int j=0;j<dims[1];j++)
335 for(int jfact=0;jfact<fact1;jfact++)
337 for(int i=0;i<dims[0];i++)
339 double *loc(outPtr+(kk+i)*nbCompo);
340 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
342 if(jfact!=0 || ifact!=0)
343 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
345 std::copy(inPtr,inPtr+nbCompo,loc);
355 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]);
356 for(int k=0;k<dims[2];k++)
358 for(int kfact=0;kfact<fact2;kfact++)
360 for(int j=0;j<dims[1];j++)
362 for(int jfact=0;jfact<fact1;jfact++)
364 for(int i=0;i<dims[0];i++)
366 double *loc(outPtr+(kk+i+j*coarseSt[0])*nbCompo);
367 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
369 if(kfact!=0 || jfact!=0 || ifact!=0)
370 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
372 std::copy(inPtr,inPtr+nbCompo,loc);
378 kk+=coarseSt[0]*coarseSt[1];
383 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : only dimensions 1, 2 and 3 supported !");
388 * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
389 * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
390 * 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.
392 * \param [in] coarseSt The cell structure of coarse mesh.
393 * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
394 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
395 * \param [in] facts The refinement coefficient per axis.
396 * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
397 * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
399 * \sa CondenseFineToCoarse,SpreadCoarseToFineGhost
401 void MEDCouplingIMesh::CondenseFineToCoarseGhost(const std::vector<int>& coarseSt, const DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, DataArrayDouble *coarseDA, int ghostSize)
404 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : ghost level has to be >= 0 !");
405 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
406 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : All input vectors (dimension) must have the same size !");
407 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
408 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the parameters 1 or 3 are NULL or not allocated !");
409 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
410 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
411 int nbCompo(fineDA->getNumberOfComponents());
412 if(coarseDA->getNumberOfComponents()!=nbCompo)
413 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the number of components of fine DA and coarse one mismatches !");
414 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
415 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
416 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
418 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples in coarse DataArray having " << coarseDA->getNumberOfTuples() << " !";
419 throw INTERP_KERNEL::Exception(oss.str().c_str());
422 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
423 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
424 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
425 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
426 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
428 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
429 throw INTERP_KERNEL::Exception(oss.str().c_str());
432 double *outPtr(coarseDA->getPointer());
433 const double *inPtr(fineDA->begin());
435 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
440 int offset(fineLocInCoarse[0].first+ghostSize),fact0(facts[0]);
441 inPtr+=ghostSize*nbCompo;
442 for(int i=0;i<dims[0];i++)
444 double *loc(outPtr+(offset+i)*nbCompo);
445 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
448 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
450 std::copy(inPtr,inPtr+nbCompo,loc);
457 int nxwg(coarseSt[0]+2*ghostSize);
458 int kk(fineLocInCoarse[0].first+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize)),fact1(facts[1]),fact0(facts[0]);
459 inPtr+=(dims[0]*fact0+2*ghostSize)*nbCompo;
460 for(int j=0;j<dims[1];j++)
462 for(int jfact=0;jfact<fact1;jfact++)
464 inPtr+=ghostSize*nbCompo;
465 for(int i=0;i<dims[0];i++)
467 double *loc(outPtr+(kk+i)*nbCompo);
468 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
470 if(jfact!=0 || ifact!=0)
471 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
473 std::copy(inPtr,inPtr+nbCompo,loc);
476 inPtr+=ghostSize*nbCompo;
483 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : only dimensions 1, 2 supported !");
488 * This method spreads the values of coarse data \a coarseDA into \a fineDA.
490 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
491 * \param [in] coarseSt The cell structure of coarse mesh.
492 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
493 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
494 * \param [in] facts The refinement coefficient per axis.
495 * \sa SpreadCoarseToFineGhost, CondenseFineToCoarse
497 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)
499 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
500 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : All input vectors (dimension) must have the same size !");
501 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
502 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the parameters 1 or 3 are NULL or not allocated !");
503 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
504 int nbCompo(fineDA->getNumberOfComponents());
505 if(coarseDA->getNumberOfComponents()!=nbCompo)
506 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the number of components of fine DA and coarse one mismatches !");
507 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
508 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) !");
509 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
511 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
512 throw INTERP_KERNEL::Exception(oss.str().c_str());
514 int nbTuplesFine(fineDA->getNumberOfTuples());
515 if(nbTuplesFine%nbOfTuplesInFineExp!=0)
516 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : Invalid nb of tuples in fine DataArray regarding its structure !");
517 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
518 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
520 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
521 throw INTERP_KERNEL::Exception(oss.str().c_str());
523 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
524 double *outPtr(fineDA->getPointer());
525 const double *inPtr(coarseDA->begin());
527 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
532 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
533 for(int i=0;i<dims[0];i++)
535 const double *loc(inPtr+(offset+i)*nbCompo);
536 for(int ifact=0;ifact<fact0;ifact++)
537 outPtr=std::copy(loc,loc+nbCompo,outPtr);
543 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact0(facts[0]),fact1(facts[1]);
544 for(int j=0;j<dims[1];j++)
546 for(int jfact=0;jfact<fact1;jfact++)
548 for(int i=0;i<dims[0];i++)
550 const double *loc(inPtr+(kk+i)*nbCompo);
551 for(int ifact=0;ifact<fact0;ifact++)
552 outPtr=std::copy(loc,loc+nbCompo,outPtr);
561 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]);
562 for(int k=0;k<dims[2];k++)
564 for(int kfact=0;kfact<fact2;kfact++)
566 for(int j=0;j<dims[1];j++)
568 for(int jfact=0;jfact<fact1;jfact++)
570 for(int i=0;i<dims[0];i++)
572 const double *loc(inPtr+(kk+i+j*coarseSt[0])*nbCompo);
573 for(int ifact=0;ifact<fact0;ifact++)
574 outPtr=std::copy(loc,loc+nbCompo,outPtr);
579 kk+=coarseSt[0]*coarseSt[1];
584 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : only dimensions 1, 2 and 3 supported !");
589 * This method spreads the values of coarse data \a coarseDA into \a fineDA.
591 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
592 * \param [in] coarseSt The cell structure of coarse mesh.
593 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
594 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
595 * \param [in] facts The refinement coefficient per axis.
596 * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
598 * \sa CondenseFineToCoarse, SpreadCoarseToFineGhostZone
600 void MEDCouplingIMesh::SpreadCoarseToFineGhost(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
603 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : ghost level has to be >= 0 !");
604 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
605 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : All input vectors (dimension) must have the same size !");
606 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
607 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the parameters 1 or 3 are NULL or not allocated !");
608 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
609 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
610 int nbCompo(fineDA->getNumberOfComponents());
611 if(coarseDA->getNumberOfComponents()!=nbCompo)
612 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the number of components of fine DA and coarse one mismatches !");
613 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
614 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
615 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
617 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
618 throw INTERP_KERNEL::Exception(oss.str().c_str());
621 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
622 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
623 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
624 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
625 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
627 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
628 throw INTERP_KERNEL::Exception(oss.str().c_str());
631 double *outPtr(fineDA->getPointer());
632 const double *inPtr(coarseDA->begin());
634 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
639 int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
640 for(int i=0;i<ghostSize;i++)
641 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
642 offset=fineLocInCoarse[0].first+ghostSize;
643 for(int i=0;i<dims[0];i++)
645 const double *loc(inPtr+(offset+i)*nbCompo);
646 for(int ifact=0;ifact<fact0;ifact++)
647 outPtr=std::copy(loc,loc+nbCompo,outPtr);
649 offset=fineLocInCoarse[0].second+ghostSize;
650 for(int i=0;i<ghostSize;i++)
651 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
656 int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
657 int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
658 for(int jg=0;jg<ghostSize;jg++)
660 for(int ig=0;ig<ghostSize;ig++)
661 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
663 for(int ig=0;ig<dims[0];ig++,kk0++)
664 for(int ifact=0;ifact<fact0;ifact++)
665 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
666 for(int ik=0;ik<ghostSize;ik++)
667 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
669 for(int j=0;j<dims[1];j++)
671 kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
672 for(int jfact=0;jfact<fact1;jfact++)
674 for(int ig=0;ig<ghostSize;ig++)
675 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
676 int kk0(kk+1);//1 not ghost. We make the hypothesis that factors is >= ghostlev
677 for(int i=0;i<dims[0];i++,kk0++)
679 const double *loc(inPtr+kk0*nbCompo);
680 for(int ifact=0;ifact<fact0;ifact++)
681 outPtr=std::copy(loc,loc+nbCompo,outPtr);
683 for(int ig=0;ig<ghostSize;ig++)
684 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
687 kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
688 for(int jg=0;jg<ghostSize;jg++)
690 for(int ig=0;ig<ghostSize;ig++)
691 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
693 for(int ig=0;ig<dims[0];ig++,kk0++)
694 for(int ifact=0;ifact<fact0;ifact++)
695 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
696 for(int ik=0;ik<ghostSize;ik++)
697 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
702 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : only dimensions 1, 2 supported !");
707 * This method spreads the values of coarse data \a coarseDA into \a fineDA \b ONLY \b in \b the \b ghost \b zone (contrary to SpreadCoarseToFineGhost that spread the values everywhere).
709 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
710 * \param [in] coarseSt The cell structure of coarse mesh.
711 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
712 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
713 * \param [in] facts The refinement coefficient per axis.
714 * \param [in] ghostSize - The size of the ghost zone. The ghost zone is expected to be the same for all axis and both for coarse and fine meshes.
716 * \sa SpreadCoarseToFineGhost
718 void MEDCouplingIMesh::SpreadCoarseToFineGhostZone(const DataArrayDouble *coarseDA, const std::vector<int>& coarseSt, DataArrayDouble *fineDA, const std::vector< std::pair<int,int> >& fineLocInCoarse, const std::vector<int>& facts, int ghostSize)
721 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : ghost level has to be >= 0 !");
722 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
723 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : All input vectors (dimension) must have the same size !");
724 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
725 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the parameters 1 or 3 are NULL or not allocated !");
726 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
727 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
728 int nbCompo(fineDA->getNumberOfComponents());
729 if(coarseDA->getNumberOfComponents()!=nbCompo)
730 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the number of components of fine DA and coarse one mismatches !");
731 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
732 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the size of fineLocInCoarse (4th param) and facts (5th param) must be equal to the sier of coarseSt (2nd param) !");
733 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
735 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
736 throw INTERP_KERNEL::Exception(oss.str().c_str());
739 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
740 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
741 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
742 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
743 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
745 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
746 throw INTERP_KERNEL::Exception(oss.str().c_str());
749 double *outPtr(fineDA->getPointer());
750 const double *inPtr(coarseDA->begin());
752 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
757 int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
758 for(int i=0;i<ghostSize;i++)
759 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
760 outPtr+=nbCompo*fact0*dims[0];
761 offset=fineLocInCoarse[0].second+ghostSize;
762 for(int i=0;i<ghostSize;i++)
763 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
768 int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
769 int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
770 for(int jg=0;jg<ghostSize;jg++)
772 for(int ig=0;ig<ghostSize;ig++)
773 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
775 for(int ig=0;ig<dims[0];ig++,kk0++)
776 for(int ifact=0;ifact<fact0;ifact++)
777 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
778 for(int ik=0;ik<ghostSize;ik++)
779 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
781 for(int j=0;j<dims[1];j++)
783 kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
784 for(int jfact=0;jfact<fact1;jfact++)
786 for(int ig=0;ig<ghostSize;ig++)
787 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
788 int kk0(kk+1+dims[0]);//1 not ghost. We make the hypothesis that factors is >= ghostlev
789 outPtr+=fact0*nbCompo*dims[0];
790 for(int ig=0;ig<ghostSize;ig++)
791 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
794 kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
795 for(int jg=0;jg<ghostSize;jg++)
797 for(int ig=0;ig<ghostSize;ig++)
798 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
800 for(int ig=0;ig<dims[0];ig++,kk0++)
801 for(int ifact=0;ifact<fact0;ifact++)
802 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
803 for(int ik=0;ik<ghostSize;ik++)
804 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
809 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : only dimensions 1, 2 supported !");
813 void MEDCouplingIMesh::setSpaceDimension(int spaceDim)
815 if(spaceDim==_space_dim)
817 CheckSpaceDimension(spaceDim);
822 void MEDCouplingIMesh::updateTime() const
826 std::size_t MEDCouplingIMesh::getHeapMemorySizeWithoutChildren() const
828 return MEDCouplingStructuredMesh::getHeapMemorySizeWithoutChildren();
831 std::vector<const BigMemoryObject *> MEDCouplingIMesh::getDirectChildren() const
833 return std::vector<const BigMemoryObject *>();
837 * This method copyies all tiny strings from other (name and components name).
838 * @throw if other and this have not same mesh type.
840 void MEDCouplingIMesh::copyTinyStringsFrom(const MEDCouplingMesh *other)
842 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
844 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::copyTinyStringsFrom : meshes have not same type !");
845 MEDCouplingStructuredMesh::copyTinyStringsFrom(other);
849 bool MEDCouplingIMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
852 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::isEqualIfNotWhy : input other pointer is null !");
853 const MEDCouplingIMesh *otherC(dynamic_cast<const MEDCouplingIMesh *>(other));
856 reason="mesh given in input is not castable in MEDCouplingIMesh !";
859 if(!MEDCouplingStructuredMesh::isEqualIfNotWhy(other,prec,reason))
861 if(!isEqualWithoutConsideringStrInternal(otherC,prec,reason))
863 if(_axis_unit!=otherC->_axis_unit)
865 reason="The units of axis are not the same !";
871 bool MEDCouplingIMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
873 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
877 return isEqualWithoutConsideringStrInternal(other,prec,tmp);
880 bool MEDCouplingIMesh::isEqualWithoutConsideringStrInternal(const MEDCouplingMesh *other, double prec, std::string& reason) const
882 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
885 if(_space_dim!=otherC->_space_dim)
887 std::ostringstream oss;
888 oss << "The spaceDimension of this (" << _space_dim << ") is not equal to those of other (" << otherC->_space_dim << ") !";
891 checkSpaceDimension();
892 for(int i=0;i<_space_dim;i++)
894 if(fabs(_origin[i]-otherC->_origin[i])>prec)
896 std::ostringstream oss;
897 oss << "The origin of this and other differs at " << i << " !";
902 for(int i=0;i<_space_dim;i++)
904 if(fabs(_dxyz[i]-otherC->_dxyz[i])>prec)
906 std::ostringstream oss;
907 oss << "The delta of this and other differs at " << i << " !";
912 for(int i=0;i<_space_dim;i++)
914 if(_structure[i]!=otherC->_structure[i])
916 std::ostringstream oss;
917 oss << "The structure of this and other differs at " << i << " !";
925 void MEDCouplingIMesh::checkDeepEquivalWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
926 DataArrayInt *&cellCor, DataArrayInt *&nodeCor) const
928 if(!isEqualWithoutConsideringStr(other,prec))
929 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalWith : Meshes are not the same !");
933 * Nothing is done here (except to check that the other is a ParaMEDMEM::MEDCouplingIMesh instance too).
934 * The user intend that the nodes are the same, so by construction of ParaMEDMEM::MEDCouplingIMesh, \a this and \a other are the same !
936 void MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
937 DataArrayInt *&cellCor) const
939 if(!isEqualWithoutConsideringStr(other,prec))
940 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith : Meshes are not the same !");
943 void MEDCouplingIMesh::checkCoherency() const
945 checkSpaceDimension();
946 for(int i=0;i<_space_dim;i++)
949 std::ostringstream oss; oss << "MEDCouplingIMesh::checkCoherency : On axis " << i << "/" << _space_dim << ", number of nodes is equal to " << _structure[i] << " ! must be >=1 !";
950 throw INTERP_KERNEL::Exception(oss.str().c_str());
954 void MEDCouplingIMesh::checkCoherency1(double eps) const
959 void MEDCouplingIMesh::checkCoherency2(double eps) const
961 checkCoherency1(eps);
964 void MEDCouplingIMesh::getNodeGridStructure(int *res) const
966 checkSpaceDimension();
967 std::copy(_structure,_structure+_space_dim,res);
970 std::vector<int> MEDCouplingIMesh::getNodeGridStructure() const
972 checkSpaceDimension();
973 std::vector<int> ret(_structure,_structure+_space_dim);
977 MEDCouplingStructuredMesh *MEDCouplingIMesh::buildStructuredSubPart(const std::vector< std::pair<int,int> >& cellPart) const
980 int dim(getSpaceDimension());
981 if(dim!=(int)cellPart.size())
983 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : the space dimension is " << dim << " and cell part size is " << cellPart.size() << " !";
984 throw INTERP_KERNEL::Exception(oss.str().c_str());
986 double retOrigin[3]={0.,0.,0.};
987 int retStruct[3]={0,0,0};
988 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(dynamic_cast<MEDCouplingIMesh *>(deepCpy()));
989 for(int i=0;i<dim;i++)
991 int startNode(cellPart[i].first),endNode(cellPart[i].second+1);
992 int myDelta(endNode-startNode);
993 if(startNode<0 || startNode>=_structure[i])
995 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : At dimension #" << i << " the start node id is " << startNode << " it should be in [0," << _structure[i] << ") !";
996 throw INTERP_KERNEL::Exception(oss.str().c_str());
998 if(myDelta<0 || myDelta>_structure[i])
1000 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;
1001 throw INTERP_KERNEL::Exception(oss.str().c_str());
1003 retOrigin[i]=_origin[i]+startNode*_dxyz[i];
1004 retStruct[i]=myDelta;
1006 ret->setNodeStruct(retStruct,retStruct+dim);
1007 ret->setOrigin(retOrigin,retOrigin+dim);
1008 ret->checkCoherency();
1013 * Return the space dimension of \a this.
1015 int MEDCouplingIMesh::getSpaceDimension() const
1020 void MEDCouplingIMesh::getCoordinatesOfNode(int nodeId, std::vector<double>& coo) const
1023 int spaceDim(getSpaceDimension());
1024 getSplitNodeValues(tmp);
1026 GetPosFromId(nodeId,spaceDim,tmp,tmp2);
1027 for(int j=0;j<spaceDim;j++)
1028 coo.push_back(_origin[j]+_dxyz[j]*tmp2[j]);
1031 std::string MEDCouplingIMesh::simpleRepr() const
1033 std::ostringstream ret;
1034 ret << "Image grid with name : \"" << getName() << "\"\n";
1035 ret << "Description of mesh : \"" << getDescription() << "\"\n";
1037 double tt(getTime(tmpp1,tmpp2));
1038 int spaceDim(_space_dim);
1039 ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
1040 ret << "Iteration : " << tmpp1 << " Order : " << tmpp2 << "\n";
1041 ret << "Space dimension : " << spaceDim << "\n";
1042 if(spaceDim<0 || spaceDim>3)
1044 ret << "The nodal structure is : "; std::copy(_structure,_structure+spaceDim,std::ostream_iterator<int>(ret," ")); ret << "\n";
1045 ret << "The origin position is [" << _axis_unit << "]: ";
1046 std::copy(_origin,_origin+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1047 ret << "The intervals along axis are : ";
1048 std::copy(_dxyz,_dxyz+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1052 std::string MEDCouplingIMesh::advancedRepr() const
1054 return simpleRepr();
1057 void MEDCouplingIMesh::getBoundingBox(double *bbox) const
1060 int dim(getSpaceDimension());
1061 for(int idim=0; idim<dim; idim++)
1063 bbox[2*idim]=_origin[idim];
1064 int coeff(_structure[idim]);
1065 if(_structure[idim]<0)
1067 std::ostringstream oss; oss << "MEDCouplingIMesh::getBoundingBox : on axis #" << idim << " number of nodes in structure is < 0 !";
1068 throw INTERP_KERNEL::Exception(oss.str().c_str());
1070 if(_structure[idim]>1)
1071 coeff=_structure[idim]-1;
1072 bbox[2*idim+1]=_origin[idim]+_dxyz[idim]*coeff;
1077 * Returns a new MEDCouplingFieldDouble containing volumes of cells constituting \a this
1079 * For 1D cells, the returned field contains lengths.<br>
1080 * For 2D cells, the returned field contains areas.<br>
1081 * For 3D cells, the returned field contains volumes.
1082 * \param [in] isAbs - a not used parameter.
1083 * \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble on cells
1084 * and one time . The caller is to delete this field using decrRef() as it is no
1087 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureField(bool isAbs) const
1090 std::string name="MeasureOfMesh_";
1092 int nbelem(getNumberOfCells());
1093 MEDCouplingFieldDouble *field(MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME));
1094 field->setName(name);
1095 DataArrayDouble* array(DataArrayDouble::New());
1096 array->alloc(nbelem,1);
1097 array->fillWithValue(getMeasureOfAnyCell());
1098 field->setArray(array) ;
1100 field->setMesh(const_cast<MEDCouplingIMesh *>(this));
1101 field->synchronizeTimeWithMesh();
1106 * not implemented yet !
1108 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureFieldOnNode(bool isAbs) const
1110 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::getMeasureFieldOnNode : not implemented yet !");
1114 int MEDCouplingIMesh::getCellContainingPoint(const double *pos, double eps) const
1116 int dim(getSpaceDimension()),ret(0),coeff(1);
1117 for(int i=0;i<dim;i++)
1119 int nbOfCells(_structure[i]-1);
1121 int tmp((int)((ref-_origin[i])/_dxyz[i]));
1122 if(tmp>=0 && tmp<nbOfCells)
1133 void MEDCouplingIMesh::rotate(const double *center, const double *vector, double angle)
1135 throw INTERP_KERNEL::Exception("No rotation available on IMesh : Traduce it to unstructured mesh to apply it !");
1139 * Translates all nodes of \a this mesh by a given vector. Actually, it adds each
1140 * component of the \a vector to all node coordinates of a corresponding axis.
1141 * \param [in] vector - the translation vector whose size must be not less than \a
1142 * this->getSpaceDimension().
1144 void MEDCouplingIMesh::translate(const double *vector)
1146 checkSpaceDimension();
1147 int dim(getSpaceDimension());
1148 std::transform(_origin,_origin+dim,vector,_origin,std::plus<double>());
1153 * Applies scaling transformation to all nodes of \a this mesh.
1154 * \param [in] point - coordinates of a scaling center. This array is to be of
1155 * size \a this->getSpaceDimension() at least.
1156 * \param [in] factor - a scale factor.
1158 void MEDCouplingIMesh::scale(const double *point, double factor)
1160 checkSpaceDimension();
1161 int dim(getSpaceDimension());
1162 std::transform(_origin,_origin+dim,point,_origin,std::minus<double>());
1163 std::transform(_origin,_origin+dim,_origin,std::bind2nd(std::multiplies<double>(),factor));
1164 std::transform(_dxyz,_dxyz+dim,_dxyz,std::bind2nd(std::multiplies<double>(),factor));
1165 std::transform(_origin,_origin+dim,point,_origin,std::plus<double>());
1169 MEDCouplingMesh *MEDCouplingIMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
1171 //not implemented yet !
1176 * Returns a new DataArrayDouble holding coordinates of all nodes of \a this mesh.
1177 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1178 * this->getNumberOfNodes() tuples per \a this->getSpaceDimension()
1179 * components. The caller is to delete this array using decrRef() as it is
1182 DataArrayDouble *MEDCouplingIMesh::getCoordinatesAndOwner() const
1185 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1186 int spaceDim(getSpaceDimension()),nbNodes(getNumberOfNodes());
1187 ret->alloc(nbNodes,spaceDim);
1188 double *pt(ret->getPointer());
1189 ret->setInfoOnComponents(buildInfoOnComponents());
1191 getSplitNodeValues(tmp);
1192 for(int i=0;i<nbNodes;i++)
1194 GetPosFromId(i,spaceDim,tmp,tmp2);
1195 for(int j=0;j<spaceDim;j++)
1196 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+_origin[j];
1202 * Returns a new DataArrayDouble holding barycenters of all cells. The barycenter is
1203 * computed by averaging coordinates of cell nodes.
1204 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1205 * this->getNumberOfCells() tuples per \a this->getSpaceDimension()
1206 * components. The caller is to delete this array using decrRef() as it is
1209 DataArrayDouble *MEDCouplingIMesh::getBarycenterAndOwner() const
1212 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1213 int spaceDim(getSpaceDimension()),nbCells(getNumberOfCells()),tmp[3],tmp2[3];
1214 ret->alloc(nbCells,spaceDim);
1215 double *pt(ret->getPointer()),shiftOrigin[3];
1216 std::transform(_dxyz,_dxyz+spaceDim,shiftOrigin,std::bind2nd(std::multiplies<double>(),0.5));
1217 std::transform(_origin,_origin+spaceDim,shiftOrigin,shiftOrigin,std::plus<double>());
1218 getSplitCellValues(tmp);
1219 ret->setInfoOnComponents(buildInfoOnComponents());
1220 for(int i=0;i<nbCells;i++)
1222 GetPosFromId(i,spaceDim,tmp,tmp2);
1223 for(int j=0;j<spaceDim;j++)
1224 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+shiftOrigin[j];
1229 DataArrayDouble *MEDCouplingIMesh::computeIsoBarycenterOfNodesPerCell() const
1231 return MEDCouplingIMesh::getBarycenterAndOwner();
1234 void MEDCouplingIMesh::renumberCells(const int *old2NewBg, bool check)
1236 throw INTERP_KERNEL::Exception("Functionnality of renumbering cell not available for IMesh !");
1239 void MEDCouplingIMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
1242 double time(getTime(it,order));
1245 littleStrings.clear();
1246 littleStrings.push_back(getName());
1247 littleStrings.push_back(getDescription());
1248 littleStrings.push_back(getTimeUnit());
1249 littleStrings.push_back(getAxisUnit());
1250 tinyInfo.push_back(it);
1251 tinyInfo.push_back(order);
1252 tinyInfo.push_back(_space_dim);
1253 tinyInfo.insert(tinyInfo.end(),_structure,_structure+3);
1254 tinyInfoD.push_back(time);
1255 tinyInfoD.insert(tinyInfoD.end(),_dxyz,_dxyz+3);
1256 tinyInfoD.insert(tinyInfoD.end(),_origin,_origin+3);
1259 void MEDCouplingIMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
1265 void MEDCouplingIMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
1267 a1=DataArrayInt::New();
1269 a2=DataArrayDouble::New();
1273 void MEDCouplingIMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
1274 const std::vector<std::string>& littleStrings)
1276 setName(littleStrings[0]);
1277 setDescription(littleStrings[1]);
1278 setTimeUnit(littleStrings[2]);
1279 setAxisUnit(littleStrings[3]);
1280 setTime(tinyInfoD[0],tinyInfo[0],tinyInfo[1]);
1281 _space_dim=tinyInfo[2];
1282 _structure[0]=tinyInfo[3]; _structure[1]=tinyInfo[4]; _structure[2]=tinyInfo[5];
1283 _dxyz[0]=tinyInfoD[1]; _dxyz[1]=tinyInfoD[2]; _dxyz[2]=tinyInfoD[3];
1284 _origin[0]=tinyInfoD[4]; _origin[1]=tinyInfoD[5]; _origin[2]=tinyInfoD[6];
1288 void MEDCouplingIMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
1291 std::ostringstream extent,origin,spacing;
1292 for(int i=0;i<3;i++)
1295 { extent << "0 " << _structure[i]-1 << " "; origin << _origin[i] << " "; spacing << _dxyz[i] << " "; }
1297 { extent << "0 0 "; origin << "0 "; spacing << "0 "; }
1299 ofs << " <" << getVTKDataSetType() << " WholeExtent=\"" << extent.str() << "\" Origin=\"" << origin.str() << "\" Spacing=\"" << spacing.str() << "\">\n";
1300 ofs << " <Piece Extent=\"" << extent.str() << "\">\n";
1301 ofs << " <PointData>\n" << pointData << std::endl;
1302 ofs << " </PointData>\n";
1303 ofs << " <CellData>\n" << cellData << std::endl;
1304 ofs << " </CellData>\n";
1305 ofs << " <Coordinates>\n";
1306 ofs << " </Coordinates>\n";
1307 ofs << " </Piece>\n";
1308 ofs << " </" << getVTKDataSetType() << ">\n";
1311 void MEDCouplingIMesh::reprQuickOverview(std::ostream& stream) const
1313 stream << "MEDCouplingIMesh C++ instance at " << this << ". Name : \"" << getName() << "\". Space dimension : " << _space_dim << ".";
1314 if(_space_dim<0 || _space_dim>3)
1317 std::ostringstream stream0,stream1;
1318 int nbNodes(1),nbCells(0);
1320 for(int i=0;i<_space_dim;i++)
1323 int tmpNodes(_structure[i]);
1324 stream1 << "- Axis " << tmp << " : " << tmpNodes << " nodes (orig=" << _origin[i] << ", inter=" << _dxyz[i] << ").";
1326 stream1 << std::endl;
1332 nbCells=nbCells==0?tmpNodes-1:nbCells*(tmpNodes-1);
1336 stream0 << "Number of cells : " << nbCells << ", Number of nodes : " << nbNodes;
1337 stream << stream0.str();
1339 stream << std::endl;
1341 stream << stream1.str();
1344 std::string MEDCouplingIMesh::getVTKDataSetType() const
1346 return std::string("ImageData");
1349 std::vector<std::string> MEDCouplingIMesh::buildInfoOnComponents() const
1351 checkSpaceDimension();
1352 int dim(getSpaceDimension());
1353 std::vector<std::string> ret(dim);
1354 for(int i=0;i<dim;i++)
1356 std::ostringstream oss;
1357 char tmp('X'+i); oss << tmp;
1358 ret[i]=DataArray::BuildInfoFromVarAndUnit(oss.str(),_axis_unit);
1363 void MEDCouplingIMesh::checkSpaceDimension() const
1365 CheckSpaceDimension(_space_dim);
1368 void MEDCouplingIMesh::CheckSpaceDimension(int spaceDim)
1370 if(spaceDim<0 || spaceDim>3)
1371 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CheckSpaceDimension : input spaceDim must be in [0,1,2,3] !");
1374 int MEDCouplingIMesh::FindIntRoot(int val, int order)
1379 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : input val is < 0 ! Not possible to compute a root !");
1382 if(order!=2 && order!=3)
1383 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the order available are 0,1,2 or 3 !");
1384 double valf((double)val);
1387 double retf(sqrt(valf));
1390 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect square root !");
1395 double retf(std::pow(val,0.3333333333333333));
1396 int ret((int)retf),ret2(ret+1);
1397 if(ret*ret*ret!=val && ret2*ret2*ret2!=val)
1398 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect cublic root !");
1399 if(ret*ret*ret==val)