Zoltan2
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 #include <Zoltan2_MachineRepresentation.hpp>
00061 #include <Zoltan2_TaskMapping.hpp>
00062 
00063 #include <unistd.h>
00064 
00065 #ifdef HAVE_ZOLTAN2_OVIS
00066 #include <ovis.h>
00067 #endif
00068 
00069 namespace Zoltan2{
00070 
00094 template<typename Adapter>
00095 class PartitioningProblem : public Problem<Adapter>
00096 {
00097 public:
00098 
00099   typedef typename Adapter::scalar_t scalar_t;
00100   typedef typename Adapter::gid_t gid_t;
00101   typedef typename Adapter::gno_t gno_t;
00102   typedef typename Adapter::lno_t lno_t;
00103   typedef typename Adapter::part_t part_t;
00104   typedef typename Adapter::user_t user_t;
00105   typedef typename Adapter::base_adapter_t base_adapter_t;
00106 
00107 #ifdef HAVE_ZOLTAN2_MPI
00108   typedef Teuchos::OpaqueWrapper<MPI_Comm> mpiWrapper_t;
00109 #endif
00110 
00113   ~PartitioningProblem() {};
00114 
00115 #ifdef HAVE_ZOLTAN2_MPI
00116 
00119   PartitioningProblem(Adapter *A, Teuchos::ParameterList *p, MPI_Comm comm); 
00120 
00121 #endif
00122 
00124   PartitioningProblem(Adapter *A, Teuchos::ParameterList *p);
00125 
00127   //
00128   //    \param updateInputData   If true this indicates that either
00129   //          this is the first attempt at solution, or that we
00130   //          are computing a new solution and the input data has
00131   //          changed since the previous solution was computed.  
00132   //          By input data we mean coordinates, topology, or weights.
00133   //          If false, this indicates that we are computing a
00134   //          new solution using the same input data was used for
00135   //          the previous solution, even though the parameters
00136   //          may have been changed.
00137   //
00138   //  For the sake of performance, we ask the caller to set \c updateInputData
00139   //  to false if he/she is computing a new solution using the same input data,
00140   //  but different problem parameters, than that which was used to compute 
00141   //  the most recent solution.
00142 
00143   void solve(bool updateInputData=true );
00144  
00146   //          when a point lies on a part boundary, the lowest part
00147   //          number on that boundary is returned.
00148   //          Note that not all partitioning algorithms will support
00149   //          this method.
00150   //
00151   //   \param dim : the number of dimensions specified for the point in space
00152   //   \param point : the coordinates of the point in space; array of size dim
00153   //   \return the part number of a part overlapping the given point
00154   //
00155   //   TODO:  This method might be more appropriate in a partitioning solution
00156   //   TODO:  But then the solution would need to know about the algorithm.
00157   //   TODO:  Consider moving the algorithm to the partitioning solution.
00158   part_t pointAssign(int dim, scalar_t *point) 
00159   {
00160     part_t p;
00161     try {
00162       if (this->algorithm_ == Teuchos::null)
00163         throw std::logic_error("no partitioning algorithm has been run yet");
00164 
00165       p = this->algorithm_->pointAssign(dim, point); 
00166     }
00167     Z2_FORWARD_EXCEPTIONS
00168     return p;
00169   }
00170 
00172   //   This method allocates memory for the return argument, but does not
00173   //   control that memory.  The user is responsible for freeing the 
00174   //   memory.
00175   //
00176   //   \param dim : (in) the number of dimensions specified for the box
00177   //   \param lower : (in) the coordinates of the lower corner of the box; 
00178   //                   array of size dim
00179   //   \param upper : (in) the coordinates of the upper corner of the box; 
00180   //                   array of size dim
00181   //   \param nPartsFound : (out) the number of parts overlapping the box
00182   //   \param partsFound :  (out) array of parts overlapping the box
00183   //   TODO:  This method might be more appropriate in a partitioning solution
00184   //   TODO:  But then the solution would need to know about the algorithm.
00185   //   TODO:  Consider moving the algorithm to the partitioning solution.
00186   void boxAssign(int dim, scalar_t *lower, scalar_t *upper,
00187                  size_t &nPartsFound, part_t **partsFound) 
00188   {
00189     try {
00190       if (this->algorithm_ == Teuchos::null)
00191         throw std::logic_error("no partitioning algorithm has been run yet");
00192 
00193       this->algorithm_->boxAssign(dim, lower, upper, nPartsFound, partsFound); 
00194     }
00195     Z2_FORWARD_EXCEPTIONS
00196   }
00197 
00199   //
00200   //   \return  a reference to the solution to the most recent solve().
00201 
00202   const PartitioningSolution<Adapter> &getSolution() {
00203     return *(solution_.getRawPtr());
00204   };
00205 
00214   const scalar_t getWeightImbalance(int idx=0) const {
00215     scalar_t imb = 0;
00216     if (!metrics_.is_null())
00217       metrics_->getWeightImbalance(imb, idx);
00218 
00219     return imb;
00220   }
00221 
00226   ArrayRCP<const MetricValues<scalar_t> > getMetrics() const {
00227    if (metrics_.is_null()){
00228       ArrayRCP<const MetricValues<scalar_t> > emptyMetrics;
00229       return emptyMetrics;
00230     }
00231     else
00232       return metrics_->getMetrics();
00233   }
00234 
00240   void printMetrics(std::ostream &os) const {
00241     if (metrics_.is_null())
00242       os << "No metrics available." << endl;
00243     else
00244       metrics_->printMetrics(os);
00245   };
00246 
00285   void setPartSizes(int len, part_t *partIds, scalar_t *partSizes, 
00286     bool makeCopy=true) 
00287   { 
00288     setPartSizesForCriteria(0, len, partIds, partSizes, makeCopy);
00289   }
00290 
00325   void setPartSizesForCriteria(int criteria, int len, part_t *partIds,
00326     scalar_t *partSizes, bool makeCopy=true) ;
00327 /*
00328   void setMachine(MachineRepresentation<typename Adapter::base_adapter_t::scalar_t> *machine);
00329 */
00332   void resetParameters(ParameterList *params)
00333   {
00334     Problem<Adapter>::resetParameters(params);  // creates new environment
00335     if (timer_.getRawPtr() != NULL)
00336       this->env_->setTimer(timer_);
00337   }
00338 
00343   const RCP<const Environment> & getEnvironment() const 
00344   {
00345     return this->envConst_;
00346   }
00347 
00348 private:
00349   void initializeProblem();
00350 
00351   void createPartitioningProblem(bool newData);
00352 
00353   RCP<PartitioningSolution<Adapter> > solution_;
00354 
00355   RCP<MachineRepresentation <typename Adapter::base_adapter_t::scalar_t>  > machine_;
00356 
00357   RCP<Comm<int> > problemComm_;
00358   RCP<const Comm<int> > problemCommConst_;
00359 
00360 #ifdef HAVE_ZOLTAN2_MPI
00361   MPI_Comm mpiComm_;
00362 #endif
00363 
00364   BaseAdapterType inputType_;
00365 
00366   //ModelType modelType_;
00367   bool modelAvail_[MAX_NUM_MODEL_TYPES];
00368 
00369   modelFlag_t graphFlags_;
00370   modelFlag_t idFlags_;
00371   modelFlag_t coordFlags_;
00372   std::string algName_;
00373 
00374   int numberOfWeights_;
00375 
00376   // Suppose Array<part_t> partIds = partIds_[w].  If partIds.size() > 0
00377   // then the user supplied part sizes for weight index "w", and the sizes
00378   // corresponding to the Ids in partIds are partSizes[w].
00379   //
00380   // If numberOfWeights_ >= 0, then there is an Id and Sizes array for
00381   // for each weight.  Otherwise the user did not supply object weights,
00382   // but they can still specify part sizes. 
00383   // So numberOfCriteria_ is numberOfWeights_ or one, whichever is greater.
00384 
00385   ArrayRCP<ArrayRCP<part_t> > partIds_;
00386   ArrayRCP<ArrayRCP<scalar_t> > partSizes_;
00387   int numberOfCriteria_;
00388 
00389   // Number of parts to be computed at each level in hierarchical partitioning.
00390   
00391   ArrayRCP<int> levelNumberParts_;
00392   bool hierarchical_;
00393 
00394   // Create a Timer if the user asked for timing stats.
00395 
00396   RCP<TimerManager> timer_;
00397 
00398   // Did the user request metrics?
00399 
00400   bool metricsRequested_;
00401   RCP<const PartitioningSolutionQuality<Adapter> > metrics_;
00402 };
00404 
00405 #ifdef HAVE_ZOLTAN2_MPI
00406 template <typename Adapter>
00407   PartitioningProblem<Adapter>::PartitioningProblem(Adapter *A, 
00408     ParameterList *p, MPI_Comm comm):
00409       Problem<Adapter>(A,p,comm), solution_(),
00410       problemComm_(), problemCommConst_(),
00411       inputType_(InvalidAdapterType), 
00412       graphFlags_(), idFlags_(), coordFlags_(), algName_(),
00413       numberOfWeights_(), partIds_(), partSizes_(), 
00414       numberOfCriteria_(), levelNumberParts_(), hierarchical_(false), 
00415       timer_(), metricsRequested_(false), metrics_()
00416 {
00417 
00418   for(int i=0;i<MAX_NUM_MODEL_TYPES;i++)
00419   {
00420     modelAvail_[i]=false;
00421   }
00422 
00423   initializeProblem();
00424 
00425 }
00426 #endif
00427 /*
00428 template <typename Adapter>
00429 void PartitioningProblem<Adapter>::setMachine(MachineRepresentation<typename Adapter::base_adapter_t::scalar_t> *machine){
00430   this->machine_ = RCP<MachineRepresentation<typename Adapter::base_adapter_t::scalar_t> > (machine, false);
00431 }
00432 */
00433 template <typename Adapter>
00434   PartitioningProblem<Adapter>::PartitioningProblem(Adapter *A, 
00435     ParameterList *p):
00436       Problem<Adapter>(A,p), solution_(),
00437       problemComm_(), problemCommConst_(),
00438       inputType_(InvalidAdapterType), 
00439       graphFlags_(), idFlags_(), coordFlags_(), algName_(),
00440       numberOfWeights_(), 
00441       partIds_(), partSizes_(), numberOfCriteria_(), 
00442       levelNumberParts_(), hierarchical_(false), timer_(),
00443       metricsRequested_(false), metrics_()
00444 {
00445   for(int i=0;i<MAX_NUM_MODEL_TYPES;i++)
00446   {
00447     modelAvail_[i]=false;
00448   }
00449 
00450   initializeProblem();
00451 }
00452 
00453 template <typename Adapter>
00454   void PartitioningProblem<Adapter>::initializeProblem()
00455 {
00456   HELLO;
00457 
00458   this->env_->debug(DETAILED_STATUS, "PartitioningProblem::initializeProblem");
00459 
00460   if (getenv("DEBUGME")){
00461     std::cout << getpid() << std::endl;
00462     sleep(15);
00463   }
00464 
00465 #ifdef HAVE_ZOLTAN2_OVIS
00466   ovis_enabled(this->comm_->getRank());
00467 #endif
00468 
00469   // Create a copy of the user's communicator.  
00470 
00471   problemComm_ = this->comm_->duplicate();
00472   problemCommConst_ = rcp_const_cast<const Comm<int> > (problemComm_);
00473 
00474   machine_ = RCP <Zoltan2::MachineRepresentation<typename Adapter::scalar_t> >(new Zoltan2::MachineRepresentation<typename Adapter::scalar_t>(problemComm_));
00475 #ifdef HAVE_ZOLTAN2_MPI
00476 
00477   // TPLs may want an MPI communicator
00478 
00479   mpiComm_ = Teuchos2MPI(problemComm_);
00480 
00481 #endif
00482 
00483   // Number of criteria is number of user supplied weights if non-zero.
00484   // Otherwise it is 1 and uniform weight is implied.
00485 
00486   numberOfWeights_ = this->inputAdapter_->getNumWeightsPerID();
00487 
00488   numberOfCriteria_ = (numberOfWeights_ > 1) ? numberOfWeights_ : 1;
00489 
00490   inputType_ = this->inputAdapter_->adapterType();
00491 
00492   // The Caller can specify part sizes in setPartSizes().  If he/she
00493   // does not, the part size arrays are empty.
00494 
00495   ArrayRCP<part_t> *noIds = new ArrayRCP<part_t> [numberOfCriteria_];
00496   ArrayRCP<scalar_t> *noSizes = new ArrayRCP<scalar_t> [numberOfCriteria_];
00497 
00498   partIds_ = arcp(noIds, 0, numberOfCriteria_, true);
00499   partSizes_ = arcp(noSizes, 0, numberOfCriteria_, true);
00500 
00501   if (this->env_->getDebugLevel() >= DETAILED_STATUS){
00502     std::ostringstream msg;
00503     msg << problemComm_->getSize() << " procs,"
00504       << numberOfWeights_ << " user-defined weights\n";
00505     this->env_->debug(DETAILED_STATUS, msg.str());
00506   }
00507 
00508   this->env_->memory("After initializeProblem");
00509 }
00510 
00511 template <typename Adapter>
00512   void PartitioningProblem<Adapter>::setPartSizesForCriteria(
00513     int criteria, int len, part_t *partIds, scalar_t *partSizes, bool makeCopy) 
00514 {
00515   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid length", 
00516     len>= 0, BASIC_ASSERTION);
00517 
00518   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid criteria", 
00519     criteria >= 0 && criteria < numberOfCriteria_, BASIC_ASSERTION);
00520 
00521   if (len == 0){
00522     partIds_[criteria] = ArrayRCP<part_t>();
00523     partSizes_[criteria] = ArrayRCP<scalar_t>();
00524     return;
00525   }
00526 
00527   this->env_->localInputAssertion(__FILE__, __LINE__, "invalid arrays", 
00528     partIds && partSizes, BASIC_ASSERTION);
00529 
00530   // The global validity of the partIds and partSizes arrays is performed
00531   // by the PartitioningSolution, which computes global part distribution and
00532   // part sizes.
00533 
00534   part_t *z2_partIds = NULL;
00535   scalar_t *z2_partSizes = NULL;
00536   bool own_memory = false;
00537 
00538   if (makeCopy){
00539     z2_partIds = new part_t [len];
00540     z2_partSizes = new scalar_t [len];
00541     this->env_->localMemoryAssertion(__FILE__, __LINE__, len, z2_partSizes);
00542     memcpy(z2_partIds, partIds, len * sizeof(part_t));
00543     memcpy(z2_partSizes, partSizes, len * sizeof(scalar_t));
00544     own_memory=true;
00545   }
00546   else{
00547     z2_partIds = partIds;
00548     z2_partSizes = partSizes;
00549   }
00550 
00551   partIds_[criteria] = arcp(z2_partIds, 0, len, own_memory);
00552   partSizes_[criteria] = arcp(z2_partSizes, 0, len, own_memory);
00553 }
00554 
00555 template <typename Adapter>
00556 void PartitioningProblem<Adapter>::solve(bool updateInputData)
00557 {
00558   HELLO;
00559   this->env_->debug(DETAILED_STATUS, "Entering solve");
00560 
00561   // Create the computational model.
00562 
00563   this->env_->timerStart(MACRO_TIMERS, "create problem");
00564 
00565   createPartitioningProblem(updateInputData);
00566 
00567   this->env_->timerStop(MACRO_TIMERS, "create problem");
00568 
00569   // TODO: If hierarchical_
00570 
00571   // Create the solution. The algorithm will query the Solution
00572   //   for part and weight information. The algorithm will
00573   //   update the solution with part assignments and quality
00574   //   metrics.  The Solution object itself will convert our internal
00575   //   global numbers back to application global Ids if needed.
00576 
00577   RCP<const IdentifierMap<user_t> > idMap = 
00578     this->baseModel_->getIdentifierMap();
00579 
00580   PartitioningSolution<Adapter> *soln = NULL;
00581 
00582   this->env_->timerStart(MACRO_TIMERS, "create solution");
00583 
00584   try{
00585     soln = new PartitioningSolution<Adapter>( 
00586       this->envConst_, problemCommConst_, idMap, numberOfWeights_, 
00587       partIds_.view(0, numberOfCriteria_), 
00588       partSizes_.view(0, numberOfCriteria_));
00589   }
00590   Z2_FORWARD_EXCEPTIONS;
00591 
00592   solution_ = rcp(soln);
00593 
00594   this->env_->timerStop(MACRO_TIMERS, "create solution");
00595 
00596   this->env_->memory("After creating Solution");
00597 
00598   // Call the algorithm
00599 
00600   try {
00601     if (algName_ == std::string("scotch")) {
00602 
00603       this->algorithm_ = rcp(new AlgPTScotch<Adapter>(this->envConst_,
00604                                             problemComm_,
00605 #ifdef HAVE_ZOLTAN2_MPI
00606                                             mpiComm_,
00607 #endif
00608                                             this->graphModel_));
00609       this->algorithm_->partition(*solution_);
00610     }
00611 
00612     else if (algName_ == std::string("block")) {
00613 
00614       this->algorithm_ = rcp(new AlgBlock<Adapter>(this->envConst_,
00615                                          problemComm_, this->identifierModel_));
00616       this->algorithm_->partition(*solution_);
00617     }
00618 
00619     else if (algName_ == std::string("rcb")) {
00620 
00621       this->algorithm_ = rcp(new AlgRCB<Adapter>(this->envConst_, problemComm_,
00622                                                  this->coordinateModel_));
00623       this->algorithm_->partition(*solution_);
00624     }
00625 
00626     else if (algName_ == std::string("multijagged")) {
00627 
00628       this->algorithm_ = rcp(new Zoltan2_AlgMJ<Adapter>(this->envConst_,
00629                                               problemComm_,
00630                                               this->coordinateModel_));
00631       this->algorithm_->partition(*solution_);
00632     }
00633 
00634     else if (algName_ == std::string("wolf")) {
00635 
00636       this->algorithm_ = rcp(new AlgWolf<Adapter>(this->envConst_,
00637                                         problemComm_,this->graphModel_,
00638                                         this->coordinateModel_));
00639 
00640       // need to add coordModel, make sure this is built
00641       this->algorithm_->partition(*solution_);
00642     }
00643     else {
00644       throw std::logic_error("partitioning algorithm not supported");
00645     }
00646   }
00647   Z2_FORWARD_EXCEPTIONS;
00648 
00649   //if mapping is requested
00650   const Teuchos::ParameterEntry *pe = this->envConst_->getParameters().getEntryPtr("mapping_type");
00651   int mapping_type = -1;
00652   if (pe){
00653     mapping_type = pe->getValue(&mapping_type);
00654   }
00655   //if mapping is 0 -- coordinate mapping
00656 
00657   if (mapping_type == 0){
00658 
00659     //part_t *task_communication_xadj = NULL, *task_communication_adj = NULL;
00660 
00661     Zoltan2::CoordinateTaskMapper <Adapter, part_t> *ctm=
00662                   new Zoltan2::CoordinateTaskMapper<Adapter,part_t>(
00663                           problemComm_.getRawPtr(),
00664                           machine_.getRawPtr(),
00665                           this->coordinateModel_.getRawPtr(),
00666                           solution_.getRawPtr(),
00667                           this->envConst_.getRawPtr()
00668                           //,task_communication_xadj,
00669                           //task_communication_adj
00670                           );
00671     //for now just delete the object.
00672     delete ctm;
00673   }
00674   else if (mapping_type == 1){
00675     //if mapping is 1 -- graph mapping
00676   }
00677 
00678 #ifdef KDDKDD_SHOULD_NEVER_CHANGE_PROBLEMCOMM
00679 #ifdef HAVE_ZOLTAN2_MPI
00680 
00681   // The algorithm may have changed the communicator.  Change it back.
00682   // KDD:  Why would the algorithm change the communicator? TODO
00683   // KDD:  Should we allow such a side effect? TODO
00684 
00685   RCP<const mpiWrapper_t > wrappedComm = rcp(new mpiWrapper_t(mpiComm_));
00686   problemComm_ = rcp(new Teuchos::MpiComm<int>(wrappedComm));
00687   problemCommConst_ = rcp_const_cast<const Comm<int> > (problemComm_);
00688 
00689 #endif
00690 #endif
00691 
00692   if (metricsRequested_){
00693     typedef PartitioningSolution<Adapter> ps_t;
00694     typedef PartitioningSolutionQuality<Adapter> psq_t;
00695 
00696     psq_t *quality = NULL;
00697     RCP<const ps_t> solutionConst = rcp_const_cast<const ps_t>(solution_);
00698     RCP<const Adapter> adapter = rcp(this->inputAdapter_, false);
00699 
00700     try{
00701       quality = new psq_t(this->envConst_, problemCommConst_, adapter, 
00702         solutionConst);
00703     }
00704     Z2_FORWARD_EXCEPTIONS
00705 
00706     metrics_ = rcp(quality);
00707   }
00708 
00709   this->env_->debug(DETAILED_STATUS, "Exiting solve");
00710 }
00711 
00712 template <typename Adapter>
00713 void PartitioningProblem<Adapter>::createPartitioningProblem(bool newData)
00714 {
00715   HELLO;
00716   this->env_->debug(DETAILED_STATUS, 
00717     "PartitioningProblem::createPartitioningProblem");
00718 
00719   using std::string;
00720   using Teuchos::ParameterList;
00721 
00722   // A Problem object may be reused.  The input data may have changed and
00723   // new parameters or part sizes may have been set.
00724   //
00725   // Save these values in order to determine if we need to create a new model.
00726 
00727   //ModelType previousModel = modelType_;
00728   bool prevModelAvail[MAX_NUM_MODEL_TYPES];
00729   for(int i=0;i<MAX_NUM_MODEL_TYPES;i++)
00730   {
00731     prevModelAvail[i] = modelAvail_[i];
00732   }
00733 
00734 
00735   modelFlag_t previousGraphModelFlags = graphFlags_;
00736   modelFlag_t previousIdentifierModelFlags = idFlags_;
00737   modelFlag_t previousCoordinateModelFlags = coordFlags_;
00738 
00739   //modelType_ = InvalidModel;
00740   for(int i=0;i<MAX_NUM_MODEL_TYPES;i++)
00741   {
00742     modelAvail_[i] = false;
00743   }
00744 
00745   graphFlags_.reset();
00746   idFlags_.reset();
00747   coordFlags_.reset();
00748 
00750   // It's possible at this point that the Problem may want to
00751   // add problem parameters to the parameter list in the Environment. 
00752   //
00753   // Since the parameters in the Environment have already been
00754   // validated in its constructor, a new Environment must be created:
00756   // Teuchos::RCP<const Teuchos::Comm<int> > oldComm = this->env_->comm_;
00757   // const ParameterList &oldParams = this->env_->getUnvalidatedParameters();
00758   // 
00759   // ParameterList newParams = oldParams;
00760   // newParams.set("new_parameter", "new_value");
00761   // 
00762   // ParameterList &newPartParams = newParams.sublist("partitioning");
00763   // newPartParams.set("new_partitioning_parameter", "its_value");
00764   // 
00765   // this->env_ = rcp(new Environment(newParams, oldComm));
00767 
00768   this->env_->debug(DETAILED_STATUS, "    parameters");
00769   Environment &env = *(this->env_);
00770   ParameterList &pl = env.getParametersNonConst();
00771 
00772   std::string defString("default");
00773 
00774   // Did the user ask for computation of quality metrics?
00775 
00776   int yesNo=0;
00777   const Teuchos::ParameterEntry *pe = pl.getEntryPtr("compute_metrics");
00778   if (pe){
00779     yesNo = pe->getValue<int>(&yesNo);
00780     metricsRequested_ = true;
00781   }
00782 
00783   // Did the user specify a computational model?
00784 
00785   std::string model(defString);
00786   pe = pl.getEntryPtr("model");
00787   if (pe)
00788     model = pe->getValue<std::string>(&model);
00789 
00790   // Did the user specify an algorithm?
00791 
00792   std::string algorithm(defString);
00793   pe = pl.getEntryPtr("algorithm");
00794   if (pe)
00795     algorithm = pe->getValue<std::string>(&algorithm);
00796 
00797   // Possible algorithm requirements that must be conveyed to the model:
00798 
00799   bool needConsecutiveGlobalIds = false;
00800   bool removeSelfEdges= false;
00801 
00803   // Determine algorithm, model, and algorithm requirements.  This
00804   // is a first pass.  Feel free to change this and add to it.
00805 
00806   if (algorithm != defString)
00807   {
00808 
00809     // Figure out the model required by the algorithm
00810     if (algorithm == std::string("block") ||
00811         algorithm == std::string("random") ||
00812         algorithm == std::string("cyclic") ){
00813 
00814       //modelType_ = IdentifierModelType;
00815       modelAvail_[IdentifierModelType] = true;
00816 
00817       algName_ = algorithm;
00818       needConsecutiveGlobalIds = true;
00819     }
00820     else if (algorithm == std::string("rcb") ||
00821              algorithm == std::string("rib") ||
00822              algorithm == std::string("multijagged") ||
00823              algorithm == std::string("hsfc"))
00824     {
00825       //modelType_ = CoordinateModelType;
00826       modelAvail_[CoordinateModelType]=true;
00827     
00828       algName_ = algorithm;
00829     }
00830     else if (algorithm == std::string("metis") ||
00831              algorithm == std::string("parmetis") ||
00832              algorithm == std::string("scotch") ||
00833              algorithm == std::string("ptscotch"))
00834     {
00835 
00836       //modelType_ = GraphModelType;
00837       modelAvail_[GraphModelType]=true;
00838 
00839       algName_ = algorithm;
00840       removeSelfEdges = true;
00841       needConsecutiveGlobalIds = true;
00842     }
00843     else if (algorithm == std::string("patoh") ||
00844              algorithm == std::string("phg"))
00845     {
00846       // if ((modelType_ != GraphModelType) &&
00847       //     (modelType_ != HypergraphModelType) )
00848       if ((modelAvail_[GraphModelType]==false) &&
00849           (modelAvail_[HypergraphModelType]==false) )
00850       {
00851         //modelType_ = HypergraphModelType;
00852         modelAvail_[HypergraphModelType]=true;
00853       }
00854       algName_ = algorithm;
00855       needConsecutiveGlobalIds = true;
00856     }
00857     else if (algorithm == std::string("wolf"))
00858     {
00859       modelAvail_[GraphModelType]=true;
00860       modelAvail_[CoordinateModelType]=true;
00861       algName_ = algorithm;
00862     }
00863     else
00864     {
00865       // Parameter list should ensure this does not happen.
00866       throw std::logic_error("parameter list algorithm is invalid");
00867     }
00868   }
00869   else if (model != defString)
00870   {
00871     // Figure out the algorithm suggested by the model.
00872     if (model == std::string("hypergraph"))
00873     {      
00874       //modelType_ = HypergraphModelType;
00875       modelAvail_[HypergraphModelType]=true;
00876 
00877       if (problemComm_->getSize() > 1)
00878         algName_ = std::string("phg"); 
00879       else
00880         algName_ = std::string("patoh"); 
00881       needConsecutiveGlobalIds = true;
00882     }
00883     else if (model == std::string("graph"))
00884     {
00885       //modelType_ = GraphModelType;
00886       modelAvail_[GraphModelType]=true;
00887 
00888 #ifdef HAVE_ZOLTAN2_SCOTCH
00889       if (problemComm_->getSize() > 1)
00890         algName_ = std::string("ptscotch"); 
00891       else
00892         algName_ = std::string("scotch"); 
00893       removeSelfEdges = true;
00894       needConsecutiveGlobalIds = true;
00895 #else
00896 #ifdef HAVE_ZOLTAN2_PARMETIS
00897       if (problemComm_->getSize() > 1)
00898         algName_ = std::string("parmetis"); 
00899       else
00900         algName_ = std::string("metis"); 
00901       removeSelfEdges = true;
00902       needConsecutiveGlobalIds = true;
00903 #else
00904       if (problemComm_->getSize() > 1)
00905         algName_ = std::string("phg"); 
00906       else
00907         algName_ = std::string("patoh"); 
00908       removeSelfEdges = true;
00909       needConsecutiveGlobalIds = true;
00910 #endif
00911 #endif
00912     }
00913     else if (model == std::string("geometry"))
00914     {
00915       //modelType_ = CoordinateModelType;
00916       modelAvail_[CoordinateModelType]=true;
00917 
00918       algName_ = std::string("rcb");
00919     }
00920     else if (model == std::string("ids"))
00921     {
00922       //modelType_ = IdentifierModelType;
00923       modelAvail_[IdentifierModelType]=true;
00924 
00925       algName_ = std::string("block");
00926       needConsecutiveGlobalIds = true;
00927     }
00928     else
00929     {
00930       // Parameter list should ensure this does not happen.
00931       env.localBugAssertion(__FILE__, __LINE__, 
00932         "parameter list model type is invalid", 1, BASIC_ASSERTION);
00933     }
00934   }
00935   else
00936   {   
00937     // Determine an algorithm and model suggested by the input type.
00938     //   TODO: this is a good time to use the time vs. quality parameter
00939     //     in choosing an algorithm, and setting some parameters
00940 
00941     if (inputType_ == MatrixAdapterType)
00942     {
00943       //modelType_ = HypergraphModelType;
00944       modelAvail_[HypergraphModelType]=true;
00945       
00946       if (problemComm_->getSize() > 1)
00947         algName_ = std::string("phg"); 
00948       else
00949         algName_ = std::string("patoh"); 
00950     }
00951     else if (inputType_ == GraphAdapterType ||
00952         inputType_ == MeshAdapterType)
00953     {
00954       //modelType_ = GraphModelType;
00955       modelAvail_[GraphModelType]=true;
00956 
00957       if (problemComm_->getSize() > 1)
00958         algName_ = std::string("phg"); 
00959       else
00960         algName_ = std::string("patoh"); 
00961     }
00962     else if (inputType_ == CoordinateAdapterType)
00963     {
00964       //modelType_ = CoordinateModelType;
00965       modelAvail_[CoordinateModelType]=true;
00966 
00967       if(algName_ != std::string("multijagged"))
00968       algName_ = std::string("rcb");
00969     }
00970     else if (inputType_ == VectorAdapterType ||
00971              inputType_ == IdentifierAdapterType)
00972     {
00973       //modelType_ = IdentifierModelType;
00974       modelAvail_[IdentifierModelType]=true;
00975 
00976       algName_ = std::string("block");
00977     }
00978     else{
00979       // This should never happen
00980       throw std::logic_error("input type is invalid");
00981     }
00982   }
00983 
00984   // Hierarchical partitioning?
00985 
00986   Array<int> valueList;
00987   pe = pl.getEntryPtr("topology");
00988 
00989   if (pe){
00990     valueList = pe->getValue<Array<int> >(&valueList);
00991 
00992     if (!Zoltan2::noValuesAreInRangeList<int>(valueList)){
00993       int *n = new int [valueList.size() + 1];
00994       levelNumberParts_ = arcp(n, 0, valueList.size() + 1, true);
00995       int procsPerNode = 1;
00996       for (int i=0; i < valueList.size(); i++){
00997         levelNumberParts_[i+1] = valueList[i];
00998         procsPerNode *= valueList[i];
00999       }
01000       // Number of parts in the first level
01001       levelNumberParts_[0] = env.numProcs_ / procsPerNode;
01002 
01003       if (env.numProcs_ % procsPerNode > 0)
01004         levelNumberParts_[0]++;
01005     }
01006   }
01007   else{
01008     levelNumberParts_.clear();
01009   }
01010 
01011   hierarchical_ = levelNumberParts_.size() > 0;
01012 
01013   // Object to be partitioned? (rows, columns, etc)
01014 
01015   std::string objectOfInterest(defString);
01016   pe = pl.getEntryPtr("objects_to_partition");
01017   if (pe)
01018     objectOfInterest = pe->getValue<std::string>(&objectOfInterest);
01019 
01021   // Set model creation flags, if any.
01022 
01023   this->env_->debug(DETAILED_STATUS, "    models");
01024   //  if (modelType_ == GraphModelType)
01025   if (modelAvail_[GraphModelType]==true)
01026   {
01027 
01028     // Any parameters in the graph sublist?
01029 
01030     std::string symParameter(defString);
01031     pe = pl.getEntryPtr("symmetrize_graph");
01032     if (pe){
01033       symParameter = pe->getValue<std::string>(&symParameter);
01034       if (symParameter == std::string("transpose"))
01035         graphFlags_.set(SYMMETRIZE_INPUT_TRANSPOSE);
01036       else if (symParameter == std::string("bipartite"))
01037         graphFlags_.set(SYMMETRIZE_INPUT_BIPARTITE);
01038     } 
01039 
01040     int sgParameter = 0;
01041     pe = pl.getEntryPtr("subset_graph");
01042     if (pe)
01043       sgParameter = pe->getValue<int>(&sgParameter);
01044 
01045     if (sgParameter == 1)
01046         graphFlags_.set(GRAPH_IS_A_SUBSET_GRAPH);
01047 
01048     // Any special behaviors required by the algorithm?
01049     
01050     if (removeSelfEdges)
01051       graphFlags_.set(SELF_EDGES_MUST_BE_REMOVED);
01052 
01053     if (needConsecutiveGlobalIds)
01054       graphFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
01055 
01056     // How does user input map to vertices and edges?
01057 
01058     if (inputType_ == MatrixAdapterType){
01059       if (objectOfInterest == defString ||
01060           objectOfInterest == std::string("matrix_rows") )
01061         graphFlags_.set(VERTICES_ARE_MATRIX_ROWS);
01062       else if (objectOfInterest == std::string("matrix_columns"))
01063         graphFlags_.set(VERTICES_ARE_MATRIX_COLUMNS);
01064       else if (objectOfInterest == std::string("matrix_nonzeros"))
01065         graphFlags_.set(VERTICES_ARE_MATRIX_NONZEROS);
01066     }
01067 
01068     else if (inputType_ == MeshAdapterType){
01069       if (objectOfInterest == defString ||
01070           objectOfInterest == std::string("mesh_nodes") )
01071         graphFlags_.set(VERTICES_ARE_MESH_NODES);
01072       else if (objectOfInterest == std::string("mesh_elements"))
01073         graphFlags_.set(VERTICES_ARE_MESH_ELEMENTS);
01074     } 
01075   }
01076   //MMW is it ok to remove else?
01077   //  else if (modelType_ == IdentifierModelType)
01078   if (modelAvail_[IdentifierModelType]==true)
01079   {
01080 
01081     // Any special behaviors required by the algorithm?
01082     
01083     if (needConsecutiveGlobalIds)
01084       idFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
01085   }
01086   //  else if (modelType_ == CoordinateModelType)
01087   if (modelAvail_[CoordinateModelType]==true)
01088   {
01089 
01090     // Any special behaviors required by the algorithm?
01091     
01092     if (needConsecutiveGlobalIds)
01093       coordFlags_.set(IDS_MUST_BE_GLOBALLY_CONSECUTIVE);
01094   }
01095 
01096 
01097   if ( newData ||
01098        (modelAvail_[GraphModelType]!=prevModelAvail[GraphModelType]) ||
01099        (modelAvail_[HypergraphModelType]!=prevModelAvail[HypergraphModelType]) ||
01100        (modelAvail_[CoordinateModelType]!=prevModelAvail[CoordinateModelType]) ||
01101        (modelAvail_[IdentifierModelType]!=prevModelAvail[IdentifierModelType]) ||
01102   //       (modelType_ != previousModel) ||
01103        (graphFlags_ != previousGraphModelFlags) ||
01104        (coordFlags_ != previousCoordinateModelFlags) ||
01105        (idFlags_ != previousIdentifierModelFlags) ) 
01106   {
01107 
01108     // Create the computational model.
01109     // Models are instantiated for base input adapter types (mesh,
01110     // matrix, graph, and so on).  We pass a pointer to the input
01111     // adapter, cast as the base input type.
01112 
01113     //KDD Not sure why this shadow declaration is needed
01114     //KDD Comment out for now; revisit later if problems.
01115     //KDD const Teuchos::ParameterList pl = this->envConst_->getParameters();
01116     //bool exceptionThrow = true;
01117 
01118     if(modelAvail_[GraphModelType]==false && modelAvail_[HypergraphModelType]==false &&
01119        modelAvail_[CoordinateModelType]==false && modelAvail_[IdentifierModelType]==false)
01120     {
01121       cout << __func__ << " Invalid model"  << endl;
01122     }
01123     else
01124     {
01125       if(modelAvail_[GraphModelType]==true)
01126       {
01127         this->env_->debug(DETAILED_STATUS, "    building graph model");
01128         this->graphModel_ = rcp(new GraphModel<base_adapter_t>(
01129           this->baseInputAdapter_, this->envConst_, problemComm_, graphFlags_));
01130 
01131         this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(this->graphModel_);
01132       }
01133       if(modelAvail_[HypergraphModelType]==true)
01134       {
01135   std::cout << "Hypergraph model not implemented yet..." << std::endl;
01136       }
01137 
01138       if(modelAvail_[CoordinateModelType]==true)
01139       {
01140         this->env_->debug(DETAILED_STATUS, "    building coordinate model");
01141         this->coordinateModel_ = rcp(new CoordinateModel<base_adapter_t>(
01142                    this->baseInputAdapter_, this->envConst_, problemComm_, coordFlags_));
01143 
01144         this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(this->coordinateModel_);
01145       }
01146 
01147       if(modelAvail_[IdentifierModelType]==true)
01148       {
01149         this->env_->debug(DETAILED_STATUS, "    building identifier model");
01150         this->identifierModel_ = rcp(new IdentifierModel<base_adapter_t>(
01151                                      this->baseInputAdapter_, this->envConst_, problemComm_, idFlags_));
01152 
01153         this->baseModel_ = rcp_implicit_cast<const Model<base_adapter_t> >(this->identifierModel_);
01154       }
01155   
01156 
01157     }
01158 
01159 
01160 
01161     this->env_->memory("After creating Model");
01162     this->env_->debug(DETAILED_STATUS, "createPartitioningProblem done");
01163   }
01164 
01165 }
01166 
01167 }  // namespace Zoltan2
01168 #endif