[sdpb] 121/233: Added parameters to SDPSolver object; got rid of SDPSolver.initialize; few other small fixes
Tobias Hansen
thansen at moszumanska.debian.org
Thu Mar 9 04:06:28 UTC 2017
This is an automated email from the git hooks/post-receive script.
thansen pushed a commit to branch master
in repository sdpb.
commit e7f620f53c788b0befca53cf85250e1e694d2484
Author: David Simmons-Duffin <davidsd at gmail.com>
Date: Fri Jan 9 12:49:04 2015 +1300
Added parameters to SDPSolver object; got rid of SDPSolver.initialize; few other small fixes
---
src/SDPSolver.cpp | 32 +++++++++++++++-----------------
src/SDPSolver.h | 37 +++++++++++++++++++++++++------------
src/main.cpp | 7 ++-----
3 files changed, 42 insertions(+), 34 deletions(-)
diff --git a/src/SDPSolver.cpp b/src/SDPSolver.cpp
index 242e0bf..88e2b4f 100644
--- a/src/SDPSolver.cpp
+++ b/src/SDPSolver.cpp
@@ -18,8 +18,9 @@ using boost::filesystem::path;
using boost::timer::nanosecond_type;
using std::cout;
-SDPSolver::SDPSolver(const SDP &sdp):
+SDPSolver::SDPSolver(const SDP &sdp, const SDPSolverParameters ¶meters):
sdp(sdp),
+ parameters(parameters),
x(sdp.primalObjective.size(), 0),
X(sdp.psdMatrixBlockDims()),
y(sdp.dualObjective.size(), 0),
@@ -53,6 +54,11 @@ SDPSolver::SDPSolver(const SDP &sdp):
eigenvaluesWorkspace.push_back(Vector(X.blocks[b].rows));
QRWorkspace.push_back(Vector(3*X.blocks[b].rows - 1));
}
+
+ // X = \Omega_p I
+ X.addDiagonal(parameters.initialMatrixScalePrimal);
+ // Y = \Omega_d I
+ Y.addDiagonal(parameters.initialMatrixScaleDual);
}
// result = b'^T a b', where b' = b \otimes 1
@@ -603,26 +609,21 @@ void SDPSolver::computeSearchDirection(const Real &beta,
dY *= -1;
}
-void SDPSolver::initialize(const SDPSolverParameters ¶meters) {
- fillVector(x, 0);
- X.setZero();
- X.addDiagonal(parameters.initialMatrixScalePrimal);
- fillVector(y, 0);
- Y.setZero();
- Y.addDiagonal(parameters.initialMatrixScaleDual);
-}
-
-SDPSolverTerminateReason SDPSolver::run(const SDPSolverParameters ¶meters,
- const path checkpointFile) {
+SDPSolverTerminateReason SDPSolver::run(const path checkpointFile) {
Real primalStepLength;
Real dualStepLength;
+ printSolverHeader();
+
for (int iteration = 1;; iteration++) {
if (timers["Last checkpoint"].elapsed().wall >= parameters.checkpointInterval * 1000000000LL)
saveCheckpoint(checkpointFile);
if (timers["Solver runtime"].elapsed().wall >= parameters.maxRuntime * 1000000000LL)
return MaxRuntimeExceeded;
+ primalObjective = sdp.objectiveConst + dotProduct(sdp.primalObjective, x);
+ dualObjective = sdp.objectiveConst + dotProduct(sdp.dualObjective, y);
+
timers["cholesky"].resume();
choleskyDecomposition(X, XCholesky);
choleskyDecomposition(Y, YCholesky);
@@ -637,14 +638,11 @@ SDPSolverTerminateReason SDPSolver::run(const SDPSolverParameters ¶meters,
timers["bilinear pairings"].stop();
// d_k = c_k - Tr(F_k Y) - (D y)_k
computeDualResidues(sdp, y, BilinearPairingsY, dualResidues);
+ dualError = maxAbsVector(dualResidues);
// PrimalResidues = sum_p F_p x_p - X - F_0 (F_0 is zero for now)
computePrimalResidues(sdp, x, X, PrimalResidues);
-
- primalError = PrimalResidues.maxAbs();
- dualError = maxAbsVector(dualResidues);
- primalObjective = sdp.objectiveConst + dotProduct(sdp.primalObjective, x);
- dualObjective = sdp.objectiveConst + dotProduct(sdp.dualObjective, y);
+ primalError = PrimalResidues.maxAbs();
const bool isPrimalFeasible = primalError < parameters.primalErrorThreshold;
const bool isDualFeasible = dualError < parameters.dualErrorThreshold;
diff --git a/src/SDPSolver.h b/src/SDPSolver.h
index cb18afa..3c9aef2 100644
--- a/src/SDPSolver.h
+++ b/src/SDPSolver.h
@@ -95,13 +95,17 @@ public:
// SDP to solve.
SDP sdp;
+ // parameters for initialization and iteration
+ SDPSolverParameters parameters;
+
/********************************************/
// Current point
// a Vector of length P = sdp.primalObjective.size()
Vector x;
- // a BlockDiagonalMatrix with structure given by sdp.psdMatrixBlockDims()
+ // a BlockDiagonalMatrix with block sizes given by
+ // sdp.psdMatrixBlockDims()
BlockDiagonalMatrix X;
// a Vector of length N = sdp.dualObjective.size()
@@ -136,7 +140,7 @@ public:
// PrimalResidues = \sum_p A_p x_p - X
//
BlockDiagonalMatrix PrimalResidues;
- Real primalError; // maxAbs(PrimalResidues)
+ Real primalError; // maxAbs(PrimalResidues)
// Discrepancy in the dual equality constraints, a Vector of length
// P, called 'd' in the manual:
@@ -149,18 +153,18 @@ public:
/********************************************/
// Intermediate computations.
- // the Cholesky decompositions of X and Y, each
- // BlockDiagonalMatrices with the same structure as X and Y
+ // the Cholesky decompositions of X and Y, each lower-triangular
+ // BlockDiagonalMatrices with the same block sizes as X and Y
BlockDiagonalMatrix XCholesky;
BlockDiagonalMatrix YCholesky;
// Z = X^{-1} (PrimalResidues Y - R), a BlockDiagonalMatrix with the
- // same structure as X and Y
+ // same block sizes as X and Y
BlockDiagonalMatrix Z;
// R = mu I - X Y for the predictor step
// R = mu I - X Y - dX dY for the corrector step
- // R is a BlockDiagonalMatrix with the same structure as X and Y.
+ // R has the same block sizes as X and Y.
BlockDiagonalMatrix R;
// Bilinear pairings needed for computing the Schur complement
@@ -194,11 +198,11 @@ public:
//
BlockDiagonalMatrix SchurBlocks;
- // SchurBlocksCholesky = L': the Cholesky decomposition of the
+ // SchurBlocksCholesky = L', the Cholesky decomposition of the
// stabilized Schur complement matrix S' = S + U U^T.
BlockDiagonalMatrix SchurBlocksCholesky;
- // SchurOffDiagonal = L'^{-1} FreeVarMatrix: needed in solving the
+ // SchurOffDiagonal = L'^{-1} FreeVarMatrix, needed in solving the
// Schur complement equation.
Matrix SchurOffDiagonal;
@@ -272,20 +276,29 @@ public:
// semidefiniteness.
BlockDiagonalMatrix StepMatrixWorkspace;
- //
+ // Needed during the computation of BilinearPairingsY, BilinearPairingsXInv
vector<Matrix> bilinearPairingsWorkspace;
+
+ // Needed during the step-length computation, where we must compute
+ // eigenvalues of a BlockDiagonalMatrix
vector<Vector> eigenvaluesWorkspace;
+
+ // Needed during the step-length computation, where we must compute
+ // the QR decomposition of a BlockDiagonalMatrix
vector<Vector> QRWorkspace;
- SDPSolver(const SDP &sdp);
+ /********************************************/
+ // Methods
+
+ // Create a new solver for a given SDP
+ SDPSolver(const SDP &sdp, const SDPSolverParameters ¶meters);
Real dualityGap() const {
return abs(primalObjective - dualObjective) /
max(Real(abs(primalObjective) + abs(dualObjective)), Real(1));
}
- void initialize(const SDPSolverParameters ¶meters);
- SDPSolverTerminateReason run(const SDPSolverParameters ¶meters, const path checkpointFile);
+ SDPSolverTerminateReason run(const path checkpointFile);
void initializeSchurComplementSolver(const BlockDiagonalMatrix &BilinearPairingsXInv,
const BlockDiagonalMatrix &BilinearPairingsY,
const Real &choleskyStabilizeThreshold);
diff --git a/src/main.cpp b/src/main.cpp
index 774d85a..81fde43 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -54,17 +54,14 @@ int solveSDP(const path &sdpFile,
cout << parameters << endl;
const SDP sdp = readBootstrapSDP(sdpFile);
- SDPSolver solver(sdp);
+ SDPSolver solver(sdp, parameters);
if (exists(checkpointFile))
solver.loadCheckpoint(checkpointFile);
- else
- solver.initialize(parameters);
timers["Solver runtime"].start();
timers["Last checkpoint"].start();
- printSolverHeader();
- SDPSolverTerminateReason reason = solver.run(parameters, checkpointFile);
+ SDPSolverTerminateReason reason = solver.run(checkpointFile);
timers["Solver runtime"].stop();
cout << "-----" << setfill('-') << setw(116) << std::left << reason << endl;
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/sdpb.git
More information about the debian-science-commits
mailing list