Dynare++ sylvester equation solver: remove the SylvMemory class, not actually used

parent 3ce051d8
......@@ -15,7 +15,7 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols,
const ConstVector &dc, const ConstVector &dd,
const SylvParams &ps)
: pars(ps),
mem_driver(pars, 1, m, n, ord), order(ord), a(Vector{da}, n),
order(ord), a(Vector{da}, n),
b(Vector{db}, n, n-zero_cols), c(Vector{dc}, m), d(Vector{dd}, n, power(m, order)),
solved(false)
{
......@@ -27,7 +27,7 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols,
const ConstVector &dc, Vector &dd,
const SylvParams &ps)
: pars(ps),
mem_driver(pars, 0, m, n, ord), order(ord), a(Vector{da}, n),
order(ord), a(Vector{da}, n),
b(Vector{db}, n, n-zero_cols), c(Vector{dc}, m), d(dd, n, power(m, order)),
solved(false)
{
......@@ -39,7 +39,7 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols,
const ConstVector &dc, const ConstVector &dd,
bool alloc_for_check)
: pars(alloc_for_check),
mem_driver(pars, 1, m, n, ord), order(ord), a(Vector{da}, n),
order(ord), a(Vector{da}, n),
b(Vector{db}, n, n-zero_cols), c(Vector{dc}, m), d(Vector{dd}, n, power(m, order)),
solved(false)
{
......@@ -51,7 +51,7 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols,
const ConstVector &dc, Vector &dd,
bool alloc_for_check)
: pars(alloc_for_check),
mem_driver(pars, 0, m, n, ord), order(ord), a(Vector{da}, n),
order(ord), a(Vector{da}, n),
b(Vector{db}, n, n-zero_cols), c(Vector{dc}, m), d(dd, n, power(m, order)),
solved(false)
{
......@@ -83,8 +83,6 @@ GeneralSylvester::solve()
if (solved)
throw SYLV_MES_EXCEPTION("Attempt to run solve() more than once.");
mem_driver.setStackMode(true);
clock_t start = clock();
// multiply d
d.multLeftITrans(bdecomp->getQ());
......@@ -99,8 +97,6 @@ GeneralSylvester::solve()
clock_t end = clock();
pars.cpu_time = ((double) (end-start))/CLOCKS_PER_SEC;
mem_driver.setStackMode(false);
solved = true;
}
......@@ -110,8 +106,6 @@ GeneralSylvester::check(const ConstVector &ds)
if (!solved)
throw SYLV_MES_EXCEPTION("Cannot run check on system, which is not solved yet.");
mem_driver.setStackMode(true);
// calculate xcheck = AX+BXC^i-D
SylvMatrix dcheck(d.numRows(), d.numCols());
dcheck.multLeft(b.numRows()-b.numCols(), b, d);
......@@ -124,6 +118,4 @@ GeneralSylvester::check(const ConstVector &ds)
pars.mat_errF = dcheck.getData().getNorm()/d.getData().getNorm();
pars.vec_err1 = dcheck.getData().getNorm1()/d.getData().getNorm1();
pars.vec_errI = dcheck.getData().getMax()/d.getData().getMax();
mem_driver.setStackMode(false);
}
......@@ -6,7 +6,6 @@
#define GENERAL_SYLVESTER_H
#include "SylvMatrix.hh"
#include "SylvMemory.hh"
#include "SimilarityDecomp.hh"
#include "SylvesterSolver.hh"
......@@ -15,7 +14,6 @@
class GeneralSylvester
{
SylvParams pars;
SylvMemoryDriver mem_driver;
int order;
const SqSylvMatrix a;
const SylvMatrix b;
......
......@@ -30,8 +30,6 @@ libsylv_a_SOURCES = \
SylvException.hh \
SylvMatrix.cc \
SylvMatrix.hh \
SylvMemory.cc \
SylvMemory.hh \
SylvParams.cc \
SylvParams.hh \
SylvesterSolver.hh \
......
......@@ -7,9 +7,7 @@
#include <string>
#include "SylvMemory.hh"
class SylvException : public MallocAllocator
class SylvException
{
protected:
std::string file;
......
/* $Header: /var/lib/cvs/dynare_cpp/sylv/cc/SylvMemory.cpp,v 1.1.1.1 2004/06/04 13:00:49 kamenik Exp $ */
/* Tag $Name: $ */
#include "SylvMemory.hh"
#include "SylvException.hh"
#include "KronVector.hh"
#if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
# include <dynmex.h>
#endif
#include <cmath>
#include <cstdlib>
/**********************************************************/
/* SylvMemoryPool */
/**********************************************************/
SylvMemoryPool memory_pool;
void
SylvMemoryPool::init(size_t size)
{
#ifdef USE_MEMORY_POOL
length = size;
# if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
if (base)
throw SYLV_MES_EXCEPTION("Attempt to use matlab memory pool twice.");
base = (char *) mxMalloc(length);
# else
base = (char *) malloc(length);
# endif
#else
throw SYLV_MES_EXCEPTION("SylvMemoryPool::init() called for non memory pool code.");
#endif
}
void *
SylvMemoryPool::allocate(size_t size)
{
#ifdef USE_MEMORY_POOL
if (allocated + size < length)
{
char *res = base + allocated;
allocated += size;
return res;
}
else
{
throw SYLV_MES_EXCEPTION("Run out of memory space");
}
#else
throw SYLV_MES_EXCEPTION("SylvMemoryPool::allocate() called for non memory pool code.");
#endif
}
void
SylvMemoryPool::free(void *p)
{
#ifdef USE_MEMORY_POOL
int offset = ((char *) p) - base;
# ifdef DEBUG
if (offset < 0)
throw SYLV_MES_EXCEPTION("SylvMemoryPool::free() frees wrong address < begin.");
if (offset >= (int) length)
throw SYLV_MES_EXCEPTION("SylvMemoryPool::free() frees wrong address > end.");
# endif
if (stack_mode && offset >= 0 && offset < (int) allocated)
allocated = offset;
#else
throw SYLV_MES_EXCEPTION("SylvMemoryPool::free() called for non memory pool code.");
#endif
}
void
SylvMemoryPool::setStackMode(bool mode)
{
stack_mode = mode;
}
SylvMemoryPool::~SylvMemoryPool()
{
reset();
}
void
SylvMemoryPool::reset()
{
#if !defined(MATLAB_MEX_FILE) && !defined(OCTAVE_MEX_FILE)
delete [] base;
base = 0;
allocated = 0;
length = 0;
stack_mode = false;
#endif
}
/**********************************************************/
/* global new and delete */
/**********************************************************/
#ifdef USE_MEMORY_POOL
void *
operator new(size_t size)
{
return memory_pool.allocate(size);
}
void *
operator new[](size_t size)
{
return memory_pool.allocate(size);
}
void
operator delete(void *p)
{
memory_pool.free(p);
}
void
operator delete[](void *p)
{
memory_pool.free(p);
}
#endif
/**********************************************************/
/* saved version of global new and delete */
/**********************************************************/
#ifdef USE_MEMORY_POOL
void *
MallocAllocator::operator new(size_t size)
{
# if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
throw SYLV_MES_EXCEPTION("Attempt to call wrong memory allocator.");
# else
void *res = malloc(size);
if (!res)
throw SYLV_MES_EXCEPTION("Malloc unable to allocate memory.");
return res;
# endif
}
void *
MallocAllocator::operator new[](size_t size)
{
# if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
throw SYLV_MES_EXCEPTION("Attempt to call wrong memory allocator.");
# else
void *res = malloc(size);
if (!res)
throw SYLV_MES_EXCEPTION("Malloc unable allocate memory.");
return res;
# endif
}
void
MallocAllocator::operator delete(void *p)
{
# if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
throw SYLV_MES_EXCEPTION("Attempt to call wrong memory destructor.");
# else
free(p);
# endif
}
void
MallocAllocator::operator delete[](void *p)
{
# if defined(MATLAB_MEX_FILE) || defined(OCTAVE_MEX_FILE)
throw SYLV_MES_EXCEPTION("Attempt to call wrong memory destructor.");
# else
free(p);
# endif
}
#endif
/**********************************************************/
/* SylvMemoryDriver */
/**********************************************************/
void
SylvMemoryDriver::allocate(int num_d, int m, int n, int order)
{
#ifdef USE_MEMORY_POOL
int x_cols = power(m, order);
int total = num_d*x_cols*n; // storage for big matrices
total += x_cols; // storage for one extra row of a big matrix
int dig_vectors = (int) ceil(((double) (power(m, order)-1))/(m-1));
total += 8*n*dig_vectors; // storage for kron vectors instantiated during solv
total += 50*(m*m+n*n); // some storage for small square matrices
total *= sizeof(double); // everything in doubles
memory_pool.init(total);
#endif
}
SylvMemoryDriver::SylvMemoryDriver(int num_d, int m, int n, int order)
{
allocate(num_d, m, n, order);
}
SylvMemoryDriver::SylvMemoryDriver(const SylvParams &pars, int num_d,
int m, int n, int order)
{
if (*(pars.method) == SylvParams::iter)
num_d++;
if (*(pars.want_check))
num_d++;
allocate(num_d, m, n, order);
}
SylvMemoryDriver::~SylvMemoryDriver()
{
memory_pool.reset();
}
void
SylvMemoryDriver::setStackMode(bool mode)
{
memory_pool.setStackMode(mode);
}
/* $Header: /var/lib/cvs/dynare_cpp/sylv/cc/SylvMemory.h,v 1.1.1.1 2004/06/04 13:00:49 kamenik Exp $ */
/* Tag $Name: $ */
#ifndef SYLV_MEMORY_H
#define SYLV_MEMORY_H
#include "SylvParams.hh"
#include <new>
class MallocAllocator
{
#ifdef USE_MEMORY_POOL
public:
void *operator new(size_t size);
void *operator new[](size_t size);
void operator delete(void *p);
void operator delete[](void *p);
#endif
};
#ifdef USE_MEMORY_POOL
void *operator new(size_t size);
void *operator new[](size_t size);
void operator delete(void *p);
void operator delete[](void *p);
#endif
class SylvMemoryPool
{
char *base{nullptr};
size_t length{0};
size_t allocated{0};
bool stack_mode{false};
public:
SylvMemoryPool() = default;
SylvMemoryPool(const SylvMemoryPool &) = delete;
SylvMemoryPool &operator=(const SylvMemoryPool &) = delete;
~SylvMemoryPool();
void init(size_t size);
void *allocate(size_t size);
void free(void *p);
void reset();
void setStackMode(bool);
};
class SylvMemoryDriver
{
public:
SylvMemoryDriver(int num_d, int m, int n, int order);
SylvMemoryDriver(const SylvParams &pars, int num_d, int m, int n, int order);
SylvMemoryDriver(const SylvMemoryDriver &) = delete;
SylvMemoryDriver &operator=(const SylvMemoryDriver &) = delete;
static void setStackMode(bool);
~SylvMemoryDriver();
protected:
void allocate(int num_d, int m, int n, int order);
};
#endif /* SYLV_MEMORY_H */
......@@ -6,13 +6,12 @@
#define MM_MATRIX_H
#include "GeneralMatrix.hh"
#include "SylvMemory.hh"
#include <string>
#include <utility>
#include <memory>
class MMException : public MallocAllocator
class MMException
{
std::string message;
public:
......@@ -26,7 +25,7 @@ public:
}
};
class MMMatrixIn : public MallocAllocator
class MMMatrixIn
{
std::shared_ptr<double> data;
int rows;
......@@ -56,7 +55,7 @@ public:
}
};
class MMMatrixOut : public MallocAllocator
class MMMatrixOut
{
public:
static void write(const std::string &fname, const GeneralMatrix &m);
......
......@@ -10,7 +10,6 @@
#include "KronUtils.hh"
#include "TriangularSylvester.hh"
#include "GeneralSylvester.hh"
#include "SylvMemory.hh"
#include "SchurDecompEig.hh"
#include "SimilarityDecomp.hh"
#include "IterativeSylvester.hh"
......@@ -26,7 +25,7 @@
#include <iomanip>
#include <memory>
class TestRunnable : public MallocAllocator
class TestRunnable
{
public:
const std::string name;
......@@ -89,7 +88,6 @@ TestRunnable::quasi_solve(bool trans, const std::string &mname, const std::strin
MMMatrixIn mmt(mname);
MMMatrixIn mmv(vname);
SylvMemoryDriver memdriver(1, mmt.row(), mmt.row(), 1);
std::unique_ptr<QuasiTriangular> t;
std::unique_ptr<QuasiTriangular> tsave;
if (mmt.row() == mmt.col())
......@@ -147,7 +145,6 @@ TestRunnable::mult_kron(bool trans, const std::string &mname, const std::string
return false;
}
SylvMemoryDriver memdriver(1, m, n, depth);
QuasiTriangular t(mmt.getData(), mmt.row());
Vector vraw{mmv.getData()};
KronVector v(vraw, m, n, depth);
......@@ -184,7 +181,6 @@ TestRunnable::level_kron(bool trans, const std::string &mname, const std::string
return false;
}
SylvMemoryDriver memdriver(1, m, n, depth);
QuasiTriangular t(mmt.getData(), mmt.row());
Vector vraw{mmv.getData()};
ConstKronVector v(vraw, m, n, depth);
......@@ -224,7 +220,6 @@ TestRunnable::kron_power(const std::string &m1name, const std::string &m2name, c
return false;
}
SylvMemoryDriver memdriver(2, m, n, depth);
QuasiTriangular t1(mmt1.getData(), mmt1.row());
QuasiTriangular t2(mmt2.getData(), mmt2.row());
Vector vraw{mmv.getData()};
......@@ -232,9 +227,7 @@ TestRunnable::kron_power(const std::string &m1name, const std::string &m2name, c
Vector craw{mmc.getData()};
KronVector c(craw, m, n, depth);
KronVector x(v);
memdriver.setStackMode(true);
KronUtils::multKron(t1, t2, x);
memdriver.setStackMode(false);
x.add(-1, c);
double norm = x.getNorm();
std::cout << "\terror norm = " << norm << std::endl;
......@@ -265,7 +258,6 @@ TestRunnable::lin_eval(const std::string &m1name, const std::string &m2name, con
return false;
}
SylvMemoryDriver memdriver(1, m, n, depth);
QuasiTriangular t1(mmt1.getData(), mmt1.row());
QuasiTriangular t2(mmt2.getData(), mmt2.row());
TriangularSylvester ts(t2, t1);
......@@ -279,9 +271,7 @@ TestRunnable::lin_eval(const std::string &m1name, const std::string &m2name, con
ConstKronVector c2(craw2, m, n, depth);
KronVector x1(m, n, depth);
KronVector x2(m, n, depth);
memdriver.setStackMode(true);
ts.linEval(alpha, beta1, beta2, x1, x2, v1, v2);
memdriver.setStackMode(false);
x1.add(-1, c1);
x2.add(-1, c2);
double norm1 = x1.getNorm();
......@@ -315,7 +305,6 @@ TestRunnable::qua_eval(const std::string &m1name, const std::string &m2name, con
return false;
}
SylvMemoryDriver memdriver(3, m, n, depth);
QuasiTriangular t1(mmt1.getData(), mmt1.row());
QuasiTriangular t2(mmt2.getData(), mmt2.row());
TriangularSylvester ts(t2, t1);
......@@ -329,9 +318,7 @@ TestRunnable::qua_eval(const std::string &m1name, const std::string &m2name, con
ConstKronVector c2(craw2, m, n, depth);
KronVector x1(m, n, depth);
KronVector x2(m, n, depth);
memdriver.setStackMode(true);
ts.quaEval(alpha, betas, gamma, delta1, delta2, x1, x2, v1, v2);
memdriver.setStackMode(false);
x1.add(-1, c1);
x2.add(-1, c2);
double norm1 = x1.getNorm();
......@@ -361,8 +348,6 @@ TestRunnable::tri_sylv(const std::string &m1name, const std::string &m2name, con
return false;
}
SylvMemoryDriver memdriver(4, m, n, depth); // need extra 2 for checks done via KronUtils::multKron
memdriver.setStackMode(true);
QuasiTriangular t1(mmt1.getData(), mmt1.row());
QuasiTriangular t2(mmt2.getData(), mmt2.row());
TriangularSylvester ts(t2, t1);
......@@ -382,7 +367,6 @@ TestRunnable::tri_sylv(const std::string &m1name, const std::string &m2name, con
double max = dcheck.getMax();
double xmax = v.getMax();
std::cout << "\trel. error max = " << max/xmax << std::endl;
memdriver.setStackMode(false);
return (norm < xnorm*eps_norm);
}
......@@ -430,7 +414,6 @@ TestRunnable::eig_bubble(const std::string &aname, int from, int to)
}
int n = mma.row();
SylvMemoryDriver memdriver(3, n, n, 2);
QuasiTriangular orig(mma.getData(), n);
SchurDecompEig dec((const QuasiTriangular &)orig);
QuasiTriangular::diag_iter itf = dec.getT().diag_begin();
......@@ -466,7 +449,6 @@ TestRunnable::block_diag(const std::string &aname, double log10norm)
}
int n = mma.row();
SylvMemoryDriver memdriver(3, n, n, 2);
SqSylvMatrix orig(mma.getData(), n);
SimilarityDecomp dec(orig.getData(), orig.numRows(), log10norm);
dec.getB().printInfo();
......@@ -515,8 +497,6 @@ TestRunnable::iter_sylv(const std::string &m1name, const std::string &m2name, co
return false;
}
SylvMemoryDriver memdriver(4, m, n, depth); // need extra 2 for checks done via KronUtils::multKron
memdriver.setStackMode(true);
QuasiTriangular t1(mmt1.getData(), mmt1.row());
QuasiTriangular t2(mmt2.getData(), mmt2.row());
IterativeSylvester is(t2, t1);
......@@ -537,7 +517,6 @@ TestRunnable::iter_sylv(const std::string &m1name, const std::string &m2name, co
double max = dcheck.getMax();
double xmax = v.getMax();
std::cout << "\trel. error max = " << max/xmax << std::endl;
memdriver.setStackMode(false);
return (cnorm < xnorm*eps_norm);
}
......
......@@ -59,8 +59,6 @@ SYLV_SRCS = \
$(TOPDIR)/sylv/cc/SylvException.hh \
$(TOPDIR)/sylv/cc/SylvMatrix.cc \
$(TOPDIR)/sylv/cc/SylvMatrix.hh \
$(TOPDIR)/sylv/cc/SylvMemory.cc \
$(TOPDIR)/sylv/cc/SylvMemory.hh \
$(TOPDIR)/sylv/cc/SylvParams.cc \
$(TOPDIR)/sylv/cc/SylvParams.hh \
$(TOPDIR)/sylv/cc/SylvesterSolver.hh \
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment