Albany: a Trilinos-based PDE code: MOR_EpetraUtils.cpp Source File

Go to the documentation of this file.
00001 //*****************************************************************//
00002 //    Albany 2.0:  Copyright 2012 Sandia Corporation               //
00003 //    This Software is released under the BSD license detailed     //
00004 //    in the file "license.txt" in the top-level Albany directory  //
00005 //*****************************************************************//
00006 
00007 #include "MOR_EpetraUtils.hpp"
00008 
00009 #include "Epetra_Map.h"
00010 
00011 #include <algorithm>
00012 #include <iterator>
00013 #include <functional>
00014 
00015 namespace MOR {
00016 
00017 Teuchos::Array<EpetraGlobalIndex> getMyLIDs(
00018     const Epetra_BlockMap &map,
00019     const Teuchos::ArrayView<const EpetraGlobalIndex> &selectedGIDs)
00020 {
00021   Teuchos::Array<EpetraGlobalIndex> sortedMyGIDs(map.MyGlobalElements(), map.MyGlobalElements() + map.NumMyElements());
00022   std::sort(sortedMyGIDs.begin(), sortedMyGIDs.end());
00023 
00024   Teuchos::Array<EpetraGlobalIndex> sortedSelectedGIDs(selectedGIDs);
00025   std::sort(sortedSelectedGIDs.begin(), sortedSelectedGIDs.end());
00026 
00027   Teuchos::Array<EpetraGlobalIndex> mySelectedGIDs;
00028   std::set_intersection(sortedMyGIDs.begin(), sortedMyGIDs.end(),
00029                         sortedSelectedGIDs.begin(), sortedSelectedGIDs.end(),
00030                         std::back_inserter(mySelectedGIDs));
00031 
00032   Teuchos::Array<EpetraGlobalIndex> result;
00033   result.reserve(mySelectedGIDs.size());
00034 
00035   std::transform(
00036       mySelectedGIDs.begin(), mySelectedGIDs.end(),
00037       std::back_inserter(result),
00038       std::bind1st(std::mem_fun_ref(static_cast<int(Epetra_BlockMap::*)(EpetraGlobalIndex) const>(&Epetra_BlockMap::LID)), map));
00039 
00040   return result;
00041 }
00042 
00043 Teuchos::RCP<Epetra_Map> mapDowncast(const Epetra_BlockMap &in)
00044 {
00045   if (in.ConstantElementSize() && in.ElementSize() == 1) {
00046     return Teuchos::rcp(new Epetra_Map(static_cast<const Epetra_Map &>(in)));
00047   }
00048   return Teuchos::null;
00049 }
00050 
00051 
00052 namespace Detail {
00053 
00054 Teuchos::RCP<Epetra_Vector> memberViewImpl(const Teuchos::RCP<const Epetra_MultiVector> &mv, int i)
00055 {
00056   if (Teuchos::nonnull(mv)) {
00057     return Teuchos::rcpWithEmbeddedObjPostDestroy(new Epetra_Vector(View, *mv, i), mv);
00058   } else {
00059     return Teuchos::null;
00060   }
00061 }
00062 
00063 Teuchos::RCP<Epetra_Vector> headViewImpl(const Teuchos::RCP<const Epetra_MultiVector> &mv)
00064 {
00065   return memberViewImpl(mv, 0);
00066 }
00067 
00068 Teuchos::RCP<Epetra_MultiVector> rangeViewImpl(
00069     const Teuchos::RCP<const Epetra_MultiVector> &mv,
00070     int firstVectorRank,
00071     int vectorCount)
00072 {
00073   if (vectorCount > 0) {
00074     return Teuchos::rcpWithEmbeddedObjPostDestroy(new Epetra_MultiVector(View, *mv, firstVectorRank, vectorCount), mv);
00075   } else {
00076     return Teuchos::null;
00077   }
00078 }
00079 
00080 Teuchos::RCP<Epetra_MultiVector> tailViewImpl(const Teuchos::RCP<const Epetra_MultiVector> &mv)
00081 {
00082   if (Teuchos::nonnull(mv)) {
00083     const int remainderVecCount = mv->NumVectors() - 1;
00084     return rangeViewImpl(mv, 0, remainderVecCount);
00085   } else {
00086     return Teuchos::null;
00087   }
00088 }
00089 
00090 Teuchos::RCP<Epetra_MultiVector> truncatedViewImpl(const Teuchos::RCP<const Epetra_MultiVector> &mv, int vectorCountMax)
00091 {
00092   if (Teuchos::nonnull(mv)) {
00093     const int vectorCount = std::min(mv->NumVectors(), vectorCountMax);
00094     return rangeViewImpl(mv, 0, vectorCount);
00095   } else {
00096     return Teuchos::null;
00097   }
00098 }
00099 
00100 } // end namespace Detail
00101 
00102 
00103 Teuchos::RCP<const Epetra_Vector> headView(const Teuchos::RCP<const Epetra_MultiVector> &mv)
00104 {
00105   return Detail::headViewImpl(mv);
00106 }
00107 
00108 Teuchos::RCP<Epetra_Vector> nonConstHeadView(const Teuchos::RCP<Epetra_MultiVector> &mv)
00109 {
00110   return Detail::headViewImpl(mv);
00111 }
00112 
00113 Teuchos::RCP<Epetra_MultiVector> nonConstTailView(const Teuchos::RCP<Epetra_MultiVector> &mv)
00114 {
00115   return Detail::tailViewImpl(mv);
00116 }
00117 
00118 Teuchos::RCP<const Epetra_MultiVector> tailView(const Teuchos::RCP<const Epetra_MultiVector> &mv)
00119 {
00120   return Detail::tailViewImpl(mv);
00121 }
00122 
00123 Teuchos::RCP<const Epetra_MultiVector> truncatedView(const Teuchos::RCP<const Epetra_MultiVector> &mv, int vectorCountMax)
00124 {
00125   return Detail::truncatedViewImpl(mv, vectorCountMax);
00126 }
00127 
00128 Teuchos::RCP<Epetra_MultiVector> nonConstTruncatedView(const Teuchos::RCP<Epetra_MultiVector> &mv, int vectorCountMax)
00129 {
00130   return Detail::truncatedViewImpl(mv, vectorCountMax);
00131 }
00132 
00133 Teuchos::RCP<const Epetra_Vector> memberView(const Teuchos::RCP<const Epetra_MultiVector> &mv, int i)
00134 {
00135   return Detail::memberViewImpl(mv, i);
00136 }
00137 
00138 Teuchos::RCP<Epetra_Vector> nonConstMemberView(const Teuchos::RCP<Epetra_MultiVector> &mv, int i)
00139 {
00140   return Detail::memberViewImpl(mv, i);
00141 }
00142 
00143 
00144 void normalize(Epetra_Vector &v)
00145 {
00146   double norm2[1];
00147   v.Norm2(&norm2[0]);
00148   v.Scale(1.0 / norm2[0]);
00149 }
00150 
00151 } // namespace MOR