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(nbTuplesFine%nbOfTuplesInFineExp!=0)
298 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : Invalid nb of tuples in fine DataArray regarding its structure !");
299 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
300 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
302 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarse : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
303 throw INTERP_KERNEL::Exception(oss.str().c_str());
305 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
306 double *outPtr(coarseDA->getPointer());
307 const double *inPtr(fineDA->begin());
309 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
314 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
315 for(int i=0;i<dims[0];i++)
317 double *loc(outPtr+(offset+i)*nbCompo);
318 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
321 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
323 std::copy(inPtr,inPtr+nbCompo,loc);
330 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact1(facts[1]),fact0(facts[0]);
331 for(int j=0;j<dims[1];j++)
333 for(int jfact=0;jfact<fact1;jfact++)
335 for(int i=0;i<dims[0];i++)
337 double *loc(outPtr+(kk+i)*nbCompo);
338 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
340 if(jfact!=0 || ifact!=0)
341 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
343 std::copy(inPtr,inPtr+nbCompo,loc);
353 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]);
354 for(int k=0;k<dims[2];k++)
356 for(int kfact=0;kfact<fact2;kfact++)
358 for(int j=0;j<dims[1];j++)
360 for(int jfact=0;jfact<fact1;jfact++)
362 for(int i=0;i<dims[0];i++)
364 double *loc(outPtr+(kk+i+j*coarseSt[0])*nbCompo);
365 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
367 if(kfact!=0 || jfact!=0 || ifact!=0)
368 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
370 std::copy(inPtr,inPtr+nbCompo,loc);
376 kk+=coarseSt[0]*coarseSt[1];
381 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarse : only dimensions 1, 2 and 3 supported !");
386 * This static method is useful to condense field on cells of a MEDCouplingIMesh instance coming from a refinement ( MEDCouplingIMesh::refineWithFactor for example)
387 * to a coarse MEDCouplingIMesh instance. So this method can be seen as a specialization in P0P0 conservative interpolation non overlaping from fine image mesh
388 * 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.
390 * \param [in] coarseSt The cell structure of coarse mesh.
391 * \param [in] fineDA The DataArray containing the cell field on uniformly refined mesh
392 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
393 * \param [in] facts The refinement coefficient per axis.
394 * \param [in,out] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
395 * \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.
397 * \sa CondenseFineToCoarse,SpreadCoarseToFineGhost
399 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)
402 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : ghost level has to be >= 0 !");
403 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
404 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : All input vectors (dimension) must have the same size !");
405 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
406 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the parameters 1 or 3 are NULL or not allocated !");
407 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
408 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
409 int nbCompo(fineDA->getNumberOfComponents());
410 if(coarseDA->getNumberOfComponents()!=nbCompo)
411 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : the number of components of fine DA and coarse one mismatches !");
412 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
413 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) !");
414 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
416 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples in coarse DataArray having " << coarseDA->getNumberOfTuples() << " !";
417 throw INTERP_KERNEL::Exception(oss.str().c_str());
420 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
421 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
422 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
423 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
424 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
426 std::ostringstream oss; oss << "MEDCouplingIMesh::CondenseFineToCoarseGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
427 throw INTERP_KERNEL::Exception(oss.str().c_str());
430 double *outPtr(coarseDA->getPointer());
431 const double *inPtr(fineDA->begin());
433 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
438 int offset(fineLocInCoarse[0].first+ghostSize),fact0(facts[0]);
439 inPtr+=ghostSize*nbCompo;
440 for(int i=0;i<dims[0];i++)
442 double *loc(outPtr+(offset+i)*nbCompo);
443 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
446 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
448 std::copy(inPtr,inPtr+nbCompo,loc);
455 int nxwg(coarseSt[0]+2*ghostSize);
456 int kk(fineLocInCoarse[0].first+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize)),fact1(facts[1]),fact0(facts[0]);
457 inPtr+=(dims[0]*fact0+2*ghostSize)*nbCompo;
458 for(int j=0;j<dims[1];j++)
460 for(int jfact=0;jfact<fact1;jfact++)
462 inPtr+=ghostSize*nbCompo;
463 for(int i=0;i<dims[0];i++)
465 double *loc(outPtr+(kk+i)*nbCompo);
466 for(int ifact=0;ifact<fact0;ifact++,inPtr+=nbCompo)
468 if(jfact!=0 || ifact!=0)
469 std::transform(inPtr,inPtr+nbCompo,loc,loc,std::plus<double>());
471 std::copy(inPtr,inPtr+nbCompo,loc);
474 inPtr+=ghostSize*nbCompo;
481 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CondenseFineToCoarseGhost : only dimensions 1, 2 supported !");
486 * This method spreads the values of coarse data \a coarseDA into \a fineDA.
488 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
489 * \param [in] coarseSt The cell structure of coarse mesh.
490 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
491 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
492 * \param [in] facts The refinement coefficient per axis.
493 * \sa SpreadCoarseToFineGhost, CondenseFineToCoarse
495 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)
497 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
498 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : All input vectors (dimension) must have the same size !");
499 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
500 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the parameters 1 or 3 are NULL or not allocated !");
501 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseSt)),nbOfTuplesInFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(fineLocInCoarse));
502 int nbCompo(fineDA->getNumberOfComponents());
503 if(coarseDA->getNumberOfComponents()!=nbCompo)
504 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : the number of components of fine DA and coarse one mismatches !");
505 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
506 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) !");
507 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
509 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
510 throw INTERP_KERNEL::Exception(oss.str().c_str());
512 int nbTuplesFine(fineDA->getNumberOfTuples());
513 if(nbTuplesFine%nbOfTuplesInFineExp!=0)
514 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : Invalid nb of tuples in fine DataArray regarding its structure !");
515 int fact(std::accumulate(facts.begin(),facts.end(),1,std::multiplies<int>()));
516 if(nbTuplesFine!=fact*nbOfTuplesInFineExp)
518 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFine : Invalid number of tuples (" << nbTuplesFine << ") of fine dataarray is invalid ! Must be " << fact*nbOfTuplesInFineExp << "!";
519 throw INTERP_KERNEL::Exception(oss.str().c_str());
521 // to improve use jump-iterator. Factorizes with SwitchOnIdsFrom BuildExplicitIdsFrom
522 double *outPtr(fineDA->getPointer());
523 const double *inPtr(coarseDA->begin());
525 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
530 int offset(fineLocInCoarse[0].first),fact0(facts[0]);
531 for(int i=0;i<dims[0];i++)
533 const double *loc(inPtr+(offset+i)*nbCompo);
534 for(int ifact=0;ifact<fact0;ifact++)
535 outPtr=std::copy(loc,loc+nbCompo,outPtr);
541 int kk(fineLocInCoarse[0].first+coarseSt[0]*fineLocInCoarse[1].first),fact0(facts[0]),fact1(facts[1]);
542 for(int j=0;j<dims[1];j++)
544 for(int jfact=0;jfact<fact1;jfact++)
546 for(int i=0;i<dims[0];i++)
548 const double *loc(inPtr+(kk+i)*nbCompo);
549 for(int ifact=0;ifact<fact0;ifact++)
550 outPtr=std::copy(loc,loc+nbCompo,outPtr);
559 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]);
560 for(int k=0;k<dims[2];k++)
562 for(int kfact=0;kfact<fact2;kfact++)
564 for(int j=0;j<dims[1];j++)
566 for(int jfact=0;jfact<fact1;jfact++)
568 for(int i=0;i<dims[0];i++)
570 const double *loc(inPtr+(kk+i+j*coarseSt[0])*nbCompo);
571 for(int ifact=0;ifact<fact0;ifact++)
572 outPtr=std::copy(loc,loc+nbCompo,outPtr);
577 kk+=coarseSt[0]*coarseSt[1];
582 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFine : only dimensions 1, 2 and 3 supported !");
587 * This method spreads the values of coarse data \a coarseDA into \a fineDA.
589 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
590 * \param [in] coarseSt The cell structure of coarse mesh.
591 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
592 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
593 * \param [in] facts The refinement coefficient per axis.
594 * \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.
596 * \sa CondenseFineToCoarse,SpreadCoarseToFineGhostZone
598 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)
601 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : ghost level has to be >= 0 !");
602 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
603 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : All input vectors (dimension) must have the same size !");
604 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
605 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the parameters 1 or 3 are NULL or not allocated !");
606 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
607 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
608 int nbCompo(fineDA->getNumberOfComponents());
609 if(coarseDA->getNumberOfComponents()!=nbCompo)
610 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : the number of components of fine DA and coarse one mismatches !");
611 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
612 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) !");
613 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
615 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
616 throw INTERP_KERNEL::Exception(oss.str().c_str());
619 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
620 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
621 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
622 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
623 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
625 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhost : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
626 throw INTERP_KERNEL::Exception(oss.str().c_str());
629 double *outPtr(fineDA->getPointer());
630 const double *inPtr(coarseDA->begin());
632 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
637 int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
638 for(int i=0;i<ghostSize;i++)
639 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
640 offset=fineLocInCoarse[0].first+ghostSize;
641 for(int i=0;i<dims[0];i++)
643 const double *loc(inPtr+(offset+i)*nbCompo);
644 for(int ifact=0;ifact<fact0;ifact++)
645 outPtr=std::copy(loc,loc+nbCompo,outPtr);
647 offset=fineLocInCoarse[0].second+ghostSize;
648 for(int i=0;i<ghostSize;i++)
649 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
654 int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
655 int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
656 for(int jg=0;jg<ghostSize;jg++)
658 for(int ig=0;ig<ghostSize;ig++)
659 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
661 for(int ig=0;ig<dims[0];ig++,kk0++)
662 for(int ifact=0;ifact<fact0;ifact++)
663 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
664 for(int ik=0;ik<ghostSize;ik++)
665 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
667 for(int j=0;j<dims[1];j++)
669 kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
670 for(int jfact=0;jfact<fact1;jfact++)
672 for(int ig=0;ig<ghostSize;ig++)
673 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
674 int kk0(kk+1);//1 not ghost. We make the hypothesis that factors is >= ghostlev
675 for(int i=0;i<dims[0];i++,kk0++)
677 const double *loc(inPtr+kk0*nbCompo);
678 for(int ifact=0;ifact<fact0;ifact++)
679 outPtr=std::copy(loc,loc+nbCompo,outPtr);
681 for(int ig=0;ig<ghostSize;ig++)
682 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
685 kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
686 for(int jg=0;jg<ghostSize;jg++)
688 for(int ig=0;ig<ghostSize;ig++)
689 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
691 for(int ig=0;ig<dims[0];ig++,kk0++)
692 for(int ifact=0;ifact<fact0;ifact++)
693 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
694 for(int ik=0;ik<ghostSize;ik++)
695 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
700 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhost : only dimensions 1, 2 supported !");
705 * 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).
707 * \param [in] coarseDA The DataArrayDouble corresponding to the a cell field of a coarse mesh whose cell structure is defined by \a coarseSt.
708 * \param [in] coarseSt The cell structure of coarse mesh.
709 * \param [in,out] fineDA The DataArray containing the cell field on uniformly refined mesh
710 * \param [in] fineLocInCoarse The cell localization of refined mesh into the coarse one.
711 * \param [in] facts The refinement coefficient per axis.
712 * \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.
714 * \sa SpreadCoarseToFineGhost
716 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)
719 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : ghost level has to be >= 0 !");
720 if(coarseSt.size()!=fineLocInCoarse.size() || coarseSt.size()!=facts.size())
721 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : All input vectors (dimension) must have the same size !");
722 if(!coarseDA || !coarseDA->isAllocated() || !fineDA || !fineDA->isAllocated())
723 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the parameters 1 or 3 are NULL or not allocated !");
724 std::vector<int> coarseStG(coarseSt.size()); std::transform(coarseSt.begin(),coarseSt.end(),coarseStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
725 int meshDim((int)coarseSt.size()),nbOfTuplesInCoarseExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(coarseStG));
726 int nbCompo(fineDA->getNumberOfComponents());
727 if(coarseDA->getNumberOfComponents()!=nbCompo)
728 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : the number of components of fine DA and coarse one mismatches !");
729 if(meshDim!=(int)fineLocInCoarse.size() || meshDim!=(int)facts.size())
730 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) !");
731 if(coarseDA->getNumberOfTuples()!=nbOfTuplesInCoarseExp)
733 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbOfTuplesInCoarseExp << " tuples having " << coarseDA->getNumberOfTuples() << " !";
734 throw INTERP_KERNEL::Exception(oss.str().c_str());
737 std::vector<int> fineStG(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
738 std::transform(fineStG.begin(),fineStG.end(),facts.begin(),fineStG.begin(),std::multiplies<int>());
739 std::transform(fineStG.begin(),fineStG.end(),fineStG.begin(),std::bind2nd(std::plus<int>(),2*ghostSize));
740 int nbTuplesFine(fineDA->getNumberOfTuples()),nbTuplesFineExp(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(fineStG));
741 if(fineDA->getNumberOfTuples()!=nbTuplesFineExp)
743 std::ostringstream oss; oss << "MEDCouplingIMesh::SpreadCoarseToFineGhostZone : Expecting " << nbTuplesFineExp << " tuples in fine DataArray having " << nbTuplesFine << " !";
744 throw INTERP_KERNEL::Exception(oss.str().c_str());
747 double *outPtr(fineDA->getPointer());
748 const double *inPtr(coarseDA->begin());
750 std::vector<int> dims(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(fineLocInCoarse));
755 int offset(fineLocInCoarse[0].first+ghostSize-1),fact0(facts[0]);//offset is always >=0 thanks to the fact that ghostSize>=1 !
756 for(int i=0;i<ghostSize;i++)
757 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
758 outPtr+=nbCompo*fact0*dims[0];
759 offset=fineLocInCoarse[0].second+ghostSize;
760 for(int i=0;i<ghostSize;i++)
761 outPtr=std::copy(inPtr+offset*nbCompo,inPtr+(offset+1)*nbCompo,outPtr);
766 int nxwg(coarseSt[0]+2*ghostSize),fact0(facts[0]),fact1(facts[1]);
767 int kk(fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].first+ghostSize-1));//kk is always >=0 thanks to the fact that ghostSize>=1 !
768 for(int jg=0;jg<ghostSize;jg++)
770 for(int ig=0;ig<ghostSize;ig++)
771 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
773 for(int ig=0;ig<dims[0];ig++,kk0++)
774 for(int ifact=0;ifact<fact0;ifact++)
775 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
776 for(int ik=0;ik<ghostSize;ik++)
777 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
779 for(int j=0;j<dims[1];j++)
781 kk=fineLocInCoarse[0].first-1+ghostSize+nxwg*(fineLocInCoarse[1].first+ghostSize+j);
782 for(int jfact=0;jfact<fact1;jfact++)
784 for(int ig=0;ig<ghostSize;ig++)
785 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
786 int kk0(kk+1+dims[0]);//1 not ghost. We make the hypothesis that factors is >= ghostlev
787 outPtr+=fact0*nbCompo*dims[0];
788 for(int ig=0;ig<ghostSize;ig++)
789 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
792 kk=fineLocInCoarse[0].first+ghostSize-1+nxwg*(fineLocInCoarse[1].second+ghostSize);
793 for(int jg=0;jg<ghostSize;jg++)
795 for(int ig=0;ig<ghostSize;ig++)
796 outPtr=std::copy(inPtr+kk*nbCompo,inPtr+(kk+1)*nbCompo,outPtr);
798 for(int ig=0;ig<dims[0];ig++,kk0++)
799 for(int ifact=0;ifact<fact0;ifact++)
800 outPtr=std::copy(inPtr+(kk0)*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
801 for(int ik=0;ik<ghostSize;ik++)
802 outPtr=std::copy(inPtr+kk0*nbCompo,inPtr+(kk0+1)*nbCompo,outPtr);
807 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::SpreadCoarseToFineGhostZone : only dimensions 1, 2 supported !");
811 void MEDCouplingIMesh::setSpaceDimension(int spaceDim)
813 if(spaceDim==_space_dim)
815 CheckSpaceDimension(spaceDim);
820 void MEDCouplingIMesh::updateTime() const
824 std::size_t MEDCouplingIMesh::getHeapMemorySizeWithoutChildren() const
826 return MEDCouplingStructuredMesh::getHeapMemorySizeWithoutChildren();
829 std::vector<const BigMemoryObject *> MEDCouplingIMesh::getDirectChildren() const
831 return std::vector<const BigMemoryObject *>();
835 * This method copyies all tiny strings from other (name and components name).
836 * @throw if other and this have not same mesh type.
838 void MEDCouplingIMesh::copyTinyStringsFrom(const MEDCouplingMesh *other)
840 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
842 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::copyTinyStringsFrom : meshes have not same type !");
843 MEDCouplingStructuredMesh::copyTinyStringsFrom(other);
847 bool MEDCouplingIMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
850 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::isEqualIfNotWhy : input other pointer is null !");
851 const MEDCouplingIMesh *otherC(dynamic_cast<const MEDCouplingIMesh *>(other));
854 reason="mesh given in input is not castable in MEDCouplingIMesh !";
857 if(!MEDCouplingStructuredMesh::isEqualIfNotWhy(other,prec,reason))
859 if(!isEqualWithoutConsideringStrInternal(otherC,prec,reason))
861 if(_axis_unit!=otherC->_axis_unit)
863 reason="The units of axis are not the same !";
869 bool MEDCouplingIMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
871 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
875 return isEqualWithoutConsideringStrInternal(other,prec,tmp);
878 bool MEDCouplingIMesh::isEqualWithoutConsideringStrInternal(const MEDCouplingMesh *other, double prec, std::string& reason) const
880 const MEDCouplingIMesh *otherC=dynamic_cast<const MEDCouplingIMesh *>(other);
883 if(_space_dim!=otherC->_space_dim)
885 std::ostringstream oss;
886 oss << "The spaceDimension of this (" << _space_dim << ") is not equal to those of other (" << otherC->_space_dim << ") !";
889 checkSpaceDimension();
890 for(int i=0;i<_space_dim;i++)
892 if(fabs(_origin[i]-otherC->_origin[i])>prec)
894 std::ostringstream oss;
895 oss << "The origin of this and other differs at " << i << " !";
900 for(int i=0;i<_space_dim;i++)
902 if(fabs(_dxyz[i]-otherC->_dxyz[i])>prec)
904 std::ostringstream oss;
905 oss << "The delta of this and other differs at " << i << " !";
910 for(int i=0;i<_space_dim;i++)
912 if(_structure[i]!=otherC->_structure[i])
914 std::ostringstream oss;
915 oss << "The structure of this and other differs at " << i << " !";
923 void MEDCouplingIMesh::checkDeepEquivalWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
924 DataArrayInt *&cellCor, DataArrayInt *&nodeCor) const
926 if(!isEqualWithoutConsideringStr(other,prec))
927 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalWith : Meshes are not the same !");
931 * Nothing is done here (except to check that the other is a ParaMEDMEM::MEDCouplingIMesh instance too).
932 * The user intend that the nodes are the same, so by construction of ParaMEDMEM::MEDCouplingIMesh, \a this and \a other are the same !
934 void MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith(const MEDCouplingMesh *other, int cellCompPol, double prec,
935 DataArrayInt *&cellCor) const
937 if(!isEqualWithoutConsideringStr(other,prec))
938 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::checkDeepEquivalOnSameNodesWith : Meshes are not the same !");
941 void MEDCouplingIMesh::checkCoherency() const
943 checkSpaceDimension();
944 for(int i=0;i<_space_dim;i++)
947 std::ostringstream oss; oss << "MEDCouplingIMesh::checkCoherency : On axis " << i << "/" << _space_dim << ", number of nodes is equal to " << _structure[i] << " ! must be >=1 !";
948 throw INTERP_KERNEL::Exception(oss.str().c_str());
952 void MEDCouplingIMesh::checkCoherency1(double eps) const
957 void MEDCouplingIMesh::checkCoherency2(double eps) const
959 checkCoherency1(eps);
962 void MEDCouplingIMesh::getNodeGridStructure(int *res) const
964 checkSpaceDimension();
965 std::copy(_structure,_structure+_space_dim,res);
968 std::vector<int> MEDCouplingIMesh::getNodeGridStructure() const
970 checkSpaceDimension();
971 std::vector<int> ret(_structure,_structure+_space_dim);
975 MEDCouplingStructuredMesh *MEDCouplingIMesh::buildStructuredSubPart(const std::vector< std::pair<int,int> >& cellPart) const
978 int dim(getSpaceDimension());
979 if(dim!=(int)cellPart.size())
981 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : the space dimension is " << dim << " and cell part size is " << cellPart.size() << " !";
982 throw INTERP_KERNEL::Exception(oss.str().c_str());
984 double retOrigin[3]={0.,0.,0.};
985 int retStruct[3]={0,0,0};
986 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> ret(dynamic_cast<MEDCouplingIMesh *>(deepCpy()));
987 for(int i=0;i<dim;i++)
989 int startNode(cellPart[i].first),endNode(cellPart[i].second+1);
990 int myDelta(endNode-startNode);
991 if(startNode<0 || startNode>=_structure[i])
993 std::ostringstream oss; oss << "MEDCouplingIMesh::buildStructuredSubPart : At dimension #" << i << " the start node id is " << startNode << " it should be in [0," << _structure[i] << ") !";
994 throw INTERP_KERNEL::Exception(oss.str().c_str());
996 if(myDelta<0 || myDelta>_structure[i])
998 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;
999 throw INTERP_KERNEL::Exception(oss.str().c_str());
1001 retOrigin[i]=_origin[i]+startNode*_dxyz[i];
1002 retStruct[i]=myDelta;
1004 ret->setNodeStruct(retStruct,retStruct+dim);
1005 ret->setOrigin(retOrigin,retOrigin+dim);
1006 ret->checkCoherency();
1011 * Return the space dimension of \a this.
1013 int MEDCouplingIMesh::getSpaceDimension() const
1018 void MEDCouplingIMesh::getCoordinatesOfNode(int nodeId, std::vector<double>& coo) const
1021 int spaceDim(getSpaceDimension());
1022 getSplitNodeValues(tmp);
1024 GetPosFromId(nodeId,spaceDim,tmp,tmp2);
1025 for(int j=0;j<spaceDim;j++)
1026 coo.push_back(_origin[j]+_dxyz[j]*tmp2[j]);
1029 std::string MEDCouplingIMesh::simpleRepr() const
1031 std::ostringstream ret;
1032 ret << "Image grid with name : \"" << getName() << "\"\n";
1033 ret << "Description of mesh : \"" << getDescription() << "\"\n";
1035 double tt(getTime(tmpp1,tmpp2));
1036 int spaceDim(_space_dim);
1037 ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
1038 ret << "Iteration : " << tmpp1 << " Order : " << tmpp2 << "\n";
1039 ret << "Space dimension : " << spaceDim << "\n";
1040 if(spaceDim<0 || spaceDim>3)
1042 ret << "The nodal structure is : "; std::copy(_structure,_structure+spaceDim,std::ostream_iterator<int>(ret," ")); ret << "\n";
1043 ret << "The origin position is [" << _axis_unit << "]: ";
1044 std::copy(_origin,_origin+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1045 ret << "The intervals along axis are : ";
1046 std::copy(_dxyz,_dxyz+spaceDim,std::ostream_iterator<double>(ret," ")); ret << "\n";
1050 std::string MEDCouplingIMesh::advancedRepr() const
1052 return simpleRepr();
1055 void MEDCouplingIMesh::getBoundingBox(double *bbox) const
1058 int dim(getSpaceDimension());
1059 for(int idim=0; idim<dim; idim++)
1061 bbox[2*idim]=_origin[idim];
1062 int coeff(_structure[idim]);
1063 if(_structure[idim]<0)
1065 std::ostringstream oss; oss << "MEDCouplingIMesh::getBoundingBox : on axis #" << idim << " number of nodes in structure is < 0 !";
1066 throw INTERP_KERNEL::Exception(oss.str().c_str());
1068 if(_structure[idim]>1)
1069 coeff=_structure[idim]-1;
1070 bbox[2*idim+1]=_origin[idim]+_dxyz[idim]*coeff;
1075 * Returns a new MEDCouplingFieldDouble containing volumes of cells constituting \a this
1077 * For 1D cells, the returned field contains lengths.<br>
1078 * For 2D cells, the returned field contains areas.<br>
1079 * For 3D cells, the returned field contains volumes.
1080 * \param [in] isAbs - a not used parameter.
1081 * \return MEDCouplingFieldDouble * - a new instance of MEDCouplingFieldDouble on cells
1082 * and one time . The caller is to delete this field using decrRef() as it is no
1085 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureField(bool isAbs) const
1088 std::string name="MeasureOfMesh_";
1090 int nbelem(getNumberOfCells());
1091 MEDCouplingFieldDouble *field(MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME));
1092 field->setName(name);
1093 DataArrayDouble* array(DataArrayDouble::New());
1094 array->alloc(nbelem,1);
1095 array->fillWithValue(getMeasureOfAnyCell());
1096 field->setArray(array) ;
1098 field->setMesh(const_cast<MEDCouplingIMesh *>(this));
1099 field->synchronizeTimeWithMesh();
1104 * not implemented yet !
1106 MEDCouplingFieldDouble *MEDCouplingIMesh::getMeasureFieldOnNode(bool isAbs) const
1108 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::getMeasureFieldOnNode : not implemented yet !");
1112 int MEDCouplingIMesh::getCellContainingPoint(const double *pos, double eps) const
1114 int dim(getSpaceDimension()),ret(0),coeff(1);
1115 for(int i=0;i<dim;i++)
1117 int nbOfCells(_structure[i]-1);
1119 int tmp((int)((ref-_origin[i])/_dxyz[i]));
1120 if(tmp>=0 && tmp<nbOfCells)
1131 void MEDCouplingIMesh::rotate(const double *center, const double *vector, double angle)
1133 throw INTERP_KERNEL::Exception("No rotation available on IMesh : Traduce it to unstructured mesh to apply it !");
1137 * Translates all nodes of \a this mesh by a given vector. Actually, it adds each
1138 * component of the \a vector to all node coordinates of a corresponding axis.
1139 * \param [in] vector - the translation vector whose size must be not less than \a
1140 * this->getSpaceDimension().
1142 void MEDCouplingIMesh::translate(const double *vector)
1144 checkSpaceDimension();
1145 int dim(getSpaceDimension());
1146 std::transform(_origin,_origin+dim,vector,_origin,std::plus<double>());
1151 * Applies scaling transformation to all nodes of \a this mesh.
1152 * \param [in] point - coordinates of a scaling center. This array is to be of
1153 * size \a this->getSpaceDimension() at least.
1154 * \param [in] factor - a scale factor.
1156 void MEDCouplingIMesh::scale(const double *point, double factor)
1158 checkSpaceDimension();
1159 int dim(getSpaceDimension());
1160 std::transform(_origin,_origin+dim,point,_origin,std::minus<double>());
1161 std::transform(_origin,_origin+dim,_origin,std::bind2nd(std::multiplies<double>(),factor));
1162 std::transform(_dxyz,_dxyz+dim,_dxyz,std::bind2nd(std::multiplies<double>(),factor));
1163 std::transform(_origin,_origin+dim,point,_origin,std::plus<double>());
1167 MEDCouplingMesh *MEDCouplingIMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
1169 //not implemented yet !
1174 * Returns a new DataArrayDouble holding coordinates of all nodes of \a this mesh.
1175 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1176 * this->getNumberOfNodes() tuples per \a this->getSpaceDimension()
1177 * components. The caller is to delete this array using decrRef() as it is
1180 DataArrayDouble *MEDCouplingIMesh::getCoordinatesAndOwner() const
1183 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1184 int spaceDim(getSpaceDimension()),nbNodes(getNumberOfNodes());
1185 ret->alloc(nbNodes,spaceDim);
1186 double *pt(ret->getPointer());
1187 ret->setInfoOnComponents(buildInfoOnComponents());
1189 getSplitNodeValues(tmp);
1190 for(int i=0;i<nbNodes;i++)
1192 GetPosFromId(i,spaceDim,tmp,tmp2);
1193 for(int j=0;j<spaceDim;j++)
1194 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+_origin[j];
1200 * Returns a new DataArrayDouble holding barycenters of all cells. The barycenter is
1201 * computed by averaging coordinates of cell nodes.
1202 * \return DataArrayDouble * - a new instance of DataArrayDouble, of size \a
1203 * this->getNumberOfCells() tuples per \a this->getSpaceDimension()
1204 * components. The caller is to delete this array using decrRef() as it is
1207 DataArrayDouble *MEDCouplingIMesh::getBarycenterAndOwner() const
1210 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
1211 int spaceDim(getSpaceDimension()),nbCells(getNumberOfCells()),tmp[3],tmp2[3];
1212 ret->alloc(nbCells,spaceDim);
1213 double *pt(ret->getPointer()),shiftOrigin[3];
1214 std::transform(_dxyz,_dxyz+spaceDim,shiftOrigin,std::bind2nd(std::multiplies<double>(),0.5));
1215 std::transform(_origin,_origin+spaceDim,shiftOrigin,shiftOrigin,std::plus<double>());
1216 getSplitCellValues(tmp);
1217 ret->setInfoOnComponents(buildInfoOnComponents());
1218 for(int i=0;i<nbCells;i++)
1220 GetPosFromId(i,spaceDim,tmp,tmp2);
1221 for(int j=0;j<spaceDim;j++)
1222 pt[i*spaceDim+j]=_dxyz[j]*tmp2[j]+shiftOrigin[j];
1227 DataArrayDouble *MEDCouplingIMesh::computeIsoBarycenterOfNodesPerCell() const
1229 return MEDCouplingIMesh::getBarycenterAndOwner();
1232 void MEDCouplingIMesh::renumberCells(const int *old2NewBg, bool check)
1234 throw INTERP_KERNEL::Exception("Functionnality of renumbering cell not available for IMesh !");
1237 void MEDCouplingIMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
1240 double time(getTime(it,order));
1243 littleStrings.clear();
1244 littleStrings.push_back(getName());
1245 littleStrings.push_back(getDescription());
1246 littleStrings.push_back(getTimeUnit());
1247 littleStrings.push_back(getAxisUnit());
1248 tinyInfo.push_back(it);
1249 tinyInfo.push_back(order);
1250 tinyInfo.push_back(_space_dim);
1251 tinyInfo.insert(tinyInfo.end(),_structure,_structure+3);
1252 tinyInfoD.push_back(time);
1253 tinyInfoD.insert(tinyInfoD.end(),_dxyz,_dxyz+3);
1254 tinyInfoD.insert(tinyInfoD.end(),_origin,_origin+3);
1257 void MEDCouplingIMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
1263 void MEDCouplingIMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
1265 a1=DataArrayInt::New();
1267 a2=DataArrayDouble::New();
1271 void MEDCouplingIMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
1272 const std::vector<std::string>& littleStrings)
1274 setName(littleStrings[0]);
1275 setDescription(littleStrings[1]);
1276 setTimeUnit(littleStrings[2]);
1277 setAxisUnit(littleStrings[3]);
1278 setTime(tinyInfoD[0],tinyInfo[0],tinyInfo[1]);
1279 _space_dim=tinyInfo[2];
1280 _structure[0]=tinyInfo[3]; _structure[1]=tinyInfo[4]; _structure[2]=tinyInfo[5];
1281 _dxyz[0]=tinyInfoD[1]; _dxyz[1]=tinyInfoD[2]; _dxyz[2]=tinyInfoD[3];
1282 _origin[0]=tinyInfoD[4]; _origin[1]=tinyInfoD[5]; _origin[2]=tinyInfoD[6];
1286 void MEDCouplingIMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const
1289 std::ostringstream extent,origin,spacing;
1290 for(int i=0;i<3;i++)
1293 { extent << "0 " << _structure[i]-1 << " "; origin << _origin[i] << " "; spacing << _dxyz[i] << " "; }
1295 { extent << "0 0 "; origin << "0 "; spacing << "0 "; }
1297 ofs << " <" << getVTKDataSetType() << " WholeExtent=\"" << extent.str() << "\" Origin=\"" << origin.str() << "\" Spacing=\"" << spacing.str() << "\">\n";
1298 ofs << " <Piece Extent=\"" << extent.str() << "\">\n";
1299 ofs << " <PointData>\n" << pointData << std::endl;
1300 ofs << " </PointData>\n";
1301 ofs << " <CellData>\n" << cellData << std::endl;
1302 ofs << " </CellData>\n";
1303 ofs << " <Coordinates>\n";
1304 ofs << " </Coordinates>\n";
1305 ofs << " </Piece>\n";
1306 ofs << " </" << getVTKDataSetType() << ">\n";
1309 void MEDCouplingIMesh::reprQuickOverview(std::ostream& stream) const
1311 stream << "MEDCouplingIMesh C++ instance at " << this << ". Name : \"" << getName() << "\". Space dimension : " << _space_dim << ".";
1312 if(_space_dim<0 || _space_dim>3)
1315 std::ostringstream stream0,stream1;
1316 int nbNodes(1),nbCells(0);
1318 for(int i=0;i<_space_dim;i++)
1321 int tmpNodes(_structure[i]);
1322 stream1 << "- Axis " << tmp << " : " << tmpNodes << " nodes (orig=" << _origin[i] << ", inter=" << _dxyz[i] << ").";
1324 stream1 << std::endl;
1330 nbCells=nbCells==0?tmpNodes-1:nbCells*(tmpNodes-1);
1334 stream0 << "Number of cells : " << nbCells << ", Number of nodes : " << nbNodes;
1335 stream << stream0.str();
1337 stream << std::endl;
1339 stream << stream1.str();
1342 std::string MEDCouplingIMesh::getVTKDataSetType() const
1344 return std::string("ImageData");
1347 std::vector<std::string> MEDCouplingIMesh::buildInfoOnComponents() const
1349 checkSpaceDimension();
1350 int dim(getSpaceDimension());
1351 std::vector<std::string> ret(dim);
1352 for(int i=0;i<dim;i++)
1354 std::ostringstream oss;
1355 char tmp('X'+i); oss << tmp;
1356 ret[i]=DataArray::BuildInfoFromVarAndUnit(oss.str(),_axis_unit);
1361 void MEDCouplingIMesh::checkSpaceDimension() const
1363 CheckSpaceDimension(_space_dim);
1366 void MEDCouplingIMesh::CheckSpaceDimension(int spaceDim)
1368 if(spaceDim<0 || spaceDim>3)
1369 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::CheckSpaceDimension : input spaceDim must be in [0,1,2,3] !");
1372 int MEDCouplingIMesh::FindIntRoot(int val, int order)
1377 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : input val is < 0 ! Not possible to compute a root !");
1380 if(order!=2 && order!=3)
1381 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the order available are 0,1,2 or 3 !");
1382 double valf((double)val);
1385 double retf(sqrt(valf));
1388 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect square root !");
1393 double retf(std::pow(val,0.3333333333333333));
1394 int ret((int)retf),ret2(ret+1);
1395 if(ret*ret*ret!=val && ret2*ret2*ret2!=val)
1396 throw INTERP_KERNEL::Exception("MEDCouplingIMesh::FindIntRoot : the input val is not a perfect cublic root !");
1397 if(ret*ret*ret==val)