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 #ifndef ERGO_MAT_ACC_EXTRAPOLATE_HEADER
00029 #define ERGO_MAT_ACC_EXTRAPOLATE_HEADER
00030
00031
00032 #include <vector>
00033
00034
00035 #include "matrix_utilities.h"
00036
00037
00038
00039 template<class Treal, class Tworker>
00040 class MatAccInvestigator
00041 {
00042 public:
00043 explicit MatAccInvestigator(mat::SizesAndBlocks const & matrix_size_block_info_);
00044 void Scan(const Tworker & worker,
00045 Treal firstParam,
00046 Treal stepFactor,
00047 int nSteps);
00048 void GetScanResult(Treal* threshList_,
00049 Treal* errorList_frob_,
00050 Treal* errorList_eucl_,
00051 Treal* errorList_maxe_,
00052 Treal* timeList_);
00053 private:
00054 mat::SizesAndBlocks matrix_size_block_info;
00055 int nScanSteps;
00056 Treal baseThresh;
00057 std::vector<Treal> threshList;
00058 std::vector<Treal> errorList_frob;
00059 std::vector<Treal> errorList_eucl;
00060 std::vector<Treal> errorList_maxe;
00061 std::vector<Treal> timeList;
00062 };
00063
00064
00065 template<class Treal, class Tworker>
00066 MatAccInvestigator<Treal, Tworker>::MatAccInvestigator(mat::SizesAndBlocks const & matrix_size_block_info_)
00067 : matrix_size_block_info(matrix_size_block_info_)
00068 {}
00069
00070
00071 template<class Treal, class Tworker>
00072 void MatAccInvestigator<Treal, Tworker>
00073 ::Scan(const Tworker & worker,
00074 Treal firstParam,
00075 Treal stepFactor,
00076 int nSteps)
00077 {
00078 nScanSteps = nSteps;
00079 baseThresh = firstParam;
00080 threshList.resize(nSteps);
00081 errorList_frob.resize(nSteps);
00082 errorList_eucl.resize(nSteps);
00083 errorList_maxe.resize(nSteps);
00084 timeList.resize(nSteps);
00085
00086
00087 symmMatrix accurateMatrix;
00088 accurateMatrix.resetSizesAndBlocks(matrix_size_block_info,
00089 matrix_size_block_info);
00090 symmMatrix otherMatrix;
00091 otherMatrix.resetSizesAndBlocks(matrix_size_block_info,
00092 matrix_size_block_info);
00093 symmMatrix errorMatrix;
00094 errorMatrix.resetSizesAndBlocks(matrix_size_block_info,
00095 matrix_size_block_info);
00096
00097
00098 worker.ComputeMatrix(firstParam, accurateMatrix);
00099
00100 Treal currParam = firstParam;
00101 for(int i = 0; i < nSteps; i++)
00102 {
00103 currParam *= stepFactor;
00104 time_t startTime, endTime;
00105 time(&startTime);
00106 worker.ComputeMatrix(currParam, otherMatrix);
00107 time(&endTime);
00108 timeList[i] = endTime - startTime;
00109 threshList[i] = currParam;
00110
00111 errorMatrix = otherMatrix;
00112 errorMatrix += (ergo_real)(-1) * accurateMatrix;
00113
00114
00115 errorList_frob[i] = errorMatrix.frob();
00116
00117 Treal euclAcc = 1e-11;
00118 errorList_eucl[i] = errorMatrix.eucl(euclAcc);
00119
00120 errorList_maxe[i] = compute_maxabs_sparse(errorMatrix);
00121 }
00122
00123 }
00124
00125
00126 template<class Treal, class Tworker>
00127 void MatAccInvestigator<Treal, Tworker>
00128 ::GetScanResult(Treal* threshList_,
00129 Treal* errorList_frob_,
00130 Treal* errorList_eucl_,
00131 Treal* errorList_maxe_,
00132 Treal* timeList_)
00133 {
00134 for(int i = 0; i < nScanSteps; i++)
00135 {
00136 threshList_[i] = threshList[i];
00137 errorList_frob_[i] = errorList_frob[i];
00138 errorList_eucl_[i] = errorList_eucl[i];
00139 errorList_maxe_[i] = errorList_maxe[i];
00140 timeList_ [i] = timeList [i];
00141 }
00142 }
00143
00144
00145
00146
00147 #endif