Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

runScal.C

Go to the documentation of this file.
00001 
00002 #include "mpi.h"
00003 #include "petsc.h"
00004 #include "sys.h"
00005 #include "octUtils.h"
00006 #include "TreeNode.h"
00007 #include "parUtils.h"
00008 #include <cstdlib>
00009 #include <cstring>
00010 #include "externVars.h"
00011 #include "dendro.h"
00012 
00013 #ifdef MPI_WTIME_IS_GLOBAL
00014 #undef MPI_WTIME_IS_GLOBAL
00015 #endif
00016 
00017 int main(int argc, char ** argv ) {     
00018   int size, rank;
00019   double startTime, endTime, minTime;
00020   unsigned int incInt = 1;
00021   bool incCorner = 1;  
00022   unsigned int numPts;
00023   unsigned int writePOut = 0;
00024   unsigned int writeBOut = 0;
00025   unsigned int readPtsFile = 0;
00026   unsigned int readOctFile =0;
00027   unsigned int runOpt =2;//2 for both, 1 for pt alone and 0 for bal alone.
00028   char Kstr[20];
00029   char inpFileName[50], p2nIn[50], p2nOut[50], n2oOut[50], balOut[50];
00030   double gSize[3];
00031 
00032   PetscInitialize(&argc,&argv,"options",NULL);
00033   ot::RegisterEvents();
00034 
00035 #ifdef PETSC_USE_LOG
00036   int stages[4];
00037   PetscLogStageRegister(&stages[0],"Prepare Input.");
00038   PetscLogStageRegister(&stages[1],"P2O");
00039   PetscLogStageRegister(&stages[2],"N2O");
00040   PetscLogStageRegister(&stages[3],"Bal");
00041 #endif
00042 
00043   MPI_Comm_size(MPI_COMM_WORLD,&size);
00044   MPI_Comm_rank(MPI_COMM_WORLD,&rank);
00045   if(argc < 3) {
00046     std::cerr << "Usage: " << argv[0] << " inpfile complx incInt[1] runOpt[2]"<<
00047       "readPtsFile[0] readOctFile[0] maxDepth[30] writePOut[0] writeBOut[0]"<<
00048       " dim[3] maxNumPtsPerOctant[1] incCorner[0] "  << std::endl;
00049     return -1;
00050   }
00051   unsigned int maxNumPts= 1;
00052   int complx = atoi(argv[2]);   
00053   unsigned int dim=3;
00054   unsigned int maxDepth=30;
00055   if(argc > 3) {
00056     incInt = atoi(argv[3]);
00057   }
00058   if(argc > 4) {
00059     runOpt = atoi(argv[4]);
00060   }
00061   if(argc > 5) {
00062     readPtsFile = atoi(argv[5]);
00063   }
00064   if(argc > 6) {
00065     readOctFile = atoi(argv[6]);
00066   }
00067   if(argc > 7) {
00068     maxDepth = atoi(argv[7]);
00069   }
00070   if(argc > 8) {
00071     writePOut = atoi(argv[8]);
00072   }
00073   if(argc > 9) {
00074     writeBOut = atoi(argv[9]);
00075   }
00076   if(argc > 10) {
00077     dim = atoi(argv[10]);
00078   }
00079   if(argc > 11) {
00080     maxNumPts = atoi(argv[11]);
00081   }
00082   if(argc > 12) { incCorner = (bool)(atoi(argv[12]));}
00083 
00084   strcpy(inpFileName, argv[1]);
00085   strcpy(balOut,inpFileName);
00086   ot::int2str(rank,Kstr);
00087   strcat(balOut,Kstr);
00088   strcat(balOut,"_\0");
00089   ot::int2str(size,Kstr);
00090   strcat(balOut,Kstr);
00091   strcpy(n2oOut,balOut);
00092   strcpy(p2nIn,balOut);
00093   strcpy(p2nOut,balOut);
00094   strcat(balOut,"_Bal.ot\0");
00095   strcat(n2oOut,"_Con.ot\0");
00096   strcat(p2nIn,".pts\0");
00097   strcat(p2nOut,"_Out.pts\0");
00098   strcat(inpFileName,".inp\0"); 
00099   std::vector<ot::TreeNode> nodes;
00100   std::vector<double> pts;
00101   unsigned int ptsLen;
00102 
00103   // print out the format ...
00104   /*
00105   if (!rank) {
00106     cout << "============================================================================================" << endl;
00107     cout << "#Points p2n_time p2n_imbal #cellsConstr #cellsBal bal_time bal_imbal maxLevBef minLevBef avgLevBef maxLevAft minLevAft avgLevAft" << endl;
00108     cout << "============================================================================================" << endl << endl;
00109   }
00110   */
00111 #ifdef PETSC_USE_LOG
00112   PetscLogStagePush(stages[0]);
00113 #endif
00114   if(runOpt == 1 || runOpt == 2)  {
00115     if(readPtsFile) {
00116       if(!rank){
00117         std::cout << " reading  "<<p2nIn<<std::endl; // Point size
00118       }
00119       ot::readPtsFromFile(p2nIn, pts);
00120       if(!rank){
00121         std::cout << " finished reading  "<<p2nIn<<std::endl; // Point size
00122       }
00123       ptsLen = pts.size();
00124       std::vector<ot::TreeNode> tmpNodes;
00125       for(int i = 0; i < ptsLen; i += 3) {
00126         if( (pts[i] > 0.0) &&
00127                  (pts[i+1] > 0.0)  
00128                 && (pts[i+2] > 0.0) &&
00129               ( ((unsigned int)(pts[i]*((double)(1u<<maxDepth)))) < (1u<<maxDepth))  &&
00130               ( ((unsigned int)(pts[i+1]*((double)(1u<<maxDepth)))) < (1u<<maxDepth))  &&
00131               ( ((unsigned int)(pts[i+2]*((double)(1u<<maxDepth)))) < (1u<<maxDepth)) ) {
00132         tmpNodes.push_back( ot::TreeNode((unsigned int)(pts[i]*(double)(1u<<maxDepth)),
00133             (unsigned int)(pts[i+1]*(double)(1u<<maxDepth)),
00134             (unsigned int)(pts[i+2]*(double)(1u<<maxDepth)),
00135             maxDepth,dim,maxDepth) );
00136         }
00137       }//end for i
00138       pts.clear();
00139     //  std::cout<<"Before Removing Duplicates, the number of pts: "<<tmpNodes.size()<<std::endl; 
00140       par::removeDuplicates<ot::TreeNode>(tmpNodes,false,MPI_COMM_WORLD);       
00141      // std::cout<<"After Removing Duplicates, the number of pts: "<<tmpNodes.size()<<std::endl; 
00142       nodes = tmpNodes;
00143       tmpNodes.clear();
00144       par::partitionW<ot::TreeNode>(nodes, NULL,MPI_COMM_WORLD);
00145       // reduce and only print the total ...
00146       DendroIntL locSz, totalSz;
00147       locSz = nodes.size();
00148 
00149       par::Mpi_Reduce<DendroIntL>(&locSz, &totalSz, 1, MPI_SUM, 0, MPI_COMM_WORLD);
00150       
00151       if(rank==0) {
00152       std::cout<<"rank= "<<rank<<" size= " << size << " total pts= " << totalSz<<std::endl;
00153       }
00154 
00155       pts.resize(3*(nodes.size()));
00156       ptsLen = (3*(nodes.size()));
00157       for(int i=0;i<nodes.size();i++) {
00158         pts[3*i] = (((double)(nodes[i].getX())) + 0.5)/((double)(1u<<maxDepth));
00159         pts[(3*i)+1] = (((double)(nodes[i].getY())) +0.5)/((double)(1u<<maxDepth));
00160         pts[(3*i)+2] = (((double)(nodes[i].getZ())) +0.5)/((double)(1u<<maxDepth));
00161       }//end for i
00162       nodes.clear();
00163     }//end if readPts
00164     gSize[0] = 1.0;
00165     gSize[1] = 1.0;
00166    gSize[2] = 1.0;
00167     if(writePOut) { 
00168       if(!rank) {
00169         std::cout<<"Writing pts to: "<<p2nOut<<std::endl;
00170       }
00171       ot::writePtsToFile(p2nOut,pts);
00172     }
00173  #ifdef PETSC_USE_LOG
00174   PetscLogStagePop();
00175 #endif
00176     MPI_Barrier(MPI_COMM_WORLD);        
00177 #ifdef PETSC_USE_LOG
00178   PetscLogStagePush(stages[1]);
00179 #endif
00180     startTime = MPI_Wtime();
00181     ot::points2Octree(pts, gSize, nodes, dim, maxDepth, maxNumPts, MPI_COMM_WORLD);
00182     endTime = MPI_Wtime();
00183 #ifdef PETSC_USE_LOG
00184   PetscLogStagePop();
00185 #endif
00186     double locTime, totalTime, minTime;
00187     locTime = endTime - startTime;
00188     par::Mpi_Reduce<double>(&locTime, &totalTime, 1, MPI_MAX, 0, MPI_COMM_WORLD);
00189     par::Mpi_Reduce<double>(&locTime, &minTime, 1, MPI_MIN, 0, MPI_COMM_WORLD);
00190     if(!rank){
00191       std::cout <<"P2n Time: "<<totalTime << " " << "secs imbalance: "<< totalTime/minTime << std::endl;
00192     }
00193     pts.clear();
00194   }
00195   //TestBuildOCt from here on:
00196   if(runOpt == 0 || runOpt == 2) {
00197     if(readOctFile) {
00198       ot::readNodesFromFile(n2oOut,nodes);
00199     }
00200     std::vector<ot::TreeNode > linOct;  
00201 #ifdef PETSC_USE_LOG
00202   PetscLogStagePush(stages[2]);
00203 #endif
00204     DendroIntL inputNodes = nodes.size();
00205     DendroIntL totInp;
00206     par::Mpi_Reduce<DendroIntL>(&inputNodes, &totInp, 1, MPI_SUM, 0, MPI_COMM_WORLD);
00207     if(!rank) {
00208       std::cout << "Nodes after p2o: "<< totInp  << std::endl;       
00209     }
00210     ot::completeOctree(nodes,linOct,dim,maxDepth,false, false, false, MPI_COMM_WORLD);
00211 #ifdef PETSC_USE_LOG
00212   PetscLogStagePop();
00213 #endif
00214 
00215     assert(!linOct.empty());
00216     if(writeBOut) { 
00217       if(!rank) {
00218         std::cout<<"Writing octree after p2o to: "<<n2oOut<<std::endl; 
00219       }
00220       ot::writeNodesToFile(n2oOut,linOct);
00221     }
00222     inputNodes = linOct.size();
00223     nodes.clear();
00224     std::vector<ot::TreeNode > balOct;
00225 
00226     par::Mpi_Reduce<DendroIntL>(&inputNodes, &totInp, 1, MPI_SUM, 0, MPI_COMM_WORLD);
00227     if(!rank) {
00228       std::cout << "Nodes before balancing: "<< totInp  << std::endl;       
00229     }
00230     MPI_Barrier(MPI_COMM_WORLD);        
00231 #ifdef PETSC_USE_LOG
00232   PetscLogStagePush(stages[3]);
00233 #endif
00234     startTime = MPI_Wtime();
00235 
00236     ot::balanceOctree (linOct, balOct, dim, maxDepth, incCorner, MPI_COMM_WORLD, NULL, NULL);
00237     endTime = MPI_Wtime();
00238 #ifdef PETSC_USE_LOG
00239   PetscLogStagePop();
00240 #endif
00241     assert(!balOct.empty());
00242     if(writeBOut) { 
00243       if(!rank) {
00244         std::cout<<"Writing octree after balancing to: "<<balOut<<std::endl; 
00245       }
00246       ot::writeNodesToFile(balOut,balOct);
00247     }
00248     // compute total inp size and output size
00249     DendroIntL locBalNodes = balOct.size();
00250     DendroIntL totBal;
00251     double balTime, totTime;
00252     balTime = endTime - startTime;
00253     par::Mpi_Reduce<DendroIntL>(&locBalNodes, &totBal, 1, MPI_SUM, 0, MPI_COMM_WORLD);
00254     par::Mpi_Reduce<double>(&balTime, &totTime, 1, MPI_MAX, 0, MPI_COMM_WORLD);
00255     par::Mpi_Reduce<double>(&balTime, &minTime, 1, MPI_MIN, 0, MPI_COMM_WORLD);
00256 
00257     if(!rank) {
00258       std::cout << "Nodes after balancing: "<< totBal << std::endl;       
00259       std::cout << "bal Time: "<<totTime << " " << "secs imbalance: "<< totTime/minTime << std::endl;
00260     }
00261     linOct.clear();
00262     balOct.clear();
00263   }
00264 
00265   PetscFinalize();
00266 
00267 }
00268 

Generated on Tue Mar 24 16:14:00 2009 for DENDRO by  doxygen 1.3.9.1