include/deal.II/lac/shifted_matrix.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: shifted_matrix.h 25345 2012-03-31 08:37:04Z bangerth @f$
00003 //
00004 //    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 #ifndef __deal2__shifted_matrix_h
00013 #define __deal2__shifted_matrix_h
00014 
00015 
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/smartpointer.h>
00018 
00019 DEAL_II_NAMESPACE_OPEN
00020 
00034 template<class MATRIX>
00035 class ShiftedMatrix
00036 {
00037   public:
00042     ShiftedMatrix (const MATRIX& A, const double sigma);
00043 
00047     void shift (const double sigma);
00048 
00052     double shift () const;
00053 
00057     template <class VECTOR>
00058     void vmult (VECTOR& dst, const VECTOR& src) const;
00059 
00063     template <class VECTOR>
00064     double residual (VECTOR& dst, const VECTOR& src, const VECTOR& rhs) const;
00065 
00066   private:
00070     SmartPointer<const MATRIX,ShiftedMatrix<MATRIX> > A;
00071 
00075     //    VECTOR aux;
00079     double sigma;
00080 };
00081 
00082 
00083 
00093 template<class MATRIX, class MASSMATRIX, class VECTOR>
00094 class ShiftedMatrixGeneralized
00095 {
00096   public:
00101     ShiftedMatrixGeneralized (const MATRIX& A,
00102                               const MASSMATRIX& M,
00103                               const double sigma);
00104 
00108     void shift (const double sigma);
00109 
00113     double shift () const;
00114 
00118     void vmult (VECTOR& dst, const VECTOR& src) const;
00119 
00123     double residual (VECTOR& dst, const VECTOR& src, const VECTOR& rhs) const;
00124 
00125   private:
00129     SmartPointer<const MATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VECTOR> > A;
00133     SmartPointer<const MASSMATRIX,ShiftedMatrixGeneralized<MATRIX,MASSMATRIX,VECTOR> > M;
00134 
00138     VECTOR aux;
00139 
00143     double sigma;
00144 };
00145 
00146 
00148 //---------------------------------------------------------------------------
00149 
00150 template <class MATRIX>
00151 inline
00152 ShiftedMatrix<MATRIX>::ShiftedMatrix (const MATRIX& A, const double sigma)
00153                 :
00154                 A(&A), sigma(sigma)
00155 {}
00156 
00157 
00158 
00159 template <class MATRIX>
00160 inline void
00161 ShiftedMatrix<MATRIX>::shift (const double s)
00162 {
00163   sigma = s;
00164 }
00165 
00166 
00167 template <class MATRIX>
00168 inline double
00169 ShiftedMatrix<MATRIX>::shift () const
00170 {
00171   return sigma;
00172 }
00173 
00174 
00175 
00176 template <class MATRIX>
00177 template <class VECTOR>
00178 inline void
00179 ShiftedMatrix<MATRIX>::vmult (VECTOR& dst, const VECTOR& src) const
00180 {
00181   A->vmult(dst, src);
00182   if (sigma != 0.)
00183     dst.add(sigma, src);
00184 }
00185 
00186 
00187 template <class MATRIX>
00188 template <class VECTOR>
00189 inline double
00190 ShiftedMatrix<MATRIX>::residual (VECTOR& dst,
00191                                  const VECTOR& src,
00192                                  const VECTOR& rhs) const
00193 {
00194   A->vmult(dst, src);
00195   if (sigma != 0.)
00196     dst.add(sigma, src);
00197   dst.sadd(-1.,1.,rhs);
00198   return dst.l2_norm ();
00199 }
00200 
00201 
00202 //---------------------------------------------------------------------------
00203 template <class MATRIX, class MASSMATRIX, class VECTOR>
00204 inline
00205 ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>
00206 ::ShiftedMatrixGeneralized (const MATRIX& A,
00207                             const MASSMATRIX& M,
00208                             const double sigma)
00209                 :
00210                 A(&A), M(&M), sigma(sigma)
00211 {}
00212 
00213 
00214 template <class MATRIX, class MASSMATRIX, class VECTOR>
00215 inline void
00216 ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::shift (const double s)
00217 {
00218   sigma = s;
00219 }
00220 
00221 template <class MATRIX, class MASSMATRIX, class VECTOR>
00222 inline double
00223 ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::shift () const
00224 {
00225   return sigma;
00226 }
00227 
00228 
00229 template <class MATRIX, class MASSMATRIX, class VECTOR>
00230 inline void
00231 ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::vmult (VECTOR& dst,
00232                                                      const VECTOR& src) const
00233 {
00234   A->vmult(dst, src);
00235   if (sigma != 0.)
00236     {
00237       aux.reinit(dst);
00238       M.vmult(aux, src);
00239       dst.add(sigma, aux);
00240     }
00241 }
00242 
00243 
00244 template <class MATRIX, class MASSMATRIX, class VECTOR>
00245 inline double
00246 ShiftedMatrixGeneralized<MATRIX, MASSMATRIX, VECTOR>::residual (VECTOR& dst,
00247                                                         const VECTOR& src,
00248                                                         const VECTOR& rhs) const
00249 {
00250   A->vmult(dst, src);
00251   if (sigma != 0.)
00252     {
00253       aux.reinit(dst);
00254       M.vmult(aux, src);
00255       dst.add(sigma, aux);
00256     }
00257   dst.sadd(-1.,1.,rhs);
00258   return dst.l2_norm ();
00259 }
00260 
00261 
00262 DEAL_II_NAMESPACE_CLOSE
00263 
00264 #endif
00265 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Wed May 23 2012 06:07:39 by doxygen 1.7.3