PartitionBaseT.cpp

Go to the documentation of this file.
00001 // This is a -*- C++ -*-  file
00002 //------------------------------------------------------------------------------------//
00003 //                                                                                    //
00004 //   Copyright (c) 2001                                                               //
00005 //      INRIA                                                                         //
00006 //      54600 VILLERS LES NANCY                                                       //
00007 //      France                                                                        //
00008 //                                                                                    //
00009 //------------------------------------------------------------------------------------//
00010 //                                                                                    //
00011 //                 * NOTICE OF PROPRIETARY INFORMATION *                              //
00012 //                                                                                    //
00013 // The information contained in this file is considered proprietary and the exclusive //
00014 // property of  INRIA. This information may not be disclosed, duplicated              //
00015 // or used, in whole or in part, for  any purpose  whatsoever without express written //
00016 // authorization from INRIA                                                           //
00017 //     
00018 #include <stdlib.h>
00019 #include <iomanip>
00020 #include <iostream>
00021 
00022 //
00023 #include "PartitionBaseT.h"
00024 //
00025 #include "tools.h"
00026 #include "Atoms.h"
00027 #include "pdbstruct.h"
00028 
00029 
00030 template <class TYPEOFATOMS>
00031 PartitionBase<TYPEOFATOMS>::PartitionBase() : _partitionType("QC_Standard_Part"),_numberOfPartition(0),_charge(NULL),_setOfAtoms(NULL),
00032                                               _setOfFragments(NULL),
00033                                               _subDomains(NULL),_atomsInZone(NULL),_Rbuffer1(3.0),_Rbuffer2(1.5),
00034                                               _useFrag(false)
00035 {}
00036 template <class TYPEOFATOMS>
00037 PartitionBase<TYPEOFATOMS>::PartitionBase(const int size, const double R1, const double R2) : 
00038                                                _partitionType("QC_Standard_Part"),_numberOfPartition(size),_charge(NULL),_setOfAtoms(NULL),
00039                                                _setOfFragments(NULL),_subDomains(new subDomain[size]),
00040                                                _atomsInZone(NULL),_Rbuffer1(R1),_Rbuffer2(R2), 
00041                                                _useFrag(false)
00042 {}
00043 
00044 template <class TYPEOFATOMS>
00045 PartitionBase<TYPEOFATOMS>::~PartitionBase()
00046 {
00047   _setOfAtoms = NULL;
00048   if(_subDomains != NULL)
00049     {delete [] _subDomains ; _numberOfPartition = 0 ;}
00050   if(_charge != NULL)
00051     {delete [] _charge ; _charge = NULL ;}
00052 }
00053 //
00054 template <class TYPEOFATOMS>
00055 void PartitionBase<TYPEOFATOMS>::resetData(){
00056   //
00057   std::cout << "  ---> resetData "  <<std::endl;
00058   for (int nsd = 0 ; nsd < _numberOfPartition ; ++nsd){
00059     _subDomains[nsd].resetData();
00060   }
00061   //  Le clear() de vector ne marche pas
00062   delete [] _atomsInZone ;
00063   _atomsInZone =  new std::vector<AtomIn>[this->_setOfAtoms->numberOfAtoms()];
00064 }
00065 template <class TYPEOFATOMS>
00066 bool PartitionBase<TYPEOFATOMS>::attachAtoms( TYPEOFATOMS *atoms )
00067 { 
00068   if( atoms == NULL){
00069     std::cerr <<"PartitionBase::attachAtoms :No Atoms to attach. " <<std::endl ;
00070     exit(1) ;
00071   }
00072   _useFrag    = false ;
00073   _setOfAtoms = atoms ;
00074   int nbAtoms = _setOfAtoms ->numberOfAtoms() ;
00075   //  std::cout <<" attachAtoms " << nbAtoms << std::endl ;
00076   _charge     = new int[_setOfAtoms ->numberOfAtoms() ] ;
00077   int num ;
00078   for ( int i = 0 ; i< _setOfAtoms ->numberOfAtoms() ; ++i){
00079     num        = giveNumberOfOA(_setOfAtoms->getAtomName(i)) ;
00080     if ( num != 0) {
00081       _charge[i] = num ; 
00082     }
00083     //    else{
00084 //       std::cout << "ERREUR for atomName : " << _setOfAtoms->getAtomName(i)
00085 //              << " the number of Orbital is   " << num << std::endl;
00086 //       std::cout << "     giveNumberOfOA(AtomName)" << giveNumberOfOA(_setOfAtoms->getAtomName(i)) <<std::endl;
00087       //      exit(EXIT_FAILURE) ;}
00088   }
00089   _atomsInZone =  new std::vector<AtomIn>[nbAtoms];
00090   //
00091   return true ;
00092 }
00093 template <class TYPEOFATOMS>
00094 bool PartitionBase<TYPEOFATOMS>::attachAtoms( setOfFragments *fragments )
00095 { 
00096   if( fragments == NULL){
00097     std::cerr <<"PartitionBase::attachAtoms :No Atoms to attach. " <<std::endl ;
00098     exit(1) ;
00099   }
00100     std::cerr <<"PartitionBase::attachAtoms :No Atoms to attach. " <<std::endl ;
00101     exit(1) ;
00102   return true ;
00103 }
00104 template <class TYPEOFATOMS>
00105 bool PartitionBase<TYPEOFATOMS>::build_fragmentStructure( ){
00106   //
00107   _nextATofFrag.resize(this->_setOfAtoms->numberOfAtoms())  ;
00108   _firstATofFrag.resize(this->_setOfAtoms->numberOfFragments()) ;
00109   //
00110   for(unsigned int i = 0 ; i < _nextATofFrag.size() ; ++i){
00111     _nextATofFrag[i] = -1  ; 
00112   }
00113  for(unsigned int i = 0 ; i < _firstATofFrag.size() ; ++i){
00114     _firstATofFrag[i] = -1  ;
00115   }
00116  //
00117  // On remplit la structure
00118  //  
00119  int fragNum, j ; 
00120  for (int i = 0 ; i < this->_setOfAtoms->numberOfAtoms() ; ++i){
00121    fragNum = _setOfAtoms->getElement(i).getFragmentNumber() -1; 
00122    //   std::cout <<" i " << i <<" fragNum "<<fragNum <<std::endl;
00123    if( _firstATofFrag[fragNum] == -1){
00124      _firstATofFrag[fragNum] = i ;}
00125    else{
00126      j = _firstATofFrag[fragNum]  ;
00127      //     std::cout << " i " << i << "  " <<_firstATofFrag[fragNum] << " " << _nextATofFrag[j] << " ";
00128      while ( _nextATofFrag[j] != -1){   
00129        j = _nextATofFrag[j]  ;
00130        //       std::cout << _nextATofFrag[j] <<"  " ;
00131      }  
00132      //     std::cout << std::endl;
00133      _nextATofFrag[j]   = i ;
00134   }
00135  }
00136  //
00137  //
00138 //  int  jj ;
00139 //  for (unsigned int i = 0 ; i < _firstATofFrag.size() ; ++i){
00140 //    std::cout << "Frag num : " << i << " Atomes :  " <<  _firstATofFrag[i];
00141 //    jj =  _firstATofFrag[i] ;
00142 //    while( _nextATofFrag[jj] != -1 ){
00143 //      std::cout << "  "<<_nextATofFrag[jj] ;
00144 //      jj = _nextATofFrag[jj] ;
00145 //    }
00146 //    std::cout <<  std::endl;
00147 //  }
00148  return true;
00149 }
00150 
00151 //
00152 template <class TYPEOFATOMS>
00153 bool PartitionBase<TYPEOFATOMS>::fragmentsToAtoms()
00154 { 
00155   std::cout << "BUGS BUGS  _atomsInZone[] not allocated " <<std::endl;
00156   if(!_useFrag)
00157     {return  true ;}
00158   //
00159   int nbAtoms = _setOfFragments->totalNumberOfAtoms()  ;
00160   //
00161   std::vector<AtomIn> * atomsInZone ;
00162   atomsInZone  =  new std::vector<AtomIn>[nbAtoms];
00163   //
00164   int *posInAtoms, *currentLocalNumber;
00165   _setOfAtoms        = new TYPEOFATOMS(nbAtoms) ;
00166   posInAtoms         = new int[_setOfFragments-> numberOfFraments() ] ;
00167   currentLocalNumber = new int[_numberOfPartition ];
00168   //
00169   int i, j, nsd, numberOfAtomsInFragment, currentAtoms;
00170   unsigned int k ;
00171   posInAtoms[0] = 0 ;
00172   for ( i = 1 ; i< _setOfFragments->numberOfFraments() ; ++i){
00173      posInAtoms[i] = posInAtoms[0] + _setOfFragments->getElement(i-1).numberOfAtoms();
00174   }
00175   //  std::cout << "Fin allocation " << nbAtoms << std::endl;
00176   for ( nsd = 1 ; nsd< _numberOfPartition ; ++nsd){
00177      currentLocalNumber[nsd] = 0; 
00178   }
00179   //
00180   // il faut reconstruire le tableau des correspondances
00181   //
00182   currentAtoms = 0 ;
00183   //  std::cout << " tableau des correspondances " << std::endl;
00184   for ( i = 0 ; i< _setOfFragments->numberOfFraments() ; ++i){
00185     FragmentElement & frag  = _setOfFragments->getElement(i);
00186     numberOfAtomsInFragment = frag.numberOfAtoms() ;
00187     // on recopie le vecteur pour tous les atoms 
00188     //    std::cout << " "<< i << "    "<< frag.numberOfAtoms() << std::endl;
00189     for (j = 0; j < numberOfAtomsInFragment ; ++j){
00190       //      copy(_atomsInZone[i].begin(),_atomsInZone[i].end(),atomsInZone[currentAtoms].begin()) ;
00191       for( k = 0 ; k <  _atomsInZone[j].size() ; ++k){
00192         atomsInZone[currentAtoms].push_back(_atomsInZone[i][k]) ;
00193       }
00194       ++currentAtoms; 
00195     }
00196     // on corrige la numerotation local;
00197 //     for (j = 0; j <  _setOfFragments->numberOfFraments(); ++j){
00198 //       std::cout << "      --  "<< j  << "    "<< atomsInZone[j].size();
00199 //       for( k = 0 ; k <  atomsInZone[j].size() ; ++k){
00200 //      std::cout << " (  "<< (_atomsInZone[j])[k]._numDomain<< ", "
00201 //                << (_atomsInZone[j])[k]._localNum<< ", "
00202 //                << (_atomsInZone[j])[k]._typeZone << " )  " ;
00203 //       }
00204 //       std::cout << std::endl;
00205 //     } 
00206     //   atomIn._numDomain = _atomsInZone[i]
00207 //     for ( j = 0 ; j < _atomsInZone[i].size ; ++j){
00208 //       fragIn._numDomain = (_atomsInZone[i])[j]._numDomain ;
00209 //       fragIn._localNum  = (_atomsInZone[i])[j]._localNum  ;
00210 //       fragIn._typeZone  = (_atomsInZone[i])[j]._typeZone  ;
00211   }
00212   //
00213   return true ;
00214 }
00215 //
00216 template <class TYPEOFATOMS>
00217 void PartitionBase<TYPEOFATOMS>::setLoad(){
00218   //
00219   int nsd, k, loadOverlap, loadCut,loadKernel;
00220   unsigned int at;
00221   
00222   //
00223   for(nsd = 0 ; nsd <  _numberOfPartition  ; ++nsd) {
00224     //
00225     loadOverlap = loadCut = loadKernel = 0 ;
00226     std::vector<int> & indexK         = _subDomains[nsd].getAtomsIndexInKernel() ;
00227     //
00228     for (k = 0 ; k < _subDomains[nsd].getNumberOfNeighbors() ; ++k){
00229       std::vector<int> & indexO = _subDomains[nsd].getAtomsIndexInOverlap(k) ;
00230       std::vector<int> & indexC = _subDomains[nsd].getAtomsIndexInCut(k) ;
00231       for(at = 0 ; at < indexO.size() ; ++at){
00232         loadOverlap += _charge[indexO[at]] ;
00233       }
00234       for(at = 0 ; at < indexC.size() ; ++at){
00235         loadCut  += _charge[indexC[at]] ;
00236       }
00237     }
00238     for(at = 0 ; at < indexK.size() ; ++at){
00239       loadKernel  += _charge[indexK[at]] ;
00240     }
00241 //     std::cout << " Subdomain : "<< nsd<<std::endl
00242 //            << "      loadKernel  : "<< loadKernel<<std::endl
00243 //            << "      loadOverlap : "<< loadOverlap<<std::endl
00244 //            << "      loadCut     : "<< loadCut<<std::endl 
00245 //            << "       Total      : "<< loadCut + loadOverlap + loadKernel <<std::endl;
00246     //
00247     _subDomains[nsd].loadInKernel()  = loadKernel ;
00248     _subDomains[nsd].loadInOverlap() = loadOverlap ;
00249     _subDomains[nsd].loadInCut()     = loadCut ;
00250   }
00251 }
00252 //
00253 template <class TYPEOFATOMS>
00254 void PartitionBase<TYPEOFATOMS>::fixAtomsInZone(){
00255   //
00256   bool toErase ;
00257   for( int i = 0 ; i< this->_setOfAtoms->numberOfAtoms() ; ++i ){
00258     if(_atomsInZone[i].size() <= 2){
00259       continue ;}
00260     toErase = false ;
00261     if((_atomsInZone[i])[0]._localNum == -1 ){
00262       //      std::cout << "  Erase atome " << i <<std::endl;
00263       toErase = true ; 
00264     }
00265     if(toErase)
00266       {(_atomsInZone[i]).erase(_atomsInZone[i].begin() ) ;}
00267   }
00268 }
00269 
00270 template <class TYPEOFATOMS>
00271 void PartitionBase<TYPEOFATOMS>::writeInFile(std::ofstream &out) 
00272 {
00273   //  std::cout <<" void PartitionBase::writeInFile(std::ofstream &out)" << std::endl;
00274   //
00275   int i , nsd ;
00276   //
00277   int j ;
00278    std::vector< std::vector<int> >  numOfSubDom( _numberOfPartition);
00279    std::vector< int > load(_numberOfPartition);
00280    int loadAtoms = 0 , numSubDom;
00281    for(i=0; i <  _numberOfPartition ; ++i){
00282      numOfSubDom[i].resize( _numberOfPartition) ;
00283      load[i] = 0 ;
00284      for(j = 0; j <  _numberOfPartition ; ++j){
00285        (numOfSubDom[i])[j] = 0 ;
00286      }
00287    }
00288    int numSubDom1  ; 
00289    for(i = 0 ; i < _setOfAtoms->numberOfAtoms()  ; ++i) {
00290      //
00291      loadAtoms               = giveNumberOfOA(_setOfAtoms->getAtomName(i)) ;
00292      std::vector<AtomIn> & v = _atomsInZone[i] ;
00293      for (unsigned int k = 0 ; k < v.size() ; ++k ){
00294        numSubDom        = v[k]._numDomain ; 
00295        load[numSubDom] += loadAtoms ;
00296        for (unsigned int l = 0 ; l < v.size() ; ++l ){
00297         numSubDom1                            = v[l]._numDomain ;
00298         if(v[l]._typeZone != Cut_Zone ||  numSubDom1 == numSubDom){
00299           (numOfSubDom[numSubDom])[numSubDom1] += 1  ;
00300         }
00301        }
00302      }
00303      //
00304    }
00305 //  
00306   out << "# " << std::endl 
00307       << "# Partition Format version 1.0   " << std::endl 
00308       << "#   GENERATED BY PartitionBase<TYPEOFATOMS>::writeInFile"<< std::endl 
00309       << "#"<< std::endl 
00310       << "# Zones of a subdomain and associated zone identifiers :" << std::endl 
00311       << "#   _____________________________________________ "<< std::endl 
00312       << "#  |        |       |           |       |        |"<< std::endl 
00313       << "#  |Cut Zone|Overlap|Kernel Zone|Overlap|Cut Zone|"<< std::endl 
00314       << "#  |   2    |Zone 1 |     0     |Zone 1 |   2    |"<< std::endl 
00315       << "#  |________|_______|___________|_______|________|"<< std::endl 
00316       << "# " << std::endl 
00317       << "# " << std::endl 
00318       << "# " << std::endl 
00319       <<"NUMBER OF SUBDOMAINS "<< std::endl << _numberOfPartition           << std::endl 
00320       <<"RADIUS BUFFER1 "      << std::endl << _Rbuffer1                    << std::endl 
00321       <<"RADIUS BUFFER2 "      << std::endl << _Rbuffer2                    << std::endl 
00322       <<"PARTITION TYPE"       << std::endl << _partitionType               << std::endl
00323       <<"NUMBER OF ATOMS"      << std::endl << _setOfAtoms->numberOfAtoms() << std::endl 
00324       << "# " << std::endl ;
00325   //
00326   out << "# num1 = numero of subdomain" << std::endl
00327       << "# num2 = number of atom(s) in the subdomain num1" << std::endl
00328       << "# num3 = number of AO in the subdomain num1" << std::endl
00329       << "# num4 = number of subdomain(s) overlapping num1" << std::endl
00330       << "# num5 = numero of a subdomain overlapping num1" << std::endl
00331       << "# num6 = number of atom(s) from num1 overlapped by num5"<< std::endl
00332       << "# The doublet (num5, num6) is repeated 'num4-1' time."<< std::endl
00333       << "#num1   num2    num3    num4    num5    num6 " << std::endl
00334       << "PARTITION DESCRIPTION" << std::endl ;
00335   for( nsd = 0 ; nsd < _numberOfPartition  ; ++nsd) {
00336 //      out << nsd << "  \t"<< _subDomains[nsd].getNumberOfAtoms() << "  \t"<<_subDomains[nsd].load() 
00337 //      <<"  \t" << load[nsd] <<"  \t" << _subDomains[nsd].getNumberOfNeighbors() << "   ";
00338     numSubDom = 0;
00339     for ( i = 0 ; i < _numberOfPartition ; ++i){
00340        if( (numOfSubDom[nsd])[i] != 0){
00341         numSubDom++;
00342        }
00343      }
00344     numSubDom--;
00345      out << nsd << "  \t"<< _subDomains[nsd].getNumberOfAtoms() 
00346               <<"  \t" << load[nsd] <<"  \t" << numSubDom << "   ";
00347     //
00348      for(  i = 0 ; i < _numberOfPartition ; ++i){
00349        if ( (numOfSubDom[nsd])[i] != 0 && i!=nsd){
00350         out << i << "  "<< (numOfSubDom[nsd])[i] << "  ";
00351        }
00352      } 
00353 //     const std::vector<int> & Neighboors  = _subDomains[nsd].getNeighboorsIndex() ;
00354 //     for(  i = 0 ; i < Neighboors.size() ; ++i){
00355 //       out << Neighboors[i] << "  "<< _subDomains[nsd].getNumberOfAtomsInOverlap(i ) 
00356 //      +  _subDomains[nsd].getNumberOfAtomsInCut(i ) << "  ";
00357 //     }
00358     out  <<  std::endl ;
00359   }
00360   out << "#"<< std::endl
00361       << "# num1 = index of atom"<< std::endl
00362       << "# num2 = the atom num1 is in num2 subdomain(s)"<< std::endl
00363       << "# num3 = numero of a subdomain including num1"<< std::endl
00364       << "# num4 = identifier of zone in num3 where is located num1"<< std::endl
00365       << "# num5 = local numero of num1 in this subdomain num3"<< std::endl
00366       << "# The tripplet (num3, num4, num5) is repeated 'num2-1' time"<< std::endl
00367       << "#"<< std::endl
00368       << "#num1   num2    num3    num4    num5"<< std::endl
00369       << "PARTITION"<< std::endl ;
00370   for(i = 0 ; i < _setOfAtoms->numberOfAtoms()  ; ++i) {
00371     out << "  "<< i << "  "<< _atomsInZone[i].size() ;
00372     std::vector<AtomIn> & v = _atomsInZone[i] ;
00373     for (unsigned int j = 0 ; j < v.size() ; ++j ){
00374       out << "  "<<  v[j]._numDomain <<"   "<< v[j]._typeZone  <<"   "<< v[j]._localNum ;
00375     }
00376     out << std::endl ;
00377   }
00378 }
00379 //
00380 // ==================================================================================
00381 //
00382 template <class TYPEOFATOMS>
00383 void PartitionBase<TYPEOFATOMS>::readFromFile(std::ifstream & data) 
00384 {
00385   std::string type,line("") ;
00386   std::stringstream ss ;
00387   int pos, loc =-1,version ;
00388   //
00389   pos = data.tellg() ;
00390   while ( loc < 0) {
00391      data.seekg(pos);
00392     std::getline(data, line);
00393     loc = line.find_last_of('.'); 
00394   }
00395   line.erase(line.begin(), line.begin()+loc-1);
00396   ss <<line ;  ss >> version ;
00397   data.seekg(pos); ;
00398 
00399   if(version == 2) {
00400     std::cout << "newReadFromFile"<<std::endl;
00401     this->newReadFromFile(data);
00402   }
00403   else if(version == 1) {
00404     std::cout << "oldReadFromFile 1"<<std::endl;
00405     this->oldReadFromFile(data);
00406   }
00407   else {
00408     std::cerr << " Wrong format of partition file " << std::endl;
00409     exit(EXIT_FAILURE) ;
00410   }
00411 }
00412 template <class TYPEOFATOMS>
00413 void PartitionBase<TYPEOFATOMS>::oldReadFromFile(std::ifstream & data) {
00414   std::cout << "BEGIN aoldReadFromFile  "<<std::endl;
00415   int nbAtoms;
00416   std::string line,comment("#") ;
00417   std::stringstream temp ;
00418   //
00419   int pos=0 ;
00420   data.seekg(std::ios::beg) ;
00421   //
00422   // Read The header
00423   
00424   std::getline(data,line) ;
00425   while(line[0] == comment[0]) {
00426     std::getline(data,line) ;
00427   }
00428   std::getline(data,line) ;
00429   temp.clear() ;
00430   // NUMBER OF SUBDOMAINS
00431   temp  << line ;  temp >> _numberOfPartition ; temp.clear();
00432   std::cout << _numberOfPartition <<std::endl;
00433   //
00434   std::getline(data,line) ;
00435   std::getline(data,line) ;
00436   temp  << line ;  temp >> _Rbuffer1 ; temp.clear();
00437   //
00438   std::getline(data,line) ;
00439   std::getline(data,line) ;
00440   temp  << line ;  temp >> _Rbuffer2 ; temp.clear();
00441   //
00442   std::getline(data,line) ;
00443   std::getline(data,line) ;
00444   temp  << line ;  temp >> _partitionType ; temp.clear();
00445   //
00446   std::getline(data,line) ;
00447   std::getline(data,line) ;
00448   temp  << line ;  temp >> nbAtoms ; temp.clear();
00449   //
00450   //  Allocations
00451   //
00452   _subDomains = new subDomain[_numberOfPartition] ;
00453   //
00454   std::getline(data,line) ;
00455   while(line[0] == comment[0]) {
00456     std::cout << line <<std::endl;
00457     std::getline(data,line) ;
00458   }
00459   std::cout << line <<std::endl;
00460   std::getline(data,line) ;
00461   std::cout << line <<std::endl;
00462   //
00463   // PARTITION DESCRIPTION
00464   //
00465   int nSD, nbAtnSD, nbAonSD, nbNeigh, numNeigh, nbATinBuffer1 ;
00466   std::vector< std::vector<int> > voisins(_numberOfPartition) ;
00467   for( int nsd = 0 ; nsd < _numberOfPartition ; ++nsd ){
00468   temp  << line ;  
00469   temp >> nSD >> nbAtnSD  >> nbAonSD >> nbNeigh ; 
00470   voisins[nSD].resize(nbNeigh) ;
00471   std::cout  << nSD << "  "<< nbAtnSD << "  "<<nbAonSD << "  "<< nbNeigh ; 
00472   _subDomains[nsd].initNumberOfOverlapDomain(nbNeigh);
00473   for( int n = 0 ; n < nbNeigh ; ++n ){
00474     temp >> numNeigh >> nbATinBuffer1 ; 
00475     voisins[nSD][n] = numNeigh ;
00476     _subDomains[nsd].setNumberOfNeighbor(n,  numNeigh) ;
00477     std::cout << "  "<<numNeigh << "  "<<  nbATinBuffer1 ;
00478   }
00479   std::cout  << std::endl;
00480   temp.clear();
00481   pos = data.tellg() ;
00482   std::getline(data,line) ;
00483   }
00484   while(line[0] == comment[0]) {
00485     pos = data.tellg() ;
00486     std::getline(data,line) ;
00487   }
00488   data.seekg(pos) ;
00489   std::cout << "Attention lecture des atomes " <<nbAtoms << std::endl ;
00490   int a, aInNBSD, isInSD, typeZone, localNum ;
00491   for( int at = 0 ; at < nbAtoms ; ++at ){
00492     data >> a >> aInNBSD ;
00493     std::cout << a << " "<< aInNBSD << " - ";
00494     for ( int i = 0 ; i < aInNBSD; ++i) {
00495       data  >> isInSD>> typeZone >> localNum ;
00496       if( typeZone==0){
00497         _subDomains[isInSD].addIndexInKernel(localNum, at) ;}
00498       else if( typeZone==1){
00499         //      _subDomains[isInSD].addIndexInOverlap(n,localNum) ;
00500       }
00501       else {
00502         //_subDomains[isInSD].addIndexInCut(n,nlocalNum) ;
00503       }
00504       std::cout  << "       "<<  isInSD << " "<<  typeZone  << " "<<  localNum ;
00505     }
00506     std::cout << std::endl;
00507   }
00508   std::cerr << "oldFromFile not yet implemented" <<std::endl;
00509   //  exit(0) ;
00510 
00511 }
00512 template <class TYPEOFATOMS>
00513 void PartitionBase<TYPEOFATOMS>::newReadFromFile(std::ifstream & data) 
00514 {
00515   std::stringstream temp ;
00516   int loc ;
00517   std::string line, bid ,comment("#");
00518   // Skip comment lines
00519   std::getline(data, line); 
00520   while ( line[0] == comment[0]) {
00521     std::getline(data, line); 
00522   }
00523   // READ HEADER
00524   loc = line.find(":") ; 
00525   temp << line.substr(loc+1,line.size() ) ; temp >> _numberOfPartition; temp.clear();
00526   std::getline(data,line) ;
00527   loc = line.find(":") ; 
00528   temp << line.substr(loc+1,line.size() ) ; temp >> _Rbuffer1; temp.clear();
00529   std::getline(data,line) ;
00530   loc = line.find(":") ; 
00531   temp << line.substr(loc+1,line.size() ) ; temp >> _Rbuffer2; temp.clear();
00532   std::getline(data,line) ;
00533   loc = line.find(":") ; 
00534   temp << line.substr(loc+1,line.size() ) ; temp >> _partitionType; temp.clear();
00535   //
00536   // ALLOCATE SUBDOMAIN
00537   _subDomains = new subDomain[_numberOfPartition]  ;
00538   //
00539   int num ;
00540   for(int i = 0; i <  _numberOfPartition ; ++i){
00541     //
00542     bool cont = true; 
00543     do {
00544       std::getline(data,line) ;
00545       loc = line.find("Domain");
00546       if( line.size() !=0 && line.find("Domain") != 2)
00547         { cont = false;}
00548     }   while(cont);
00549     loc = line.find(":");
00550     temp << line.substr(loc+1,line.size() ) ; temp >> num; temp.clear();
00551     //  lecture des partitions 
00552     // CORE 
00553     int nbVoisins, atomInCore, atomInOverlap, atomInCut, numV ;
00554     data >>  nbVoisins >> atomInCore>> atomInOverlap >> atomInCut ;
00555     int numAt ;
00556     std::getline(data,line) ;
00557     std::getline(data,line) ;
00558     for(int j = 0; j <  atomInCore ; ++j){
00559       data >> numAt ; 
00560       _subDomains[num].addIndexInKernel(numAt);
00561     }
00562     std::getline(data,line) ;
00563     std::getline(data,line) ;
00564     // 
00565     // ALLOCATE 
00566     _subDomains[num].initNumberOfOverlapDomain(nbVoisins);
00567     for ( int n = 0 ; n < nbVoisins ; ++n){
00568       data >> numV >>atomInOverlap;  
00569       _subDomains[num].setNumberOfNeighbor(n,  numV) ;
00570       for( int j = 0; j <  atomInOverlap ; ++j){
00571         data >> numAt ; 
00572         _subDomains[num].addIndexInOverlap(n,numAt);
00573       }
00574       std::getline(data,line) ;
00575     }
00576     std::getline(data,line) ;
00577     for ( int n = 0 ; n < nbVoisins ; ++n){
00578       data >> numV >>atomInCut;
00579       for( int j = 0; j <  atomInCut ; ++j){
00580         data >> numAt ; 
00581         _subDomains[num].addIndexInCut(n,numAt);
00582       }
00583       std::getline(data,line) ;
00584     }
00585 
00586   }
00587   
00588 }
00589 
00590 template <class TYPEOFATOMS>
00591 void PartitionBase<TYPEOFATOMS>::newWriteInFile(std::ofstream &out) 
00592 {
00593   //  std::cout <<" void PartitionBase::newWriteInFile(std::ofstream &out)" << std::endl;
00594   //
00595   //  
00596   unsigned int num, numOverlap, numCut  ;
00597   //
00598   out <<"# Partition Format Version 2.0  BUG BUG " << std::endl 
00599       << "# GENERATED BY  PartitionBase<TYPEOFATOMS>::newWriteInFile"<< std::endl ;
00600   out << "NUMBER OF PARTITION    :  " << _numberOfPartition << std::endl ;
00601   out << "Radius of Buffer1      :  " << _Rbuffer1 << std::endl
00602       << "Radius of Buffer2      :  " << _Rbuffer2 << std::endl 
00603       << "Type                   :  " << _partitionType<< std::endl ;
00604 
00605   //
00606   for(int i = 0; i <  _numberOfPartition ; ++i){
00607     numOverlap = numCut = 0 ;
00608     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00609       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ;
00610       numOverlap += indexOVERLAP.size() ;
00611       std::vector<int> & indexCUT     = _subDomains[i].getAtomsIndexInCut(n) ;
00612       numCut += indexCUT.size() ;
00613     }
00614 
00615     num = 0 ;
00616     std::vector<int> & indexCORE    = _subDomains[i].getAtomsIndexInKernel() ;
00617     out << "Sub Domain number : " << i <<std::endl ;
00618     out << " "<< _subDomains[i].getNumberOfNeighbors() <<"    "
00619         << _subDomains[i].getAtomsIndexInKernel().size() <<"  "<< numOverlap << "  " << numCut<<std::endl ;
00620     out << "   CORE ZONE    : " <<std::endl <<"      ";
00621     if(indexCORE.size() != 0 ){
00622       for(unsigned int j = 0; j <  indexCORE.size() ; ++j){
00623         out <<  indexCORE[j] << " " ;
00624       }
00625       out << std::endl ;
00626     }
00627     num += indexCORE.size() ;
00628     out << "   OVERLAP ZONE  : " <<std::endl ;
00629     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00630       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ;
00631       if  (indexOVERLAP.size() ==0) {continue ;}
00632       out << "      " << _subDomains[i].getNumberOfNeighbor(n) <<"  "<< indexOVERLAP.size()<< "      ";
00633       for(unsigned int j = 0; j <  indexOVERLAP.size() ; ++j){
00634         out << indexOVERLAP[j] << "  ";
00635       }
00636       out << std::endl;
00637       num += indexOVERLAP.size() ;
00638     }
00639     out << "   CUT ZONE      : "  <<std::endl ;
00640     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00641       std::vector<int> & indexCUT     = _subDomains[i].getAtomsIndexInCut(n) ;
00642       if  (indexCUT.size() ==0) {continue ;}
00643       out << "      " << _subDomains[i].getNumberOfNeighbor(n) << "  "<< indexCUT.size() <<"      ";
00644       for(unsigned int j = 0; j <  indexCUT.size() ; ++j){
00645         out << indexCUT[j] << "   ";
00646       }
00647       out << std::endl;
00648       num += indexCUT.size() ;
00649     }
00650     out   << std::endl << std::endl ;
00651   }
00652 }
00653 template <class TYPEOFATOMS>
00654 void PartitionBase<TYPEOFATOMS>::writeInVTKFormat(std::ofstream &out) 
00655 {
00656   //  
00657   // Affectation des couleurs
00658   //
00659   int nbAtoms = this->_setOfAtoms->numberOfAtoms() ;
00660   std::vector<int>  scalar(nbAtoms,-1) ;
00661   for(int i = 0; i <  _numberOfPartition ; ++i){
00662     std::vector<int> & indexCORE    = _subDomains[i].getAtomsIndexInKernel() ;
00663     for(unsigned int j = 0; j <  indexCORE.size() ; ++j){
00664       scalar[indexCORE[j]] = i ;
00665     }
00666     if(_partitionType == "QC_Dixon_Part"){continue;}
00667     //  Notre cas il faut gerer l'overlap
00668     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00669       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ;
00670       if(  i > _subDomains[i].getNumberOfNeighbor(n) ){
00671         for(unsigned int j = 0; j <  indexOVERLAP.size() ; ++j){
00672           scalar[indexOVERLAP[j]] = _numberOfPartition ;
00673         }
00674       }
00675     }
00676   }
00677   //
00678   out << "<VTKFile type=\"PolyData\" version=\"0.1\" byte_order=\"LittleEndian\">  "<< std::endl 
00679       <<"<PolyData>" << std::endl <<"<Piece NumberOfPoints=\"" << nbAtoms <<"\">"<< std::endl ;
00680   out << "<Points>" <<std::endl ;
00681   out << "<DataArray type=\"Float64\" Name=\"coord_atoms\"  NumberOfComponents=\"3\" format=\"ascii\">" <<std::endl ;
00682   double * xyz ;
00683   for ( int i = 0 ; i < nbAtoms ; ++i){
00684     xyz =  this->_setOfAtoms->getElement(i).pos();
00685     out << xyz[0] << " " <<xyz[1] << " " <<xyz[2] << " " <<std::endl;
00686   }
00687   out << "</DataArray>" <<std::endl ;
00688   out << "</Points>" <<std::endl ;
00689   
00690   out << "<PointData  Scalars=\"PartitionNumber\"> "<<std::endl ;
00691   out << "<DataArray type=\"Int32\" Name=\"PartitionNumber\"  format=\"ascii\"> " <<std::endl;
00692   //
00693   for ( unsigned int i = 0 ; i < scalar.size() ; ++i){
00694     out << scalar[i] << " " ;
00695   }
00696   out <<std::endl ;
00697   //
00698   out << "</DataArray>" <<std::endl ;
00699   out << "</PointData>" <<std::endl ;
00700   out << "</Piece>" <<std::endl ;
00701   out << "</PolyData>" <<std::endl ;
00702   out << "</VTKFile>" <<std::endl ;
00703   //
00704 }
00705 
00706 
00707 template <class TYPEOFATOMS>
00708 void PartitionBase<TYPEOFATOMS>::newWriteInFileOld(std::ofstream &out) 
00709 {
00710   //  std::cout <<" void PartitionBase::newWriteInFile(std::ofstream &out)" << std::endl;
00711   //
00712   //  
00713   unsigned int num, numOverlap, numCut  ;
00714   //
00715   out << "NEW FORMAT" << std::endl ;
00716   out << "NUMBER OF PARTITION  : " << _numberOfPartition << std::endl ;
00717   out << "Size of Buffer1      :  " << _Rbuffer1 << std::endl
00718       << "Size of Buffer2      :  " << _Rbuffer2 << std::endl 
00719       << "Type                 :  " << _partitionType<< std::endl ;
00720 
00721   //
00722   for(int i = 0; i <  _numberOfPartition ; ++i){
00723     numOverlap = numCut = 0 ;
00724     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00725       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ;
00726       numOverlap += indexOVERLAP.size() ;
00727       std::vector<int> & indexCUT     = _subDomains[i].getAtomsIndexInCut(n) ;
00728       numCut += indexCUT.size() ;
00729     }
00730 
00731     num = 0 ;
00732     std::vector<int> & indexCORE    = _subDomains[i].getAtomsIndexInKernel() ;
00733     out << "Sub Domain number : " << i <<std::endl ;
00734     out << " "<< _subDomains[i].getAtomsIndexInKernel().size() <<"  "<< numOverlap << "  " << numCut<<std::endl ;
00735     out << "   CORE ZONE    : " ;
00736     for(unsigned int j = 0; j <  indexCORE.size() ; ++j){
00737       out << "  { " << j << " , "<< indexCORE[j] << " } ";
00738     }
00739     out << std::endl ;
00740     num += indexCORE.size() ;
00741     out << "   OVERLAP ZONE  : " ;
00742     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00743       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ;
00744       for(unsigned int j = 0; j <  indexOVERLAP.size() ; ++j){
00745         out << "  { " << num + j << " , "<< indexOVERLAP[j] << " } ";
00746       }
00747       num += indexOVERLAP.size() ;
00748     }
00749     out << std::endl ;
00750     out << "   CUT ZONE      : " ;
00751     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00752       std::vector<int> & indexCUT     = _subDomains[i].getAtomsIndexInCut(n) ;
00753       for(unsigned int j = 0; j <  indexCUT.size() ; ++j){
00754         out << "  { " << num + j << " , "<< indexCUT[j] << " } ";
00755       }
00756       num += indexCUT.size() ;
00757     }
00758     out << std::endl << std::endl ;
00759     out <<"======================================================"<<std::endl
00760       << std::endl << std::endl ;
00761   }
00762 }
00763 template <class TYPEOFATOMS>
00764 void  PartitionBase<TYPEOFATOMS>::getAtomsIndexInPartition(const int& part, std::vector<int>& index){
00765   std::cout << _subDomains[part].getNumberOfAtoms() <<std::endl;
00766   index.resize(_subDomains[part].getNumberOfAtoms() );
00767   std::vector<int> & indexCORE    = _subDomains[part].getAtomsIndexInKernel() ;
00768   for(unsigned int j = 0; j <  indexCORE.size() ; ++j){
00769       index[j] = indexCORE[j] ;
00770     }
00771   int num =  indexCORE.size();
00772   for ( int n = 0 ; n < _subDomains[part].getNumberOfNeighbors() ; ++n){
00773     std::vector<int> & indexOVERLAP = _subDomains[part].getAtomsIndexInOverlap(n) ;
00774     for(unsigned int j = 0; j <  indexOVERLAP.size() ; ++j){
00775       index[num + j] =  indexOVERLAP[j] ;
00776       }
00777       num += indexOVERLAP.size() ;
00778     }
00779   for ( int n = 0 ; n < _subDomains[part].getNumberOfNeighbors() ; ++n){
00780     std::vector<int> & indexCUT = _subDomains[part].getAtomsIndexInCut(n) ;
00781     for(unsigned int j = 0; j <  indexCUT.size() ; ++j){
00782       index[num + j] =  indexCUT[j] ;
00783       }
00784       num += indexCUT.size() ;
00785     }
00786 }
00787 
00788 //
00789 // ==================================================================================
00790 //
00791 template <class TYPEOFATOMS>
00792 void PartitionBase<TYPEOFATOMS>::buildAtomsInZone() {
00793   //
00794   //
00795   AtomIn  atomNSD ;
00796   int     num ;
00797   //
00798   for (unsigned int i = 0 ; i < _atomsInZone->size() ; ++i){
00799     _atomsInZone[i].clear() ;
00800   }
00801   // numerotation locale du sous domaine i
00802   //  std::cout << "Number Of Partiotions : "<< _numberOfPartition <<std::endl;
00803   for(int i = 0; i <  _numberOfPartition ; ++i){
00804 //     std::cout << "Sous domaine " << i << "  nombre de voisins  " 
00805 //            <<_subDomains[i].getNumberOfNeighbors()<< std::endl;
00806     num                = 0 ;
00807     atomNSD._numDomain = i ;  
00808     //
00809     //  Core
00810     //
00811     std::vector<int> & indexCORE    = _subDomains[i].getAtomsIndexInKernel() ;
00812     atomNSD._typeZone  = Core_Zone ;
00813     for(unsigned int j = 0; j <  indexCORE.size() ; ++j){
00814       atomNSD._localNum  = j ;
00815       _atomsInZone[indexCORE[j]].push_back(atomNSD);
00816     }
00817     //
00818     //  Overlap 
00819     //
00820     num               += indexCORE.size() ;
00821     atomNSD._typeZone  = Overlap_Zone ;
00822     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00823       //      std::cout << "       BUFFER1  " <<  n << std::endl; 
00824       std::vector<int> & indexOVERLAP = _subDomains[i].getAtomsIndexInOverlap(n) ; 
00825       for(unsigned int j = 0; j <  indexOVERLAP.size() ; ++j){
00826         atomNSD._localNum  = num + j ;
00827         _atomsInZone[indexOVERLAP[j]].push_back(atomNSD);
00828       }
00829       num += indexOVERLAP.size() ;
00830     }
00831     //
00832     //  Cut
00833     //
00834     atomNSD._typeZone  = Cut_Zone ;
00835     for ( int n = 0 ; n < _subDomains[i].getNumberOfNeighbors() ; ++n){
00836 //       std::cout << "       BUFFER2 "<<  n << std::endl; 
00837 //       std::cout << "               " ;
00838       std::vector<int> & indexCUT     = _subDomains[i].getAtomsIndexInCut(n) ;
00839       for(unsigned int j = 0; j <  indexCUT.size() ; ++j){
00840         atomNSD._localNum  = num + j ;
00841         //      std::cout <<   indexCUT[j] << " " ;
00842         _atomsInZone[indexCUT[j]].push_back(atomNSD);
00843       }
00844       //      std::cout << std::endl;
00845       num += indexCUT.size() ;
00846     }
00847     //    std::cout << "END subdomain "<<  i << std::endl; 
00848   }
00849 }
00850 template <class TYPEOFATOMS>
00851 bool PartitionBase<TYPEOFATOMS>::setAtomsInOverlapAndCut(){
00852   //
00853   bool retVal = false; 
00854   // Ajout des atomes des zones overlap et cut :
00855   int k, nsd, nsd_k = 0,  atomJ, fragNum ;
00856   double dist=0., distMin ;
00857   double R = _Rbuffer1 + _Rbuffer2, R1 =  _Rbuffer1* _Rbuffer1 ;
00858   //
00859   std::vector<bool> alreadyUsed ;
00860   //
00861   R = R*R ;
00862   //
00863   std::vector<int> numFragInSD,numFragInSD_k ;
00864   //
00865   for(nsd = _numberOfPartition -1 ; nsd > -1; --nsd) {
00866     //  for(nsd = 0 ; nsd <  _numberOfPartition  ; ++nsd) {
00867     //
00868     std::vector<int>& index      = _subDomains[nsd].getAtomsIndexInKernel();
00869     //
00870     // std::cout << " index du core " <<nsd <<std::endl; 
00871     // for (unsigned int i = 0 ; i < index.size() ; ++i){    std::cout <<"  "<<index[i] ;} std::cout <<std::endl; 
00872     bool             isIn        ;
00873     if( _useFrag){
00874       numFragInSD.clear() ;
00875       for ( unsigned int i = 0 ; i < index.size() ; ++i){
00876         isIn       = false ;
00877         fragNum    = _setOfAtoms->getElement(index[i]).getFragmentNumber() -1; 
00878         for ( unsigned int k = 0 ; k < numFragInSD.size() ; ++k){
00879           if(fragNum == numFragInSD[k]){ isIn = true ; break ;}
00880         }
00881         if(!isIn)  { 
00882           numFragInSD.push_back(fragNum) ; 
00883         }
00884       }
00885     }
00886     //
00887     for (k = 0 ; k < _subDomains[nsd].getNumberOfNeighbors() ; ++k){ 
00888       nsd_k = _subDomains[nsd].getNumberOfNeighbor(k) ;
00889       //       std::cout << "    nsd_k : " << nsd_k << std::endl;
00890       std::vector<int>& indexZONE1 = _subDomains[nsd].getAtomsIndexInOverlap(k);
00891       std::vector<int>& indexZONE2 = _subDomains[nsd].getAtomsIndexInCut(k);
00892       indexZONE1.clear(); 
00893       indexZONE2.clear(); 
00894       //
00895       std::vector<int>& index_k    = _subDomains[nsd_k].getAtomsIndexInKernel();
00896       //      std::cout << "     index du core " <<nsd_k <<"  :   "; 
00897       //      for (unsigned int i = 0 ; i < index_k.size() ; ++i){    std::cout <<"  "<<index_k[i] ;} std::cout <<std::endl; 
00898       //
00899       // Boucle sur les atomes
00900       //
00901       if( !_useFrag){
00902         for ( unsigned int j = 0 ; j < index_k.size() ; ++j){
00903           atomJ               = index_k[j] ;
00904           distMin             = 9999. ;
00905           for ( unsigned int i = 0 ; i < index.size() ; ++i){
00906             dist    = this->_setOfAtoms->getElement(atomJ).dist2To(this->_setOfAtoms->getElement(index[i]).pos() );
00907             if (distMin > dist){
00908               distMin = dist ; 
00909             }
00910           } // for i
00911           if( distMin > R){
00912             continue ; }
00913           else {
00914             if (distMin < R1 ){
00915               indexZONE1.push_back(atomJ) ;
00916             }
00917             else {
00918               indexZONE2.push_back( atomJ) ;
00919             }
00920           }
00921         } // for j
00922       }
00923       else {
00924         std::cout << "Construction avec des fragments " <<std::endl;
00925         alreadyUsed.resize(index_k.size()) ;
00926         numFragInSD_k.clear() ;
00927         for ( unsigned int j = 0 ; j < index_k.size() ; ++j){
00928           alreadyUsed[j] = false ;
00929           isIn           = false ;
00930           fragNum        = _setOfAtoms->getElement(index_k[j]).getFragmentNumber() -1; 
00931           for ( unsigned int k = 0 ; k < numFragInSD_k.size() ; ++k){
00932             if(fragNum == numFragInSD_k[k]){ isIn = true ; break ;}
00933           }
00934           if(!isIn)  { 
00935             numFragInSD_k.push_back(fragNum) ; 
00936           }
00937         }
00938         //
00939         //
00940 //      std::cout << " nsd   " << nsd << " nbFrag   "<<  numFragInSD.size() <<"  : ";
00941 //      for ( unsigned int frag = 0 ; frag <  numFragInSD.size()  ; ++frag){
00942 //        std::cout << numFragInSD[frag] << " " ;}
00943 //      std::cout << std::endl;
00944 //      std::cout   <<  " nsd_k " << nsd_k << " nbFrag_k "<<  numFragInSD_k.size() << "  : ";
00945 //      for ( unsigned int frag = 0 ; frag <  numFragInSD_k.size()  ; ++frag){
00946 //        std::cout << numFragInSD_k[frag] << " " ;}
00947 //      std::cout << std::endl<< std::endl;
00948         //
00949         //
00950         //
00951         int atIdx, atIdx_k, fragNum,fragNum_k,atIdxToAdd ;
00952         bool isInShell1, isInShell2 ;
00953         for ( unsigned int frag = 0 ; frag <  numFragInSD.size()  ; ++frag){
00954           fragNum = numFragInSD[frag] ;
00955           atIdx   = _firstATofFrag[fragNum]  ;
00956           //      std::cout << " numFragInSD_k.size()  : " << numFragInSD_k.size() <<std::endl;
00957           for ( unsigned int frag_k = 0 ; frag_k <  numFragInSD_k.size()  ; ++frag_k){
00958           isInShell1 = isInShell2 = false ;
00959           fragNum_k  = numFragInSD_k[frag_k] ;
00960           //
00961           while ( atIdx != -1 && !isInShell1){   
00962             atIdx_k   = _firstATofFrag[fragNum_k]  ;
00963 //          std::cout << " fragNum   " << fragNum   << " atIdx : " << atIdx 
00964 //                    << " fragNum_k " << fragNum_k << " :  " ;
00965             while ( atIdx_k != -1 && !isInShell1){   
00966               dist    = this->_setOfAtoms->getElement(atIdx).dist2To(this->_setOfAtoms->getElement(atIdx_k).pos() );
00967               if( dist < R){
00968                 if( dist < R1){
00969                   isInShell1 = true ;
00970                 }
00971                 else {
00972                   isInShell2 = true ;
00973                 }
00974               }
00975 //            std::cout << " atIdx : " << atIdx << " atIdx_k : " << atIdx_k << "  dist : "<< dist<<" isInShell1 "
00976 //                      << isInShell1<<" isInShell2 "<< isInShell2 << "  "<<true <<std::endl;
00977               //
00978               atIdx_k = _nextATofFrag[atIdx_k]  ;
00979             } //while atoms in frag_k
00980             //
00981             atIdx     = _nextATofFrag[atIdx]  ;
00982           } //while atoms in frag
00983           atIdx   = _firstATofFrag[fragNum]  ;
00984           if(isInShell1){ 
00985             atIdxToAdd = _firstATofFrag[fragNum_k] ;
00986             while ( atIdxToAdd != -1){   
00987               indexZONE1.push_back(atIdxToAdd) ;
00988               atIdxToAdd = _nextATofFrag[atIdxToAdd]  ;
00989             }
00990             std::cout << "On ajoute les atomes du fragment "<< frag_k<< " dans le shell 1 du fragment "<<frag <<std::endl; 
00991             break ;
00992           }
00993           if(isInShell2){ 
00994             atIdxToAdd = _firstATofFrag[fragNum_k] ;
00995             while ( atIdxToAdd != -1){   
00996               indexZONE2.push_back(atIdxToAdd) ;
00997               atIdxToAdd = _nextATofFrag[atIdxToAdd]  ;
00998             }
00999             std::cout << "On ajoute les atomes du fragment "<< frag_k<< " dans le shell 2  du fragment "<<frag <<std::endl;
01000             break ;
01001           }
01002           }// frag_k
01003         }  // frag
01004       } //else
01005     }  //for nsd_k
01006   } //for nsd
01007   //
01008   return retVal;
01009 }
01010 
01011 //
01012 //
01013 template class PartitionBase<PDBstructure> ;
01014 template class PartitionBase<setOfAtoms> ;
01015 
01016 

Generated on Sat Jan 28 21:07:25 2006 for QC++ by  doxygen 1.4.4