QDP_error_exit("too many CG iterations: count = %d", n_count);
QDPIO::cerr << "Nonconvergence Warning" << endl;
void InvCG2_a(const LinearOperator< multi1d<T> >& M,	      const multi1d<T> & chi,	      multi1d<T>& psi,	      const Real& RsdCG, 	      int MaxCG, 	      int& n_count){  const int N = psi.size();  const OrderedSubset& s = M.subset();  Real chi_sq =  Real(norm2(chi,s));  QDPIO::cout << "chi_norm = " << sqrt(chi_sq) << endl;  Real rsd_sq = (RsdCG * RsdCG) * chi_sq;  //                                            +  //  r[0]  :=  Chi - A . Psi[0]    where  A = M  . M      //                      +  //  r  :=  [ Chi  -  M(u)  . M(u) . psi ]  multi1d<T> r(N), mp(N), mmp(N);  M(mp, psi, PLUS);  M(mmp, mp, MINUS);  for(int n=0; n < N; ++n)    r[n][s] = chi[n] - mmp[n];#ifdef PRINT_5D_RESID   for(int n=0; n < N; n++) {    Double norm_r = norm2(r[n],s);    if( toBool( norm_r > Double(1.0e-20)) ) {      QDPIO::cout << "Iteration 0  r[" << n << "] = " << norm_r << endl;    }  }#endif  //  p[1]  :=  r[0]  multi1d<T> p(N);  for(int n=0; n < N; ++n)    p[n][s] = r[n];    //  Cp = |r[0]|^2  Double cp = norm2(r, s);   	       	   /* 2 Nc Ns  flops */  QDPIO::cout << "InvCG: k = 0  cp = " << cp << "  rsd_sq = " << rsd_sq << endl;  //  IF |r[0]| <= RsdCG |Chi| THEN RETURN;  if ( toBool(cp  <=  rsd_sq) )  {    n_count = 0;    return;  }  //  //  FOR k FROM 1 TO MaxCG DO  //  Real a, b;  Double c, d;    for(int k = 1; k <= MaxCG; ++k)  {    //  c  =  | r[k-1] |**2    c = cp;    //  a[k] := | r[k-1] |**2 / < p[k], Ap[k] > ;    //      	       	       	       	       	  +    //  First compute  d  =  < p, A.p >  =  < p, M . M . p >  =  < M.p, M.p >    //  Mp = M(u) * p    M(mp, p, PLUS);    //  d = | mp | ** 2    d = norm2(mp, s);	/* 2 Nc Ns  flops */    a = Real(c)/Real(d);    //  Psi[k] += a[k] p[k]    for(int n=0; n < N; ++n)      psi[n][s] += a * p[n];	/* 2 Nc Ns  flops */    //  r[k] -= a[k] A . p[k] ;    //      	       +            +    //  r  =  r  -  M(u)  . Mp  =  M  . M . p  =  A . p    M(mmp, mp, MINUS);    for(int n=0; n < N; ++n)      r[n][s] -= a * mmp[n];#ifdef PRINT_5D_RESID    for(int n=0; n < N; n++) {      Double norm_r = norm2(r[n],s);      if( toBool( norm_r > Double(1.0e-20)) )	QDPIO::cout << "Iteration " << k << " r[" << n << "] = " << norm_r << endl;    }#endif    //  IF |r[k]| <= RsdCG |Chi| THEN RETURN;    //  cp  =  | r[k] |**2    cp = norm2(r, s);	                /* 2 Nc Ns  flops */    QDPIO::cout << "InvCG: k = " << k << "  cp = " << cp << endl;    if ( toBool(cp  <=  rsd_sq) )    {      n_count = k;      return;    }    //  b[k+1] := |r[k]|**2 / |r[k-1]|**2    b = Real(cp) / Real(c);    //  p[k+1] := r[k] + b[k+1] p[k]    for(int n=0; n < N; ++n)      p[n][s] = r[n] + b*p[n];	/* Nc Ns  flops */  }  n_count = MaxCG;  QDP_error_exit("too many CG iterations: count = %d", n_count);}