]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Commit this to have it in the CVS, but I will delete it immediately after again since...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 May 1999 13:30:47 +0000 (13:30 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 May 1999 13:30:47 +0000 (13:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@1332 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/solver_gmres.h

index 4d47b06bb37c1758652ea8b874a67e2e62054ffd..c5cd615dba3aa141b23163500e72f9c053c7d475 100644 (file)
@@ -136,136 +136,6 @@ class SolverGMRES : public Solver<Matrix, Vector>
 
 
 
-template <class Vector >
-void 
-Update(::Vector<double> &x, int k, FullMatrix<double> &h,
-       ::Vector<double> &s, vector<Vector> &v)
-{
-  ::Vector<double> y(s);
-
-  // Backsolve:  
-  for (int i = k; i >= 0; i--) {
-    y(i) /= h(i,i);
-    for (int j = i - 1; j >= 0; j--)
-      y(j) -= h(j,i) * y(i);
-  }
-
-  for (int j = 0; j <= k; j++)
-    x.add (y(j), v[j]);
-}
-
-
-
-
-template < class Matrix, class Vector, class Preconditioner>
-int 
-solve2(const Matrix &A, Vector &x, const Vector &b,
-       const Preconditioner &M, int &m, int &max_iter,
-       double &tol)
-{
-  FullMatrix<double> H(m+1,m);
-  
-  double resid;
-  int i, j = 1, k;
-  ::Vector<double> s(m+1), cs(m+1), sn(m+1);
-  Vector w(x.size());
-  
-  Vector tmp (b.size());
-  M (tmp, b);    
-  double normb = tmp.l2_norm();
-  Vector r(b.size());
-  A.residual (tmp, x, b);
-  M(r, tmp);
-  double beta = r.l2_norm();
-  
-  if (normb == 0.0)
-    normb = 1;
-  
-  if ((resid = r.l2_norm() / normb) <= tol) {
-    tol = resid;
-    max_iter = 0;
-    return 0;
-  }
-
-  vector<Vector> v(m+1);
-
-  while (j <= max_iter) {
-    v[0] = r;
-    v[0].scale (1.0 / beta);    // ??? r / beta
-    s = 0.0;
-    s(0) = beta;
-
-    for (i = 0; i < m && j <= max_iter; i++, j++) {
-      A.vmult (tmp, v[i]);
-      M(w, tmp);
-      for (k = 0; k <= i; k++) {
-        H(k, i) = w * v[k];
-        w.sadd (-H(k, i), v[k]);
-      }
-      H(i+1, i) = w.l2_norm();
-      v[i+1] = w;
-      v[i+1].scale (1.0 / H(i+1, i)); // ??? w / H(i+1, i)
-
-      for (k = 0; k < i; k++)
-        ApplyPlaneRotation(H(k,i), H(k+1,i), cs(k), sn(k));
-      
-      GeneratePlaneRotation(H(i,i), H(i+1,i), cs(i), sn(i));
-      ApplyPlaneRotation(H(i,i), H(i+1,i), cs(i), sn(i));
-      ApplyPlaneRotation(s(i), s(i+1), cs(i), sn(i));
-      
-      if ((resid = fabs(s(i+1)) / normb) < tol) {
-        Update(x, i, H, s, v);
-        tol = resid;
-        max_iter = j;
-        return 0;
-      }
-    }
-    Update(x, i - 1, H, s, v);
-    A.residual (tmp, x, b);
-    M(r, tmp);
-    beta = r.l2_norm();
-    if ((resid = beta / normb) < tol) {
-      tol = resid;
-      max_iter = j;
-      return 0;
-    }
-  }
-  
-  tol = resid;
-  return 1;
-}
-
-
-template<class Real> 
-void GeneratePlaneRotation(Real &dx, Real &dy, Real &cs, Real &sn)
-{
-  if (dy == 0.0) {
-    cs = 1.0;
-    sn = 0.0;
-  } else if (fabs(dy) > fabs(dx)) {
-    Real temp = dx / dy;
-    sn = 1.0 / sqrt( 1.0 + temp*temp );
-    cs = temp * sn;
-  } else {
-    Real temp = dy / dx;
-    cs = 1.0 / sqrt( 1.0 + temp*temp );
-    sn = temp * cs;
-  }
-}
-
-
-template<class Real> 
-void ApplyPlaneRotation(Real &dx, Real &dy, Real &cs, Real &sn)
-{
-  Real temp  =  cs * dx + sn * dy;
-  dy = -sn * dx + cs * dy;
-  dx = temp;
-}
-
-
-
-
-
 /* --------------------- Inline and template functions ------------------- */
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.