FEI Version of the Day
fei_Trilinos_Helpers.cpp
00001 /*
00002 // @HEADER
00003 // ************************************************************************
00004 //             FEI: Finite Element Interface to Linear Solvers
00005 //                  Copyright (2005) Sandia Corporation.
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
00008 // U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Alan Williams (william@sandia.gov) 
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 */
00042 
00043 #include <fei_macros.hpp>
00044 
00045 #include <fei_Include_Trilinos.hpp>
00046 
00047 #ifdef HAVE_FEI_EPETRA
00048 #include <fei_VectorTraits_Epetra.hpp>
00049 #include <fei_MatrixTraits_Epetra.hpp>
00050 #include <fei_LinProbMgr_EpetraBasic.hpp>
00051 #endif
00052 
00053 #include <fei_Trilinos_Helpers.hpp>
00054 #include <fei_ParameterSet.hpp>
00055 #include <fei_Vector_Impl.hpp>
00056 #include <fei_VectorReducer.hpp>
00057 #include <fei_Matrix_Impl.hpp>
00058 #include <fei_MatrixReducer.hpp>
00059 //#include <EpetraExt_BlockMapOut.h>
00060 
00061 namespace Trilinos_Helpers {
00062 
00063 #ifdef HAVE_FEI_EPETRA
00064 
00065 Epetra_Map
00066 create_Epetra_Map(MPI_Comm comm,
00067                   const std::vector<int>& local_eqns)
00068 {
00069 #ifndef FEI_SER
00070   Epetra_MpiComm EComm(comm);
00071 #else
00072   Epetra_SerialComm EComm;
00073 #endif
00074 
00075   int localSize = local_eqns.size();
00076   int globalSize = 0;
00077   EComm.SumAll(&localSize, &globalSize, 1);
00078 
00079   if (localSize < 0 || globalSize < 0) {
00080     throw std::runtime_error("Trilinos_Helpers::create_Epetra_Map: negative local or global size.");
00081   }
00082 
00083   Epetra_Map EMap(globalSize, localSize, 0, EComm);
00084   return(EMap);
00085 }
00086 
00087 Epetra_BlockMap
00088 create_Epetra_BlockMap(const fei::SharedPtr<fei::VectorSpace>& vecspace)
00089 {
00090   if (vecspace.get() == 0) {
00091     throw std::runtime_error("create_Epetra_Map needs non-null fei::VectorSpace");
00092   }
00093 
00094 #ifndef FEI_SER
00095   MPI_Comm comm = vecspace->getCommunicator();
00096   Epetra_MpiComm EComm(comm);
00097 #else
00098   Epetra_SerialComm EComm;
00099 #endif
00100 
00101   int localSizeBlk = vecspace->getNumBlkIndices_Owned();
00102   int globalSizeBlk = vecspace->getGlobalNumBlkIndices();
00103 
00104   if (localSizeBlk < 0 || globalSizeBlk < 0) {
00105     throw std::runtime_error("Trilinos_Helpers::create_Epetra_BlockMap: fei::VectorSpace has negative local or global size.");
00106   }
00107 
00108   std::vector<int> blkEqns(localSizeBlk*2);
00109   int* blkEqnsPtr = &(blkEqns[0]);
00110 
00111   int chkNum = 0;
00112   int errcode = vecspace->getBlkIndices_Owned(localSizeBlk,
00113                 blkEqnsPtr, blkEqnsPtr+localSizeBlk,
00114                 chkNum);
00115   if (errcode != 0) {
00116     FEI_OSTRINGSTREAM osstr;
00117     osstr << "create_Epetra_BlockMap ERROR, nonzero errcode="<<errcode
00118     << " returned by vecspace->getBlkIndices_Owned.";
00119     throw std::runtime_error(osstr.str());
00120   }
00121 
00122   Epetra_BlockMap EBMap(globalSizeBlk, localSizeBlk,
00123       blkEqnsPtr, blkEqnsPtr+localSizeBlk, 0, EComm);
00124 
00125   return(EBMap);
00126 }
00127 
00128 Epetra_CrsGraph
00129 create_Epetra_CrsGraph(const fei::SharedPtr<fei::MatrixGraph>& matgraph,
00130            bool blockEntries, bool orderRowsWithLocalColsFirst)
00131 {
00132   fei::SharedPtr<fei::SparseRowGraph> localSRGraph =
00133     matgraph->createGraph(blockEntries);
00134   if (localSRGraph.get() == NULL) {
00135     throw std::runtime_error("create_Epetra_CrsGraph ERROR in fei::MatrixGraph::createGraph");
00136   }
00137 
00138   int numLocallyOwnedRows = localSRGraph->rowNumbers.size();
00139   int* rowNumbers = &(localSRGraph->rowNumbers[0]);
00140   int* rowOffsets = &(localSRGraph->rowOffsets[0]);
00141   int* packedColumnIndices = &(localSRGraph->packedColumnIndices[0]);
00142 
00143   fei::SharedPtr<fei::VectorSpace> vecspace = matgraph->getRowSpace();
00144   MPI_Comm comm = vecspace->getCommunicator();
00145   std::vector<int>& local_eqns = localSRGraph->rowNumbers;
00146 
00147   Epetra_BlockMap emap = blockEntries ?
00148     create_Epetra_BlockMap(vecspace) : create_Epetra_Map(comm, local_eqns);
00149 
00150   if (orderRowsWithLocalColsFirst == true &&
00151       emap.Comm().NumProc() > 2 && !blockEntries) {
00152     bool* used_row = new bool[local_eqns.size()];
00153     for(unsigned ii=0; ii<local_eqns.size(); ++ii) used_row[ii] = false;
00154 
00155     int offset = 0;
00156     std::vector<int> ordered_local_eqns(local_eqns.size());
00157     for(unsigned ii=0; ii<local_eqns.size(); ++ii) {
00158       bool row_has_off_proc_cols = false;
00159       for(int j=rowOffsets[ii]; j<rowOffsets[ii+1]; ++j) {
00160         if (emap.MyGID(packedColumnIndices[j]) == false) {
00161           row_has_off_proc_cols = true;
00162           break;
00163         }
00164       }
00165 
00166       if (row_has_off_proc_cols == false) {
00167         ordered_local_eqns[offset++] = rowNumbers[ii];
00168         used_row[ii] = true;
00169       }
00170     }
00171 
00172     for(unsigned ii=0; ii<local_eqns.size(); ++ii) {
00173       if (used_row[ii] == true) continue;
00174       ordered_local_eqns[offset++] = rowNumbers[ii];
00175     }
00176 
00177     emap = create_Epetra_Map(comm, ordered_local_eqns);
00178     delete [] used_row;
00179   }
00180 
00181 //  EpetraExt::BlockMapToMatrixMarketFile("EBMap.np12.mm",emap,"AriaTest");
00182 
00183   std::vector<int> rowLengths; rowLengths.reserve(numLocallyOwnedRows);
00184   for(int ii=0; ii<numLocallyOwnedRows; ++ii) {
00185     rowLengths.push_back(rowOffsets[ii+1]-rowOffsets[ii]);
00186   }
00187 
00188   bool staticProfile = true;
00189   Epetra_CrsGraph egraph(Copy, emap, &(rowLengths[0]), staticProfile);
00190 
00191   const Epetra_Comm& ecomm = emap.Comm();
00192   int localProc = ecomm.MyPID();
00193 
00194   int firstLocalEqn = numLocallyOwnedRows > 0 ? rowNumbers[0] : -1;
00195 
00196   int offset = 0;
00197   for(int i=0; i<numLocallyOwnedRows; ++i) {
00198     int err = egraph.InsertGlobalIndices(firstLocalEqn+i,
00199            rowLengths[i],
00200            &(packedColumnIndices[offset]));
00201     if (err != 0) {
00202       fei::console_out() << "proc " << localProc << " err-return " << err
00203                << " inserting row " << firstLocalEqn+i<<", cols ";
00204       for(int ii=0; ii<rowLengths[i]; ++ii) {
00205         fei::console_out() << packedColumnIndices[offset+ii]<<",";
00206       }
00207       fei::console_out() << FEI_ENDL;
00208       throw std::runtime_error("... occurred in create_Epetra_CrsGraph");
00209     }
00210 
00211     offset += rowLengths[i];
00212   }
00213 
00214   //Epetra_BlockMap* domainmap = const_cast<Epetra_BlockMap*>(&(epetraGraph_->DomainMap()));
00215   //Epetra_BlockMap* rangemap = const_cast<Epetra_BlockMap*>(&(epetraGraph_->RangeMap()));
00216   egraph.FillComplete();
00217 
00218   //only optimize-storage for graph if we're not in block-matrix mode.
00219   //(Epetra_VbrMatrix can't be constructed with an optimize-storage'd graph.)
00220   if (!blockEntries) egraph.OptimizeStorage();
00221 
00222   return(egraph);
00223 }
00224 
00225 fei::SharedPtr<fei::Matrix>
00226 create_from_Epetra_Matrix(fei::SharedPtr<fei::MatrixGraph> matrixGraph,
00227                           bool blockEntryMatrix,
00228                           fei::SharedPtr<fei::Reducer> reducer,
00229                           bool orderRowsWithLocalColsFirst)
00230 {
00231   fei::SharedPtr<fei::VectorSpace> vecSpace = matrixGraph->getRowSpace();
00232   //localSize is in point-equations, because we will only use it for constructing
00233   //the fei::Matrix_Impl, and those always deal in point-equations.
00234   int localSize = vecSpace->getNumIndices_Owned();
00235 
00236   fei::SharedPtr<fei::Matrix> tmpmat;
00237   fei::SharedPtr<fei::Matrix> feimat;
00238   try {
00239     Epetra_CrsGraph egraph =
00240       Trilinos_Helpers::create_Epetra_CrsGraph(matrixGraph, blockEntryMatrix,
00241                                                orderRowsWithLocalColsFirst);
00242 
00243     if (blockEntryMatrix) {
00244       fei::SharedPtr<Epetra_VbrMatrix>
00245         epetraMatrix(new Epetra_VbrMatrix(Copy, egraph));
00246 
00247       tmpmat.reset(new fei::Matrix_Impl<Epetra_VbrMatrix>(epetraMatrix,
00248                                                    matrixGraph, localSize));
00249       zero_Epetra_VbrMatrix(epetraMatrix.get());
00250     }
00251     else {
00252       fei::SharedPtr<Epetra_CrsMatrix>
00253         epetraMatrix(new Epetra_CrsMatrix(Copy, egraph));
00254 
00255       tmpmat.reset(new fei::Matrix_Impl<Epetra_CrsMatrix>(epetraMatrix,
00256                                           matrixGraph, localSize));
00257     }
00258   }
00259   catch(std::runtime_error& exc) {
00260     fei::console_out() << "Trilinos_Helpers::create_from_Epetra_Matrix ERROR, "
00261            << "caught exception: '" << exc.what() << "', rethrowing..."
00262            << FEI_ENDL;
00263     throw exc;
00264   }
00265 
00266   if (reducer.get() != NULL) {
00267     feimat.reset(new fei::MatrixReducer(reducer, tmpmat));
00268   }
00269   else {
00270     feimat = tmpmat;
00271   }
00272 
00273   return(feimat);
00274 }
00275 
00276 fei::SharedPtr<fei::Matrix>
00277 create_from_LPM_EpetraBasic(fei::SharedPtr<fei::MatrixGraph> matrixGraph,
00278                             bool blockEntryMatrix,
00279                             fei::SharedPtr<fei::Reducer> reducer,
00280                             fei::SharedPtr<fei::LinearProblemManager>
00281                               lpm_epetrabasic)
00282 {
00283   fei::SharedPtr<fei::SparseRowGraph> srgraph =
00284     matrixGraph->createGraph(blockEntryMatrix);
00285   if (srgraph.get() == NULL) {
00286     throw std::runtime_error("create_LPM_EpetraBasic ERROR in fei::MatrixGraph::createGraph");
00287   }
00288 
00289   fei::SharedPtr<fei::SparseRowGraph> sharedsrg(srgraph);
00290   int localSize;
00291   if (reducer.get() != NULL) {
00292     std::vector<int>& reduced_eqns = reducer->getLocalReducedEqns();
00293     lpm_epetrabasic->setRowDistribution(reduced_eqns);
00294     lpm_epetrabasic->setMatrixGraph(sharedsrg);
00295     localSize = reduced_eqns.size();
00296   }
00297   else {
00298     fei::SharedPtr<fei::VectorSpace> vecSpace = matrixGraph->getRowSpace();
00299     int err = 0,chkNum;
00300     std::vector<int> indices;
00301     if (blockEntryMatrix) {
00302       localSize = vecSpace->getNumBlkIndices_Owned();
00303       indices.resize(localSize*2);
00304       err = vecSpace->getBlkIndices_Owned(localSize, &indices[0],
00305                                   &indices[localSize], chkNum);
00306     }
00307     else {
00308       localSize = vecSpace->getNumIndices_Owned();
00309       indices.resize(localSize);
00310       err = vecSpace->getIndices_Owned(indices);
00311     }
00312     if (err != 0) {
00313       throw std::runtime_error("Factory_Trilinos: createMatrix: error in vecSpace->getIndices_Owned");
00314     }
00315 
00316     lpm_epetrabasic->setRowDistribution(indices);
00317     lpm_epetrabasic->setMatrixGraph(sharedsrg);
00318   }
00319 
00320   fei::SharedPtr<fei::Matrix> tmpmat(new fei::Matrix_Impl<fei::LinearProblemManager>(lpm_epetrabasic, matrixGraph, localSize));
00321 
00322   fei::SharedPtr<fei::Matrix> feimat;
00323 
00324   if (reducer.get() != NULL) {
00325     feimat.reset(new fei::MatrixReducer(reducer, tmpmat));
00326   }
00327   else {
00328     feimat = tmpmat;
00329   }
00330 
00331   return(feimat);
00332 }
00333 #endif //HAVE_FEI_EPETRA
00334 
00335 void copy_parameterset(const fei::ParameterSet& paramset,
00336                        Teuchos::ParameterList& paramlist)
00337 {
00338   fei::ParameterSet::const_iterator
00339     iter = paramset.begin(),
00340     iter_end = paramset.end();
00341 
00342   for(; iter != iter_end; ++iter) {
00343     const fei::Param& param = *iter;
00344     fei::Param::ParamType ptype = param.getType();
00345     switch(ptype) {
00346     case fei::Param::STRING:
00347       paramlist.set(param.getName(), param.getStringValue());
00348       break;
00349     case fei::Param::DOUBLE:
00350       paramlist.set(param.getName(), param.getDoubleValue());
00351       break;
00352     case fei::Param::INT:
00353       paramlist.set(param.getName(), param.getIntValue());
00354       break;
00355     case fei::Param::BOOL:
00356       paramlist.set(param.getName(), param.getBoolValue());
00357       break;
00358     default:
00359       break;
00360     }
00361   }
00362 }
00363 
00364 void copy_parameterlist(const Teuchos::ParameterList& paramlist,
00365                         fei::ParameterSet& paramset)
00366 {
00367   Teuchos::ParameterList::ConstIterator
00368     iter = paramlist.begin(),
00369     iter_end = paramlist.end();
00370 
00371   for(; iter != iter_end; ++iter) {
00372     const Teuchos::ParameterEntry& param = paramlist.entry(iter);
00373     if (param.isType<std::string>()) {
00374       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<std::string>(param).c_str()));
00375     }
00376     else if (param.isType<double>()) {
00377       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<double>(param)));
00378     }
00379     else if (param.isType<int>()) {
00380       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<int>(param)));
00381     }
00382     else if (param.isType<bool>()) {
00383       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<bool>(param)));
00384     }
00385   }
00386 }
00387 
00388 #ifdef HAVE_FEI_EPETRA
00389 
00390 Epetra_MultiVector*
00391 get_Epetra_MultiVector(fei::Vector* feivec, bool soln_vec)
00392 {
00393   fei::Vector* vecptr = feivec;
00394   fei::VectorReducer* feireducer = dynamic_cast<fei::VectorReducer*>(feivec);
00395   if (feireducer != NULL) vecptr = feireducer->getTargetVector().get();
00396 
00397   fei::Vector_Impl<Epetra_MultiVector>* fei_epetra_vec =
00398     dynamic_cast<fei::Vector_Impl<Epetra_MultiVector>*>(vecptr);
00399   fei::Vector_Impl<fei::LinearProblemManager>* fei_lpm =
00400     dynamic_cast<fei::Vector_Impl<fei::LinearProblemManager>*>(vecptr);
00401 
00402   if (fei_epetra_vec == NULL && fei_lpm == NULL) {
00403     throw std::runtime_error("failed to obtain Epetra_MultiVector from fei::Vector.");
00404   }
00405 
00406   if (fei_epetra_vec != NULL) {
00407     return( fei_epetra_vec->getUnderlyingVector());
00408   }
00409 
00410   LinProbMgr_EpetraBasic* lpm_epetrabasic =
00411     dynamic_cast<LinProbMgr_EpetraBasic*>(fei_lpm->getUnderlyingVector());
00412   if (lpm_epetrabasic == 0) {
00413     throw std::runtime_error("fei Trilinos_Helpers: ERROR getting LinProbMgr_EpetraBasic");
00414   }
00415 
00416   if (soln_vec) {
00417     return(lpm_epetrabasic->get_solution_vector().get());
00418   }
00419 
00420   return(lpm_epetrabasic->get_rhs_vector().get());
00421 }
00422 
00423 Epetra_VbrMatrix* get_Epetra_VbrMatrix(fei::Matrix* feimat)
00424 {
00425   fei::Matrix_Impl<Epetra_VbrMatrix>* fei_epetra_vbr =
00426     dynamic_cast<fei::Matrix_Impl<Epetra_VbrMatrix>*>(feimat);
00427   fei::MatrixReducer* feireducer =
00428     dynamic_cast<fei::MatrixReducer*>(feimat);
00429 
00430   if (feireducer != NULL) {
00431     fei::SharedPtr<fei::Matrix> feimat2 = feireducer->getTargetMatrix();
00432     fei_epetra_vbr =
00433       dynamic_cast<fei::Matrix_Impl<Epetra_VbrMatrix>*>(feimat2.get());
00434   }
00435 
00436   if (fei_epetra_vbr == NULL) {
00437     throw std::runtime_error("failed to obtain Epetra_VbrMatrix from fei::Matrix.");
00438   }
00439 
00440   return(fei_epetra_vbr->getMatrix().get());
00441 }
00442 
00443 Epetra_CrsMatrix* get_Epetra_CrsMatrix(fei::Matrix* feimat)
00444 {
00445   fei::Matrix_Impl<Epetra_CrsMatrix>* fei_epetra_crs =
00446     dynamic_cast<fei::Matrix_Impl<Epetra_CrsMatrix>*>(feimat);
00447   fei::Matrix_Impl<fei::LinearProblemManager>* fei_lpm =
00448     dynamic_cast<fei::Matrix_Impl<fei::LinearProblemManager>*>(feimat);
00449   fei::MatrixReducer* feireducer =
00450     dynamic_cast<fei::MatrixReducer*>(feimat);
00451 
00452   if (feireducer != NULL) {
00453     fei::SharedPtr<fei::Matrix> feimat2 = feireducer->getTargetMatrix();
00454     fei_epetra_crs =
00455       dynamic_cast<fei::Matrix_Impl<Epetra_CrsMatrix>*>(feimat2.get());
00456     fei_lpm =
00457       dynamic_cast<fei::Matrix_Impl<fei::LinearProblemManager>*>(feimat2.get());
00458   }
00459 
00460   if (fei_epetra_crs == NULL && fei_lpm == NULL) {
00461     throw std::runtime_error("failed to obtain Epetra_CrsMatrix from fei::Matrix.");
00462   }
00463 
00464   if (fei_epetra_crs != NULL) {
00465     return(fei_epetra_crs->getMatrix().get());
00466   }
00467 
00468   LinProbMgr_EpetraBasic* lpm_epetrabasic =
00469     dynamic_cast<LinProbMgr_EpetraBasic*>(fei_lpm->getMatrix().get());
00470   if (lpm_epetrabasic == 0) {
00471     throw std::runtime_error("fei Trilinos_Helpers ERROR getting LinProbMgr_EpetraBasic");
00472   }
00473 
00474   return(lpm_epetrabasic->get_A_matrix().get());
00475 }
00476 
00477 
00478 void get_Epetra_pointers(fei::SharedPtr<fei::Matrix> feiA,
00479                          fei::SharedPtr<fei::Vector> feix,
00480                          fei::SharedPtr<fei::Vector> feib,
00481                          Epetra_CrsMatrix*& crsA,
00482                          Epetra_Operator*& opA,
00483                          Epetra_MultiVector*& x,
00484                          Epetra_MultiVector*& b)
00485 {
00486   x = get_Epetra_MultiVector(feix.get(), true);
00487   b = get_Epetra_MultiVector(feib.get(), false);
00488 
00489   const char* matname = feiA->typeName();
00490   if (!strcmp(matname, "Epetra_VbrMatrix")) {
00491     Epetra_VbrMatrix* A = get_Epetra_VbrMatrix(feiA.get());
00492     opA = A;
00493   }
00494   else {
00495     crsA = get_Epetra_CrsMatrix(feiA.get());
00496     opA = crsA;
00497   }
00498 }
00499 
00500 int zero_Epetra_VbrMatrix(Epetra_VbrMatrix* mat)
00501 {
00502   const Epetra_CrsGraph& crsgraph = mat->Graph();
00503   const Epetra_BlockMap& rowmap = crsgraph.RowMap();
00504   const Epetra_BlockMap& colmap = crsgraph.ColMap();
00505   mat->RowMatrixRowMap();//generate point objects
00506   int maxBlkRowSize = mat->GlobalMaxRowDim();
00507   int maxBlkColSize = mat->GlobalMaxColDim();
00508   std::vector<double> zeros(maxBlkRowSize*maxBlkColSize, 0);
00509   int numMyRows = rowmap.NumMyElements();
00510   int* myRows = rowmap.MyGlobalElements();
00511   for(int i=0; i<numMyRows; ++i) {
00512     int row = myRows[i];
00513     int rowlength = 0;
00514     int* colindicesView = NULL;
00515     int localrow = rowmap.LID(row);
00516     int err = crsgraph.ExtractMyRowView(localrow, rowlength, colindicesView);
00517     if (err != 0) {
00518       return err;
00519     }
00520     err = mat->BeginReplaceMyValues(localrow, rowlength, colindicesView);
00521     if (err != 0) {
00522       return err;
00523     }
00524     int blkRowSize = rowmap.ElementSize(localrow);
00525     for(int j=0; j<rowlength; ++j) {
00526       int blkColSize = colmap.ElementSize(colindicesView[j]);
00527       err = mat->SubmitBlockEntry(&zeros[0], maxBlkRowSize,
00528                             blkRowSize, blkColSize);
00529       if (err != 0) {
00530         return err;
00531       }
00532     }
00533     err = mat->EndSubmitEntries();
00534     if (err != 0) {
00535       return err;
00536     }
00537   }
00538 
00539   return 0;
00540 }
00541 
00542 #endif //HAVE_FEI_EPETRA
00543 
00544 }//namespace Trilinos_Helpers
00545 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends