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

Generated on 9 Apr 2008 with Doxygen (ver 1.5.2) [logfile].