benchmark.cxx 9.87 KB
Newer Older
1
// Xerus - A General Purpose Tensor Library
2
// Copyright (C) 2014-2017 Benjamin Huber and Sebastian Wolf. 
3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
// 
// Xerus is free software: you can redistribute it and/or modify
// it under the terms of the GNU Affero General Public License as published
// by the Free Software Foundation, either version 3 of the License,
// or (at your option) any later version.
// 
// Xerus is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Affero General Public License for more details.
// 
// You should have received a copy of the GNU Affero General Public License
// along with Xerus. If not, see <http://www.gnu.org/licenses/>.
//
// For further information on Xerus visit https://libXerus.org 
// or contact us at contact@libXerus.org.

#include <fstream>
#include <string>
#include <iostream>
#include <stdio.h>
#include <vector>
#include <set>

Benjamin Huber's avatar
Benjamin Huber committed
27 28
#include <boost/filesystem.hpp>

29 30
#include "include/xerus.h"

31
std::mt19937_64 rnd = xerus::misc::randomEngine;
32 33 34 35
std::normal_distribution<double> normalDist(0,1);

using namespace xerus;

Benjamin Huber's avatar
Benjamin Huber committed
36 37 38 39 40
// ---------------------------------------------------------------------------------------------------------------------------------
// ------------------------------------------------- general settings --------------------------------------------------------------

const value_t HISTOGRAM_BASE_CONVERGENCE_RATES = 1.2;
const value_t HISTOGRAM_BASE_END_RESIDUAL = 1.7;
41
const size_t NUM_SOLVES_PER_RUN = 10;
Benjamin Huber's avatar
Benjamin Huber committed
42 43 44 45 46 47 48

// ---------------------------------------------------------------------------------------------------------------------------------
// ------------------------------------------------- benchmark problems ------------------------------------------------------------

using LeastSquaresSolver = std::pair<std::string, std::function<double(const TTOperator&, TTTensor&, const TTTensor&, PerformanceData&)>>;

struct LeastSquaresProblem {
49 50 51 52
	std::string name;
	std::vector<size_t> dimensions;
	std::vector<size_t> x_ranks;
	std::vector<size_t> b_ranks;
Benjamin Huber's avatar
Benjamin Huber committed
53 54 55
	std::vector<LeastSquaresSolver> solver;
	LeastSquaresProblem(const std::string &_name, const std::vector<LeastSquaresSolver> &_solver)
		: name(_name), solver(_solver) {};
56 57 58 59 60 61 62
	
	virtual TTOperator get_a() const {
		std::vector<size_t> dim(dimensions);
		dim.insert(dim.end(), dimensions.begin(), dimensions.end());
		return TTOperator::identity(dim);
	}
	virtual TTTensor get_x() const {
63
		TTTensor x = TTTensor::random(dimensions, x_ranks, normalDist);
64 65
		x /= frob_norm(x);
		return x;
66 67
	};
	virtual TTTensor get_b() const {
68
		TTTensor b = TTTensor::random(dimensions, b_ranks, normalDist);
69 70
		b /= frob_norm(b);
		return b;
71 72 73 74
	};
};

namespace ls {
Benjamin Huber's avatar
Benjamin Huber committed
75 76 77
	struct approximation : public LeastSquaresProblem {
		approximation(size_t _n, size_t _d, size_t _rankB, size_t _rankX, const std::vector<LeastSquaresSolver> &_solver)
			: LeastSquaresProblem("approximation", _solver)
78 79 80 81 82 83 84
		{
			dimensions = std::vector<size_t>(_d, _n);
			x_ranks = std::vector<size_t>(_d-1, _rankX);
			b_ranks = std::vector<size_t>(_d-1, _rankB);
		};
	};
	
Benjamin Huber's avatar
Benjamin Huber committed
85
	struct random : public LeastSquaresProblem {
86 87
		std::vector<size_t> a_ranks;
		
Benjamin Huber's avatar
Benjamin Huber committed
88 89
		random(size_t _n, size_t _d, size_t _rankA, size_t _rankB, size_t _rankX, const std::vector<LeastSquaresSolver> &_solver)
			: LeastSquaresProblem("random", _solver)
90 91 92 93 94 95 96 97 98 99
		{
			dimensions = std::vector<size_t>(_d, _n);
			a_ranks = std::vector<size_t>(_d-1, _rankA);
			x_ranks = std::vector<size_t>(_d-1, _rankX);
			b_ranks = std::vector<size_t>(_d-1, _rankB);
		};
		
		TTOperator get_a() const override {
			std::vector<size_t> dim(dimensions);
			dim.insert(dim.end(), dimensions.begin(), dimensions.end());
100
			TTOperator A = TTOperator::random(dim, a_ranks, normalDist);
101 102
			A /= frob_norm(A);
			return A;
103 104 105
		}
	};
	
Benjamin Huber's avatar
Benjamin Huber committed
106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
	struct symmetric_posdef_random : public LeastSquaresProblem {
		std::vector<size_t> a_ranks;
		
		symmetric_posdef_random(size_t _n, size_t _d, size_t _rankA, size_t _rankB, size_t _rankX, const std::vector<LeastSquaresSolver> &_solver)
			: LeastSquaresProblem("symmetric_posdef_random", _solver)
		{
			dimensions = std::vector<size_t>(_d, _n);
			a_ranks = std::vector<size_t>(_d-1, _rankA);
			x_ranks = std::vector<size_t>(_d-1, _rankX);
			b_ranks = std::vector<size_t>(_d-1, _rankB);
		};
		
		TTOperator get_a() const override {
			std::vector<size_t> dim(dimensions);
			dim.insert(dim.end(), dimensions.begin(), dimensions.end());
121
			TTOperator A = TTOperator::random(dim, a_ranks, normalDist);
Benjamin Huber's avatar
Benjamin Huber committed
122 123
			Index i,j,k;
			A(i,j) = A(i,k) * A(j,k);
124
			A /= frob_norm(A);
Benjamin Huber's avatar
Benjamin Huber committed
125 126 127
			return A;
		}
	};
128 129 130
}


Benjamin Huber's avatar
Benjamin Huber committed
131
std::vector<LeastSquaresSolver> leastSquaresAlgorithms{
132
	{"ALS", ALSVariant(1, 0, 1e-8, ALSVariant::lapack_solver, true)}, 
Benjamin Huber's avatar
Benjamin Huber committed
133
	{"CG", GeometricCGVariant(0, 0, 1e-8, false, SubmanifoldRetractionI, ProjectiveVectorTransport)}, 
134
	{"SteepestDescent_submanifold", SteepestDescentVariant(0, 1e-8, false, SubmanifoldRetractionII)},
Benjamin Huber's avatar
Benjamin Huber committed
135 136 137
	{"SteepestDescent_als", SteepestDescentVariant(0, 1e-8, false, ALSRetractionII)},
	{"SteepestDescent_hosvd", SteepestDescentVariant(0, 1e-8, false, HOSVDRetraction(3))}, //TODO
};
138

Benjamin Huber's avatar
Benjamin Huber committed
139
std::vector<LeastSquaresSolver> leastSquaresAlgorithmsSPD{
140
	{"ALS", ALSVariant(1, 0, 1e-8, ALSVariant::lapack_solver, true)}, 
Benjamin Huber's avatar
Benjamin Huber committed
141
	{"CG", GeometricCGVariant(0, 0, 1e-8, true, SubmanifoldRetractionI, ProjectiveVectorTransport)}, 
142
	{"SteepestDescent_submanifold", SteepestDescentVariant(0, 1e-8, true, SubmanifoldRetractionII)},
Benjamin Huber's avatar
Benjamin Huber committed
143 144
	{"SteepestDescent_als", SteepestDescentVariant(0, 1e-8, true, ALSRetractionII)},
	{"SteepestDescent_hosvd", SteepestDescentVariant(0, 1e-8, true, HOSVDRetraction(3))}, //TODO
145 146
};

Benjamin Huber's avatar
Benjamin Huber committed
147 148 149 150 151 152 153 154 155 156 157
struct Approximation_Variant {
	std::function<double(TTTensor&, const TTTensor&, PerformanceData&)> solver;
	Approximation_Variant(std::function<double(TTTensor&, const TTTensor&, PerformanceData&)> _solver)
		: solver(_solver) {}
	double operator()(const TTOperator&, TTTensor &_x, const TTTensor &_b, PerformanceData &_pd) const {
		return solver(_x, _b, _pd);
	}
};

std::vector<LeastSquaresProblem> leastSquaresProblems{
	ls::approximation(2, 10, 4, 2, std::vector<LeastSquaresSolver>{
158
		{"ALS", Approximation_Variant(ALSVariant(1, 0, 1e-8, ALSVariant::lapack_solver, true))}, 
Benjamin Huber's avatar
Benjamin Huber committed
159
		{"CG", Approximation_Variant(GeometricCGVariant(0, 0, 1e-8, true, SubmanifoldRetractionI, ProjectiveVectorTransport))}, 
160
		{"SteepestDescent_submanifold", Approximation_Variant(SteepestDescentVariant(0, 1e-8, true, SubmanifoldRetractionII))},
Benjamin Huber's avatar
Benjamin Huber committed
161 162 163 164 165 166 167 168 169 170 171 172 173 174 175
		{"SteepestDescent_als", Approximation_Variant(SteepestDescentVariant(0, 1e-8, true, ALSRetractionII))},
		{"SteepestDescent_hosvd", Approximation_Variant(SteepestDescentVariant(0, 1e-8, true, HOSVDRetraction(2)))}, //TODO
	}),
	ls::random(2, 10, 3, 3, 3, leastSquaresAlgorithms),
	ls::symmetric_posdef_random(2, 10, 2, 3, 3, leastSquaresAlgorithmsSPD)
};




// ---------------------------------------------------------------------------------------------------------------------------------
// ------------------------------------------------- benchmark routines ------------------------------------------------------------


std::string generate_profile_name() {
176
	std::string profileName;
177
#ifdef XERUS_TEST_COVERAGE
178 179
	static_assert(false, "test coverage checking nonsensical with benchmark run");
#endif
Benjamin Huber's avatar
Benjamin Huber committed
180 181 182 183 184 185 186
	
#ifdef XERUS_VERSION
	profileName += STRINGIFY(XERUS_VERSION);
#else
	profileName += "unknownVersion";
#endif
	
187
#ifdef LOW_OPTIMIZATION
Benjamin Huber's avatar
Benjamin Huber committed
188
	profileName += "_lowOpt";
189
#elif defined(HIGH_OPTIMIZATION)
Benjamin Huber's avatar
Benjamin Huber committed
190
	profileName += "_highOpt";
191
#elif defined(DANGEROUS_OPTIMIZATION)
Benjamin Huber's avatar
Benjamin Huber committed
192
	profileName += "_dangerousOpt";
193
#elif defined(RIDICULOUS_OPTIMIZATION)
Benjamin Huber's avatar
Benjamin Huber committed
194
	profileName += "_ridiculousOpt";
195
#else
Benjamin Huber's avatar
Benjamin Huber committed
196
	profileName += "_noOpt";
197 198 199 200 201
#endif
	
#ifdef USE_LTO
	profileName += "_lto";
#endif
202
#ifdef XERUS_DISABLE_RUNTIME_CHECKS
203 204
	profileName += "_noChecks";
#endif
205
#ifdef XERUS_REPLACE_ALLOCATOR
206 207
	profileName += "_replaceAlloc";
#endif
208
#ifdef XERUS_PERFORMANCE_ANALYSIS
209 210
	profileName += "_perfAnalysis";
#endif
Benjamin Huber's avatar
Benjamin Huber committed
211 212 213 214 215 216
	return profileName;
}


int main() {
	std::string profileName(generate_profile_name());
217
	LOG(benchmark, "running profile " << profileName);
218 219 220 221 222
	while (true) {
		for (const LeastSquaresProblem &prob : leastSquaresProblems) {
			std::vector<TTOperator> A;
			std::vector<TTTensor> X;
			std::vector<TTTensor> B;
Benjamin Huber's avatar
Benjamin Huber committed
223
			for (size_t i=0; i< NUM_SOLVES_PER_RUN; ++i) {
224 225 226
				A.emplace_back(prob.get_a());
				X.emplace_back(prob.get_x());
				B.emplace_back(prob.get_b());
Benjamin Huber's avatar
Benjamin Huber committed
227
			}
228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263
			for (const LeastSquaresSolver &solver : prob.solver) {
				LOG(benchmark, "solving " << prob.name << " with " << solver.first);
				PerformanceData perfData;
				misc::LogHistogram speedHist(HISTOGRAM_BASE_CONVERGENCE_RATES);
				misc::LogHistogram residualHist(HISTOGRAM_BASE_END_RESIDUAL);
				
				for (size_t i=0; i< NUM_SOLVES_PER_RUN; ++i) {
					perfData.reset();
					// solving the system
					TTTensor xCpy(X[i]);
					solver.second(A[i], xCpy, B[i], perfData);
					
					// generate histograms of this run
					speedHist += perfData.get_histogram(HISTOGRAM_BASE_CONVERGENCE_RATES, true);
					residualHist.add(perfData.data.back().residual);
				}
				
				// merge histograms with data on disk
				std::string fileName = std::string("benchmark/")+profileName+"/"+prob.name+"/"+solver.first;
				if (boost::filesystem::exists(fileName+"_speed.tsv")) {
					misc::LogHistogram speedHistFile = misc::LogHistogram::read_from_file(fileName+"_speed.tsv");
					speedHist += speedHistFile;
				} else {
					// will fail silently if the directories already exist
					boost::filesystem::create_directories(std::string("benchmark/")+profileName+"/"+prob.name);
				}
				speedHist.dump_to_file(fileName+"_speed.tsv");
				
				if (boost::filesystem::exists(fileName+"_residual.tsv")) {
					misc::LogHistogram residualHistFile = misc::LogHistogram::read_from_file(fileName+"_residual.tsv");
					residualHist += residualHistFile;
				} else {
					// will fail silently if the directories already exist
					boost::filesystem::create_directories(std::string("benchmark/")+profileName+"/"+prob.name);
				}
				residualHist.dump_to_file(fileName+"_residual.tsv");
Benjamin Huber's avatar
Benjamin Huber committed
264 265 266
			}
		}
	}
267
}