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       bool zeroSharedRows = false;
00248       tmpmat.reset(new fei::Matrix_Impl<Epetra_VbrMatrix>(epetraMatrix,
00249                                      matrixGraph, localSize, zeroSharedRows));
00250       zero_Epetra_VbrMatrix(epetraMatrix.get());
00251     }
00252     else {
00253       fei::SharedPtr<Epetra_CrsMatrix>
00254         epetraMatrix(new Epetra_CrsMatrix(Copy, egraph));
00255 
00256       tmpmat.reset(new fei::Matrix_Impl<Epetra_CrsMatrix>(epetraMatrix,
00257                                           matrixGraph, localSize));
00258     }
00259   }
00260   catch(std::runtime_error& exc) {
00261     fei::console_out() << "Trilinos_Helpers::create_from_Epetra_Matrix ERROR, "
00262            << "caught exception: '" << exc.what() << "', rethrowing..."
00263            << FEI_ENDL;
00264     throw exc;
00265   }
00266 
00267   if (reducer.get() != NULL) {
00268     feimat.reset(new fei::MatrixReducer(reducer, tmpmat));
00269   }
00270   else {
00271     feimat = tmpmat;
00272   }
00273 
00274   return(feimat);
00275 }
00276 
00277 fei::SharedPtr<fei::Matrix>
00278 create_from_LPM_EpetraBasic(fei::SharedPtr<fei::MatrixGraph> matrixGraph,
00279                             bool blockEntryMatrix,
00280                             fei::SharedPtr<fei::Reducer> reducer,
00281                             fei::SharedPtr<fei::LinearProblemManager>
00282                               lpm_epetrabasic)
00283 {
00284   fei::SharedPtr<fei::SparseRowGraph> srgraph =
00285     matrixGraph->createGraph(blockEntryMatrix);
00286   if (srgraph.get() == NULL) {
00287     throw std::runtime_error("create_LPM_EpetraBasic ERROR in fei::MatrixGraph::createGraph");
00288   }
00289 
00290   fei::SharedPtr<fei::SparseRowGraph> sharedsrg(srgraph);
00291   int localSize;
00292   if (reducer.get() != NULL) {
00293     std::vector<int>& reduced_eqns = reducer->getLocalReducedEqns();
00294     lpm_epetrabasic->setRowDistribution(reduced_eqns);
00295     lpm_epetrabasic->setMatrixGraph(sharedsrg);
00296     localSize = reduced_eqns.size();
00297   }
00298   else {
00299     fei::SharedPtr<fei::VectorSpace> vecSpace = matrixGraph->getRowSpace();
00300     int err = 0,chkNum;
00301     std::vector<int> indices;
00302     if (blockEntryMatrix) {
00303       localSize = vecSpace->getNumBlkIndices_Owned();
00304       indices.resize(localSize*2);
00305       err = vecSpace->getBlkIndices_Owned(localSize, &indices[0],
00306                                   &indices[localSize], chkNum);
00307     }
00308     else {
00309       localSize = vecSpace->getNumIndices_Owned();
00310       indices.resize(localSize);
00311       err = vecSpace->getIndices_Owned(indices);
00312     }
00313     if (err != 0) {
00314       throw std::runtime_error("Factory_Trilinos: createMatrix: error in vecSpace->getIndices_Owned");
00315     }
00316 
00317     lpm_epetrabasic->setRowDistribution(indices);
00318     lpm_epetrabasic->setMatrixGraph(sharedsrg);
00319   }
00320 
00321   fei::SharedPtr<fei::Matrix> tmpmat(new fei::Matrix_Impl<fei::LinearProblemManager>(lpm_epetrabasic, matrixGraph, localSize));
00322 
00323   fei::SharedPtr<fei::Matrix> feimat;
00324 
00325   if (reducer.get() != NULL) {
00326     feimat.reset(new fei::MatrixReducer(reducer, tmpmat));
00327   }
00328   else {
00329     feimat = tmpmat;
00330   }
00331 
00332   return(feimat);
00333 }
00334 #endif //HAVE_FEI_EPETRA
00335 
00336 void copy_parameterset(const fei::ParameterSet& paramset,
00337                        Teuchos::ParameterList& paramlist)
00338 {
00339   fei::ParameterSet::const_iterator
00340     iter = paramset.begin(),
00341     iter_end = paramset.end();
00342 
00343   for(; iter != iter_end; ++iter) {
00344     const fei::Param& param = *iter;
00345     fei::Param::ParamType ptype = param.getType();
00346     switch(ptype) {
00347     case fei::Param::STRING:
00348       paramlist.set(param.getName(), param.getStringValue());
00349       break;
00350     case fei::Param::DOUBLE:
00351       paramlist.set(param.getName(), param.getDoubleValue());
00352       break;
00353     case fei::Param::INT:
00354       paramlist.set(param.getName(), param.getIntValue());
00355       break;
00356     case fei::Param::BOOL:
00357       paramlist.set(param.getName(), param.getBoolValue());
00358       break;
00359     default:
00360       break;
00361     }
00362   }
00363 }
00364 
00365 void copy_parameterlist(const Teuchos::ParameterList& paramlist,
00366                         fei::ParameterSet& paramset)
00367 {
00368   Teuchos::ParameterList::ConstIterator
00369     iter = paramlist.begin(),
00370     iter_end = paramlist.end();
00371 
00372   for(; iter != iter_end; ++iter) {
00373     const Teuchos::ParameterEntry& param = paramlist.entry(iter);
00374     if (param.isType<std::string>()) {
00375       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<std::string>(param).c_str()));
00376     }
00377     else if (param.isType<double>()) {
00378       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<double>(param)));
00379     }
00380     else if (param.isType<int>()) {
00381       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<int>(param)));
00382     }
00383     else if (param.isType<bool>()) {
00384       paramset.add(fei::Param(paramlist.name(iter).c_str(), Teuchos::getValue<bool>(param)));
00385     }
00386   }
00387 }
00388 
00389 #ifdef HAVE_FEI_EPETRA
00390 
00391 Epetra_MultiVector*
00392 get_Epetra_MultiVector(fei::Vector* feivec, bool soln_vec)
00393 {
00394   fei::Vector* vecptr = feivec;
00395   fei::VectorReducer* feireducer = dynamic_cast<fei::VectorReducer*>(feivec);
00396   if (feireducer != NULL) vecptr = feireducer->getTargetVector().get();
00397 
00398   fei::Vector_Impl<Epetra_MultiVector>* fei_epetra_vec =
00399     dynamic_cast<fei::Vector_Impl<Epetra_MultiVector>*>(vecptr);
00400   fei::Vector_Impl<fei::LinearProblemManager>* fei_lpm =
00401     dynamic_cast<fei::Vector_Impl<fei::LinearProblemManager>*>(vecptr);
00402 
00403   if (fei_epetra_vec == NULL && fei_lpm == NULL) {
00404     throw std::runtime_error("failed to obtain Epetra_MultiVector from fei::Vector.");
00405   }
00406 
00407   if (fei_epetra_vec != NULL) {
00408     return( fei_epetra_vec->getUnderlyingVector());
00409   }
00410 
00411   LinProbMgr_EpetraBasic* lpm_epetrabasic =
00412     dynamic_cast<LinProbMgr_EpetraBasic*>(fei_lpm->getUnderlyingVector());
00413   if (lpm_epetrabasic == 0) {
00414     throw std::runtime_error("fei Trilinos_Helpers: ERROR getting LinProbMgr_EpetraBasic");
00415   }
00416 
00417   if (soln_vec) {
00418     return(lpm_epetrabasic->get_solution_vector().get());
00419   }
00420 
00421   return(lpm_epetrabasic->get_rhs_vector().get());
00422 }
00423 
00424 Epetra_VbrMatrix* get_Epetra_VbrMatrix(fei::Matrix* feimat)
00425 {
00426   fei::Matrix_Impl<Epetra_VbrMatrix>* fei_epetra_vbr =
00427     dynamic_cast<fei::Matrix_Impl<Epetra_VbrMatrix>*>(feimat);
00428   fei::MatrixReducer* feireducer =
00429     dynamic_cast<fei::MatrixReducer*>(feimat);
00430 
00431   if (feireducer != NULL) {
00432     fei::SharedPtr<fei::Matrix> feimat2 = feireducer->getTargetMatrix();
00433     fei_epetra_vbr =
00434       dynamic_cast<fei::Matrix_Impl<Epetra_VbrMatrix>*>(feimat2.get());
00435   }
00436 
00437   if (fei_epetra_vbr == NULL) {
00438     throw std::runtime_error("failed to obtain Epetra_VbrMatrix from fei::Matrix.");
00439   }
00440 
00441   return(fei_epetra_vbr->getMatrix().get());
00442 }
00443 
00444 Epetra_CrsMatrix* get_Epetra_CrsMatrix(fei::Matrix* feimat)
00445 {
00446   fei::Matrix_Impl<Epetra_CrsMatrix>* fei_epetra_crs =
00447     dynamic_cast<fei::Matrix_Impl<Epetra_CrsMatrix>*>(feimat);
00448   fei::Matrix_Impl<fei::LinearProblemManager>* fei_lpm =
00449     dynamic_cast<fei::Matrix_Impl<fei::LinearProblemManager>*>(feimat);
00450   fei::MatrixReducer* feireducer =
00451     dynamic_cast<fei::MatrixReducer*>(feimat);
00452 
00453   if (feireducer != NULL) {
00454     fei::SharedPtr<fei::Matrix> feimat2 = feireducer->getTargetMatrix();
00455     fei_epetra_crs =
00456       dynamic_cast<fei::Matrix_Impl<Epetra_CrsMatrix>*>(feimat2.get());
00457     fei_lpm =
00458       dynamic_cast<fei::Matrix_Impl<fei::LinearProblemManager>*>(feimat2.get());
00459   }
00460 
00461   if (fei_epetra_crs == NULL && fei_lpm == NULL) {
00462     throw std::runtime_error("failed to obtain Epetra_CrsMatrix from fei::Matrix.");
00463   }
00464 
00465   if (fei_epetra_crs != NULL) {
00466     return(fei_epetra_crs->getMatrix().get());
00467   }
00468 
00469   LinProbMgr_EpetraBasic* lpm_epetrabasic =
00470     dynamic_cast<LinProbMgr_EpetraBasic*>(fei_lpm->getMatrix().get());
00471   if (lpm_epetrabasic == 0) {
00472     throw std::runtime_error("fei Trilinos_Helpers ERROR getting LinProbMgr_EpetraBasic");
00473   }
00474 
00475   return(lpm_epetrabasic->get_A_matrix().get());
00476 }
00477 
00478 
00479 void get_Epetra_pointers(fei::SharedPtr<fei::Matrix> feiA,
00480                          fei::SharedPtr<fei::Vector> feix,
00481                          fei::SharedPtr<fei::Vector> feib,
00482                          Epetra_CrsMatrix*& crsA,
00483                          Epetra_Operator*& opA,
00484                          Epetra_MultiVector*& x,
00485                          Epetra_MultiVector*& b)
00486 {
00487   x = get_Epetra_MultiVector(feix.get(), true);
00488   b = get_Epetra_MultiVector(feib.get(), false);
00489 
00490   const char* matname = feiA->typeName();
00491   if (!strcmp(matname, "Epetra_VbrMatrix")) {
00492     Epetra_VbrMatrix* A = get_Epetra_VbrMatrix(feiA.get());
00493     opA = A;
00494   }
00495   else {
00496     crsA = get_Epetra_CrsMatrix(feiA.get());
00497     opA = crsA;
00498   }
00499 }
00500 
00501 int zero_Epetra_VbrMatrix(Epetra_VbrMatrix* mat)
00502 {
00503   const Epetra_CrsGraph& crsgraph = mat->Graph();
00504   const Epetra_BlockMap& rowmap = crsgraph.RowMap();
00505   const Epetra_BlockMap& colmap = crsgraph.ColMap();
00506   mat->RowMatrixRowMap();//generate point objects
00507   int maxBlkRowSize = mat->GlobalMaxRowDim();
00508   int maxBlkColSize = mat->GlobalMaxColDim();
00509   std::vector<double> zeros(maxBlkRowSize*maxBlkColSize, 0);
00510   int numMyRows = rowmap.NumMyElements();
00511   int* myRows = rowmap.MyGlobalElements();
00512   for(int i=0; i<numMyRows; ++i) {
00513     int row = myRows[i];
00514     int rowlength = 0;
00515     int* colindicesView = NULL;
00516     int localrow = rowmap.LID(row);
00517     int err = crsgraph.ExtractMyRowView(localrow, rowlength, colindicesView);
00518     if (err != 0) {
00519       return err;
00520     }
00521     err = mat->BeginReplaceMyValues(localrow, rowlength, colindicesView);
00522     if (err != 0) {
00523       return err;
00524     }
00525     int blkRowSize = rowmap.ElementSize(localrow);
00526     for(int j=0; j<rowlength; ++j) {
00527       int blkColSize = colmap.ElementSize(colindicesView[j]);
00528       err = mat->SubmitBlockEntry(&zeros[0], maxBlkRowSize,
00529                             blkRowSize, blkColSize);
00530       if (err != 0) {
00531         return err;
00532       }
00533     }
00534     err = mat->EndSubmitEntries();
00535     if (err != 0) {
00536       return err;
00537     }
00538   }
00539 
00540   return 0;
00541 }
00542 
00543 #endif //HAVE_FEI_EPETRA
00544 
00545 }//namespace Trilinos_Helpers
00546 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends