1
0
mirror of https://github.com/paboyle/Grid.git synced 2024-11-15 02:05:37 +00:00

update convergence check

This commit is contained in:
Yong-Chull Jang 2018-03-05 18:37:51 -05:00
parent 11219a8f7a
commit 386b4fcb04

View File

@ -440,41 +440,41 @@ public:
// Convergence test // Convergence test
Glog <<" #Convergence test: "<<std::endl; Glog <<" #Convergence test: "<<std::endl;
Nconv = 0;
for(int k = 0; k<Nr; ++k) B[k]=0.0; for(int k = 0; k<Nr; ++k) B[k]=0.0;
for(int j = 0; j<Nr; j+=Nconv_test_interval){ for(int j = 0; j<Nr; j+=Nconv_test_interval){
Glog <<" #rotation for evec" if ( j/Nconv_test_interval == Nconv ) {
Glog <<" #rotation for next check point evec"
<< std::setw(4)<< std::setiosflags(std::ios_base::right) << std::setw(4)<< std::setiosflags(std::ios_base::right)
<< "["<< j <<"]" <<std::endl; << "["<< j <<"]" <<std::endl;
for(int k = 0; k<Nr; ++k){ for(int k = 0; k<Nr; ++k){
B[j].checkerboard = evec[k].checkerboard; B[j].checkerboard = evec[k].checkerboard;
B[j] += evec[k]*Qt(k,j); B[j] += evec[k]*Qt(k,j);
} }
}
Nconv = 0; _Linop.HermOp(B[j],v);
for(int i=0; i<Nr; i+=Nconv_test_interval){ RealD vnum = real(innerProduct(B[j],v)); // HermOp.
RealD vden = norm2(B[j]);
_Linop.HermOp(B[i],v); eval2[j] = vnum/vden;
RealD vnum = real(innerProduct(B[i],v)); // HermOp. v -= eval2[j]*B[j];
RealD vden = norm2(B[i]);
eval2[i] = vnum/vden;
v -= eval2[i]*B[i];
RealD vv = norm2(v); RealD vv = norm2(v);
resid[i] = vv; resid[j] = vv;
std::cout.precision(13); std::cout.precision(13);
std::cout << "[" << std::setw(4)<< std::setiosflags(std::ios_base::right) <<i<<"] "; std::cout << "[" << std::setw(4)<< std::setiosflags(std::ios_base::right) <<j<<"] ";
std::cout << "eval = "<<std::setw(20)<< std::setiosflags(std::ios_base::left)<< eval2[i]; std::cout << "eval = "<<std::setw(20)<< std::setiosflags(std::ios_base::left)<< eval2[j];
std::cout << " resid^2 = "<< std::setw(20)<< std::setiosflags(std::ios_base::right)<< vv<< std::endl; std::cout << " resid^2 = "<< std::setw(20)<< std::setiosflags(std::ios_base::right)<< vv<< std::endl;
// change the criteria as evals are supposed to be sorted, all evals smaller(larger) than Nstop should have converged // change the criteria as evals are supposed to be sorted, all evals smaller(larger) than Nstop should have converged
//if( (vv<eresid*eresid) && (i == Nconv) ){ //if( (vv<eresid*eresid) && (i == Nconv) ){
if (vv<eresid*eresid) { if (vv<eresid*eresid) {
Iconv[Nconv] = i; Iconv[Nconv] = j;
++Nconv; ++Nconv;
} }
} else {
} // i-loop end break;
}
} // j-loop end
Glog <<" #modes converged: "<<Nconv<<std::endl; Glog <<" #modes converged: "<<Nconv<<std::endl;
for(int i=0; i<Nconv; ++i){ for(int i=0; i<Nconv; ++i){