src/MEDCoupling/MEDCouplingSkyLineArray.cxx

   1 // Copyright (C) 2007-2016  CEA/DEN, EDF R&D
   2 //
   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.
   7 //
   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.
  12 //
  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
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "MEDCouplingSkyLineArray.hxx"
  21
  22 #include <sstream>
  23
  24 using namespace MEDCoupling;
  25
  26 MEDCouplingSkyLineArray::MEDCouplingSkyLineArray():
  27   _super_index( DataArrayInt::New() ), _index( DataArrayInt::New() ), _values( DataArrayInt::New() )
  28 {
  29 }
  30
  31 MEDCouplingSkyLineArray::~MEDCouplingSkyLineArray()
  32 {
  33 }
  34
  35 MEDCouplingSkyLineArray* MEDCouplingSkyLineArray::New()
  36 {
  37   return new MEDCouplingSkyLineArray();
  38 }
  39
  40 MEDCouplingSkyLineArray* MEDCouplingSkyLineArray::New( const std::vector<int>& index,
  41                                                        const std::vector<int>& value )
  42 {
  43   MEDCouplingSkyLineArray * ret = new MEDCouplingSkyLineArray();
  44   ret->_index->reserve( index.size() );
  45   ret->_index->insertAtTheEnd( index.begin(), index.end() );
  46   ret->_values->reserve( value.size() );
  47   ret->_values->insertAtTheEnd( value.begin(), value.end() );
  48   return ret;
  49 }
  50
  51 MEDCouplingSkyLineArray* MEDCouplingSkyLineArray::New( DataArrayInt* index, DataArrayInt* value )
  52 {
  53   MEDCouplingSkyLineArray* ret = new MEDCouplingSkyLineArray();
  54   ret->set(index, value);
  55   return ret;
  56 }
  57
  58 MEDCouplingSkyLineArray* MEDCouplingSkyLineArray::New( const MEDCouplingSkyLineArray & other )
  59 {
  60   MEDCouplingSkyLineArray* ret = new MEDCouplingSkyLineArray();
  61   ret->_super_index = other._super_index;
  62   ret->_index = other._index;
  63   ret->_values = other._values;
  64   return ret;
  65 }
  66
  67 /**! Build a three level SkyLine array from the dynamic connectivity of a dynamic mesh (i.e. containing only
  68  * polyhedrons or polygons).
  69  * The input arrays are deep copied, contrary to the other ctors.
  70  */
  71 MEDCouplingSkyLineArray * MEDCouplingSkyLineArray::BuildFromPolyhedronConn( const DataArrayInt* c, const DataArrayInt* cI )
  72 {
  73   using namespace std;
  74
  75   MEDCouplingSkyLineArray* ret = new MEDCouplingSkyLineArray();
  76
  77   const int * cP(c->begin()), * cIP(cI->begin());
  78   int prev = -1;
  79   if ((int)c->getNbOfElems() != *(cI->end()-1))
  80     throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::BuildFromDynamicConn: misformatted connectivity (wrong nb of tuples)!");
  81   for (std::size_t i=0; i < cI->getNbOfElems(); i++)
  82     {
  83       int j = cIP[i];
  84       if (cIP[i] < prev)
  85         throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::BuildFromDynamicConn: misformatted connectivity (indices not monotonic ascending)!");
  86       prev = cIP[i];
  87       if (i!=cI->getNbOfElems()-1)
  88         if (cP[j] != INTERP_KERNEL::NORM_POLYHED)
  89           throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::BuildFromDynamicConn: connectivity containing other types than POLYHED!");
  90     }
  91
  92   vector<int> superIdx, idx, vals;
  93   int cnt = 0, cnt2 = 0;
  94   superIdx.reserve(cI->getNbOfElems());
  95   superIdx.push_back(0);
  96   idx.push_back(0);
  97   vals.resize(c->getNbOfElems()); // too much because of the type and the -1, but still better than push_back().
  98   for (std::size_t i=0; i < cI->getNbOfElems()-1; i++)
  99     {
 100       int start = cIP[i]+1, end = cIP[i+1];
 101       int * work = vals.data() + cnt;
 102       const int * w = cP+start;
 103       const int * w2 = find(w, cP+end, -1);
 104       while (w2 != cP+end)
 105         {
 106           copy(w, w2, work);
 107           int d = distance(w, w2);
 108           cnt += d; work +=d;
 109           idx.push_back(cnt); cnt2++;
 110           w = w2+1;  // skip the -1
 111           w2 = find(w, cP+end, -1);
 112         }
 113       copy(w, cP+end, work);
 114       cnt += distance(w, cP+end);
 115       idx.push_back(cnt); cnt2++;
 116       superIdx.push_back(cnt2);
 117     }
 118   ret->_super_index->alloc(superIdx.size(),1);
 119   copy(superIdx.begin(), superIdx.end(), ret->_super_index->getPointer());
 120   ret->_index->alloc(idx.size(),1);
 121   copy(idx.begin(), idx.end(), ret->_index->getPointer());
 122   ret->_values->alloc(cnt,1);
 123   copy(vals.begin(), vals.begin()+cnt, ret->_values->getPointer());
 124
 125   return ret;
 126 }
 127
 128 /**
 129  * Convert a three-level SkyLineArray into a polyhedral connectivity.
 130  * The super-packs are interpreted as cell description, and the packs represent the face connectivity.
 131  */
 132 void MEDCouplingSkyLineArray::convertToPolyhedronConn( MCAuto<DataArrayInt>& c,  MCAuto<DataArrayInt>& cI) const
 133 {
 134   // TODO: in this case an iterator would be nice
 135   using namespace std;
 136
 137   checkSuperIndex("convertToPolyhedronConn");
 138
 139   const int * siP(_super_index->begin()), * iP(_index->begin()), *vP(_values->begin());
 140   int cnt = 0;
 141   cI->alloc(_super_index->getNbOfElems(),1);  // same number of super packs as number of cells
 142   int * cIVecP(cI->getPointer());
 143   MCAuto <DataArrayInt> dsi = _index->deltaShiftIndex();
 144   int sz = dsi->accumulate(0) + dsi->getNbOfElems();  // think about it: one slot for the type, -1 at the end of each face of the cell
 145   c->alloc(sz, 1);
 146   int * cVecP(c->getPointer());
 147
 148   for ( std::size_t i=0; i < _super_index->getNbOfElems()-1; i++)
 149      {
 150        cIVecP[i]= cnt;
 151        int endId = siP[i+1];
 152        cVecP[cnt++] = INTERP_KERNEL::NORM_POLYHED;
 153        for (int j=siP[i]; j < endId; j++)
 154          {
 155            int startId2 = iP[j], endId2 = iP[j+1];
 156            copy(vP+startId2, vP+endId2, cVecP+cnt);
 157            cnt += endId2-startId2;
 158            if(j != endId-1)
 159              cVecP[cnt++] = -1;
 160          }
 161      }
 162   cIVecP[_super_index->getNbOfElems()-1] = cnt;
 163 }
 164
 165 std::size_t MEDCouplingSkyLineArray::getHeapMemorySizeWithoutChildren() const
 166 {
 167   return _index->getHeapMemorySizeWithoutChildren()+_values->getHeapMemorySizeWithoutChildren()+_super_index->getHeapMemorySizeWithoutChildren();
 168 }
 169
 170 std::vector<const BigMemoryObject *> MEDCouplingSkyLineArray::getDirectChildrenWithNull() const
 171 {
 172   std::vector<const BigMemoryObject *> ret;
 173   ret.push_back(_super_index);
 174   ret.push_back(_index);
 175   ret.push_back(_values);
 176   return ret;
 177 }
 178
 179
 180 void MEDCouplingSkyLineArray::set( DataArrayInt* index, DataArrayInt* value )
 181 {
 182   _index=index;
 183   _values=value;
 184   if ( (DataArrayInt*)_index ) _index->incrRef();
 185   else                         _index = DataArrayInt::New();
 186   if ( (DataArrayInt*)_values ) _values->incrRef();
 187   else                         _values = DataArrayInt::New();
 188 }
 189
 190 void MEDCouplingSkyLineArray::set3( DataArrayInt* superIndex, DataArrayInt* index, DataArrayInt* value )
 191 {
 192   _super_index=superIndex;
 193   if ( (DataArrayInt*)_super_index ) _super_index->incrRef();
 194   else                         _super_index = DataArrayInt::New();
 195   set(index, value);
 196 }
 197
 198 DataArrayInt* MEDCouplingSkyLineArray::getSuperIndexArray() const
 199 {
 200   return const_cast<MEDCouplingSkyLineArray*>(this)->_super_index;
 201 }
 202
 203
 204 DataArrayInt* MEDCouplingSkyLineArray::getIndexArray() const
 205 {
 206   return const_cast<MEDCouplingSkyLineArray*>(this)->_index;
 207 }
 208
 209 DataArrayInt* MEDCouplingSkyLineArray::getValuesArray() const
 210 {
 211   return const_cast<MEDCouplingSkyLineArray*>(this)->_values;
 212 }
 213
 214 void MEDCouplingSkyLineArray::checkSuperIndex(const std::string& func) const
 215 {
 216   if (!_super_index->getNbOfElems())
 217     {
 218       std::ostringstream oss;
 219       oss << "MEDCouplingSkyLineArray::"<< func << ": not a three level SkyLineArray! Method is not available for two-level SkyLineArray.";
 220       throw INTERP_KERNEL::Exception(oss.str());
 221     }
 222 }
 223
 224 void MEDCouplingSkyLineArray::validSuperIndex(const std::string& func, int superIndex) const
 225 {
 226   if(superIndex < 0 || superIndex >= (int)_super_index->getNbOfElems())
 227     {
 228       std::ostringstream oss;
 229       oss << "MEDCouplingSkyLineArray::" << func <<  ": invalid super index!";
 230       throw INTERP_KERNEL::Exception(oss.str());
 231     }
 232 }
 233
 234 void MEDCouplingSkyLineArray::validIndex(const std::string& func, int idx) const
 235 {
 236   if(idx < 0 || idx >= (int)_index->getNbOfElems())
 237     {
 238       std::ostringstream oss;
 239       oss << "MEDCouplingSkyLineArray::" << func <<  ": invalid index!";
 240       throw INTERP_KERNEL::Exception(oss.str());
 241     }
 242 }
 243
 244 void MEDCouplingSkyLineArray::validSuperIndexAndIndex(const std::string& func, int superIndex, int index) const
 245 {
 246   validSuperIndex(func, superIndex);
 247   int idx = _super_index->begin()[superIndex] + index;
 248   if(idx < 0 || idx >= (int)_index->getNbOfElems())
 249     {
 250       std::ostringstream oss;
 251       oss << "MEDCouplingSkyLineArray::" << func <<  ": invalid index!";
 252       throw INTERP_KERNEL::Exception(oss.str());
 253     }
 254 }
 255
 256 std::string MEDCouplingSkyLineArray::simpleRepr() const
 257 {
 258   std::ostringstream oss;
 259   oss << "MEDCouplingSkyLineArray (" << this << ")" << std::endl;
 260   MCAuto<DataArrayInt> super_index = _super_index->deepCopy();
 261   if (_super_index->getNbOfElems())
 262     oss << "   Nb of super-packs: " << getSuperNumberOf() << std::endl;
 263   else
 264     {
 265       super_index->alloc(2,1);
 266       super_index->setIJSilent(0,0,0);
 267       super_index->setIJSilent(1,0,_index->getNbOfElems()-1);
 268     }
 269   oss << "   Nb of packs: " << getNumberOf() << std::endl;
 270   oss << "   Nb of values: " << getLength() << std::endl;
 271
 272   if (_super_index->getNbOfElems())
 273     {
 274       oss << "   Super-indices:" << std::endl;
 275       oss << "   ";
 276       const int * i = _super_index->begin();
 277       for ( ; i != _super_index->end(); ++i )
 278         oss << *i << " ";
 279       oss << std::endl;
 280     }
 281
 282   oss << "   Indices:" << std::endl;
 283   oss << "   ";
 284   const int * i = _index->begin();
 285   for ( ; i != _index->end(); ++i )
 286     oss << *i << " ";
 287   oss << std::endl;
 288   oss << "   Values:" << std::endl;
 289   oss << "     ";
 290   const int * v = _values->begin();
 291   int cnt = 0, cntI = 0;
 292   i = _index->begin();
 293   for ( const int * si = super_index->begin()+1; v != _values->end(); ++v, ++cnt )
 294     {
 295       if ( cnt == *i )
 296         {
 297           if ( cntI == *si && cnt != 0)
 298             {
 299               oss << std::endl << "     ";
 300               ++si;
 301             }
 302
 303           oss << "| ";
 304           ++i; ++cntI;
 305         }
 306       oss << *v << " ";
 307     }
 308   oss << std::endl;
 309
 310   return oss.str();
 311 }
 312
 313 /**
 314  * For a 2- or 3-level SkyLine array, return a copy of the absolute pack with given identifier.
 315  */
 316 void MEDCouplingSkyLineArray::getSimplePackSafe(const int absolutePackId, std::vector<int> & pack) const
 317 {
 318   if(absolutePackId < 0 || absolutePackId >= (int)_index->getNbOfElems())
 319     throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::getPackSafe: invalid index!");
 320   const int * iP(_index->begin()), *vP(_values->begin());
 321   int sz = iP[absolutePackId+1]-iP[absolutePackId];
 322   pack.resize(sz);
 323   std::copy(vP+iP[absolutePackId], vP+iP[absolutePackId+1],pack.begin());
 324 }
 325
 326 /**
 327  * Same as getPackSafe, but directly returns a pointer to the internal data with the size of the pack.
 328  */
 329 const int * MEDCouplingSkyLineArray::getSimplePackSafePtr(const int absolutePackId, int & packSize) const
 330 {
 331   if(absolutePackId < 0 || absolutePackId >= (int)_index->getNbOfElems())
 332     throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::getPackSafe: invalid index!");
 333   const int * iP(_index->begin()), *vP(_values->begin());
 334   packSize = iP[absolutePackId+1]-iP[absolutePackId];
 335   return vP+iP[absolutePackId];
 336 }
 337
 338
 339 /**!
 340  * For each given super-pack ID, provide the sub-index of the first matching pack. If no matching pack is found for the
 341  * given super-pack -1 is returned.
 342  * \param[in] superPackIndices the list of super-packs that should be inspected
 343  * \param[in] packBg the pack that the function is looking for in each of the provided super-pack
 344  * \param[in] packEnd the pack that the function is looking for in each of the provided super-pack
 345  * \param[out] a vector of int, having the same size as superPackIndices and containing for each inspected super-pack
 346  * the index of the first matching pack, or -1 if none found.
 347  */
 348 void MEDCouplingSkyLineArray::findPackIds(const std::vector<int> & superPackIndices,
 349                                           const int *packBg, const int *packEnd,
 350                                           std::vector<int>& out) const
 351 {
 352   using namespace std;
 353
 354   checkSuperIndex("findPackIds");
 355
 356   int packSz = std::distance(packBg, packEnd);
 357   if (!packSz)
 358     throw INTERP_KERNEL::Exception("MEDCouplingSkyLineArray::findPackIds: void pack!");
 359
 360   out.resize(superPackIndices.size());
 361   int i = 0;
 362   const int * siP(_super_index->begin()), * iP(_index->begin()), *vP(_values->begin());
 363   for(vector<int>::const_iterator it=superPackIndices.begin(); it!=superPackIndices.end(); ++it, i++)
 364     {
 365       out[i] = -1;
 366       const int sPackIdx = *it;
 367       // for each pack
 368       for (int idx=siP[sPackIdx], j=0; idx < siP[sPackIdx+1]; idx++, j++)
 369         {
 370           if (packSz == (iP[idx+1] - iP[idx]))
 371             if (equal(&vP[iP[idx]], &vP[iP[idx+1]], packBg))
 372               {
 373                 out[i] = j;
 374                 break;
 375               }
 376         }
 377     }
 378 }
 379
 380 /**!
 381  * Delete pack number 'idx' in super-pack number 'superIdx'.
 382  * \param[in] superIdx is the super-pack number
 383  * \param[in] idx is the pack index inside the super-pack 'superIdx'.
 384  */
 385 void MEDCouplingSkyLineArray::deletePack(const int superIdx, const int idx)
 386 {
 387   using namespace std;
 388
 389   checkSuperIndex("deletePack");
 390   validSuperIndexAndIndex("deletePack", superIdx, idx);
 391
 392   int * vP = _values->getPointer();
 393   int * siP(_super_index->getPointer()), *iP(_index->getPointer());
 394   const int start = iP[siP[superIdx]+idx], end = iP[siP[superIdx]+idx+1];
 395   // _values
 396   copy(vP+end, vP+_values->getNbOfElems(), vP+start);
 397   _values->reAlloc(_values->getNbOfElems() - (end-start));
 398
 399   // _index
 400   int nt = _index->getNbOfElems();
 401   copy(iP+siP[superIdx]+idx+1, iP+nt, iP+siP[superIdx]+idx);
 402   _index->reAlloc(nt-1); iP = _index->getPointer();  // better not forget this ...
 403   for(int ii = siP[superIdx]+idx; ii < nt-1; ii++)
 404     iP[ii] -= (end-start);
 405
 406   // _super_index
 407   for(int ii = superIdx+1; ii < (int)_super_index->getNbOfElems(); ii++)
 408     (siP[ii])--;
 409 }
 410
 411 /**!
 412  * Insert a new pack in super-pack at index 'superIdx'. The pack is inserted at the end of the pack list of the chosen super-pack.
 413  */
 414 void MEDCouplingSkyLineArray::pushBackPack(const int superIdx, const int * packBg, const int * packEnd)
 415 {
 416   using namespace std;
 417
 418   checkSuperIndex("pushBackPack");
 419   validSuperIndex("pushBackPack", superIdx);
 420
 421   int *siP(_super_index->getPointer()), *iP(_index->getPointer());
 422   const int sz(distance(packBg, packEnd));
 423
 424   // _values
 425   _values->reAlloc(_values->getNbOfElems()+sz);
 426   int * vPE(_values->getPointer()+_values->getNbOfElems());
 427   int *vP(_values->getPointer());
 428   copy(vP+iP[siP[superIdx+1]], vPE-sz, vP+iP[siP[superIdx+1]]+sz);
 429   // insert pack
 430   copy(packBg, packEnd, vP+iP[siP[superIdx+1]]);
 431
 432   // _index
 433   int nt = _index->getNbOfElems();
 434   _index->reAlloc(nt+1); iP = _index->getPointer();
 435   copy(iP+siP[superIdx+1]+1, iP+nt, iP+siP[superIdx+1]+2);
 436   iP[siP[superIdx+1]+1] = iP[siP[superIdx+1]] + sz;
 437   for(int ii = siP[superIdx+1]+2; ii < nt+1; ii++)
 438     iP[ii] += sz;
 439
 440   // _super_index
 441   for(int ii = superIdx+1; ii < (int)_super_index->getNbOfElems(); ii++)
 442     (siP[ii])++;
 443 }
 444
 445 /**
 446  * Replace pack with absolute index 'idx' with the provided new pack. Function can be used either
 447  * for 2-level SkyLine or 3-level SkyLine.
 448  */
 449 void MEDCouplingSkyLineArray::replaceSimplePack(const int idx, const int * packBg, const int * packEnd)
 450 {
 451   using namespace std;
 452
 453   validIndex("replaceSimplePack", idx);
 454
 455   int * iP(_index->getPointer());
 456   int newSz = distance(packBg, packEnd);
 457   const int start = iP[idx], end = iP[idx+1];
 458
 459   // _values
 460   int initValSz = _values->getNbOfElems();
 461   int deltaSz = newSz-(end-start);  // can be negative
 462   if (deltaSz)
 463     {
 464       if (deltaSz > 0)
 465         _values->reAlloc(initValSz+deltaSz);
 466       int *vP(_values->getPointer());
 467       copy(vP+end, vP+initValSz, vP+end+deltaSz);
 468       if (deltaSz < 0)
 469         _values->reAlloc(initValSz+deltaSz);
 470     }
 471
 472   // copy new pack
 473   copy(packBg, packEnd, _values->getPointer()+start);
 474
 475   // _index
 476   for(int ii = idx+1; ii < (int)_index->getNbOfElems(); ii++)
 477     iP[ii] += deltaSz;
 478 }
 479
 480 /**
 481  * Replace pack with super index 'superIdx' and index 'idx' with the provided new pack.
 482  * Function can be used only for 3-level SkyLine.
 483  */
 484 void MEDCouplingSkyLineArray::replacePack(const int superIdx, const int idx, const int * packBg, const int * packEnd)
 485 {
 486   using namespace std;
 487
 488   checkSuperIndex("replacePack");
 489   validSuperIndexAndIndex("replacePack", superIdx, idx);
 490
 491   int * siP(_super_index->getPointer()), *iP(_index->getPointer());
 492   int newSz = distance(packBg, packEnd);
 493   const int start = iP[siP[superIdx]+idx], end = iP[siP[superIdx]+idx+1];
 494
 495   // _values
 496   int initValSz = _values->getNbOfElems();
 497   int deltaSz = newSz-(end-start);  // can be negative
 498   if (deltaSz)
 499     {
 500       if (deltaSz > 0)
 501         _values->reAlloc(initValSz+deltaSz);
 502       int *vP(_values->getPointer());
 503       copy(vP+end, vP+initValSz, vP+end+deltaSz);
 504       if (deltaSz < 0)
 505         _values->reAlloc(initValSz+deltaSz);
 506     }
 507
 508   // copy new pack
 509   copy(packBg, packEnd, _values->getPointer()+start);
 510
 511   // _index
 512   for(int ii = siP[superIdx]+idx+1; ii < (int)_index->getNbOfElems(); ii++)
 513     iP[ii] += deltaSz;
 514 }