00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00040 #ifndef ERGO_MAT_ACC_EXTRAPOLATE_HEADER
00041 #define ERGO_MAT_ACC_EXTRAPOLATE_HEADER
00042
00043
00044 #include <vector>
00045
00046
00047 #include "matrix_utilities.h"
00048
00049
00050
00051 template<class Treal, class Tworker>
00052 class MatAccInvestigator
00053 {
00054 public:
00055 explicit MatAccInvestigator(mat::SizesAndBlocks const & matrix_size_block_info_);
00056 void Scan(const Tworker & worker,
00057 Treal firstParam,
00058 Treal stepFactor,
00059 int nSteps);
00060 void GetScanResult(Treal* threshList_,
00061 Treal* errorList_frob_,
00062 Treal* errorList_eucl_,
00063 Treal* errorList_maxe_,
00064 Treal* timeList_);
00065 private:
00066 mat::SizesAndBlocks matrix_size_block_info;
00067 int nScanSteps;
00068 Treal baseThresh;
00069 std::vector<Treal> threshList;
00070 std::vector<Treal> errorList_frob;
00071 std::vector<Treal> errorList_eucl;
00072 std::vector<Treal> errorList_maxe;
00073 std::vector<Treal> timeList;
00074 };
00075
00076
00077 template<class Treal, class Tworker>
00078 MatAccInvestigator<Treal, Tworker>::MatAccInvestigator(mat::SizesAndBlocks const & matrix_size_block_info_)
00079 : matrix_size_block_info(matrix_size_block_info_)
00080 {}
00081
00082
00083 template<class Treal, class Tworker>
00084 void MatAccInvestigator<Treal, Tworker>
00085 ::Scan(const Tworker & worker,
00086 Treal firstParam,
00087 Treal stepFactor,
00088 int nSteps)
00089 {
00090 nScanSteps = nSteps;
00091 baseThresh = firstParam;
00092 threshList.resize(nSteps);
00093 errorList_frob.resize(nSteps);
00094 errorList_eucl.resize(nSteps);
00095 errorList_maxe.resize(nSteps);
00096 timeList.resize(nSteps);
00097
00098
00099 symmMatrix accurateMatrix;
00100 accurateMatrix.resetSizesAndBlocks(matrix_size_block_info,
00101 matrix_size_block_info);
00102 symmMatrix otherMatrix;
00103 otherMatrix.resetSizesAndBlocks(matrix_size_block_info,
00104 matrix_size_block_info);
00105 symmMatrix errorMatrix;
00106 errorMatrix.resetSizesAndBlocks(matrix_size_block_info,
00107 matrix_size_block_info);
00108
00109
00110 worker.ComputeMatrix(firstParam, accurateMatrix);
00111
00112 Treal currParam = firstParam;
00113 for(int i = 0; i < nSteps; i++)
00114 {
00115 currParam *= stepFactor;
00116 time_t startTime, endTime;
00117 time(&startTime);
00118 worker.ComputeMatrix(currParam, otherMatrix);
00119 time(&endTime);
00120 timeList[i] = endTime - startTime;
00121 threshList[i] = currParam;
00122
00123 errorMatrix = otherMatrix;
00124 errorMatrix += (ergo_real)(-1) * accurateMatrix;
00125
00126
00127 errorList_frob[i] = errorMatrix.frob();
00128
00129 Treal euclAcc = 1e-11;
00130 errorList_eucl[i] = errorMatrix.eucl(euclAcc);
00131
00132 errorList_maxe[i] = compute_maxabs_sparse(errorMatrix);
00133 }
00134
00135 }
00136
00137
00138 template<class Treal, class Tworker>
00139 void MatAccInvestigator<Treal, Tworker>
00140 ::GetScanResult(Treal* threshList_,
00141 Treal* errorList_frob_,
00142 Treal* errorList_eucl_,
00143 Treal* errorList_maxe_,
00144 Treal* timeList_)
00145 {
00146 for(int i = 0; i < nScanSteps; i++)
00147 {
00148 threshList_[i] = threshList[i];
00149 errorList_frob_[i] = errorList_frob[i];
00150 errorList_eucl_[i] = errorList_eucl[i];
00151 errorList_maxe_[i] = errorList_maxe[i];
00152 timeList_ [i] = timeList [i];
00153 }
00154 }
00155
00156
00157
00158
00159 #endif