From b0a7826cb91d2533faa49898e8832ee0588a7f1a Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Tue, 15 Jan 2019 18:55:09 +0100
Subject: [PATCH] Dynare++ / sylvester equation solver: first batch of
modernizations
---
dynare++/sylv/cc/BlockDiagonal.hh | 6 ++-
dynare++/sylv/cc/GeneralMatrix.hh | 10 ++--
dynare++/sylv/cc/GeneralSylvester.cc | 15 ++----
dynare++/sylv/cc/GeneralSylvester.hh | 11 +++--
dynare++/sylv/cc/IterativeSylvester.cc | 7 +--
dynare++/sylv/cc/KronVector.cc | 41 +++++----------
dynare++/sylv/cc/KronVector.hh | 4 +-
dynare++/sylv/cc/QuasiTriangular.cc | 7 +--
dynare++/sylv/cc/QuasiTriangular.hh | 38 ++++++--------
dynare++/sylv/cc/QuasiTriangularZero.hh | 18 ++++---
dynare++/sylv/cc/SchurDecomp.cc | 28 +++--------
dynare++/sylv/cc/SchurDecomp.hh | 15 +++---
dynare++/sylv/cc/SimilarityDecomp.cc | 13 ++---
dynare++/sylv/cc/SimilarityDecomp.hh | 11 +++--
dynare++/sylv/cc/SylvException.cc | 17 +------
dynare++/sylv/cc/SylvException.hh | 7 +--
dynare++/sylv/cc/SylvMatrix.hh | 6 +--
dynare++/sylv/cc/SylvMemory.cc | 4 --
dynare++/sylv/cc/SylvMemory.hh | 2 +-
dynare++/sylv/cc/TriangularSylvester.cc | 18 ++-----
dynare++/sylv/cc/TriangularSylvester.hh | 8 +--
dynare++/sylv/cc/Vector.cc | 66 +++++++++++--------------
dynare++/sylv/cc/Vector.hh | 25 ++++------
23 files changed, 138 insertions(+), 239 deletions(-)
diff --git a/dynare++/sylv/cc/BlockDiagonal.hh b/dynare++/sylv/cc/BlockDiagonal.hh
index 3dedcaef4c..3b108bbd32 100644
--- a/dynare++/sylv/cc/BlockDiagonal.hh
+++ b/dynare++/sylv/cc/BlockDiagonal.hh
@@ -5,6 +5,8 @@
#ifndef BLOCK_DIAGONAL_H
#define BLOCK_DIAGONAL_H
+#include <memory>
+
#include "QuasiTriangular.hh"
class BlockDiagonal : public QuasiTriangular
@@ -39,10 +41,10 @@ public:
col_iter col_begin(const DiagonalBlock &b) override;
const_row_iter row_end(const DiagonalBlock &b) const override;
row_iter row_end(const DiagonalBlock &b) override;
- QuasiTriangular *
+ std::unique_ptr<QuasiTriangular>
clone() const override
{
- return new BlockDiagonal(*this);
+ return std::make_unique<BlockDiagonal>(*this);
}
private:
void setZerosToRU(diag_iter edge);
diff --git a/dynare++/sylv/cc/GeneralMatrix.hh b/dynare++/sylv/cc/GeneralMatrix.hh
index 49053aa709..e7ff5611db 100644
--- a/dynare++/sylv/cc/GeneralMatrix.hh
+++ b/dynare++/sylv/cc/GeneralMatrix.hh
@@ -27,8 +27,7 @@ public:
ConstGeneralMatrix(const GeneralMatrix &m);
ConstGeneralMatrix(const GeneralMatrix &m, int i, int j, int nrows, int ncols);
ConstGeneralMatrix(const ConstGeneralMatrix &m, int i, int j, int nrows, int ncols);
- virtual ~ConstGeneralMatrix()
- = default;
+ virtual ~ConstGeneralMatrix() = default;
const double &
get(int i, int j) const
@@ -146,11 +145,8 @@ public:
GeneralMatrix(const GeneralMatrix &a, const char *dum1,
const GeneralMatrix &b, const char *dum2);
- virtual
- ~GeneralMatrix();
- GeneralMatrix &
- operator=(const GeneralMatrix &m)
- = default;
+ virtual ~GeneralMatrix();
+ GeneralMatrix &operator=(const GeneralMatrix &m) = default;
const double &
get(int i, int j) const
diff --git a/dynare++/sylv/cc/GeneralSylvester.cc b/dynare++/sylv/cc/GeneralSylvester.cc
index fbe3e51e1b..bc5c031e1f 100644
--- a/dynare++/sylv/cc/GeneralSylvester.cc
+++ b/dynare++/sylv/cc/GeneralSylvester.cc
@@ -67,14 +67,14 @@ GeneralSylvester::init()
a.multInvLeft2(ainvb, d, rcond1, rcondinf);
pars.rcondA1 = rcond1;
pars.rcondAI = rcondinf;
- bdecomp = new SchurDecompZero(ainvb);
- cdecomp = new SimilarityDecomp(c.getData().base(), c.numRows(), *(pars.bs_norm));
+ bdecomp = std::make_unique<SchurDecompZero>(ainvb);
+ cdecomp = std::make_unique<SimilarityDecomp>(c.getData().base(), c.numRows(), *(pars.bs_norm));
cdecomp->check(pars, c);
cdecomp->infoToPars(pars);
if (*(pars.method) == SylvParams::recurse)
- sylv = new TriangularSylvester(*bdecomp, *cdecomp);
+ sylv = std::make_unique<TriangularSylvester>(*bdecomp, *cdecomp);
else
- sylv = new IterativeSylvester(*bdecomp, *cdecomp);
+ sylv = std::make_unique<IterativeSylvester>(*bdecomp, *cdecomp);
}
void
@@ -128,10 +128,3 @@ GeneralSylvester::check(const double *ds)
mem_driver.setStackMode(false);
}
-
-GeneralSylvester::~GeneralSylvester()
-{
- delete bdecomp;
- delete cdecomp;
- delete sylv;
-}
diff --git a/dynare++/sylv/cc/GeneralSylvester.hh b/dynare++/sylv/cc/GeneralSylvester.hh
index 107e70e7da..489169ae5d 100644
--- a/dynare++/sylv/cc/GeneralSylvester.hh
+++ b/dynare++/sylv/cc/GeneralSylvester.hh
@@ -10,6 +10,8 @@
#include "SimilarityDecomp.hh"
#include "SylvesterSolver.hh"
+#include <memory>
+
class GeneralSylvester
{
SylvParams pars;
@@ -20,9 +22,9 @@ class GeneralSylvester
const SqSylvMatrix c;
SylvMatrix d;
bool solved;
- SchurDecompZero *bdecomp;
- SimilarityDecomp *cdecomp;
- SylvesterSolver *sylv;
+ std::unique_ptr<SchurDecompZero> bdecomp;
+ std::unique_ptr<SimilarityDecomp> cdecomp;
+ std::unique_ptr<SylvesterSolver> sylv;
public:
/* construct with my copy of d*/
GeneralSylvester(int ord, int n, int m, int zero_cols,
@@ -42,8 +44,7 @@ public:
const double *da, const double *db,
const double *dc, double *dd,
const SylvParams &ps);
- virtual
- ~GeneralSylvester();
+ virtual ~GeneralSylvester() = default;
int
getM() const
{
diff --git a/dynare++/sylv/cc/IterativeSylvester.cc b/dynare++/sylv/cc/IterativeSylvester.cc
index e653a3a446..874829d0d9 100644
--- a/dynare++/sylv/cc/IterativeSylvester.cc
+++ b/dynare++/sylv/cc/IterativeSylvester.cc
@@ -13,8 +13,8 @@ IterativeSylvester::solve(SylvParams &pars, KronVector &x) const
double max_norm = *(pars.convergence_tol);
double norm = performFirstStep(x);
- QuasiTriangular *kpow = matrixK->clone();
- QuasiTriangular *fpow = matrixF->clone();
+ auto kpow = matrixK->clone();
+ auto fpow = matrixF->clone();
while (steps < max_steps &&norm > max_norm)
{
kpow->multRight(SqSylvMatrix(*kpow)); // be careful to make copy
@@ -23,9 +23,6 @@ IterativeSylvester::solve(SylvParams &pars, KronVector &x) const
steps++;
}
- delete fpow;
- delete kpow;
-
pars.converged = (norm <= max_norm);
pars.iter_last_norm = norm;
pars.num_iter = steps;
diff --git a/dynare++/sylv/cc/KronVector.cc b/dynare++/sylv/cc/KronVector.cc
index 7897d571e8..494a8f13b1 100644
--- a/dynare++/sylv/cc/KronVector.cc
+++ b/dynare++/sylv/cc/KronVector.cc
@@ -10,9 +10,7 @@ power(int m, int depth)
{
int p = 1;
for (int i = 0; i < depth; i++)
- {
- p *= m;
- }
+ p *= m;
return p;
}
@@ -24,11 +22,8 @@ KronVector::KronVector(int mm, int nn, int dp)
KronVector::KronVector(Vector &v, int mm, int nn, int dp)
: Vector(v), m(mm), n(nn), depth(dp)
{
- len = power(m, depth)*n;
- if (v.length() != length())
- {
- throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
- }
+ if (length() != power(m, depth)*n)
+ throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
}
KronVector::KronVector(KronVector &v, int i)
@@ -36,9 +31,7 @@ KronVector::KronVector(KronVector &v, int i)
depth(v.depth-1)
{
if (depth < 0)
- {
- throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
- }
+ throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
}
KronVector::KronVector(const ConstKronVector &v)
@@ -61,9 +54,7 @@ const KronVector &
KronVector::operator=(const Vector &v)
{
if (length() != v.length())
- {
- throw SYLV_MES_EXCEPTION("Wrong lengths for vector operator =.");
- }
+ throw SYLV_MES_EXCEPTION("Wrong lengths for vector operator =.");
Vector::operator=(v);
return *this;
}
@@ -82,21 +73,15 @@ ConstKronVector::ConstKronVector(const ConstKronVector &v)
ConstKronVector::ConstKronVector(const Vector &v, int mm, int nn, int dp)
: ConstVector(v), m(mm), n(nn), depth(dp)
{
- len = power(m, depth)*n;
- if (v.length() != length())
- {
- throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
- }
+ if (length() != power(m, depth)*n)
+ throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
}
ConstKronVector::ConstKronVector(const ConstVector &v, int mm, int nn, int dp)
: ConstVector(v), m(mm), n(nn), depth(dp)
{
- len = power(m, depth)*n;
- if (v.length() != length())
- {
- throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
- }
+ if (length() != power(m, depth)*n)
+ throw SYLV_MES_EXCEPTION("Bad conversion KronVector from Vector.");
}
ConstKronVector::ConstKronVector(const KronVector &v, int i)
@@ -105,9 +90,7 @@ ConstKronVector::ConstKronVector(const KronVector &v, int i)
m(v.getM()), n(v.getN()), depth(v.getDepth()-1)
{
if (depth < 0)
- {
- throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
- }
+ throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
}
ConstKronVector::ConstKronVector(const ConstKronVector &v, int i)
@@ -115,7 +98,5 @@ ConstKronVector::ConstKronVector(const ConstKronVector &v, int i)
m(v.getM()), n(v.getN()), depth(v.getDepth()-1)
{
if (depth < 0)
- {
- throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
- }
+ throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
}
diff --git a/dynare++/sylv/cc/KronVector.hh b/dynare++/sylv/cc/KronVector.hh
index c6e7e1bc94..8bddd29f79 100644
--- a/dynare++/sylv/cc/KronVector.hh
+++ b/dynare++/sylv/cc/KronVector.hh
@@ -23,9 +23,7 @@ public:
KronVector(Vector &v, int mm, int nn, int dp); // conversion
KronVector(KronVector &, int i); // picks i-th subvector
KronVector(const ConstKronVector &v); // new instance and copy
- KronVector &
- operator=(const KronVector &v)
- = default;
+ KronVector &operator=(const KronVector &v) = default;
const KronVector &operator=(const ConstKronVector &v);
const KronVector &operator=(const Vector &v);
int
diff --git a/dynare++/sylv/cc/QuasiTriangular.cc b/dynare++/sylv/cc/QuasiTriangular.cc
index 2b0cdeac3e..348403a72f 100644
--- a/dynare++/sylv/cc/QuasiTriangular.cc
+++ b/dynare++/sylv/cc/QuasiTriangular.cc
@@ -431,9 +431,6 @@ QuasiTriangular::QuasiTriangular(const double *d, int d_size)
{
}
-QuasiTriangular::~QuasiTriangular()
-= default;
-
QuasiTriangular::QuasiTriangular(int p, const QuasiTriangular &t)
: SqSylvMatrix(t.numRows()), diagonal(getData().base(), t.diagonal)
{
@@ -477,9 +474,7 @@ QuasiTriangular::QuasiTriangular(const SchurDecompZero &decomp)
}
}
// construct diagonal
- auto *const d = new Diagonal(getData().base(), decomp.getDim());
- diagonal = *d;
- delete d;
+ diagonal = Diagonal{getData().base(), decomp.getDim()};
}
void
diff --git a/dynare++/sylv/cc/QuasiTriangular.hh b/dynare++/sylv/cc/QuasiTriangular.hh
index edbd38d401..d4f39847a3 100644
--- a/dynare++/sylv/cc/QuasiTriangular.hh
+++ b/dynare++/sylv/cc/QuasiTriangular.hh
@@ -10,6 +10,7 @@
#include "SylvMatrix.hh"
#include <list>
+#include <memory>
using namespace std;
@@ -21,8 +22,7 @@ private:
double *a1;
double *a2;
public:
- DiagPair()
- = default;
+ DiagPair() = default;
DiagPair(double *aa1, double *aa2)
{
a1 = aa1; a2 = aa2;
@@ -31,9 +31,7 @@ public:
{
a1 = p.a1; a2 = p.a2;
}
- DiagPair &
- operator=(const DiagPair &p)
- = default;
+ DiagPair &operator=(const DiagPair &p) = default;
DiagPair &
operator=(double v)
{
@@ -70,8 +68,7 @@ private:
}
public:
- DiagonalBlock()
- = default;
+ DiagonalBlock() = default;
DiagonalBlock(int jb, bool r, double *a1, double *a2,
double *b1, double *b2)
: alpha(a1, a2)
@@ -205,8 +202,7 @@ private:
int num_real{0};
void copy(const Diagonal &);
public:
- Diagonal()
- = default;
+ Diagonal() = default;
Diagonal(double *data, int d_size);
Diagonal(double *data, const Diagonal &d);
Diagonal(const Diagonal &d)
@@ -218,8 +214,7 @@ public:
{
copy(d); return *this;
}
- virtual ~Diagonal()
- = default;
+ virtual ~Diagonal() = default;
int
getNumComplex() const
@@ -290,8 +285,7 @@ public:
{
ptr = base; d_size = ds; real = r;
}
- virtual ~_matrix_iter()
- = default;
+ virtual ~_matrix_iter() = default;
const _Self &
operator=(const _Self &it)
{
@@ -414,7 +408,7 @@ public:
QuasiTriangular(const SchurDecompZero &decomp);
QuasiTriangular(const QuasiTriangular &t);
- ~QuasiTriangular() override;
+ ~QuasiTriangular() override = default;
const Diagonal &
getDiagonal() const
{
@@ -484,25 +478,25 @@ public:
virtual row_iter row_end(const DiagonalBlock &b);
/* clone */
- virtual QuasiTriangular *
+ virtual std::unique_ptr<QuasiTriangular>
clone() const
{
- return new QuasiTriangular(*this);
+ return std::make_unique<QuasiTriangular>(*this);
}
- virtual QuasiTriangular *
+ virtual std::unique_ptr<QuasiTriangular>
clone(int p, const QuasiTriangular &t) const
{
- return new QuasiTriangular(p, t);
+ return std::make_unique<QuasiTriangular>(p, t);
}
- virtual QuasiTriangular *
+ virtual std::unique_ptr<QuasiTriangular>
clone(double r) const
{
- return new QuasiTriangular(r, *this);
+ return std::make_unique<QuasiTriangular>(r, *this);
}
- virtual QuasiTriangular *
+ virtual std::unique_ptr<QuasiTriangular>
clone(double r, double rr, const QuasiTriangular &tt) const
{
- return new QuasiTriangular(r, *this, rr, tt);
+ return std::make_unique<QuasiTriangular>(r, *this, rr, tt);
}
protected:
void setMatrix(double r, const QuasiTriangular &t);
diff --git a/dynare++/sylv/cc/QuasiTriangularZero.hh b/dynare++/sylv/cc/QuasiTriangularZero.hh
index ec408aca17..f4bb4a9e5a 100644
--- a/dynare++/sylv/cc/QuasiTriangularZero.hh
+++ b/dynare++/sylv/cc/QuasiTriangularZero.hh
@@ -8,6 +8,8 @@
#include "QuasiTriangular.hh"
#include "GeneralMatrix.hh"
+#include <memory>
+
class QuasiTriangularZero : public QuasiTriangular
{
int nz; // number of zero columns
@@ -31,25 +33,25 @@ public:
void multKronTrans(KronVector &x) const override;
void multLeftOther(GeneralMatrix &a) const override;
/* clone */
- QuasiTriangular *
+ std::unique_ptr<QuasiTriangular>
clone() const override
{
- return new QuasiTriangularZero(*this);
+ return std::make_unique<QuasiTriangularZero>(*this);
}
- QuasiTriangular *
+ std::unique_ptr<QuasiTriangular>
clone(int p, const QuasiTriangular &t) const override
{
- return new QuasiTriangularZero(p, (const QuasiTriangularZero &) t);
+ return std::make_unique<QuasiTriangularZero>(p, (const QuasiTriangularZero &) t);
}
- QuasiTriangular *
+ std::unique_ptr<QuasiTriangular>
clone(double r) const override
{
- return new QuasiTriangularZero(r, *this);
+ return std::make_unique<QuasiTriangularZero>(r, *this);
}
- QuasiTriangular *
+ std::unique_ptr<QuasiTriangular>
clone(double r, double rr, const QuasiTriangular &tt) const override
{
- return new QuasiTriangularZero(r, *this, rr, (const QuasiTriangularZero &) tt);
+ return std::make_unique<QuasiTriangularZero>(r, *this, rr, (const QuasiTriangularZero &) tt);
}
void print() const override;
};
diff --git a/dynare++/sylv/cc/SchurDecomp.cc b/dynare++/sylv/cc/SchurDecomp.cc
index 0612e176d8..1989dabdcb 100644
--- a/dynare++/sylv/cc/SchurDecomp.cc
+++ b/dynare++/sylv/cc/SchurDecomp.cc
@@ -7,10 +7,9 @@
#include <dynlapack.h>
SchurDecomp::SchurDecomp(const SqSylvMatrix &m)
- : q_destroy(true), t_destroy(true)
+ : q(m.numRows())
{
lapack_int rows = m.numRows();
- q = new SqSylvMatrix(rows);
SqSylvMatrix auxt(m);
lapack_int sdim;
auto *const wr = new double[rows];
@@ -19,36 +18,25 @@ SchurDecomp::SchurDecomp(const SqSylvMatrix &m)
auto *const work = new double[lwork];
lapack_int info;
dgees("V", "N", nullptr, &rows, auxt.base(), &rows, &sdim,
- wr, wi, q->base(), &rows,
+ wr, wi, q.base(), &rows,
work, &lwork, nullptr, &info);
delete [] work;
delete [] wi;
delete [] wr;
- t = new QuasiTriangular(auxt.base(), rows);
+ t_storage = std::make_unique<QuasiTriangular>(auxt.base(), rows);
+ t = t_storage.get();
}
SchurDecomp::SchurDecomp(const QuasiTriangular &tr)
- : q_destroy(true), t_destroy(true)
+ : q(tr.numRows()), t_storage{std::make_unique<QuasiTriangular>(tr)}, t{t_storage.get()}
{
- q = new SqSylvMatrix(tr.numRows());
- q->setUnit();
- t = new QuasiTriangular(tr);
+ q.setUnit();
}
SchurDecomp::SchurDecomp(QuasiTriangular &tr)
- : q_destroy(true), t_destroy(false)
+ : q(tr.numRows()), t{&tr}
{
- q = new SqSylvMatrix(tr.numRows());
- q->setUnit();
- t = &tr;
-}
-
-SchurDecomp::~SchurDecomp()
-{
- if (t_destroy)
- delete t;
- if (q_destroy)
- delete q;
+ q.setUnit();
}
int
diff --git a/dynare++/sylv/cc/SchurDecomp.hh b/dynare++/sylv/cc/SchurDecomp.hh
index 7a6a31b2e0..f93e5aaf6e 100644
--- a/dynare++/sylv/cc/SchurDecomp.hh
+++ b/dynare++/sylv/cc/SchurDecomp.hh
@@ -8,12 +8,14 @@
#include "SylvMatrix.hh"
#include "QuasiTriangular.hh"
+#include <memory>
+
class QuasiTriangular;
class SchurDecomp
{
- bool q_destroy;
- SqSylvMatrix *q;
- bool t_destroy;
+ SqSylvMatrix q;
+ // Stores t if is owned
+ std::unique_ptr<QuasiTriangular> t_storage;
QuasiTriangular *t;
public:
SchurDecomp(const SqSylvMatrix &m);
@@ -22,7 +24,7 @@ public:
const SqSylvMatrix &
getQ() const
{
- return *q;
+ return q;
}
const QuasiTriangular &
getT() const
@@ -32,7 +34,7 @@ public:
SqSylvMatrix &
getQ()
{
- return *q;
+ return q;
}
QuasiTriangular &
getT()
@@ -40,8 +42,7 @@ public:
return *t;
}
virtual int getDim() const;
- virtual
- ~SchurDecomp();
+ virtual ~SchurDecomp() = default;
};
class SchurDecompZero : public SchurDecomp
diff --git a/dynare++/sylv/cc/SimilarityDecomp.cc b/dynare++/sylv/cc/SimilarityDecomp.cc
index 27ac760245..3de3445001 100644
--- a/dynare++/sylv/cc/SimilarityDecomp.cc
+++ b/dynare++/sylv/cc/SimilarityDecomp.cc
@@ -14,22 +14,15 @@
SimilarityDecomp::SimilarityDecomp(const double *d, int d_size, double log10norm)
{
SchurDecomp sd(SqSylvMatrix(d, d_size));
- q = new SqSylvMatrix(sd.getQ());
- b = new BlockDiagonal(sd.getT());
- invq = new SqSylvMatrix(d_size);
+ q = std::make_unique<SqSylvMatrix>(sd.getQ());
+ b = std::make_unique<BlockDiagonal>(sd.getT());
+ invq = std::make_unique<SqSylvMatrix>(d_size);
invq->setUnit();
invq->multLeftTrans(sd.getQ());
double norm = pow(10.0, log10norm);
diagonalize(norm);
}
-SimilarityDecomp::~SimilarityDecomp()
-{
- delete invq;
- delete b;
- delete q;
-}
-
void
SimilarityDecomp::getXDim(diag_iter start, diag_iter end,
int &rows, int &cols) const
diff --git a/dynare++/sylv/cc/SimilarityDecomp.hh b/dynare++/sylv/cc/SimilarityDecomp.hh
index f546143b27..b831c821d2 100644
--- a/dynare++/sylv/cc/SimilarityDecomp.hh
+++ b/dynare++/sylv/cc/SimilarityDecomp.hh
@@ -9,16 +9,17 @@
#include "BlockDiagonal.hh"
#include "SylvParams.hh"
+#include <memory>
+
class SimilarityDecomp
{
- SqSylvMatrix *q;
- BlockDiagonal *b;
- SqSylvMatrix *invq;
+ std::unique_ptr<SqSylvMatrix> q;
+ std::unique_ptr<BlockDiagonal> b;
+ std::unique_ptr<SqSylvMatrix> invq;
using diag_iter = BlockDiagonal::diag_iter;
public:
SimilarityDecomp(const double *d, int d_size, double log10norm = 3.0);
- virtual
- ~SimilarityDecomp();
+ virtual ~SimilarityDecomp() = default;
const SqSylvMatrix &
getQ() const
{
diff --git a/dynare++/sylv/cc/SylvException.cc b/dynare++/sylv/cc/SylvException.cc
index bc8ae17a86..1065c283fc 100644
--- a/dynare++/sylv/cc/SylvException.cc
+++ b/dynare++/sylv/cc/SylvException.cc
@@ -7,19 +7,10 @@
#include <cstring>
#include <cstdio>
-SylvException::SylvException(const char *f, int l, const SylvException *s)
+SylvException::SylvException(const char *f, int l)
{
strcpy(file, f);
line = l;
- source = s;
-}
-
-SylvException::~SylvException()
-{
- if (source != nullptr)
- {
- delete source;
- }
}
void
@@ -35,10 +26,6 @@ int
SylvException::printMessage(char *str, int maxlen) const
{
int remain = maxlen;
- if (source != nullptr)
- {
- remain = source->printMessage(str, maxlen);
- }
char aux[100];
sprintf(aux, "From %s:%d\n", file, line);
int newremain = remain - strlen(aux);
@@ -53,7 +40,7 @@ SylvException::printMessage(char *str, int maxlen) const
SylvExceptionMessage::SylvExceptionMessage(const char *f, int i,
const char *mes)
- : SylvException(f, i, nullptr)
+ : SylvException(f, i)
{
strcpy(message, mes);
}
diff --git a/dynare++/sylv/cc/SylvException.hh b/dynare++/sylv/cc/SylvException.hh
index e16cd86f6e..786a2a69a0 100644
--- a/dynare++/sylv/cc/SylvException.hh
+++ b/dynare++/sylv/cc/SylvException.hh
@@ -12,11 +12,9 @@ class SylvException : public MallocAllocator
protected:
char file[50];
int line;
- const SylvException *source;
public:
- SylvException(const char *f, int l, const SylvException *s);
- virtual
- ~SylvException();
+ SylvException(const char *f, int l);
+ virtual ~SylvException() = default;
virtual int printMessage(char *str, int maxlen) const;
void printMessage() const;
};
@@ -30,7 +28,6 @@ public:
};
// define macros:
-#define SYLV_EXCEPTION(exc) (SylvException(__FILE__, __LINE__, exc))
#define SYLV_MES_EXCEPTION(mes) (SylvExceptionMessage(__FILE__, __LINE__, mes))
#endif /* SYLV_EXCEPTION_H */
diff --git a/dynare++/sylv/cc/SylvMatrix.hh b/dynare++/sylv/cc/SylvMatrix.hh
index 4478d3a9a7..9c97913956 100644
--- a/dynare++/sylv/cc/SylvMatrix.hh
+++ b/dynare++/sylv/cc/SylvMatrix.hh
@@ -74,8 +74,7 @@ public:
SqSylvMatrix(double *d, int m) : SylvMatrix(d, m, m)
{
}
- SqSylvMatrix(const SqSylvMatrix &m)
- = default;
+ SqSylvMatrix(const SqSylvMatrix &m) = default;
SqSylvMatrix(const GeneralMatrix &m, int i, int j, int nrows)
: SylvMatrix(m, i, j, nrows, nrows)
{
@@ -88,7 +87,8 @@ public:
const SqSylvMatrix &
operator=(const SqSylvMatrix &m)
{
- GeneralMatrix::operator=(m); return *this;
+ GeneralMatrix::operator=(m);
+ return *this;
}
/* x = (this \otimes this..\otimes this)*d */
void multVecKron(KronVector &x, const KronVector &d) const;
diff --git a/dynare++/sylv/cc/SylvMemory.cc b/dynare++/sylv/cc/SylvMemory.cc
index 145a3339c7..a156637704 100644
--- a/dynare++/sylv/cc/SylvMemory.cc
+++ b/dynare++/sylv/cc/SylvMemory.cc
@@ -20,10 +20,6 @@
SylvMemoryPool memory_pool;
-SylvMemoryPool::SylvMemoryPool()
-
-= default;
-
void
SylvMemoryPool::init(size_t size)
{
diff --git a/dynare++/sylv/cc/SylvMemory.hh b/dynare++/sylv/cc/SylvMemory.hh
index ca937781b8..cd7e7f4a7e 100644
--- a/dynare++/sylv/cc/SylvMemory.hh
+++ b/dynare++/sylv/cc/SylvMemory.hh
@@ -34,7 +34,7 @@ class SylvMemoryPool
size_t allocated{0};
bool stack_mode{false};
public:
- SylvMemoryPool();
+ SylvMemoryPool() = default;
SylvMemoryPool(const SylvMemoryPool &) = delete;
const SylvMemoryPool &operator=(const SylvMemoryPool &) = delete;
~SylvMemoryPool();
diff --git a/dynare++/sylv/cc/TriangularSylvester.cc b/dynare++/sylv/cc/TriangularSylvester.cc
index 621f56ed3c..e637339579 100644
--- a/dynare++/sylv/cc/TriangularSylvester.cc
+++ b/dynare++/sylv/cc/TriangularSylvester.cc
@@ -17,7 +17,7 @@ TriangularSylvester::TriangularSylvester(const QuasiTriangular &k,
const QuasiTriangular &f)
: SylvesterSolver(k, f),
matrixKK(matrixK->clone(2, *matrixK)),
- matrixFF(new QuasiTriangular(2, *matrixF))
+ matrixFF(std::make_unique<QuasiTriangular>(2, *matrixF))
{
}
@@ -25,7 +25,7 @@ TriangularSylvester::TriangularSylvester(const SchurDecompZero &kdecomp,
const SchurDecomp &fdecomp)
: SylvesterSolver(kdecomp, fdecomp),
matrixKK(matrixK->clone(2, *matrixK)),
- matrixFF(new QuasiTriangular(2, *matrixF))
+ matrixFF(std::make_unique<QuasiTriangular>(2, *matrixF))
{
}
@@ -33,16 +33,10 @@ TriangularSylvester::TriangularSylvester(const SchurDecompZero &kdecomp,
const SimilarityDecomp &fdecomp)
: SylvesterSolver(kdecomp, fdecomp),
matrixKK(matrixK->clone(2, *matrixK)),
- matrixFF(new BlockDiagonal(2, *matrixF))
+ matrixFF(std::make_unique<BlockDiagonal>(2, *matrixF))
{
}
-TriangularSylvester::~TriangularSylvester()
-{
- delete matrixKK;
- delete matrixFF;
-}
-
void
TriangularSylvester::print() const
{
@@ -65,9 +59,8 @@ TriangularSylvester::solvi(double r, KronVector &d, double &eig_min) const
{
if (d.getDepth() == 0)
{
- QuasiTriangular *t = matrixK->clone(r);
+ auto t = matrixK->clone(r);
t->solvePre(d, eig_min);
- delete t;
}
else
{
@@ -114,9 +107,8 @@ TriangularSylvester::solviip(double alpha, double betas,
if (d.getDepth() == 0)
{
double aspbs = alpha*alpha+betas;
- QuasiTriangular *t = matrixK->clone(2*alpha, aspbs, *matrixKK);
+ auto t = matrixK->clone(2*alpha, aspbs, *matrixKK);
t->solvePre(d, eig_min);
- delete t;
}
else
{
diff --git a/dynare++/sylv/cc/TriangularSylvester.hh b/dynare++/sylv/cc/TriangularSylvester.hh
index 6838da0319..e804e3343a 100644
--- a/dynare++/sylv/cc/TriangularSylvester.hh
+++ b/dynare++/sylv/cc/TriangularSylvester.hh
@@ -11,16 +11,18 @@
#include "QuasiTriangularZero.hh"
#include "SimilarityDecomp.hh"
+#include <memory>
+
class TriangularSylvester : public SylvesterSolver
{
- const QuasiTriangular *const matrixKK;
- const QuasiTriangular *const matrixFF;
+ const std::unique_ptr<const QuasiTriangular> matrixKK;
+ const std::unique_ptr<const QuasiTriangular> matrixFF;
public:
TriangularSylvester(const QuasiTriangular &k, const QuasiTriangular &f);
TriangularSylvester(const SchurDecompZero &kdecomp, const SchurDecomp &fdecomp);
TriangularSylvester(const SchurDecompZero &kdecomp, const SimilarityDecomp &fdecomp);
- ~TriangularSylvester() override;
+ ~TriangularSylvester() override = default;
void print() const;
void solve(SylvParams &pars, KronVector &d) const override;
diff --git a/dynare++/sylv/cc/Vector.cc b/dynare++/sylv/cc/Vector.cc
index 1c529113b9..c58bc089d2 100644
--- a/dynare++/sylv/cc/Vector.cc
+++ b/dynare++/sylv/cc/Vector.cc
@@ -6,12 +6,13 @@
#include <dynblas.h>
-#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <algorithm>
#include <limits>
+#include <iostream>
+#include <iomanip>
using namespace std;
@@ -29,42 +30,41 @@ Vector::Vector(const ConstVector &v)
copy(v.base(), v.skip());
}
-const Vector &
+Vector &
Vector::operator=(const Vector &v)
{
if (this == &v)
return *this;
if (v.length() != length())
- {
- throw SYLV_MES_EXCEPTION("Attempt to assign vectors with different lengths.");
- }
+ throw SYLV_MES_EXCEPTION("Attempt to assign vectors with different lengths.");
+
if (s == v.s
&& (data <= v.data && v.data < data+len*s
|| v.data <= data && data < v.data+v.len*v.s)
&& (data-v.data) % s == 0)
{
- printf("this destroy=%d, v destroy=%d, data-v.data=%lu, len=%d\n", destroy, v.destroy, (unsigned long) (data-v.data), len);
+ std::cout << "this destroy=" << destroy << ", v destroy=" << v.destroy
+ << ", data-v.data=" << (unsigned long) (data-v.data)
+ << ", len=" << len << std::endl;
throw SYLV_MES_EXCEPTION("Attempt to assign overlapping vectors.");
}
copy(v.base(), v.skip());
return *this;
}
-const Vector &
+Vector &
Vector::operator=(const ConstVector &v)
{
if (v.length() != length())
- {
- throw SYLV_MES_EXCEPTION("Attempt to assign vectors with different lengths.");
- }
+ throw SYLV_MES_EXCEPTION("Attempt to assign vectors with different lengths.");
+
if (v.skip() == 1 && skip() == 1 && (
(base() < v.base() + v.length() && base() >= v.base())
|| (base() + length() < v.base() + v.length()
&& base() + length() > v.base())))
- {
- throw SYLV_MES_EXCEPTION("Attempt to assign overlapping vectors.");
- }
+ throw SYLV_MES_EXCEPTION("Attempt to assign overlapping vectors.");
+
copy(v.base(), v.skip());
return *this;
}
@@ -175,9 +175,7 @@ Vector::infs()
Vector::~Vector()
{
if (destroy)
- {
- delete [] data;
- }
+ delete[] data;
}
void
@@ -283,28 +281,25 @@ Vector::isFinite() const
void
Vector::print() const
{
+ auto ff = std::cout.flags();
+ std::cout << std::setprecision(4);
for (int i = 0; i < length(); i++)
- {
- printf("%d\t%8.4g\n", i, operator[](i));
- }
+ std::cout << i << '\t' << std::setw(8) << operator[](i) << std::endl;
+ std::cout.flags(ff);
}
ConstVector::ConstVector(const Vector &v, int off, int l)
: BaseConstVector(l, v.skip(), v.base() + v.skip()*off)
{
if (off < 0 || off + length() > v.length())
- {
- throw SYLV_MES_EXCEPTION("Subvector not contained in supvector.");
- }
+ throw SYLV_MES_EXCEPTION("Subvector not contained in supvector.");
}
ConstVector::ConstVector(const ConstVector &v, int off, int l)
: BaseConstVector(l, v.skip(), v.base() + v.skip()*off)
{
if (off < 0 || off + length() > v.length())
- {
- throw SYLV_MES_EXCEPTION("Subvector not contained in supvector.");
- }
+ throw SYLV_MES_EXCEPTION("Subvector not contained in supvector.");
}
ConstVector::ConstVector(const double *d, int skip, int l)
@@ -353,9 +348,7 @@ ConstVector::getNorm() const
{
double s = 0;
for (int i = 0; i < length(); i++)
- {
- s += operator[](i)*operator[](i);
- }
+ s += operator[](i)*operator[](i);
return sqrt(s);
}
@@ -364,10 +357,8 @@ ConstVector::getMax() const
{
double r = 0;
for (int i = 0; i < length(); i++)
- {
- if (abs(operator[](i)) > r)
- r = abs(operator[](i));
- }
+ if (abs(operator[](i)) > r)
+ r = abs(operator[](i));
return r;
}
@@ -376,9 +367,7 @@ ConstVector::getNorm1() const
{
double norm = 0.0;
for (int i = 0; i < length(); i++)
- {
- norm += abs(operator[](i));
- }
+ norm += abs(operator[](i));
return norm;
}
@@ -405,10 +394,11 @@ ConstVector::isFinite() const
void
ConstVector::print() const
{
+ auto ff = std::cout.flags();
+ std::cout << std::setprecision(4);
for (int i = 0; i < length(); i++)
- {
- printf("%d\t%8.4g\n", i, operator[](i));
- }
+ std::cout << i << '\t' << std::setw(8) << operator[](i) << std::endl;
+ std::cout.flags(ff);
}
ZeroPad::ZeroPad()
diff --git a/dynare++/sylv/cc/Vector.hh b/dynare++/sylv/cc/Vector.hh
index 5873a0a12d..ce08790fd5 100644
--- a/dynare++/sylv/cc/Vector.hh
+++ b/dynare++/sylv/cc/Vector.hh
@@ -10,7 +10,6 @@
* members, and methods are thus duplicated */
#include <array>
-#include <cstdio>
class GeneralMatrix;
class ConstVector;
@@ -23,8 +22,7 @@ protected:
double *data{nullptr};
bool destroy{false};
public:
- Vector()
- = default;
+ Vector() = default;
Vector(int l) : len(l), data(new double[l]), destroy(true)
{
}
@@ -50,8 +48,8 @@ public:
Vector(const Vector &v, int off, int l);
Vector(GeneralMatrix &m, int col);
Vector(int row, GeneralMatrix &m);
- const Vector &operator=(const Vector &v);
- const Vector &operator=(const ConstVector &v);
+ Vector &operator=(const Vector &v);
+ Vector &operator=(const ConstVector &v);
double &
operator[](int i)
{
@@ -92,8 +90,7 @@ public:
bool operator>(const Vector &y) const;
bool operator>=(const Vector &y) const;
- virtual
- ~Vector();
+ virtual ~Vector();
void zeros();
void nans();
void infs();
@@ -150,11 +147,8 @@ public:
BaseConstVector(int l, int si, const double *d) : len(l), s(si), data(d)
{
}
- BaseConstVector(const BaseConstVector &v)
- = default;
- BaseConstVector &
- operator=(const BaseConstVector &v)
- = default;
+ BaseConstVector(const BaseConstVector &v) = default;
+ BaseConstVector &operator=(const BaseConstVector &v) = default;
const double &
operator[](int i) const
{
@@ -185,8 +179,7 @@ public:
ConstVector(const Vector &v) : BaseConstVector(v.length(), v.skip(), v.base())
{
}
- ConstVector(const ConstVector &v)
- = default;
+ ConstVector(const ConstVector &v) = default;
ConstVector(const double *d, int l) : BaseConstVector(l, 1, d)
{
}
@@ -196,8 +189,8 @@ public:
ConstVector(const ConstGeneralMatrix &m, int col);
ConstVector(int row, const ConstGeneralMatrix &m);
- virtual ~ConstVector()
- = default;
+ virtual ~ConstVector() = default;
+ ConstVector &operator=(const ConstVector &v) = default;
/** Exact equality. */
bool operator==(const ConstVector &y) const;
bool
--
GitLab