transfer.hh

Go to the documentation of this file.
00001 #ifndef DUNE_AMGTRANSFER_HH
00002 #define DUNE_AMGTRANSFER_HH
00003 
00004 #include<dune/istl/paamg/pinfo.hh>
00005 #include<dune/istl/owneroverlapcopy.hh>
00006 #include<dune/istl/paamg/aggregates.hh>
00007 #include<dune/common/exceptions.hh>
00008 #include "pinfo.hh"
00009 
00010 namespace Dune
00011 {
00012   namespace Amg
00013   {
00014     
00025     template<class V1, class V2, class T>
00026     class Transfer
00027     {
00028     
00029     public:
00030       typedef V1 Vertex;
00031       typedef V2 Vector;
00032       
00033       static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine, 
00034                              typename Vector::field_type damp);
00035 
00036       static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
00037                            T& comm);
00038     };
00039 
00040     template<class V,class B>
00041     class Transfer<V,BlockVector<B>, SequentialInformation>
00042     {
00043     public:
00044       typedef V Vertex;
00045       typedef BlockVector<B> Vector;
00046       static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine, 
00047                              typename Vector::field_type damp);
00048 
00049       static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
00050                            const SequentialInformation& comm);
00051     };
00052 
00053 #if HAVE_MPI
00054     
00055     template<class V,class B, class T>
00056     class Transfer<V,BlockVector<B>,ParallelInformation<T> > 
00057     {
00058     public:
00059       typedef V Vertex;
00060       typedef BlockVector<B> Vector;
00061       static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine, 
00062                              typename Vector::field_type damp);
00063 
00064       static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
00065                            ParallelInformation<T>& comm);
00066     };
00067 
00068     template<class V,class B, class T1, class T2>
00069     class Transfer<V,BlockVector<B>,OwnerOverlapCopyCommunication<T1,T2> > 
00070     {
00071     public:
00072       typedef V Vertex;
00073       typedef BlockVector<B> Vector;
00074       static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine, 
00075                              typename Vector::field_type damp);
00076 
00077       static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
00078                            OwnerOverlapCopyCommunication<T1,T2>& comm);
00079     };
00080 
00081 #endif
00082 
00083     template<class V, class B>
00084     inline void Transfer<V,BlockVector<B>,SequentialInformation>::prolongate(const AggregatesMap<Vertex>& aggregates,
00085                                                                              Vector& coarse, Vector& fine, 
00086                                                                              typename Vector::field_type damp)
00087     {
00088       typedef typename Vector::iterator Iterator;
00089 
00090       Iterator end = fine.end();
00091       
00092       coarse *= damp;
00093       
00094       for(Iterator block=fine.begin(); block != end; ++block){
00095         const Vertex& vertex = aggregates[block.index()];
00096         if(vertex != AggregatesMap<Vertex>::ISOLATED)
00097           *block += coarse[aggregates[block.index()]];
00098       }
00099     }
00100 
00101     template<class V, class B>
00102     inline void Transfer<V,BlockVector<B>,SequentialInformation>::restrict(const AggregatesMap<Vertex>& aggregates,
00103                                                                            Vector& coarse, 
00104                                                                            const Vector& fine,
00105                                                                            const SequentialInformation& comm)
00106     {
00107       // Set coarse vector to zero
00108       coarse=0;
00109       
00110       typedef typename Vector::const_iterator Iterator;
00111       Iterator end = fine.end();
00112       for(Iterator block=fine.begin(); block != end; ++block){
00113         const Vertex& vertex = aggregates[block.index()];
00114         if(vertex != AggregatesMap<Vertex>::ISOLATED)
00115           coarse[vertex] += *block;
00116       }
00117     }
00118 
00119 #if HAVE_MPI
00120     template<class V, class B, class T>
00121     inline void Transfer<V,BlockVector<B>,ParallelInformation<T> >::prolongate(const AggregatesMap<Vertex>& aggregates,
00122                                                                                Vector& coarse, Vector& fine, 
00123                                                                                typename Vector::field_type damp)
00124     {
00125       Transfer<V,BlockVector<B>,SequentialInformation>::prolongate(aggregates, coarse, fine, damp);
00126     }
00127     
00128     template<class V, class B, class T>
00129     inline void Transfer<V,BlockVector<B>,ParallelInformation<T> >::restrict(const AggregatesMap<Vertex>& aggregates,
00130                                                                              Vector& coarse, const Vector& fine, 
00131                                                                              ParallelInformation<T>& comm)
00132     {
00133       Transfer<V,BlockVector<B>,SequentialInformation>::restrict(aggregates, coarse, fine, SequentialInformation());
00134     }
00135 
00136     template<class V, class B, class T1, class T2>
00137     inline void Transfer<V,BlockVector<B>,OwnerOverlapCopyCommunication<T1,T2> >::prolongate(const AggregatesMap<Vertex>& aggregates,
00138                                                                                Vector& coarse, Vector& fine, 
00139                                                                                typename Vector::field_type damp)
00140     {
00141       Transfer<V,BlockVector<B>,SequentialInformation>::prolongate(aggregates, coarse, fine, damp);
00142     }
00143     
00144     template<class V, class B, class T1, class T2>
00145     inline void Transfer<V,BlockVector<B>,OwnerOverlapCopyCommunication<T1,T2> >::restrict(const AggregatesMap<Vertex>& aggregates,
00146                                                                              Vector& coarse, const Vector& fine, 
00147                                                                              OwnerOverlapCopyCommunication<T1,T2>& comm)
00148     {
00149       Transfer<V,BlockVector<B>,SequentialInformation>::restrict(aggregates, coarse, fine, SequentialInformation());
00150       //      comm.project(coarse);
00151     }
00152 #endif
00153 
00154       }// namspace Amg
00155     } // namspace Dune
00156 #endif

Generated on 12 Dec 2007 with Doxygen (ver 1.5.1)