formatting

1 files changed, 158 insertions, 156 deletions

diff --git a/unsupported/Eigen/src/IterativeSolvers/GMRES.h b/unsupported/Eigen/src/IterativeSolvers/GMRES.h
index e70a864e1..3cbc2f525 100644
--- a/unsupported/Eigen/src/IterativeSolvers/GMRES.h
+++ b/unsupported/Eigen/src/IterativeSolvers/GMRES.h
@@ -16,189 +16,191 @@ namespace Eigen {
 namespace internal {
 
 /**
- * Generalized Minimal Residual Algorithm based on the
- * Arnoldi algorithm implemented with Householder reflections.
- *
- * Parameters:
- *  \param mat       matrix of linear system of equations
- *  \param Rhs       right hand side vector of linear system of equations
- *  \param x         on input: initial guess, on output: solution
- *  \param precond   preconditioner used
- *  \param iters     on input: maximum number of iterations to perform
- *                   on output: number of iterations performed
- *  \param restart   number of iterations for a restart
- *  \param tol_error on input: relative residual tolerance
- *                   on output: residuum achieved
- *
- * \sa IterativeMethods::bicgstab()
- *
- *
- * For references, please see:
- *
- * Saad, Y. and Schultz, M. H.
- * GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems.
- * SIAM J.Sci.Stat.Comp. 7, 1986, pp. 856 - 869.
- *
- * Saad, Y.
- * Iterative Methods for Sparse Linear Systems.
- * Society for Industrial and Applied Mathematics, Philadelphia, 2003.
- *
- * Walker, H. F.
- * Implementations of the GMRES method.
- * Comput.Phys.Comm. 53, 1989, pp. 311 - 320.
- *
- * Walker, H. F.
- * Implementation of the GMRES Method using Householder Transformations.
- * SIAM J.Sci.Stat.Comp. 9, 1988, pp. 152 - 163.
- *
- */
+* Generalized Minimal Residual Algorithm based on the
+* Arnoldi algorithm implemented with Householder reflections.
+*
+* Parameters:
+*  \param mat       matrix of linear system of equations
+*  \param Rhs       right hand side vector of linear system of equations
+*  \param x         on input: initial guess, on output: solution
+*  \param precond   preconditioner used
+*  \param iters     on input: maximum number of iterations to perform
+*                   on output: number of iterations performed
+*  \param restart   number of iterations for a restart
+*  \param tol_error on input: relative residual tolerance
+*                   on output: residuum achieved
+*
+* \sa IterativeMethods::bicgstab()
+*
+*
+* For references, please see:
+*
+* Saad, Y. and Schultz, M. H.
+* GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems.
+* SIAM J.Sci.Stat.Comp. 7, 1986, pp. 856 - 869.
+*
+* Saad, Y.
+* Iterative Methods for Sparse Linear Systems.
+* Society for Industrial and Applied Mathematics, Philadelphia, 2003.
+*
+* Walker, H. F.
+* Implementations of the GMRES method.
+* Comput.Phys.Comm. 53, 1989, pp. 311 - 320.
+*
+* Walker, H. F.
+* Implementation of the GMRES Method using Householder Transformations.
+* SIAM J.Sci.Stat.Comp. 9, 1988, pp. 152 - 163.
+*
+*/
 template<typename MatrixType, typename Rhs, typename Dest, typename Preconditioner>
 bool gmres(const MatrixType & mat, const Rhs & rhs, Dest & x, const Preconditioner & precond,
-		Index &iters, const Index &restart, typename Dest::RealScalar & tol_error) {
+    Index &iters, const Index &restart, typename Dest::RealScalar & tol_error) {
 
-	using std::sqrt;
-	using std::abs;
+  using std::sqrt;
+  using std::abs;
 
-	typedef typename Dest::RealScalar RealScalar;
-	typedef typename Dest::Scalar Scalar;
-	typedef Matrix < Scalar, Dynamic, 1 > VectorType;
-	typedef Matrix < Scalar, Dynamic, Dynamic > FMatrixType;
+  typedef typename Dest::RealScalar RealScalar;
+  typedef typename Dest::Scalar Scalar;
+  typedef Matrix < Scalar, Dynamic, 1 > VectorType;
+  typedef Matrix < Scalar, Dynamic, Dynamic > FMatrixType;
 
-	RealScalar tol = tol_error;
-	const Index maxIters = iters;
-	iters = 0;
+  RealScalar tol = tol_error;
+  const Index maxIters = iters;
+  iters = 0;
 
-	const Index m = mat.rows();
+  const Index m = mat.rows();
 
-	// residual and preconditioned residual
-	VectorType p0 = rhs - mat*x;
-	VectorType r0 = precond.solve(p0);
+  // residual and preconditioned residual
+  VectorType p0 = rhs - mat*x;
+  VectorType r0 = precond.solve(p0);
   VectorType t(m), v(m), workspace(m);
 
-	const RealScalar r0Norm = r0.norm();
+  const RealScalar r0Norm = r0.norm();
 
-	// is initial guess already good enough?
-	if(r0Norm == 0) {
-		tol_error=0;
-		return true;
-	}
+  // is initial guess already good enough?
+  if(r0Norm == 0)
+  {
+    tol_error = 0;
+    return true;
+  }
 
-	// storage for Hessenberg matrix and Householder data
-	FMatrixType H = FMatrixType::Zero(m, restart + 1);
-	VectorType w = VectorType::Zero(restart + 1);
-	VectorType tau = VectorType::Zero(restart + 1);
+  // storage for Hessenberg matrix and Householder data
+  FMatrixType H   = FMatrixType::Zero(m, restart + 1);
+  VectorType w    = VectorType::Zero(restart + 1);
+  VectorType tau  = VectorType::Zero(restart + 1);
 
-	// storage for Jacobi rotations
-	std::vector < JacobiRotation < Scalar > > G(restart);
+  // storage for Jacobi rotations
+  std::vector < JacobiRotation < Scalar > > G(restart);
 
-	// generate first Householder vector
+  // generate first Householder vector
   Ref<VectorType> H0_tail = H.col(0).tail(m - 1);
-	RealScalar beta;
-	r0.makeHouseholder(H0_tail, tau.coeffRef(0), beta);
-	w(0)=(Scalar) beta;
-
-	for (Index k = 1; k <= restart; ++k) {
-
-		++iters;
+  RealScalar beta;
+  r0.makeHouseholder(H0_tail, tau.coeffRef(0), beta);
+  w(0) = Scalar(beta);
 
-		v = VectorType::Unit(m, k - 1);
+  for (Index k = 1; k <= restart; ++k)
+  {
+    ++iters;
 
-		// apply Householder reflections H_{1} ... H_{k-1} to v
-		for (Index i = k - 1; i >= 0; --i) {
-			v.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
-		}
+    v = VectorType::Unit(m, k - 1);
 
-		// apply matrix M to v:  v = mat * v;
-		t=mat*v;
-		v=precond.solve(t);
+    // apply Householder reflections H_{1} ... H_{k-1} to v
+    for (Index i = k - 1; i >= 0; --i) {
+      v.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
+    }
 
-		// apply Householder reflections H_{k-1} ... H_{1} to v
-		for (Index i = 0; i < k; ++i) {
-			v.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
-		}
+    // apply matrix M to v:  v = mat * v;
+    t=mat*v;
+    v=precond.solve(t);
 
-		if (v.tail(m - k).norm() != 0.0) {
-			if (k <= restart) {
+    // apply Householder reflections H_{k-1} ... H_{1} to v
+    for (Index i = 0; i < k; ++i) {
+      v.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
+    }
 
-				// generate new Householder vector
+    if (v.tail(m - k).norm() != 0.0)
+    {
+      if (k <= restart)
+      {
+        // generate new Householder vector
         Ref<VectorType> Hk_tail = H.col(k).tail(m - k - 1);
-        
-				v.tail(m - k).makeHouseholder(Hk_tail, tau.coeffRef(k), beta);
-
-				// apply Householder reflection H_{k} to v
-				v.tail(m - k).applyHouseholderOnTheLeft(Hk_tail, tau.coeffRef(k), workspace.data());
-
-			}
-		}
-
-		if (k > 1) {
-			for (Index i = 0; i < k - 1; ++i) {
-				// apply old Givens rotations to v
-				v.applyOnTheLeft(i, i + 1, G[i].adjoint());
-			}
-		}
-
-		if (k<m && v(k) != (Scalar) 0) {
-
-			// determine next Givens rotation
-			G[k - 1].makeGivens(v(k - 1), v(k));
-
-			// apply Givens rotation to v and w
-			v.applyOnTheLeft(k - 1, k, G[k - 1].adjoint());
-			w.applyOnTheLeft(k - 1, k, G[k - 1].adjoint());
-		}
-
-		// insert coefficients into upper matrix triangle
-		H.col(k - 1).head(k) = v.head(k);
-
-		bool stop=(k==m || abs(w(k)) < tol * r0Norm || iters == maxIters);
+        v.tail(m - k).makeHouseholder(Hk_tail, tau.coeffRef(k), beta);
 
-		if (stop || k == restart) {
-
-			// solve upper triangular system
-			VectorType y = w.head(k);
-			H.topLeftCorner(k, k).template triangularView < Eigen::Upper > ().solveInPlace(y);
-
-			// use Horner-like scheme to calculate solution vector
-			VectorType x_new = y(k - 1) * VectorType::Unit(m, k - 1);
-
-			// apply Householder reflection H_{k} to x_new
-			x_new.tail(m - k + 1).applyHouseholderOnTheLeft(H.col(k - 1).tail(m - k), tau.coeffRef(k - 1), workspace.data());
-
-			for (Index i = k - 2; i >= 0; --i) {
-				x_new += y(i) * VectorType::Unit(m, i);
-				// apply Householder reflection H_{i} to x_new
-				x_new.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
-			}
-
-			x += x_new;
-
-			if (stop) {
-				return true;
-			} else {
-				k=0;
-
-				// reset data for restart
-				p0 = rhs - mat*x;
-				r0 = precond.solve(p0);
+        // apply Householder reflection H_{k} to v
+        v.tail(m - k).applyHouseholderOnTheLeft(Hk_tail, tau.coeffRef(k), workspace.data());
+      }
+    }
 
-				// clear Hessenberg matrix and Householder data
-				H = FMatrixType::Zero(m, restart + 1);
-				w = VectorType::Zero(restart + 1);
-				tau = VectorType::Zero(restart + 1);
+    if (k > 1)
+    {
+      for (Index i = 0; i < k - 1; ++i)
+      {
+        // apply old Givens rotations to v
+        v.applyOnTheLeft(i, i + 1, G[i].adjoint());
+      }
+    }
 
-				// generate first Householder vector
-				r0.makeHouseholder(H0_tail, tau.coeffRef(0), beta);
-				w(0)=(Scalar) beta;
+    if (k<m && v(k) != (Scalar) 0)
+    {
+      // determine next Givens rotation
+      G[k - 1].makeGivens(v(k - 1), v(k));
 
-			}
+      // apply Givens rotation to v and w
+      v.applyOnTheLeft(k - 1, k, G[k - 1].adjoint());
+      w.applyOnTheLeft(k - 1, k, G[k - 1].adjoint());
+    }
 
-		}
+    // insert coefficients into upper matrix triangle
+    H.col(k - 1).head(k) = v.head(k);
 
+    bool stop = (k==m || abs(w(k)) < tol * r0Norm || iters == maxIters);
 
-	}
+    if (stop || k == restart)
+    {
+      // solve upper triangular system
+      VectorType y = w.head(k);
+      H.topLeftCorner(k, k).template triangularView <Upper>().solveInPlace(y);
+
+      // use Horner-like scheme to calculate solution vector
+      VectorType x_new = y(k - 1) * VectorType::Unit(m, k - 1);
+
+      // apply Householder reflection H_{k} to x_new
+      x_new.tail(m - k + 1).applyHouseholderOnTheLeft(H.col(k - 1).tail(m - k), tau.coeffRef(k - 1), workspace.data());
+
+      for (Index i = k - 2; i >= 0; --i)
+      {
+        x_new += y(i) * VectorType::Unit(m, i);
+        // apply Householder reflection H_{i} to x_new
+        x_new.tail(m - i).applyHouseholderOnTheLeft(H.col(i).tail(m - i - 1), tau.coeffRef(i), workspace.data());
+      }
+
+      x += x_new;
+
+      if (stop)
+      {
+        return true;
+      }
+      else
+      {
+        k=0;
+
+        // reset data for restart
+        p0 = rhs - mat*x;
+        r0 = precond.solve(p0);
+
+        // clear Hessenberg matrix and Householder data
+        H = FMatrixType::Zero(m, restart + 1);
+        w = VectorType::Zero(restart + 1);
+        tau = VectorType::Zero(restart + 1);
+
+        // generate first Householder vector
+        r0.makeHouseholder(H0_tail, tau.coeffRef(0), beta);
+        w(0)=(Scalar) beta;
+      }
+    }
+  }
 
-	return false;
+  return false;
 
 }
 
@@ -314,8 +316,8 @@ public:
         failed = true;
     }
     m_info = failed ? NumericalIssue
-           : m_error <= Base::m_tolerance ? Success
-           : NoConvergence;
+          : m_error <= Base::m_tolerance ? Success
+          : NoConvergence;
     m_isInitialized = true;
   }