* Initialize for a new matrix.
*/
void initialize (const FullMatrix<number>& A);
+
+ void operator() (const unsigned int level,
+ VECTOR &dst,
+ const VECTOR &src) const;
+ private:
/**
- * Release matrix pointer.
+ * Matrix for QR-factorization.
*/
- void clear ();
+ Householder<number> householder;
+};
+/**
+ * Coarse grid solver using singular value decomposition of LAPACK matrices.
+ *
+ * Upon initialization, the singular value decomposition of the matrix is
+ * computed. then, the operator() uses
+ *
+ * @author Guido Kanschat, 2003, 2012
+ */
+template<typename number = double, class VECTOR = Vector<number> >
+class MGCoarseGridSVD : public MGCoarseGridBase<VECTOR>
+{
+ public:
+ /**
+ * Constructor leaving an
+ * uninitialized object.
+ */
+ MGCoarseGridSVD ();
+
/**
- * Solution operator, defined in
- * the base class.
+ * Initialize for a new
+ * matrix. This resets the
+ * dimensions to the
*/
+ void initialize (const FullMatrix<number>& A, const double threshold = 1.e-12);
+
void operator() (const unsigned int level,
VECTOR &dst,
const VECTOR &src) const;
-
+
+ /**
+ * Write the singular values to #deallog.
+ */
+ void log () const;
+
private:
+
/**
- * Matrix for QR-factorization.
+ * Matrix for singular value decomposition.
*/
- Householder<number> householder;
+ LAPACKFullMatrix<number> matrix;
};
/*@}*/
template<typename number, class VECTOR>
void
-MGCoarseGridHouseholder<number, VECTOR>::clear()
+MGCoarseGridHouseholder<number, VECTOR>::operator() (
+ const unsigned int,
+ VECTOR &dst,
+ const VECTOR &src) const
+{
+ householder.least_squares(dst, src);
+}
+
+//---------------------------------------------------------------------------
+
+template<typename number, class VECTOR>
+inline
+MGCoarseGridSVD<number, VECTOR>::MGCoarseGridSVD()
{}
template<typename number, class VECTOR>
void
-MGCoarseGridHouseholder<number, VECTOR>::operator() (
+MGCoarseGridSVD<number, VECTOR>::initialize(
+ const FullMatrix<number>& A,
+ double threshold)
+{
+ matrix.reinit(A.n_rows(), A.n_cols());
+ matrix = A;
+ matrix.compute_inverse_svd(threshold);
+}
+
+
+template<typename number, class VECTOR>
+void
+MGCoarseGridSVD<number, VECTOR>::operator() (
const unsigned int,
VECTOR &dst,
const VECTOR &src) const
{
- householder.least_squares(dst, src);
+ matrix.vmult(dst, src);
}
+
+template<typename number, class VECTOR>
+void
+MGCoarseGridSVD<number, VECTOR>::log() const
+{
+ const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
+
+ for (unsigned int i=0;i<n;++i)
+ deallog << ' ' << matrix.singular_value(i);
+ deallog << std::endl;
+}
+
+
#endif // DOXYGEN
DEAL_II_NAMESPACE_CLOSE