Zoltan 2 Version 0.5
Zoltan2_PartitioningProblem.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 //
00003 // ***********************************************************************
00004 //
00005 //   Zoltan2: A package of combinatorial algorithms for scientific computing
00006 //                  Copyright 2012 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Karen Devine      (kddevin@sandia.gov)
00039 //                    Erik Boman        (egboman@sandia.gov)
00040 //                    Siva Rajamanickam (srajama@sandia.gov)
00041 //
00042 // ***********************************************************************
00043 //
00044 // @HEADER
00045 
00050 #ifndef _ZOLTAN2_PARTITIONINGPROBLEM_HPP_
00051 #define _ZOLTAN2_PARTITIONINGPROBLEM_HPP_
00052 
00053 #include <Zoltan2_Problem.hpp>
00054 #include <Zoltan2_PartitioningAlgorithms.hpp>
00055 #include <Zoltan2_PartitioningSolution.hpp>
00056 #include <Zoltan2_PartitioningSolutionQuality.hpp>
00057 #include <Zoltan2_GraphModel.hpp>
00058 #include <Zoltan2_IdentifierModel.hpp>
00059 #include <Zoltan2_IntegerRangeList.hpp>
00060 
00061 #ifdef HAVE_ZOLTAN2_OVIS
00062 #include <ovis.h>
00063 #endif
00064 
00065 namespace Zoltan2{
00066 
00086 template<typename Adapter>
00087 class PartitioningProblem : public Problem<Adapter>
00088 {
00089 public:
00090 
00091   typedef typename Adapter::scalar_t scalar_t;
00092   typedef typename Adapter::gid_t gid_t;
00093   typedef typename Adapter::gno_t gno_t;
00094   typedef typename Adapter::lno_t lno_t;
00095   typedef typename Adapter::user_t user_t;
00096   typedef typename Adapter::base_adapter_t base_adapter_t;
00097 
00098 #ifdef HAVE_ZOLTAN2_MPI
00099   typedef Teuchos::OpaqueWrapper<MPI_Comm> mpiWrapper_t;
00100 #endif
00101 
00104   ~PartitioningProblem() {};
00105 
00106 #ifdef HAVE_ZOLTAN2_MPI
00107 
00110   PartitioningProblem(Adapter *A, Teuchos::ParameterList *p, MPI_Comm comm); 
00111 
00112 #endif
00113 
00115   PartitioningProblem(Adapter *A, Teuchos::ParameterList *p) ;
00116 
00118   //
00119   //    \param updateInputData   If true this indicates that either
00120   //          this is the first attempt at solution, or that we
00121   //          are computing a new solution and the input data has
00122   //          changed since the previous solution was computed.  
00123   //          By input data we mean coordinates, topology, or weights.
00124   //          If false, this indicates that we are computing a
00125   //          new solution using the same input data was used for
00126   //          the previous solution, even though the parameters
00127   //          may have been changed.
00128   //
00129   //  For the sake of performance, we ask the caller to set \c updateInputData
00130   //  to false if he/she is computing a new solution using the same input data,
00131   //  but different problem parameters, than that which was used to compute 
00132   //  the most recent solution.
00133 
00134   void solve(bool updateInputData=true );
00135  
00137   //
00138   //   \return  a reference to the solution to the most recent solve().
00139 
00140   const PartitioningSolution<Adapter> &getSolution() {
00141     return *(solution_.getRawPtr());
00142   };
00143 
00152   const scalar_t getImbalance(int dim=0) const {
00153     scalar_t imb = 0;
00154     if (!metrics_.is_null())
00155       metrics_->getWeightImbalance(imb, dim);
00156 
00157     return imb;
00158   }
00159 
00164   ArrayRCP<const MetricValues<scalar_t> > getMetrics() const {
00165    if (metrics_.is_null()){
00166       ArrayRCP<const MetricValues<scalar_t> > emptyMetrics;
00167       return emptyMetrics;
00168     }
00169     else
00170       return metrics_->getMetrics();
00171   }
00172 
00178   void printMetrics(ostream &os) const {
00179     if (metrics_.is_null())
00180       os << "No metrics available." << endl;
00181     else
00182       metrics_->printMetrics(os);
00183   };
00184 
00223   void setPartSizes(int len, partId_t *partIds, scalar_t *partSizes, 
00224     bool makeCopy=true) 
00225   { 
00226     setPartSizesForCritiera(0, len, partIds, partSizes, makeCopy);
00227   }
00228 
00263   void setPartSizesForCritiera(int criteria, int len, partId_t *partIds, 
00264     scalar_t *partSizes, bool makeCopy=true) ;
00265 
00268   void resetParameters(ParameterList *params)
00269   {
00270     Problem<Adapter>::resetParameters(params);  // creates new environment
00271     if (timer_.getRawPtr() != NULL)
00272       this->env_->setTimer(timer_);
00273   }
00274 
00279   const RCP<const Environment> & getEnvironment() const 
00280   {
00281     return this->envConst_;
00282   }
00283 
00284 private:
00285   void initializeProblem();
00286 
00287   void createPartitioningProblem(bool newData);
00288 
00289   RCP<PartitioningSolution<Adapter> > solution_;
00290 
00291   RCP<Comm<int> > problemComm_;
00292   RCP<const Comm<int> > problemCommConst_;
00293 
00294 #ifdef HAVE_ZOLTAN2_MPI
00295   MPI_Comm mpiComm_;
00296 #endif
00297 
00298   InputAdapterType inputType_;
00299   ModelType modelType_;
00300   modelFlag_t graphFlags_;
00301   modelFlag_t idFlags_;
00302   modelFlag_t coordFlags_;
00303   std::string algorithm_;
00304 
00305   int numberOfWeights_;
00306 
00307   // Suppose Array<partId_t> partIds = partIds_[w].  If partIds.size() > 0
00308   // then the user supplied part sizes for weight index "w", and the sizes
00309   // corresponding to the Ids in partIds are partSizes[w].
00310   //
00311   // If numberOfWeights_ >= 0, then there is an Id and Sizes array for
00312   // for each weight.  Otherwise the user did not supply object weights,
00313   // but they can still specify part sizes. 
00314   // So numberOfCriteria_ is numberOfWeights_ or one, whichever is greater.
00315 
00316   ArrayRCP<ArrayRCP<partId_t> > partIds_;
00317   ArrayRCP<ArrayRCP<scalar_t> > partSizes_;
00318   int numberOfCriteria_;
00319 
00320   // Number of parts to be computed at each level in hierarchical partitioning.
00321   
00322   ArrayRCP<int> levelNumberParts_;
00323   bool hierarchical_;
00324 
00325   // Create a Timer if the user asked for timing stats.
00326 
00327   RCP<TimerManager> timer_;
00328 
00329   // Did the user request metrics?
00330 
00331   bool metricsRequested_;
00332   RCP<const PartitioningSolutionQuality<Adapter> > metrics_;
00333 };
00335 
00336 #ifdef HAVE_ZOLTAN2_MPI
00337 template <typename Adapter>
00338   PartitioningProblem<Adapter>::PartitioningProblem(Adapter *A, 
00339     ParameterList *p, MPI_Comm comm):
00340       Problem<Adapter>(A,p,comm), solution_(),
00341       problemComm_(), problemCommConst_(),
00342       inputType_(InvalidAdapterType), modelType_(InvalidModel), 
00343       graphFlags_(), idFlags_(), coordFlags_(), algorithm_(),
00344       numberOfWeights_(), partIds_(), partSizes_(), 
00345       numberOfCriteria_(), levelNumberParts_(), hierarchical_(false), 
00346       timer_(), metricsRequested_(false), metrics_()
00347 {
00348 
00349   initializeProblem();
00350 }
00351 #endif
00352 
00353 template <typename Adapter>
00354   PartitioningProblem<Adapter>::PartitioningProblem(Adapter *A, 
00355     ParameterList *p):
00356       Problem<Adapter>(A,p), solution_(),
00357       problemComm_(), problemCommConst_(),
00358       inputType_(InvalidAdapterType), modelType_(InvalidModel), 
00359       graphFlags_(), idFlags_(), coordFlags_(), algorithm_(),
00360       numberOfWeights_(), 
00361       partIds_(), partSizes_(), numberOfCriteria_(), 
00362       levelNumberParts_(), hierarchical_(false), timer_(),
00363       metricsRequested_(false), metrics_()
00364 {
00365   initializeProblem();
00366 }
00367 
00368 template <typename Adapter>
00369   void PartitioningProblem<Adapter>::initializeProblem()
00370 {
00371   HELLO;
00372 
00373   this->env_->debug(DETAILED_STATUS, "PartitioningProblem::initializeProblem");
00374 
00375   if (getenv("DEBUGME")){
00376     std::cout << getpid() << std::endl;
00377     sleep(15);
00378   }
00379 
00380 #ifdef HAVE_ZOLTAN2_OVIS
00381   ovis_enabled(this->comm_->getRank());
00382 #endif
00383 
00384   // Create a copy of the user's communicator.  
00385 
00386   problemComm_ = this->comm_->duplicate();
00387   problemCommConst_ = rcp_const_cast<const Comm<int> > (problemComm_);
00388 
00389 #ifdef HAVE_ZOLTAN2_MPI
00390 
00391   // TPLs may want an MPI communicator
00392 
00393   mpiComm_ = Teuchos2MPI(problemComm_);
00394 
00395 #endif
00396 
00397   // Number of criteria is number of user supplied weights if non-zero.
00398   // Otherwise it is 1 and uniform weight is implied.
00399 
00400   numberOfWeights_ = this->inputAdapter_->getNumberOfWeightsPerObject();
00401 
00402   numberOfCriteria_ = (numberOfWeights_ > 1) ? numberOfWeights_ : 1;
00403 
00404   inputType_ = this->inputAdapter_->inputAdapterType();
00405 
00406   // The Caller can specify part sizes in setPartSizes().  If he/she
00407   // does not, the part size arrays are empty.
00408 
00409   ArrayRCP<partId_t> *noIds = new ArrayRCP<partId_t> [numberOfCriteria_];
00410   ArrayRCP<scalar_t> *noSizes = new ArrayRCP<scalar_t> [numberOfCriteria_];
00411 
00412   partIds_ = arcp(noIds, 0, numberOfCriteria_, true);
00413   partSizes_ = arcp(noSizes, 0, numberOfCriteria_, true);
00414 
00415   if (this->env_->getDebugLevel() >= DETAILED_STATUS){
00416     ostringstream msg;
00417     msg << problemComm_->getSize() << " procs,"
00418       << numberOfWeights_ << " user-defined weights, "
00419       << this->inputAdapter_->inputAdapterName() 
00420       << " input adapter type.\n";
00421     this->env_->debug(DETAILED_STATUS, msg.str());
00422   }
00423 
00424   this->env_->memory("After initializeProblem");
00425 }
00426 
00427 template <typename Adapter>
00428   void PartitioningProblem<Adapter>::setPartSizesForCritiera(
00429     int criteria, int len, partId_t *partIds, scalar_t *partSizes, bool makeCopy) 
00430 {
00431   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid length", 
00432     len>= 0, BASIC_ASSERTION);
00433 
00434   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid criteria", 
00435     criteria >= 0 && criteria < numberOfCriteria_, BASIC_ASSERTION);
00436 
00437   if (len == 0){
00438     partIds_[criteria] = ArrayRCP<partId_t>();
00439     partSizes_[criteria] = ArrayRCP<scalar_t>();
00440     return;
00441   }
00442 
00443   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid arrays", 
00444     partIds && partSizes, BASIC_ASSERTION);
00445 
00446   // The global validity of the partIds and partSizes arrays is performed
00447   // by the PartitioningSolution, which computes global part distribution and
00448   // part sizes.
00449 
00450   partId_t *z2_partIds = NULL;
00451   scalar_t *z2_partSizes = NULL;
00452   bool own_memory = false;
00453 
00454   if (makeCopy){
00455     z2_partIds = new partId_t [len];
00456     z2_partSizes = new scalar_t [len];
00457     this->env_->localMemoryAssertion(__FILE__, __LINE__, len, z2_partSizes);
00458     memcpy(z2_partIds, partIds, len * sizeof(partId_t));
00459     memcpy(z2_partSizes, partSizes, len * sizeof(scalar_t));
00460     own_memory=true;
00461   }
00462   else{
00463     z2_partIds = partIds;
00464     z2_partSizes = partSizes;
00465   }
00466 
00467   partIds_[criteria] = arcp(z2_partIds, 0, len, own_memory);
00468   partSizes_[criteria] = arcp(z2_partSizes, 0, len, own_memory);
00469 }
00470 
00471 template <typename Adapter>
00472 void PartitioningProblem<Adapter>::solve(bool updateInputData)
00473 {
00474   HELLO;
00475   this->env_->debug(DETAILED_STATUS, "Entering solve");
00476 
00477   // Create the computational model.
00478 
00479   this->env_->timerStart(MACRO_TIMERS, "create problem");
00480 
00481   createPartitioningProblem(updateInputData);
00482 
00483   this->env_->timerStop(MACRO_TIMERS, "create problem");
00484 
00485   // TODO: If hierarchical_
00486 
00487   // Create the solution. The algorithm will query the Solution
00488   //   for part and weight information. The algorithm will
00489   //   update the solution with part assignments and quality
00490   //   metrics.  The Solution object itself will convert our internal
00491   //   global numbers back to application global Ids if needed.
00492 
00493   RCP<const IdentifierMap<user_t> > idMap = 
00494     this->baseModel_->getIdentifierMap();
00495 
00496   PartitioningSolution<Adapter> *soln = NULL;
00497 
00498   this->env_->timerStart(MACRO_TIMERS, "create solution");
00499 
00500   try{
00501     soln = new PartitioningSolution<Adapter>( 
00502       this->envConst_, problemCommConst_, idMap, numberOfWeights_, 
00503       partIds_.view(0, numberOfCriteria_), 
00504       partSizes_.view(0, numberOfCriteria_));
00505   }
00506   Z2_FORWARD_EXCEPTIONS;
00507 
00508   solution_ = rcp(soln);
00509 
00510   this->env_->timerStop(MACRO_TIMERS, "create solution");
00511 
00512   this->env_->memory("After creating Solution");
00513 
00514   // Call the algorithm
00515 
00516   try {
00517     if (algorithm_ == string("scotch")){
00518       AlgPTScotch<Adapter>(this->envConst_, problemComm_,
00519 #ifdef HAVE_ZOLTAN2_MPI
00520         mpiComm_,
00521 #endif
00522         this->graphModel_, solution_);
00523     }
00524     else if (algorithm_ == string("block")){
00525       AlgBlock<Adapter>(this->envConst_, problemComm_,
00526         this->identifierModel_, solution_);
00527     }
00528     else if (algorithm_ == string("rcb")){
00529       AlgRCB<Adapter>(this->envConst_, problemComm_,
00530         this->coordinateModel_, solution_);
00531     }
00532     else if (algorithm_ == string("multijagged")){
00533       AlgPQJagged<Adapter>(this->envConst_, problemComm_,
00534               this->coordinateModel_, solution_);
00535     }
00536     else{
00537       throw std::logic_error("partitioning algorithm not supported yet");
00538     }
00539   }
00540   Z2_FORWARD_EXCEPTIONS;
00541 
00542 #ifdef HAVE_ZOLTAN2_MPI
00543 
00544   // The algorithm may have changed the communicator.  Change it back.
00545   // KDD:  Why would the algorithm change the communicator? TODO
00546   // KDD:  Should we allow such a side effect? TODO
00547 
00548   RCP<const mpiWrapper_t > wrappedComm = rcp(new mpiWrapper_t(mpiComm_));
00549   problemComm_ = rcp(new Teuchos::MpiComm<int>(wrappedComm));
00550   problemCommConst_ = rcp_const_cast<const Comm<int> > (problemComm_);
00551 
00552 #endif
00553 
00554   if (metricsRequested_){
00555     typedef PartitioningSolution<Adapter> ps_t;
00556     typedef PartitioningSolutionQuality<Adapter> psq_t;
00557 
00558     psq_t *quality = NULL;
00559     RCP<const ps_t> solutionConst = rcp_const_cast<const ps_t>(solution_);
00560     RCP<const Adapter> adapter = rcp(this->inputAdapter_, false);
00561 
00562     try{
00563       quality = new psq_t(this->envConst_, problemCommConst_, adapter, 
00564         solutionConst);
00565     }
00566     Z2_FORWARD_EXCEPTIONS
00567 
00568     metrics_ = rcp(quality);
00569   }
00570 
00571   this->env_->debug(DETAILED_STATUS, "Exiting solve");
00572 }
00573 
00574 template <typename Adapter>
00575 void PartitioningProblem<Adapter>::createPartitioningProblem(bool newData)
00576 {
00577   HELLO;
00578   this->env_->debug(DETAILED_STATUS, 
00579     "PartitioningProblem::createPartitioningProblem");
00580 
00581   using std::string;
00582   using Teuchos::ParameterList;
00583 
00584   // A Problem object may be reused.  The input data may have changed and
00585   // new parameters or part sizes may have been set.
00586   //
00587   // Save these values in order to determine if we need to create a new model.
00588 
00589   ModelType previousModel = modelType_;
00590   modelFlag_t previousGraphModelFlags = graphFlags_;
00591   modelFlag_t previousIdentifierModelFlags = idFlags_;
00592   modelFlag_t previousCoordinateModelFlags = coordFlags_;
00593 
00594   modelType_ = InvalidModel;
00595   graphFlags_.reset();
00596   idFlags_.reset();
00597   coordFlags_.reset();
00598 
00600   // It's possible at this point that the Problem may want to
00601   // add problem parameters to the parameter list in the Environment. 
00602   //
00603   // Since the parameters in the Environment have already been
00604   // validated in its constructor, a new Environment must be created:
00606   // Teuchos::RCP<const Teuchos::Comm<int> > oldComm = this->env_->comm_;
00607   // const ParameterList &oldParams = this->env_->getUnvalidatedParameters();
00608   // 
00609   // ParameterList newParams = oldParams;
00610   // newParams.set("new_parameter", "new_value");
00611   // 
00612   // ParameterList &newPartParams = newParams.sublist("partitioning");
00613   // newPartParams.set("new_partitioning_parameter", "its_value");
00614   // 
00615   // this->env_ = rcp(new Environment(newParams, oldComm));
00617 
00618   Environment &env = *(this->env_);
00619   ParameterList &pl = env.getParametersNonConst();
00620 
00621   string defString("default");
00622 
00623   // Did the user ask for computation of quality metrics?
00624 
00625   int yesNo=0;
00626   const Teuchos::ParameterEntry *pe = pl.getEntryPtr("compute_metrics");
00627   if (pe){
00628     yesNo = pe->getValue<int>(&yesNo);
00629     metricsRequested_ = true;
00630   }
00631 
00632   // Did the user specify a computational model?
00633 
00634   string model(defString);
00635   pe = pl.getEntryPtr("model");
00636   if (pe)
00637     model = pe->getValue<string>(&model);
00638 
00639   // Did the user specify an algorithm?
00640 
00641   string algorithm(defString);
00642   pe = pl.getEntryPtr("algorithm");
00643   if (pe)
00644     algorithm = pe->getValue<string>(&algorithm);
00645 
00646   // Possible algorithm requirements that must be conveyed to the model:
00647 
00648   bool needConsecutiveGlobalIds = false;
00649   bool removeSelfEdges= false;
00650 
00652   // Determine algorithm, model, and algorithm requirements.  This
00653   // is a first pass.  Feel free to change this and add to it.
00654   
00655   if (algorithm != defString){
00656     // Figure out the model required by the algorithm
00657     if (algorithm == string("block") ||
00658         algorithm == string("random") ||
00659         algorithm == string("cyclic") ){
00660 
00661       modelType_ = IdentifierModelType;
00662       algorithm_ = algorithm;
00663       needConsecutiveGlobalIds = true;
00664     }
00665     else if (algorithm == string("rcb") ||
00666              algorithm == string("rib") ||
00667              algorithm == string("multijagged") ||
00668              algorithm == string("hsfc")){
00669 
00670       modelType_ = CoordinateModelType;
00671       algorithm_ = algorithm;
00672     }
00673     else if (algorithm == string("metis") ||
00674              algorithm == string("parmetis") ||
00675              algorithm == string("scotch") ||
00676              algorithm == string("ptscotch")){
00677 
00678       modelType_ = GraphModelType;
00679       algorithm_ = algorithm;
00680       removeSelfEdges = true;
00681       needConsecutiveGlobalIds = true;
00682     }
00683     else if (algorithm == string("patoh") ||
00684              algorithm == string("phg")){
00685 
00686       if ((modelType_ != GraphModelType) &&
00687           (modelType_ != HypergraphModelType) ){
00688         modelType_ = HypergraphModelType;
00689       }
00690       algorithm_ = algorithm;
00691       needConsecutiveGlobalIds = true;
00692     }
00693     else{
00694       // Parameter list should ensure this does not happen.
00695       throw std::logic_error("parameter list algorithm is invalid");
00696     }
00697   }
00698   else if (model != defString){
00699     // Figure out the algorithm suggested by the model.
00700     if (model == string("hypergraph")){
00701       modelType_ = HypergraphModelType;
00702       if (problemComm_->getSize() > 1)
00703         algorithm_ = string("phg"); 
00704       else
00705         algorithm_ = string("patoh"); 
00706       needConsecutiveGlobalIds = true;
00707     }
00708     else if (model == string("graph")){
00709       modelType_ = GraphModelType;
00710 #ifdef HAVE_ZOLTAN2_SCOTCH
00711       if (problemComm_->getSize() > 1)
00712         algorithm_ = string("ptscotch"); 
00713       else
00714         algorithm_ = string("scotch"); 
00715       removeSelfEdges = true;
00716       needConsecutiveGlobalIds = true;
00717 #else
00718 #ifdef HAVE_ZOLTAN2_PARMETIS
00719       if (problemComm_->getSize() > 1)
00720         algorithm_ = string("parmetis"); 
00721       else
00722         algorithm_ = string("metis"); 
00723       removeSelfEdges = true;
00724       needConsecutiveGlobalIds = true;
00725 #else
00726       if (problemComm_->getSize() > 1)
00727         algorithm_ = string("phg"); 
00728       else
00729         algorithm_ = string("patoh"); 
00730       removeSelfEdges = true;
00731       needConsecutiveGlobalIds = true;
00732 #endif
00733 #endif
00734     }
00735     else if (model == string("geometry")){
00736       modelType_ = CoordinateModelType;
00737       algorithm_ = string("rcb");
00738     }
00739     else if (model == string("ids")){
00740       modelType_ = IdentifierModelType;
00741       algorithm_ = string("block");
00742       needConsecutiveGlobalIds = true;
00743     }
00744     else{
00745       // Parameter list should ensure this does not happen.
00746       env.localBugAssertion(__FILE__, __LINE__, 
00747         "parameter list model type is invalid", 1, BASIC_ASSERTION);
00748     }
00749   }
00750   else{   
00751     // Determine an algorithm and model suggested by the input type.
00752     //   TODO: this is a good time to use the time vs. quality parameter
00753     //     in choosing an algorithm, and setting some parameters
00754 
00755     if (inputType_ == MatrixAdapterType){
00756       modelType_ = HypergraphModelType;
00757       if (problemComm_->getSize() > 1)
00758         algorithm_ = string("phg"); 
00759       else
00760         algorithm_ = string("patoh"); 
00761     }
00762     else if (inputType_ == GraphAdapterType ||
00763         inputType_ == MeshAdapterType){
00764       modelType_ = GraphModelType;
00765       if (problemComm_->getSize() > 1)
00766         algorithm_ = string("phg"); 
00767       else
00768         algorithm_ = string("patoh"); 
00769     }
00770     else if (inputType_ == CoordinateAdapterType){
00771       modelType_ = CoordinateModelType;
00772       if(algorithm_ != string("multijagged"))
00773       algorithm_ = string("rcb");
00774     }
00775     else if (inputType_ == VectorAdapterType ||
00776              inputType_ == IdentifierAdapterType){
00777       modelType_ = IdentifierModelType;
00778       algorithm_ = string("block");
00779     }
00780     else{
00781       // This should never happen
00782       throw std::logic_error("input type is invalid");
00783     }
00784   }
00785 
00786   // Hierarchical partitioning?
00787 
00788   Array<int> valueList;
00789   pe = pl.getEntryPtr("topology");
00790 
00791   if (pe){
00792     valueList = pe->getValue<Array<int> >(&valueList);
00793 
00794     if (!Zoltan2::noValuesAreInRangeList<int>(valueList)){
00795       int *n = new int [valueList.size() + 1];
00796       levelNumberParts_ = arcp(n, 0, valueList.size() + 1, true);
00797       int procsPerNode = 1;
00798       for (int i=0; i < valueList.size(); i++){
00799         levelNumberParts_[i+1] = valueList[i];
00800         procsPerNode *= valueList[i];
00801       }
00802       // Number of parts in the first level
00803       levelNumberParts_[0] = env.numProcs_ / procsPerNode;
00804 
00805       if (env.numProcs_ % procsPerNode > 0)
00806         levelNumberParts_[0]++;
00807     }
00808   }
00809   else{
00810     levelNumberParts_.clear();
00811   }
00812 
00813   hierarchical_ = levelNumberParts_.size() > 0;
00814 
00815   // Object to be partitioned? (rows, columns, etc)
00816 
00817   string objectOfInterest(defString);
00818   pe = pl.getEntryPtr("objects_to_partition");
00819   if (pe)
00820     objectOfInterest = pe->getValue<string>(&objectOfInterest);
00821 
00823   // Set model creation flags, if any.
00824 
00825   if (modelType_ == GraphModelType){
00826 
00827     // Any parameters in the graph sublist?
00828 
00829     string symParameter(defString);
00830     pe = pl.getEntryPtr("symmetrize_graph");
00831     if (pe){
00832       symParameter = pe->getValue<string>(&symParameter);
00833       if (symParameter == string("transpose"))
00834         graphFlags_.set(SYMMETRIZE_INPUT_TRANSPOSE);
00835       else if (symParameter == string("bipartite"))
00836         graphFlags_.set(SYMMETRIZE_INPUT_BIPARTITE);
00837     } 
00838 
00839     int sgParameter = 0;
00840     pe = pl.getEntryPtr("subset_graph");
00841     if (pe)
00842       sgParameter = pe->getValue<int>(&sgParameter);
00843 
00844     if (sgParameter == 1)
00845         graphFlags_.set(GRAPH_IS_A_SUBSET_GRAPH);
00846 
00847     // Any special behaviors required by the algorithm?
00848     
00849     if (removeSelfEdges)
00850       graphFlags_.set(SELF_EDGES_MUST_BE_REMOVED);
00851 
00852     if (needConsecutiveGlobalIds)
00853       graphFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
00854 
00855     // How does user input map to vertices and edges?
00856 
00857     if (inputType_ == MatrixAdapterType){
00858       if (objectOfInterest == defString ||
00859           objectOfInterest == string("matrix_rows") )
00860         graphFlags_.set(VERTICES_ARE_MATRIX_ROWS);
00861       else if (objectOfInterest == string("matrix_columns"))
00862         graphFlags_.set(VERTICES_ARE_MATRIX_COLUMNS);
00863       else if (objectOfInterest == string("matrix_nonzeros"))
00864         graphFlags_.set(VERTICES_ARE_MATRIX_NONZEROS);
00865     }
00866 
00867     else if (inputType_ == MeshAdapterType){
00868       if (objectOfInterest == defString ||
00869           objectOfInterest == string("mesh_nodes") )
00870         graphFlags_.set(VERTICES_ARE_MESH_NODES);
00871       else if (objectOfInterest == string("mesh_elements"))
00872         graphFlags_.set(VERTICES_ARE_MESH_ELEMENTS);
00873     } 
00874   }
00875   else if (modelType_ == IdentifierModelType){
00876 
00877     // Any special behaviors required by the algorithm?
00878     
00879     if (needConsecutiveGlobalIds)
00880       idFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
00881   }
00882   else if (modelType_ == CoordinateModelType){
00883 
00884     // Any special behaviors required by the algorithm?
00885     
00886     if (needConsecutiveGlobalIds)
00887       coordFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
00888   }
00889 
00890 
00891   if (  newData ||
00892        (modelType_ != previousModel) ||
00893        (graphFlags_ != previousGraphModelFlags) ||
00894        (coordFlags_ != previousCoordinateModelFlags) ||
00895        (idFlags_ != previousIdentifierModelFlags) ) {
00896 
00897     // Create the computational model.
00898     // Models are instantiated for base input adapter types (mesh,
00899     // matrix, graph, and so on).  We pass a pointer to the input
00900     // adapter, cast as the base input type.
00901 
00902     typedef typename Adapter::base_adapter_t base_adapter_t;
00903     const Teuchos::ParameterList pl = this->envConst_->getParameters();
00904     bool exceptionThrow = true;
00905 
00906     switch (modelType_) {
00907 
00908     case GraphModelType:
00909       this->graphModel_ = rcp(new GraphModel<base_adapter_t>(
00910         this->baseInputAdapter_, this->envConst_, problemComm_, graphFlags_));
00911 
00912       this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(
00913         this->graphModel_);
00914 
00915       break;
00916 
00917     case HypergraphModelType:
00918       break;
00919   
00920     case CoordinateModelType:
00921       this->coordinateModel_ = rcp(new CoordinateModel<base_adapter_t>(
00922         this->baseInputAdapter_, this->envConst_, problemComm_, coordFlags_));
00923 
00925       // It's possible at this point that the Problem may want to
00926       // add problem parameters to the parameter list in the Environment.
00927       //
00928       // Since the parameters in the Environment have already been
00929       // validated in its constructor, a new Environment must be created:
00931       // Teuchos::RCP<const Teuchos::Comm<int> > oldComm = this->env_->comm_;
00932       // const ParameterList &oldParams = this->env_->getUnvalidatedParameters();
00933       //
00934       // ParameterList newParams = oldParams;
00935       // newParams.set("new_parameter", "new_value");
00936       //
00937       // ParameterList &newPartParams = newParams.sublist("partitioning");
00938       // newPartParams.set("new_partitioning_parameter", "its_value");
00939       //
00940       // this->env_ = rcp(new Environment(newParams, oldComm));
00942 
00943       if(algorithm == string("multijagged")){
00944           //int coordinateCnt = this->coordinateModel_->getCoordinateDim();
00945           //cout << coordinateCnt << " " << pl.getPtr<Array <int> >("pqParts")->size() << endl;
00946           //exceptionThrow = coordinateCnt == pl.getPtr<Array <int> >("pqParts")->size();
00947           int arraySize = pl.getPtr<Array <int> >("pqParts")->size() - 1;
00948           exceptionThrow = arraySize > 0;
00949           this->envConst_->localInputAssertion(__FILE__, __LINE__, "invalid length of cut lines. Size of cut lines should be at least 1.",
00950                         exceptionThrow, BASIC_ASSERTION);
00951 
00952 
00953           int totalPartCount = 1;
00954           for(int i = 0; i < arraySize; ++i){
00955             //cout <<  pl.getPtr<Array <int> >("pqParts")->getRawPtr()[i] << " ";
00956             totalPartCount *= pl.getPtr<Array <int> >("pqParts")->getRawPtr()[i];
00957           }
00958           Teuchos::ParameterList newParams = pl;
00959           Teuchos::ParameterList &parParams = newParams.sublist("partitioning");
00960 
00961           parParams.set("num_global_parts", totalPartCount);
00962 
00963 
00964           //cout << endl;
00965           Teuchos::RCP<const Teuchos::Comm<int> > oldComm = this->envConst_->comm_;
00966 
00967           //this->envConst_ = rcp(new Environment(newParams, oldComm));
00968 
00969 
00970       }
00971 
00972 
00973       this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(
00974         this->coordinateModel_);
00975       break;
00976 
00977     case IdentifierModelType:
00978       this->identifierModel_ = rcp(new IdentifierModel<base_adapter_t>(
00979         this->baseInputAdapter_, this->envConst_, problemComm_, idFlags_));
00980 
00981       this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(
00982         this->identifierModel_);
00983       break;
00984 
00985     default:
00986       cout << __func__ << " Invalid model" << modelType_ << endl;
00987       break;
00988     }
00989 
00990     this->env_->memory("After creating Model");
00991   }
00992 
00993   /*
00994   Teuchos::RCP<const Teuchos::Comm<int> > oldComm = this->env_->comm_;
00995   const ParameterList &oldParams = this->env_->getUnvalidatedParameters();
00996 
00997   ParameterList newParams = oldParams;
00998   int totalPartCount = 1;
00999   const int *partNo = pl.getPtr<Array <int> >("pqParts")->getRawPtr();
01000 
01001   for (int i = 0; i < coordDim; ++i){
01002     totalPartCount *= partNo[i];
01003   }
01004   newParams.set("num_global_parts", totalPartCount);
01005 
01006 
01007   this->env_ = rcp(new Environment(newParams, oldComm));
01008 */
01009 }
01010 
01011 }  // namespace Zoltan2
01012 #endif