Changeset 8505904 in tol

Timestamp:

May 24, 2010, 2:04:31 PM (16 years ago)

Author:

Víctor de Buen <vdebuen@…>

Branches:

1941-docker-image-to-test, devel, master

Children:

b353f4f

Parents:

ab485e9

Message:

Fixes #899

git-svn-id: ​https://www.tol-project.org/svn/tolp/trunk@2193 be4846ae-f051-483a-831f-b64e7d08527f

File:

• tol/btol/bmodel/modcalc.cpp (modified) (3 diffs)

tol/btol/bmodel/modcalc.cpp

@@ -1426 +1426 @@
     PIW_ = PIW_.Sub(n-N,0,N,I);   
     BMat CI(I,1);  CI.SetAllValuesTo(0);
-    CI = MinimumResidualsSolve(PIW_,-A_, CI);
+    CI = CholeskiMinimumResidualsSolve(PIW_,-A_);
     BBool CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
     BMat err = PIW_.T()*(PIW_*CI+A_);
@@ -1432 +1432 @@
     if(!CI_isKnown ||(test>Sqrt(BDat::Tolerance())))
     {
-    //Std(BText("\nBModel::CleanInterruptions test>Sqrt(BDat::Tolerance())"));
+      Warning(BText("There are numerical problems to calculate values for "
+      "interruptions by means of Cholesky decomposition. So SVD method will "
+      "be tried and process could become slower"));
       BText method = BText("Golub_Reinsch_Mod");
       BMat U(N,I), V(I,I);
       BDiagMatrix<BDat> D(I,1.0), Dp;
       gsl_SingularValueDecomposition(PIW_,U,D,V, method);
+    //Std(BText("\nBModel::CleanInterruptions U=\n")+U.Name()+"\n");
+    //Std(BText("\nBModel::CleanInterruptions D=\n")+D.Name()+"\n");
+    //Std(BText("\nBModel::CleanInterruptions V=\n")+V.Name()+"\n");
       Dp=D.P(Sqrt(DEpsilon()));
-      CI=V*Dp*(A_.T()*U).T();
+    //Std(BText("\nBModel::CleanInterruptions Dp=\n")+Dp.Name()+"\n");
+      CI=-V*Dp*(A_.T()*U).T();
+    //Std(BText("\nBModel::CleanInterruptions CI=\n")+CI.T().Name()+"\n");
+    //Std(BText("\nBModel::CleanInterruptions PIW_*CI=\n")+(PIW_*CI).T().Name()+"\n");
       err = PIW_.T()*(PIW_*CI+A_);
+    //Std(BText("\nBModel::CleanInterruptions err=\n")+(err).T().Name()+"\n");
       test = err.MaxAbsNorm(); 
+    //Std(BText("\nBModel::CleanInterruptions test=")+test+"\n");
     }
     CI_isKnown = FrobeniusNormU(CI.Data()).IsKnown();
@@ -1460 +1470 @@
       A_  = E + PIWC;
       sigma = FrobeniusNormU(A.Data());
-    //Std(BText("\nCleanInterruptions sigma after = ") + sigma);
+      //Std(BText("\nCleanInterruptions sigma after = ") + sigma);
       if(sigma.IsKnown() && 
          !gsl_isinf(sigma.Value()) &&