stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Dynare++ tensor library: modernize the Symmetry class 

We now use a initializer list constructor for creating symmetries of the form
$y^n$, $y^n u^m$, $y^nu^m\sigma^k$.

The constructor taking a single integer is used to initialize a symmetry of a
given length.

Similar changes are made to IntSequence.

This behavior is similar to std::vector.
time-shift
Sébastien Villemot 2019-02-08 18:38:05 +01:00
parent 219d2bb31b
commit 1f7d3beddc
No known key found for this signature in database
GPG Key ID: 2CECE9350ECEBE4A
28 changed files with 206 additions and 218 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
dynare++/integ/cc/product.cc
View File

@ -10,7 +10,7 @@
prodpit::prodpit(const ProductQuadrature &q, int j0, int l) prodpit::prodpit(const ProductQuadrature &q, int j0, int l)
: prodq(q), level(l), npoints(q.uquad.numPoints(l)), : prodq(q), level(l), npoints(q.uquad.numPoints(l)),
jseq{q.dimen(), 0}, jseq(q.dimen(), 0),
end_flag(false), end_flag(false),
sig{q.dimen()}, sig{q.dimen()},
p{q.dimen()} p{q.dimen()}

2
dynare++/integ/cc/smolyak.cc
View File

@ -12,7 +12,7 @@
smolpit::smolpit(const SmolyakQuadrature &q, unsigned int isum) smolpit::smolpit(const SmolyakQuadrature &q, unsigned int isum)
: smolq(q), isummand(isum), : smolq(q), isummand(isum),
jseq{q.dimen(), 0}, jseq(q.dimen(), 0),
sig{q.dimen()}, sig{q.dimen()},
p{q.dimen()} p{q.dimen()}
{ {

4
dynare++/integ/testing/tests.cc
View File

@ -256,7 +256,7 @@ TestRunnable::smolyak_normal_moments(const GeneralMatrix &m, int imom, int level
// check against theoretical moments // check against theoretical moments
UNormalMoments moments(imom, msq); UNormalMoments moments(imom, msq);
smol_out.add(-1.0, (moments.get(Symmetry(imom)))->getData()); smol_out.add(-1.0, (moments.get(Symmetry{imom}))->getData());
std::cout << "\tError: " << std::setw(16) << std::setprecision(12) << smol_out.getMax() << std::endl; std::cout << "\tError: " << std::setw(16) << std::setprecision(12) << smol_out.getMax() << std::endl;
return smol_out.getMax() < 1.e-7; return smol_out.getMax() < 1.e-7;
} }
@ -285,7 +285,7 @@ TestRunnable::product_normal_moments(const GeneralMatrix &m, int imom, int level
// check against theoretical moments // check against theoretical moments
UNormalMoments moments(imom, msq); UNormalMoments moments(imom, msq);
prod_out.add(-1.0, (moments.get(Symmetry(imom)))->getData()); prod_out.add(-1.0, (moments.get(Symmetry{imom}))->getData());
std::cout << "\tError: " << std::setw(16) << std::setprecision(12) << prod_out.getMax() << std::endl; std::cout << "\tError: " << std::setw(16) << std::setprecision(12) << prod_out.getMax() << std::endl;
return prod_out.getMax() < 1.e-7; return prod_out.getMax() < 1.e-7;
} }

14
dynare++/kord/approximation.cc
View File

@ -80,7 +80,7 @@ Approximation::approxAtSteady()
{ {
model.calcDerivativesAtSteady(); model.calcDerivativesAtSteady();
FirstOrder fo(model.nstat(), model.npred(), model.nboth(), model.nforw(), FirstOrder fo(model.nstat(), model.npred(), model.nboth(), model.nforw(),
model.nexog(), *(model.getModelDerivatives().get(Symmetry(1))), model.nexog(), *(model.getModelDerivatives().get(Symmetry{1})),
journal, qz_criterium); journal, qz_criterium);
KORD_RAISE_IF_X(!fo.isStable(), KORD_RAISE_IF_X(!fo.isStable(),
"The model is not Blanchard-Kahn stable", "The model is not Blanchard-Kahn stable",
@ -293,21 +293,21 @@ Approximation::calcStochShift(Vector &out, double at_sigma) const
{ {
if (KOrder::is_even(d)) if (KOrder::is_even(d))
{ {
Symmetry sym(0, d, 0, 0); Symmetry sym{0, d, 0, 0};
// calculate $F_{u'^d}$ via |ZAuxContainer| // calculate $F_{u'^d}$ via |ZAuxContainer|
FGSTensor *ten = new FGSTensor(ypart.ny(), TensorDimens(sym, nvs)); FGSTensor *ten = new FGSTensor(ypart.ny(), TensorDimens(sym, nvs));
ten->zeros(); ten->zeros();
for (int l = 1; l <= d; l++) for (int l = 1; l <= d; l++)
{ {
const FSSparseTensor *f = model.getModelDerivatives().get(Symmetry(l)); const FSSparseTensor *f = model.getModelDerivatives().get(Symmetry{l});
zaux.multAndAdd(*f, *ten); zaux.multAndAdd(*f, *ten);
} }
// multiply with shocks and add to result // multiply with shocks and add to result
FGSTensor *tmp = new FGSTensor(ypart.ny(), TensorDimens(Symmetry(0, 0, 0, 0), nvs)); FGSTensor *tmp = new FGSTensor(ypart.ny(), TensorDimens(Symmetry{0, 0, 0, 0}, nvs));
tmp->zeros(); tmp->zeros();
ten->contractAndAdd(1, *tmp, *(mom.get(Symmetry(d)))); ten->contractAndAdd(1, *tmp, *(mom.get(Symmetry{d})));
out.add(pow(at_sigma, d)/dfac, tmp->getData()); out.add(pow(at_sigma, d)/dfac, tmp->getData());
delete ten; delete ten;
@ -360,8 +360,8 @@ Approximation::check(double at_sigma) const
TwoDMatrix * TwoDMatrix *
Approximation::calcYCov() const Approximation::calcYCov() const
{ {
const TwoDMatrix &gy = *(rule_ders->get(Symmetry(1, 0, 0, 0))); const TwoDMatrix &gy = *(rule_ders->get(Symmetry{1, 0, 0, 0}));
const TwoDMatrix &gu = *(rule_ders->get(Symmetry(0, 1, 0, 0))); const TwoDMatrix &gu = *(rule_ders->get(Symmetry{0, 1, 0, 0}));
TwoDMatrix G(model.numeq(), model.numeq()); TwoDMatrix G(model.numeq(), model.numeq());
G.zeros(); G.zeros();
G.place(gy, 0, model.nstat()); G.place(gy, 0, model.nstat());

10
dynare++/kord/decision_rule.hh
View File

@ -217,11 +217,11 @@ DecisionRuleImpl<t>::fillTensors(const _Tg &g, double sigma)
int j = d-i; int j = d-i;
int kfact = 1; int kfact = 1;
_Ttensor tmp(ypart.ny(), _Ttensor tmp(ypart.ny(),
TensorDimens(Symmetry(i, j), tns)); TensorDimens(Symmetry{i, j}, tns));
tmp.zeros(); tmp.zeros();
for (int k = 0; k+d <= g.getMaxDim(); k++, kfact *= k) for (int k = 0; k+d <= g.getMaxDim(); k++, kfact *= k)
{ {
Symmetry sym(i, j, 0, k); Symmetry sym{i, j, 0, k};
if (g.check(sym)) if (g.check(sym))
{ {
double mult = pow(sigma, k)/dfact/kfact; double mult = pow(sigma, k)/dfact/kfact;
@ -539,7 +539,7 @@ DRFixPoint<t>::DRFixPoint(const _Tg &g, const PartitionY &yp,
fillTensors(g, sigma); fillTensors(g, sigma);
_Tparent yspol(ypart.nstat, ypart.nys(), *this); _Tparent yspol(ypart.nstat, ypart.nys(), *this);
bigf = new _Tparent((const _Tparent &) yspol); bigf = new _Tparent((const _Tparent &) yspol);
_Ttensym *frst = bigf->get(Symmetry(1)); _Ttensym *frst = bigf->get(Symmetry{1});
for (int i = 0; i < ypart.nys(); i++) for (int i = 0; i < ypart.nys(); i++)
frst->get(i, i) = frst->get(i, i) - 1; frst->get(i, i) = frst->get(i, i) - 1;
bigfder = new _Tparent(*bigf, 0); bigfder = new _Tparent(*bigf, 0);
@ -572,9 +572,9 @@ DRFixPoint<t>::fillTensors(const _Tg &g, double sigma)
int kfact = 1; int kfact = 1;
for (int k = 0; d+k <= g.getMaxDim(); k++, kfact *= k) for (int k = 0; d+k <= g.getMaxDim(); k++, kfact *= k)
{ {
if (g.check(Symmetry(d, 0, 0, k))) if (g.check(Symmetry{d, 0, 0, k}))
{ {
const _Ttensor *ten = g.get(Symmetry(d, 0, 0, k)); const _Ttensor *ten = g.get(Symmetry{d, 0, 0, k});
double mult = pow(sigma, k)/dfact/kfact; double mult = pow(sigma, k)/dfact/kfact;
g_yd->add(mult, *ten); g_yd->add(mult, *ten);
} }

4
dynare++/kord/first_order.hh
View File

@ -72,10 +72,10 @@ public:
{ {
IntSequence nvs(4); IntSequence nvs(4);
nvs[0] = fo.ypart.nys(); nvs[1] = fo.nu; nvs[2] = fo.nu; nvs[3] = 1; nvs[0] = fo.ypart.nys(); nvs[1] = fo.nu; nvs[2] = fo.nu; nvs[3] = 1;
_Ttensor *ten = new _Ttensor(fo.ypart.ny(), TensorDimens(Symmetry(1, 0, 0, 0), nvs)); _Ttensor *ten = new _Ttensor(fo.ypart.ny(), TensorDimens(Symmetry{1, 0, 0, 0}, nvs));
ten->zeros(); ten->add(1.0, fo.gy); ten->zeros(); ten->add(1.0, fo.gy);
this->insert(ten); this->insert(ten);
ten = new _Ttensor(fo.ypart.ny(), TensorDimens(Symmetry(0, 1, 0, 0), nvs)); ten = new _Ttensor(fo.ypart.ny(), TensorDimens(Symmetry{0, 1, 0, 0}, nvs));
ten->zeros(); ten->add(1.0, fo.gu); ten->zeros(); ten->add(1.0, fo.gu);
this->insert(ten); this->insert(ten);
} }

4
dynare++/kord/global_check.cc
View File

@ -95,9 +95,9 @@ ResidFunction::setYU(const ConstVector &ys, const ConstVector &xx)
hss = new FTensorPolynomial(dr_ss, ytmp_star); hss = new FTensorPolynomial(dr_ss, ytmp_star);
ConstVector ysteady_ss(dr.c->getSteady(), model->nstat()+model->npred(), ConstVector ysteady_ss(dr.c->getSteady(), model->nstat()+model->npred(),
model->nboth()+model->nforw()); model->nboth()+model->nforw());
if (hss->check(Symmetry(0))) if (hss->check(Symmetry{0}))
{ {
hss->get(Symmetry(0))->getData().add(1.0, ysteady_ss); hss->get(Symmetry{0})->getData().add(1.0, ysteady_ss);
} }
else else
{ {

18
dynare++/kord/korder.cc
View File

@ -239,9 +239,9 @@ KOrder::KOrder(int num_stat, int num_pred, int num_both, int num_forw,
_uGstack(&_ugs, ypart.nys(), nu), _uGstack(&_ugs, ypart.nys(), nu),
_fGstack(&_fgs, ypart.nys(), nu), _fGstack(&_fgs, ypart.nys(), nu),
_um(maxk, v), _fm(_um), f(fcont), _um(maxk, v), _fm(_um), f(fcont),
matA(*(f.get(Symmetry(1))), _uZstack.getStackSizes(), gy, ypart), matA(*(f.get(Symmetry{1})), _uZstack.getStackSizes(), gy, ypart),
matS(*(f.get(Symmetry(1))), _uZstack.getStackSizes(), gy, ypart), matS(*(f.get(Symmetry{1})), _uZstack.getStackSizes(), gy, ypart),
matB(*(f.get(Symmetry(1))), _uZstack.getStackSizes()), matB(*(f.get(Symmetry{1})), _uZstack.getStackSizes()),
journal(jr) journal(jr)
{ {
KORD_RAISE_IF(gy.ncols() != ypart.nys(), KORD_RAISE_IF(gy.ncols() != ypart.nys(),
@ -265,20 +265,20 @@ KOrder::KOrder(int num_stat, int num_pred, int num_both, int num_forw,
// put $g_y$ and $g_u$ to the container // put $g_y$ and $g_u$ to the container
/* Note that $g_\sigma$ is zero by the nature and we do not insert it to /* Note that $g_\sigma$ is zero by the nature and we do not insert it to
the container. We insert a new physical copies. */ the container. We insert a new physical copies. */
UGSTensor *tgy = new UGSTensor(ny, TensorDimens(Symmetry(1, 0, 0, 0), nvs)); UGSTensor *tgy = new UGSTensor(ny, TensorDimens(Symmetry{1, 0, 0, 0}, nvs));
tgy->getData() = gy.getData(); tgy->getData() = gy.getData();
insertDerivative<unfold>(tgy); insertDerivative<unfold>(tgy);
UGSTensor *tgu = new UGSTensor(ny, TensorDimens(Symmetry(0, 1, 0, 0), nvs)); UGSTensor *tgu = new UGSTensor(ny, TensorDimens(Symmetry{0, 1, 0, 0}, nvs));
tgu->getData() = gu.getData(); tgu->getData() = gu.getData();
insertDerivative<unfold>(tgu); insertDerivative<unfold>(tgu);
// put $G_y$, $G_u$ and $G_{u'}$ to the container // put $G_y$, $G_u$ and $G_{u'}$ to the container
/* Also note that since $g_\sigma$ is zero, so $G_\sigma$. */ /* Also note that since $g_\sigma$ is zero, so $G_\sigma$. */
UGSTensor *tGy = faaDiBrunoG<unfold>(Symmetry(1, 0, 0, 0)); UGSTensor *tGy = faaDiBrunoG<unfold>(Symmetry{1, 0, 0, 0});
G<unfold>().insert(tGy); G<unfold>().insert(tGy);
UGSTensor *tGu = faaDiBrunoG<unfold>(Symmetry(0, 1, 0, 0)); UGSTensor *tGu = faaDiBrunoG<unfold>(Symmetry{0, 1, 0, 0});
G<unfold>().insert(tGu); G<unfold>().insert(tGu);
UGSTensor *tGup = faaDiBrunoG<unfold>(Symmetry(0, 0, 1, 0)); UGSTensor *tGup = faaDiBrunoG<unfold>(Symmetry{0, 0, 1, 0});
G<unfold>().insert(tGup); G<unfold>().insert(tGup);
} }
@ -306,7 +306,7 @@ KOrder::sylvesterSolve<KOrder::unfold>(ctraits<unfold>::Ttensor &der) const
{ {
KORD_RAISE_IF(!der.isFinite(), KORD_RAISE_IF(!der.isFinite(),
"RHS of Sylverster is not finite"); "RHS of Sylverster is not finite");
TwoDMatrix gs_y(*(gs<unfold>().get(Symmetry(1, 0, 0, 0)))); TwoDMatrix gs_y(*(gs<unfold>().get(Symmetry{1, 0, 0, 0})));
GeneralSylvester sylv(der.getSym()[0], ny, ypart.nys(), GeneralSylvester sylv(der.getSym()[0], ny, ypart.nys(),
ypart.nstat+ypart.npred, ypart.nstat+ypart.npred,
matA.getData(), matB.getData(), matA.getData(), matB.getData(),

36
dynare++/kord/korder.hh
View File

@ -461,7 +461,7 @@ template<int t>
void void
KOrder::recover_y(int i) KOrder::recover_y(int i)
{ {
Symmetry sym(i, 0, 0, 0); Symmetry sym{i, 0, 0, 0};
JournalRecordPair pa(journal); JournalRecordPair pa(journal);
pa << "Recovering symmetry " << sym << endrec; pa << "Recovering symmetry " << sym << endrec;
@ -475,7 +475,7 @@ KOrder::recover_y(int i)
insertDerivative<t>(g_yi); insertDerivative<t>(g_yi);
_Ttensor *gss_y = gss<t>().get(Symmetry(1, 0, 0, 0)); _Ttensor *gss_y = gss<t>().get(Symmetry{1, 0, 0, 0});
gs<t>().multAndAdd(*gss_y, *G_yi); gs<t>().multAndAdd(*gss_y, *G_yi);
_Ttensor *gss_yi = gss<t>().get(sym); _Ttensor *gss_yi = gss<t>().get(sym);
gs<t>().multAndAdd(*gss_yi, *G_yi); gs<t>().multAndAdd(*gss_yi, *G_yi);
@ -496,7 +496,7 @@ template <int t>
void void
KOrder::recover_yu(int i, int j) KOrder::recover_yu(int i, int j)
{ {
Symmetry sym(i, j, 0, 0); Symmetry sym{i, j, 0, 0};
JournalRecordPair pa(journal); JournalRecordPair pa(journal);
pa << "Recovering symmetry " << sym << endrec; pa << "Recovering symmetry " << sym << endrec;
@ -508,7 +508,7 @@ KOrder::recover_yu(int i, int j)
matA.multInv(*g_yiuj); matA.multInv(*g_yiuj);
insertDerivative<t>(g_yiuj); insertDerivative<t>(g_yiuj);
gs<t>().multAndAdd(*(gss<t>().get(Symmetry(1, 0, 0, 0))), *G_yiuj); gs<t>().multAndAdd(*(gss<t>().get(Symmetry{1, 0, 0, 0})), *G_yiuj);
} }
/* Here we solve /* Here we solve
@ -535,7 +535,7 @@ template <int t>
void void
KOrder::recover_ys(int i, int j) KOrder::recover_ys(int i, int j)
{ {
Symmetry sym(i, 0, 0, j); Symmetry sym{i, 0, 0, j};
JournalRecordPair pa(journal); JournalRecordPair pa(journal);
pa << "Recovering symmetry " << sym << endrec; pa << "Recovering symmetry " << sym << endrec;
@ -595,7 +595,7 @@ template <int t>
void void
KOrder::recover_yus(int i, int j, int k) KOrder::recover_yus(int i, int j, int k)
{ {
Symmetry sym(i, j, 0, k); Symmetry sym{i, j, 0, k};
JournalRecordPair pa(journal); JournalRecordPair pa(journal);
pa << "Recovering symmetry " << sym << endrec; pa << "Recovering symmetry " << sym << endrec;
@ -662,7 +662,7 @@ template <int t>
void void
KOrder::recover_s(int i) KOrder::recover_s(int i)
{ {
Symmetry sym(0, 0, 0, i); Symmetry sym{0, 0, 0, i};
JournalRecordPair pa(journal); JournalRecordPair pa(journal);
pa << "Recovering symmetry " << sym << endrec; pa << "Recovering symmetry " << sym << endrec;
@ -710,7 +710,7 @@ KOrder::fillG(int i, int j, int k)
{ {
if (is_even(k-m)) if (is_even(k-m))
{ {
_Ttensor *G_yiujupms = faaDiBrunoG<t>(Symmetry(i, j, m, k-m)); _Ttensor *G_yiujupms = faaDiBrunoG<t>(Symmetry{i, j, m, k-m});
G<t>().insert(G_yiujupms); G<t>().insert(G_yiujupms);
} }
} }
@ -727,12 +727,12 @@ template <int t>
_Ttensor * _Ttensor *
KOrder::calcD_ijk(int i, int j, int k) const KOrder::calcD_ijk(int i, int j, int k) const
{ {
_Ttensor *res = new _Ttensor(ny, TensorDimens(Symmetry(i, j, 0, 0), nvs)); _Ttensor *res = new _Ttensor(ny, TensorDimens(Symmetry{i, j, 0, 0}, nvs));
res->zeros(); res->zeros();
if (is_even(k)) if (is_even(k))
{ {
_Ttensor *tmp = faaDiBrunoZ<t>(Symmetry(i, j, k, 0)); _Ttensor *tmp = faaDiBrunoZ<t>(Symmetry{i, j, k, 0});
tmp->contractAndAdd(2, *res, *(m<t>().get(Symmetry(k)))); tmp->contractAndAdd(2, *res, *(m<t>().get(Symmetry{k})));
delete tmp; delete tmp;
} }
return res; return res;
@ -749,13 +749,13 @@ template <int t>
_Ttensor * _Ttensor *
KOrder::calcE_ijk(int i, int j, int k) const KOrder::calcE_ijk(int i, int j, int k) const
{ {
_Ttensor *res = new _Ttensor(ny, TensorDimens(Symmetry(i, j, 0, 0), nvs)); _Ttensor *res = new _Ttensor(ny, TensorDimens(Symmetry{i, j, 0, 0}, nvs));
res->zeros(); res->zeros();
for (int n = 2; n <= k-1; n += 2) for (int n = 2; n <= k-1; n += 2)
{ {
_Ttensor *tmp = faaDiBrunoZ<t>(Symmetry(i, j, n, k-n)); _Ttensor *tmp = faaDiBrunoZ<t>(Symmetry{i, j, n, k-n});
tmp->mult((double) (Tensor::noverk(k, n))); tmp->mult((double) (Tensor::noverk(k, n)));
tmp->contractAndAdd(2, *res, *(m<t>().get(Symmetry(n)))); tmp->contractAndAdd(2, *res, *(m<t>().get(Symmetry{n})));
delete tmp; delete tmp;
} }
return res; return res;
@ -861,7 +861,7 @@ KOrder::check(int dim) const
// check for $F_{y^iu^j}=0 // check for $F_{y^iu^j}=0
for (int i = 0; i <= dim; i++) for (int i = 0; i <= dim; i++)
{ {
Symmetry sym(dim-i, i, 0, 0); Symmetry sym{dim-i, i, 0, 0};
_Ttensor *r = faaDiBrunoZ<t>(sym); _Ttensor *r = faaDiBrunoZ<t>(sym);
double err = r->getData().getMax(); double err = r->getData().getMax();
JournalRecord(journal) << "\terror for symmetry " << sym << "\tis " << err << endrec; JournalRecord(journal) << "\terror for symmetry " << sym << "\tis " << err << endrec;
@ -879,7 +879,7 @@ KOrder::check(int dim) const
int k = (*si)[2]; int k = (*si)[2];
if (i+j > 0 && k > 0) if (i+j > 0 && k > 0)
{ {
Symmetry sym(i, j, 0, k); Symmetry sym{i, j, 0, k};
_Ttensor *r = faaDiBrunoZ<t>(sym); _Ttensor *r = faaDiBrunoZ<t>(sym);
_Ttensor *D_ijk = calcD_ijk<t>(i, j, k); _Ttensor *D_ijk = calcD_ijk<t>(i, j, k);
r->add(1.0, *D_ijk); r->add(1.0, *D_ijk);
@ -894,7 +894,7 @@ KOrder::check(int dim) const
} }
// check for $F_{\sigma^i}+D_i+E_i=0 // check for $F_{\sigma^i}+D_i+E_i=0
_Ttensor *r = faaDiBrunoZ<t>(Symmetry(0, 0, 0, dim)); _Ttensor *r = faaDiBrunoZ<t>(Symmetry{0, 0, 0, dim});
_Ttensor *D_k = calcD_k<t>(dim); _Ttensor *D_k = calcD_k<t>(dim);
r->add(1.0, *D_k); r->add(1.0, *D_k);
delete D_k; delete D_k;
@ -902,7 +902,7 @@ KOrder::check(int dim) const
r->add(1.0, *E_k); r->add(1.0, *E_k);
delete E_k; delete E_k;
double err = r->getData().getMax(); double err = r->getData().getMax();
Symmetry sym(0, 0, 0, dim); Symmetry sym{0, 0, 0, dim};
JournalRecord(journal) << "\terror for symmetry " << sym << "\tis " << err << endrec; JournalRecord(journal) << "\terror for symmetry " << sym << "\tis " << err << endrec;
if (err > maxerror) if (err > maxerror)
maxerror = err; maxerror = err;

4
dynare++/kord/korder_stoch.cc
View File

@ -35,7 +35,7 @@ KOrderStoch::KOrderStoch(const PartitionY &yp, int nu,
_uGstack(&_ugs, ypart.nys(), nu), _uGstack(&_ugs, ypart.nys(), nu),
_fGstack(&_fgs, ypart.nys(), nu), _fGstack(&_fgs, ypart.nys(), nu),
f(fcont), f(fcont),
matA(*(fcont.get(Symmetry(1))), _uZstack.getStackSizes(), *(hh.get(Symmetry(1, 0, 0, 0))), matA(*(fcont.get(Symmetry{1})), _uZstack.getStackSizes(), *(hh.get(Symmetry{1, 0, 0, 0})),
ypart) ypart)
{ {
nvs[0] = ypart.nys(); nvs[0] = ypart.nys();
@ -56,7 +56,7 @@ KOrderStoch::KOrderStoch(const PartitionY &yp, int nu,
_uGstack(&_ugs, ypart.nys(), nu), _uGstack(&_ugs, ypart.nys(), nu),
_fGstack(&_fgs, ypart.nys(), nu), _fGstack(&_fgs, ypart.nys(), nu),
f(fcont), f(fcont),
matA(*(fcont.get(Symmetry(1))), _uZstack.getStackSizes(), *(hh.get(Symmetry(1, 0, 0, 0))), matA(*(fcont.get(Symmetry{1})), _uZstack.getStackSizes(), *(hh.get(Symmetry{1, 0, 0, 0})),
ypart) ypart)
{ {
nvs[0] = ypart.nys(); nvs[0] = ypart.nys();

22
dynare++/kord/korder_stoch.hh
View File

@ -87,7 +87,7 @@ IntegDerivs<t>::IntegDerivs(int r, const IntSequence &nvs, const _Tgss &g, const
for (int i = 0; i <= d; i++) for (int i = 0; i <= d; i++)
{ {
int p = d-i; int p = d-i;
Symmetry sym(i, 0, 0, p); Symmetry sym{i, 0, 0, p};
_Ttensor *ten = new _Ttensor(r, TensorDimens(sym, nvs)); _Ttensor *ten = new _Ttensor(r, TensorDimens(sym, nvs));
// calculate derivative $h_{y^i\sigma^p}$ // calculate derivative $h_{y^i\sigma^p}$
@ -105,14 +105,14 @@ IntegDerivs<t>::IntegDerivs(int r, const IntSequence &nvs, const _Tgss &g, const
for (int m = 0; i+m+n+k <= maxd; m++, mfac *= m) for (int m = 0; i+m+n+k <= maxd; m++, mfac *= m)
{ {
double mult = (pow(at_sigma, m)*povern)/mfac; double mult = (pow(at_sigma, m)*povern)/mfac;
Symmetry sym_mn(i, m+n, 0, k); Symmetry sym_mn{i, m+n, 0, k};
if (m+n == 0 && g.check(sym_mn)) if (m+n == 0 && g.check(sym_mn))
ten->add(mult, *(g.get(sym_mn))); ten->add(mult, *(g.get(sym_mn)));
if (m+n > 0 && KOrder::is_even(m+n) && g.check(sym_mn)) if (m+n > 0 && KOrder::is_even(m+n) && g.check(sym_mn))
{ {
_Ttensor gtmp(*(g.get(sym_mn))); _Ttensor gtmp(*(g.get(sym_mn)));
gtmp.mult(mult); gtmp.mult(mult);
gtmp.contractAndAdd(1, *ten, *(mom.get(Symmetry(m+n)))); gtmp.contractAndAdd(1, *ten, *(mom.get(Symmetry{m+n})));
} }
} }
} }
@ -187,8 +187,8 @@ StochForwardDerivs<t>::StochForwardDerivs(const PartitionY &ypart, int nu,
for (int i = 0; i <= d; i++) for (int i = 0; i <= d; i++)
{ {
int k = d-i; int k = d-i;
if (g_int.check(Symmetry(i, 0, 0, k))) if (g_int.check(Symmetry{i, 0, 0, k}))
ten->addSubTensor(*(g_int.get(Symmetry(i, 0, 0, k)))); ten->addSubTensor(*(g_int.get(Symmetry{i, 0, 0, k})));
} }
g_int_sym.insert(ten); g_int_sym.insert(ten);
} }
@ -224,13 +224,13 @@ StochForwardDerivs<t>::StochForwardDerivs(const PartitionY &ypart, int nu,
true_nvs[1] = nu; true_nvs[2] = nu; true_nvs[1] = nu; true_nvs[2] = nu;
for (int d = 1; d <= maxd; d++) for (int d = 1; d <= maxd; d++)
{ {
if (g_int_cent.check(Symmetry(d))) if (g_int_cent.check(Symmetry{d}))
{ {
for (int i = 0; i <= d; i++) for (int i = 0; i <= d; i++)
{ {
Symmetry sym(i, 0, 0, d-i); Symmetry sym{i, 0, 0, d-i};
IntSequence coor(sym, pp); IntSequence coor(sym, pp);
_Ttensor *ten = new _Ttensor(*(g_int_cent.get(Symmetry(d))), ss, coor, _Ttensor *ten = new _Ttensor(*(g_int_cent.get(Symmetry{d})), ss, coor,
TensorDimens(sym, true_nvs)); TensorDimens(sym, true_nvs));
this->insert(ten); this->insert(ten);
} }
@ -274,7 +274,7 @@ GXContainer<_Ttype>::getType(int i, const Symmetry &s) const
if (i == 2) if (i == 2)
return _Stype::zero; return _Stype::zero;
if (i == 3) if (i == 3)
if (s == Symmetry(0, 0, 0, 1)) if (s == Symmetry{0, 0, 0, 1})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;
@ -319,12 +319,12 @@ ZXContainer<_Ttype>::getType(int i, const Symmetry &s) const
else else
return _Stype::matrix; return _Stype::matrix;
if (i == 2) if (i == 2)
if (s == Symmetry(1, 0, 0, 0)) if (s == Symmetry{1, 0, 0, 0})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;
if (i == 3) if (i == 3)
if (s == Symmetry(0, 1, 0, 0)) if (s == Symmetry{0, 1, 0, 0})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;

2
dynare++/kord/tests.cc
View File

@ -287,7 +287,7 @@ TestRunnable::korder_unfold_fold(int maxdim, int unfold_dim,
for (int d = 1; d <= maxdim; d++) for (int d = 1; d <= maxdim; d++)
{ {
printf("\ttensor fill for dim=%d is: %3.2f %%\n", printf("\ttensor fill for dim=%d is: %3.2f %%\n",
d, c.get(Symmetry(d))->getFillFactor()*100.0); d, c.get(Symmetry{d})->getFillFactor()*100.0);
} }
Journal jr("out.txt"); Journal jr("out.txt");
KOrder kord(nstat, npred, nboth, nforw, c, gy, gu, v, jr); KOrder kord(nstat, npred, nboth, nforw, c, gy, gu, v, jr);

2
dynare++/src/dynare3.cc
View File

@ -353,7 +353,7 @@ DynareDerEvalLoader::DynareDerEvalLoader(const ogp::FineAtoms &a,
void void
DynareDerEvalLoader::load(int i, int iord, const int *vars, double res) DynareDerEvalLoader::load(int i, int iord, const int *vars, double res)
{ {
FSSparseTensor *t = md.get(Symmetry(iord)); FSSparseTensor *t = md.get(Symmetry{iord});
IntSequence s(iord, 0); IntSequence s(iord, 0);
for (int j = 0; j < iord; j++) for (int j = 0; j < iord; j++)
s[j] = atoms.get_pos_of_all(vars[j]); s[j] = atoms.get_pos_of_all(vars[j]);

4
dynare++/tl/cc/fs_tensor.hh
View File

@ -73,7 +73,7 @@ public:
Symmetry Symmetry
getSym() const getSym() const
{ {
return Symmetry(dimen()); return Symmetry{dimen()};
} }
int getOffset(const IntSequence &v) const override; int getOffset(const IntSequence &v) const override;
@ -117,7 +117,7 @@ public:
Symmetry Symmetry
getSym() const getSym() const
{ {
return Symmetry(dimen()); return Symmetry{dimen()};
} }
int getOffset(const IntSequence &v) const override; int getOffset(const IntSequence &v) const override;

2
dynare++/tl/cc/gs_tensor.cc
View File

@ -9,7 +9,7 @@
|@<|TensorDimens| class declaration@>| for details. */ |@<|TensorDimens| class declaration@>| for details. */
TensorDimens::TensorDimens(const IntSequence &ss, const IntSequence &coor) TensorDimens::TensorDimens(const IntSequence &ss, const IntSequence &coor)
: nvs(ss), : nvs(ss),
sym(ss.size(), ""), sym(ss.size()),
nvmax(coor.size(), 0) nvmax(coor.size(), 0)
{ {
TL_RAISE_IF(!coor.isSorted(), TL_RAISE_IF(!coor.isSorted(),

2
dynare++/tl/cc/gs_tensor.hh
View File

@ -64,7 +64,7 @@ public:
{ {
} }
TensorDimens(int nvar, int dimen) TensorDimens(int nvar, int dimen)
: nvs(1), sym(dimen), nvmax(dimen, nvar) : nvs(1), sym{dimen}, nvmax(dimen, nvar)
{ {
nvs[0] = nvar; nvs[0] = nvar;
} }

36
dynare++/tl/cc/int_sequence.hh
View File

@ -27,10 +27,16 @@
#include <vector> #include <vector>
#include <memory> #include <memory>
#include <algorithm> #include <algorithm>
#include <initializer_list>
/* The implementation of |IntSequence| is straightforward. It has a pointer /* The implementation of |IntSequence| is straightforward. It has a pointer
|data|, an |offset| integer indicating the beginning of the data relatively |data|, an |offset| integer indicating the beginning of the data relatively
to the pointer and a |length| of the sequence. */ to the pointer and a |length| of the sequence.
WARNING: IntSequence(n) and IntSequence{n} are not the same. The former
initializes a sequence of length n, while the latter constructs a sequence
of a single element equal to n. This is similar to the behaviour of
std::vector. */
class Symmetry; class Symmetry;
class IntSequence class IntSequence
@ -39,42 +45,52 @@ class IntSequence
int length; int length;
int offset{0}; int offset{0};
public: public:
/* We have a constructor allocating a given length of data, constructor // Constructor allocating a given length of (uninitialized) data
allocating and then initializing all members to a given number, a copy
constructor, a conversion from |vector<int>|, a subsequence
constructor, a constructor used for calculating implied symmetry from
a more general symmetry and one equivalence class (see |Symmetry|
class). Finally we have a constructor which unfolds a sequence with
respect to a given symmetry and constructor which inserts a given
number to the ordered sequence or given number to a given position. */
explicit IntSequence(int l) explicit IntSequence(int l)
: data{new int[l], [](int *arr) { delete[] arr; }}, length{l} : data{new int[l], [](int *arr) { delete[] arr; }}, length{l}
{ {
} }
// Constructor allocating and then initializing all members to a given number
IntSequence(int l, int n) IntSequence(int l, int n)
: data{new int[l], [](int *arr) { delete[] arr; }}, length{l} : data{new int[l], [](int *arr) { delete[] arr; }}, length{l}
{ {
std::fill_n(data.get(), length, n); std::fill_n(data.get(), length, n);
} }
/* Constructor using an initializer list (gives the contents of the
IntSequence, similarly to std::vector) */
IntSequence(std::initializer_list<int> init)
: data{new int[init.size()], [](int *arr) { delete[] arr; }},
length{static_cast<int>(init.size())}
{
std::copy(init.begin(), init.end(), data.get());
}
// Copy constructor
IntSequence(const IntSequence &s) IntSequence(const IntSequence &s)
: data{new int[s.length], [](int *arr) { delete[] arr; }}, length{s.length} : data{new int[s.length], [](int *arr) { delete[] arr; }}, length{s.length}
{ {
std::copy_n(s.data.get()+s.offset, length, data.get()); std::copy_n(s.data.get()+s.offset, length, data.get());
} }
// Move constructor
IntSequence(IntSequence &&s) = default; IntSequence(IntSequence &&s) = default;
// Subsequence constructor (which shares the data pointer)
IntSequence(IntSequence &s, int i1, int i2) IntSequence(IntSequence &s, int i1, int i2)
: data{s.data}, length{i2-i1}, offset{s.offset+i1} : data{s.data}, length{i2-i1}, offset{s.offset+i1}
{ {
} }
// Subsequence constructor (without pointer sharing)
IntSequence(const IntSequence &s, int i1, int i2) IntSequence(const IntSequence &s, int i1, int i2)
: data{new int[i2-i1], [](int *arr) { delete[] arr; }}, length{i2-i1} : data{new int[i2-i1], [](int *arr) { delete[] arr; }}, length{i2-i1}
{ {
std::copy_n(s.data.get()+s.offset+i1, length, data.get()); std::copy_n(s.data.get()+s.offset+i1, length, data.get());
} }
/* Constructor used for calculating implied symmetry from a more general
symmetry and one equivalence class */
IntSequence(const Symmetry &sy, const std::vector<int> &se); IntSequence(const Symmetry &sy, const std::vector<int> &se);
// Unfolds a given integer sequence with respect to a given symmetry
IntSequence(const Symmetry &sy, const IntSequence &se); IntSequence(const Symmetry &sy, const IntSequence &se);
// Inserts an element in an ordered sequence
IntSequence(int i, const IntSequence &s); IntSequence(int i, const IntSequence &s);
// Inserts an element at a given position
IntSequence(int i, const IntSequence &s, int pos); IntSequence(int i, const IntSequence &s, int pos);
const IntSequence &operator=(const IntSequence &s); const IntSequence &operator=(const IntSequence &s);

4
dynare++/tl/cc/rfs_tensor.hh
View File

@ -73,7 +73,7 @@ public:
Symmetry Symmetry
getSym() const getSym() const
{ {
return Symmetry(dimen()); return Symmetry{dimen()};
} }
}; };
@ -113,7 +113,7 @@ public:
Symmetry Symmetry
getSym() const getSym() const
{ {
return Symmetry(dimen()); return Symmetry{dimen()};
} }
}; };

2
dynare++/tl/cc/sparse_tensor.cc
View File

@ -116,7 +116,7 @@ SparseTensor::print() const
FSSparseTensor::FSSparseTensor(int d, int nvar, int r) FSSparseTensor::FSSparseTensor(int d, int nvar, int r)
: SparseTensor(d, r, FFSTensor::calcMaxOffset(nvar, d)), : SparseTensor(d, r, FFSTensor::calcMaxOffset(nvar, d)),
nv(nvar), sym(d) nv(nvar), sym{d}
{ {
} }

18
dynare++/tl/cc/stack_container.hh
View File

@ -284,8 +284,8 @@ public:
multAndAdd(int dim, const TensorContainer<FSSparseTensor> &c, multAndAdd(int dim, const TensorContainer<FSSparseTensor> &c,
FGSTensor &out) const FGSTensor &out) const
{ {
if (c.check(Symmetry(dim))) if (c.check(Symmetry{dim}))
multAndAdd(*(c.get(Symmetry(dim))), out); multAndAdd(*(c.get(Symmetry{dim})), out);
} }
void multAndAdd(const FSSparseTensor &t, FGSTensor &out) const; void multAndAdd(const FSSparseTensor &t, FGSTensor &out) const;
void multAndAdd(int dim, const FGSContainer &c, FGSTensor &out) const; void multAndAdd(int dim, const FGSContainer &c, FGSTensor &out) const;
@ -314,8 +314,8 @@ public:
multAndAdd(int dim, const TensorContainer<FSSparseTensor> &c, multAndAdd(int dim, const TensorContainer<FSSparseTensor> &c,
UGSTensor &out) const UGSTensor &out) const
{ {
if (c.check(Symmetry(dim))) if (c.check(Symmetry{dim}))
multAndAdd(*(c.get(Symmetry(dim))), out); multAndAdd(*(c.get(Symmetry{dim})), out);
} }
void multAndAdd(const FSSparseTensor &t, UGSTensor &out) const; void multAndAdd(const FSSparseTensor &t, UGSTensor &out) const;
void multAndAdd(int dim, const UGSContainer &c, UGSTensor &out) const; void multAndAdd(int dim, const UGSContainer &c, UGSTensor &out) const;
@ -370,12 +370,12 @@ public:
else else
return _Stype::matrix; return _Stype::matrix;
if (i == 2) if (i == 2)
if (s == Symmetry(1, 0, 0, 0)) if (s == Symmetry{1, 0, 0, 0})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;
if (i == 3) if (i == 3)
if (s == Symmetry(0, 1, 0, 0)) if (s == Symmetry{0, 1, 0, 0})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;
@ -446,19 +446,19 @@ public:
getType(int i, const Symmetry &s) const override getType(int i, const Symmetry &s) const override
{ {
if (i == 0) if (i == 0)
if (s[2] > 0 || s == Symmetry(0, 0, 0, 1)) if (s[2] > 0 || s == Symmetry{0, 0, 0, 1})
return _Stype::zero; return _Stype::zero;
else else
return _Stype::matrix; return _Stype::matrix;
if (i == 1) if (i == 1)
if (s == Symmetry(0, 0, 1, 0)) if (s == Symmetry{0, 0, 1, 0})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;
if (i == 2) if (i == 2)
return _Stype::zero; return _Stype::zero;
if (i == 3) if (i == 3)
if (s == Symmetry(0, 0, 0, 1)) if (s == Symmetry{0, 0, 0, 1})
return _Stype::unit; return _Stype::unit;
else else
return _Stype::zero; return _Stype::zero;

34
dynare++/tl/cc/symmetry.cc
View File

@ -3,7 +3,7 @@
#include "symmetry.hh" #include "symmetry.hh"
#include "permutation.hh" #include "permutation.hh"
#include <cstdio> #include <iostream>
/* Construct symmetry as numbers of successively equal items in the sequence. */ /* Construct symmetry as numbers of successively equal items in the sequence. */
@ -57,28 +57,18 @@ Symmetry::isFull() const
beginning as subordinal |symiterator|. */ beginning as subordinal |symiterator|. */
symiterator::symiterator(SymmetrySet &ss) symiterator::symiterator(SymmetrySet &ss)
: s(ss), subit(nullptr), subs(nullptr), end_flag(false) : s(ss), end_flag(false)
{ {
s.sym()[0] = 0; s.sym()[0] = 0;
if (s.size() == 2) if (s.size() == 2)
{ s.sym()[1] = s.dimen();
s.sym()[1] = s.dimen();
}
else else
{ {
subs = new SymmetrySet(s, s.dimen()); subs = std::make_unique<SymmetrySet>(s, s.dimen());
subit = new symiterator(*subs); subit = std::make_unique<symiterator>(*subs);
} }
} }
symiterator::~symiterator()
{
if (subit)
delete subit;
if (subs)
delete subs;
}
/* Here we move to the next symmetry. We do so only, if we are not at /* Here we move to the next symmetry. We do so only, if we are not at
the end. If length is 2, we increase lower index and decrease upper the end. If length is 2, we increase lower index and decrease upper
index, otherwise we increase the subordinal symmetry. If we got to the index, otherwise we increase the subordinal symmetry. If we got to the
@ -101,11 +91,9 @@ symiterator::operator++()
++(*subit); ++(*subit);
if (subit->isEnd()) if (subit->isEnd())
{ {
delete subit;
delete subs;
s.sym()[0]++; s.sym()[0]++;
subs = new SymmetrySet(s, s.dimen()-s.sym()[0]); subs = std::make_unique<SymmetrySet>(s, s.dimen()-s.sym()[0]);
subit = new symiterator(*subs); subit = std::make_unique<symiterator>(*subs);
} }
} }
if (s.sym()[0] == s.dimen()+1) if (s.sym()[0] == s.dimen()+1)
@ -117,9 +105,7 @@ symiterator::operator++()
InducedSymmetries::InducedSymmetries(const Equivalence &e, const Symmetry &s) InducedSymmetries::InducedSymmetries(const Equivalence &e, const Symmetry &s)
{ {
for (const auto & i : e) for (const auto & i : e)
{ emplace_back(s, i);
push_back(Symmetry(s, i));
}
} }
// |InducedSymmetries| permuted constructor code // |InducedSymmetries| permuted constructor code
@ -129,7 +115,7 @@ InducedSymmetries::InducedSymmetries(const Equivalence &e, const Permutation &p,
for (int i = 0; i < e.numClasses(); i++) for (int i = 0; i < e.numClasses(); i++)
{ {
auto it = e.find(p.getMap()[i]); auto it = e.find(p.getMap()[i]);
push_back(Symmetry(s, *it)); emplace_back(s, *it);
} }
} }
@ -138,7 +124,7 @@ InducedSymmetries::InducedSymmetries(const Equivalence &e, const Permutation &p,
void void
InducedSymmetries::print() const InducedSymmetries::print() const
{ {
printf("Induced symmetries: %lu\n", (unsigned long) size()); std::cout << "Induced symmetries: " << size() << std::endl;
for (unsigned int i = 0; i < size(); i++) for (unsigned int i = 0; i < size(); i++)
operator[](i).print(); operator[](i).print();
} }

68
dynare++/tl/cc/symmetry.hh
View File

@ -46,63 +46,49 @@
#include <list> #include <list>
#include <vector> #include <vector>
#include <initializer_list>
#include <utility>
#include <memory>
/* Clear. The method |isFull| returns true if and only if the symmetry /* Clear. The method |isFull| returns true if and only if the symmetry
allows for any permutation of indices. */ allows for any permutation of indices.
WARNING: Symmetry(n) and Symmetry{n} are not the same. The former
initializes a symmetry of n elements, while the latter is a full symmetry of
order n. This is similar to the behaviour of std::vector. */
class Symmetry : public IntSequence class Symmetry : public IntSequence
{ {
public: public:
/* We provide three constructors for symmetries of the form $y^n$, // Constructor allocating a given length of (zero-initialized) data
$y^nu^m$, $y^nu^m\sigma^k$. Also a copy constructor, and finally a Symmetry(int len)
constructor of implied symmetry for a symmetry and an equivalence
class. It is already implemented in |IntSequence| so we only call
appropriate constructor of |IntSequence|. We also provide the
subsymmetry, which takes the given length of symmetry from the end.
The last constructor constructs a symmetry from an integer sequence
(supposed to be ordered) as a symmetry counting successively equal
items. For instance the sequence $(a,a,a,b,c,c,d,d,d,d)$ produces
symmetry $(3,1,2,4)$. */
Symmetry(int len, const char *dummy)
: IntSequence(len, 0) : IntSequence(len, 0)
{ {
} }
Symmetry(int i1) /* Constructor using an initializer list, that gives the contents of the
: IntSequence(1, i1) Symmetry. Used for symmetries of the form $y^n$, $y^n u^m$, $y^nu^m\sigma^k$ */
Symmetry(std::initializer_list<int> init)
: IntSequence(std::move(init))
{ {
} }
Symmetry(int i1, int i2) // Copy constructor
: IntSequence(2)
{
operator[](0) = i1;
operator[](1) = i2;
}
Symmetry(int i1, int i2, int i3)
: IntSequence(3)
{
operator[](0) = i1;
operator[](1) = i2;
operator[](2) = i3;
}
Symmetry(int i1, int i2, int i3, int i4)
: IntSequence(4)
{
operator[](0) = i1;
operator[](1) = i2;
operator[](2) = i3;
operator[](3) = i4;
}
Symmetry(const Symmetry &s) = default; Symmetry(const Symmetry &s) = default;
// Move constructor
Symmetry(Symmetry &&s) = default; Symmetry(Symmetry &&s) = default;
// Constructor of implied symmetry for a symmetry and an equivalence class
Symmetry(const Symmetry &s, const OrdSequence &cl) Symmetry(const Symmetry &s, const OrdSequence &cl)
: IntSequence(s, cl.getData()) : IntSequence(s, cl.getData())
{ {
} }
/* Subsymmetry, which takes the given length of symmetry from the end (shares
data pointer) */
Symmetry(Symmetry &s, int len) Symmetry(Symmetry &s, int len)
: IntSequence(s, s.size()-len, s.size()) : IntSequence(s, s.size()-len, s.size())
{ {
} }
/* Constructs a symmetry from an integer sequence (supposed to be ordered) as
a symmetry counting successively equal items. For instance the sequence
$(a,a,a,b,c,c,d,d,d,d)$ produces symmetry $(3,1,2,4)$. */
Symmetry(const IntSequence &s); Symmetry(const IntSequence &s);
int int
@ -150,7 +136,7 @@ class SymmetrySet
int dim; int dim;
public: public:
SymmetrySet(int d, int length) SymmetrySet(int d, int length)
: run(length, ""), dim(d) : run(length), dim(d)
{ {
} }
SymmetrySet(SymmetrySet &s, int d) SymmetrySet(SymmetrySet &s, int d)
@ -184,7 +170,7 @@ public:
to the |SymmetrySet| only to know dimension and for access of its to the |SymmetrySet| only to know dimension and for access of its
symmetry storage. Further we have pointers to subordinal |symiterator| symmetry storage. Further we have pointers to subordinal |symiterator|
and its |SymmetrySet|. These are pointers, since the recursion ends at and its |SymmetrySet|. These are pointers, since the recursion ends at
length equal to 2, in which case these pointers are |NULL|. length equal to 2, in which case these pointers are uninitialized.
The constructor creates the iterator which initializes to the first The constructor creates the iterator which initializes to the first
symmetry (beginning). */ symmetry (beginning). */
@ -192,12 +178,12 @@ public:
class symiterator class symiterator
{ {
SymmetrySet &s; SymmetrySet &s;
symiterator *subit; std::unique_ptr<symiterator> subit;
SymmetrySet *subs; std::unique_ptr<SymmetrySet> subs;
bool end_flag; bool end_flag;
public: public:
symiterator(SymmetrySet &ss); symiterator(SymmetrySet &ss);
~symiterator(); ~symiterator() = default;
symiterator &operator++(); symiterator &operator++();
bool bool
isEnd() const isEnd() const

54
dynare++/tl/cc/t_polynomial.hh
View File

@ -168,14 +168,14 @@ public:
const _Stype &xpow = pwp.getNext((const _Stype *) nullptr); const _Stype &xpow = pwp.getNext((const _Stype *) nullptr);
for (int j = 0; j <= tp.maxdim-i; j++) for (int j = 0; j <= tp.maxdim-i; j++)
{ {
if (tp.check(Symmetry(i+j))) if (tp.check(Symmetry{i+j}))
{ {
// initialize |ten| of dimension |j| // initialize |ten| of dimension |j|
/* The pointer |ten| is either a new tensor or got from |this| container. */ /* The pointer |ten| is either a new tensor or got from |this| container. */
_Ttype *ten; _Ttype *ten;
if (_Tparent::check(Symmetry(j))) if (_Tparent::check(Symmetry{j}))
{ {
ten = _Tparent::get(Symmetry(j)); ten = _Tparent::get(Symmetry{j});
} }
else else
{ {
@ -184,9 +184,9 @@ public:
insert(ten); insert(ten);
} }
Symmetry sym(i, j); Symmetry sym{i, j};
IntSequence coor(sym, pp); IntSequence coor(sym, pp);
_TGStype slice(*(tp.get(Symmetry(i+j))), ss, coor, TensorDimens(sym, ss)); _TGStype slice(*(tp.get(Symmetry{i+j})), ss, coor, TensorDimens(sym, ss));
slice.mult(Tensor::noverk(i+j, j)); slice.mult(Tensor::noverk(i+j, j));
_TGStype tmp(*ten); _TGStype tmp(*ten);
slice.contractAndAdd(0, tmp, xpow); slice.contractAndAdd(0, tmp, xpow);
@ -200,15 +200,15 @@ public:
simple addition. */ simple addition. */
for (int j = 0; j <= tp.maxdim; j++) for (int j = 0; j <= tp.maxdim; j++)
{ {
if (tp.check(Symmetry(j))) if (tp.check(Symmetry{j}))
{ {
// initialize |ten| of dimension |j| // initialize |ten| of dimension |j|
/* Same code as above */ /* Same code as above */
_Ttype *ten; _Ttype *ten;
if (_Tparent::check(Symmetry(j))) if (_Tparent::check(Symmetry{j}))
{ {
ten = _Tparent::get(Symmetry(j)); ten = _Tparent::get(Symmetry{j});
} }
else else
{ {
@ -217,9 +217,9 @@ public:
insert(ten); insert(ten);
} }
Symmetry sym(0, j); Symmetry sym{0, j};
IntSequence coor(sym, pp); IntSequence coor(sym, pp);
_TGStype slice(*(tp.get(Symmetry(j))), ss, coor, TensorDimens(sym, ss)); _TGStype slice(*(tp.get(Symmetry{j})), ss, coor, TensorDimens(sym, ss));
ten->add(1.0, slice); ten->add(1.0, slice);
} }
} }
@ -247,8 +247,8 @@ public:
void void
evalTrad(Vector &out, const ConstVector &v) const evalTrad(Vector &out, const ConstVector &v) const
{ {
if (_Tparent::check(Symmetry(0))) if (_Tparent::check(Symmetry{0}))
out = _Tparent::get(Symmetry(0))->getData(); out = _Tparent::get(Symmetry{0})->getData();
else else
out.zeros(); out.zeros();
@ -256,7 +256,7 @@ public:
for (int d = 1; d <= maxdim; d++) for (int d = 1; d <= maxdim; d++)
{ {
const _Stype &p = pp.getNext((const _Stype *) nullptr); const _Stype &p = pp.getNext((const _Stype *) nullptr);
Symmetry cs(d); Symmetry cs{d};
if (_Tparent::check(cs)) if (_Tparent::check(cs))
{ {
const _Ttype *t = _Tparent::get(cs); const _Ttype *t = _Tparent::get(cs);
@ -271,8 +271,8 @@ public:
void void
evalHorner(Vector &out, const ConstVector &v) const evalHorner(Vector &out, const ConstVector &v) const
{ {
if (_Tparent::check(Symmetry(0))) if (_Tparent::check(Symmetry{0}))
out = _Tparent::get(Symmetry(0))->getData(); out = _Tparent::get(Symmetry{0})->getData();
else else
out.zeros(); out.zeros();
@ -281,12 +281,12 @@ public:
_Ttype *last; _Ttype *last;
if (maxdim == 1) if (maxdim == 1)
last = new _Ttype(*(_Tparent::get(Symmetry(1)))); last = new _Ttype(*(_Tparent::get(Symmetry{1})));
else else
last = new _Ttype(*(_Tparent::get(Symmetry(maxdim))), v); last = new _Ttype(*(_Tparent::get(Symmetry{maxdim})), v);
for (int d = maxdim-1; d >= 1; d--) for (int d = maxdim-1; d >= 1; d--)
{ {
Symmetry cs(d); Symmetry cs{d};
if (_Tparent::check(cs)) if (_Tparent::check(cs))
{ {
const _Ttype *nt = _Tparent::get(cs); const _Ttype *nt = _Tparent::get(cs);
@ -337,9 +337,9 @@ public:
{ {
for (int d = 1; d <= maxdim; d++) for (int d = 1; d <= maxdim; d++)
{ {
if (_Tparent::check(Symmetry(d))) if (_Tparent::check(Symmetry{d}))
{ {
_Ttype *ten = _Tparent::get(Symmetry(d)); _Ttype *ten = _Tparent::get(Symmetry{d});
ten->mult((double) std::max((d-k), 0)); ten->mult((double) std::max((d-k), 0));
} }
} }
@ -367,14 +367,14 @@ public:
auto *res = new _Ttype(nrows(), nvars(), s); auto *res = new _Ttype(nrows(), nvars(), s);
res->zeros(); res->zeros();
if (_Tparent::check(Symmetry(s))) if (_Tparent::check(Symmetry{s}))
res->add(1.0, *(_Tparent::get(Symmetry(s)))); res->add(1.0, *(_Tparent::get(Symmetry{s})));
for (int d = s+1; d <= maxdim; d++) for (int d = s+1; d <= maxdim; d++)
{ {
if (_Tparent::check(Symmetry(d))) if (_Tparent::check(Symmetry{d}))
{ {
const _Ttype &ltmp = *(_Tparent::get(Symmetry(d))); const _Ttype &ltmp = *(_Tparent::get(Symmetry{d}));
auto *last = new _Ttype(ltmp); auto *last = new _Ttype(ltmp);
for (int j = 0; j < d - s; j++) for (int j = 0; j < d - s; j++)
{ {
@ -474,12 +474,12 @@ public:
for (Tensor::index i = dum.begin(); i != dum.end(); ++i) for (Tensor::index i = dum.begin(); i != dum.end(); ++i)
{ {
int d = i.getCoor().sum(); int d = i.getCoor().sum();
Symmetry symrun(_Ttype::dimen()-d, d); Symmetry symrun{_Ttype::dimen()-d, d};
_TGStype dumgs(0, TensorDimens(symrun, dumnvs)); _TGStype dumgs(0, TensorDimens(symrun, dumnvs));
if (pol.check(Symmetry(d))) if (pol.check(Symmetry{d}))
{ {
TwoDMatrix subt(*this, offset, dumgs.ncols()); TwoDMatrix subt(*this, offset, dumgs.ncols());
subt.add(1.0, *(pol.get(Symmetry(d)))); subt.add(1.0, *(pol.get(Symmetry{d})));
} }
offset += dumgs.ncols(); offset += dumgs.ncols();
} }

4
dynare++/tl/testing/factory.hh
View File

@ -51,7 +51,7 @@ public:
if (symnum == 1) if (symnum == 1)
{ {
// full symmetry // full symmetry
Symmetry sym(dim); Symmetry sym{dim};
auto *t = make<_Ttype>(r, sym, nvs); auto *t = make<_Ttype>(r, sym, nvs);
res->insert(t); res->insert(t);
} }
@ -60,7 +60,7 @@ public:
// general symmetry // general symmetry
for (int i = 0; i <= dim; i++) for (int i = 0; i <= dim; i++)
{ {
Symmetry sym(i, dim-i); Symmetry sym{i, dim-i};
auto *t = make<_Ttype>(r, sym, nvs); auto *t = make<_Ttype>(r, sym, nvs);
res->insert(t); res->insert(t);
} }

6
dynare++/tl/testing/monoms.cc
View File

@ -121,7 +121,7 @@ FGSTensor *
Monom1Vector::deriv(int dim) const Monom1Vector::deriv(int dim) const
{ {
FGSTensor *res FGSTensor *res
= new FGSTensor(len, TensorDimens(Symmetry(dim), IntSequence(1, nx))); = new FGSTensor(len, TensorDimens(Symmetry{dim}, IntSequence(1, nx)));
for (Tensor::index it = res->begin(); it != res->end(); ++it) for (Tensor::index it = res->begin(); it != res->end(); ++it)
{ {
Vector outcol{res->getCol(*it)}; Vector outcol{res->getCol(*it)};
@ -229,7 +229,7 @@ Monom2Vector::deriv(int maxdim) const
for (int ydim = 0; ydim <= dim; ydim++) for (int ydim = 0; ydim <= dim; ydim++)
{ {
int udim = dim - ydim; int udim = dim - ydim;
Symmetry s(ydim, udim); Symmetry s{ydim, udim};
res->insert(deriv(s)); res->insert(deriv(s));
} }
} }
@ -409,7 +409,7 @@ Monom4Vector::deriv(int dim) const
FFSTensor dummy(0, nx1+nx2+nx3+nx4, dim); FFSTensor dummy(0, nx1+nx2+nx3+nx4, dim);
for (Tensor::index run = dummy.begin(); run != dummy.end(); ++run) for (Tensor::index run = dummy.begin(); run != dummy.end(); ++run)
{ {
Symmetry ind_sym(0, 0, 0, 0); Symmetry ind_sym{0, 0, 0, 0};
IntSequence ind(run.getCoor()); IntSequence ind(run.getCoor());
for (int i = 0; i < ind.size(); i++) for (int i = 0; i < ind.size(); i++)
{ {

40
dynare++/tl/testing/tests.cc
View File

@ -226,7 +226,7 @@ TestRunnable::dense_prod(const Symmetry &bsym, const IntSequence &bnvs,
FGSContainer *cont FGSContainer *cont
= f.makeCont<FGSTensor, FGSContainer>(hnv, bnvs, bsym.dimen()-hdim+1); = f.makeCont<FGSTensor, FGSContainer>(hnv, bnvs, bsym.dimen()-hdim+1);
auto *fh auto *fh
= f.make<FGSTensor>(rows, Symmetry(hdim), IntSequence(1, hnv)); = f.make<FGSTensor>(rows, Symmetry{hdim}, IntSequence(1, hnv));
UGSTensor uh(*fh); UGSTensor uh(*fh);
FGSTensor fb(rows, TensorDimens(bsym, bnvs)); FGSTensor fb(rows, TensorDimens(bsym, bnvs));
fb.getData().zeros(); fb.getData().zeros();
@ -275,7 +275,7 @@ TestRunnable::folded_monomial(int ng, int nx, int ny, int nu, int dim)
double maxnorm = 0; double maxnorm = 0;
for (int ydim = 0; ydim <= dim; ydim++) for (int ydim = 0; ydim <= dim; ydim++)
{ {
Symmetry s(ydim, dim-ydim); Symmetry s{ydim, dim-ydim};
printf("\tSymmetry: "); s.print(); printf("\tSymmetry: "); s.print();
FGSTensor res(ng, TensorDimens(s, nvs)); FGSTensor res(ng, TensorDimens(s, nvs));
res.getData().zeros(); res.getData().zeros();
@ -314,7 +314,7 @@ TestRunnable::unfolded_monomial(int ng, int nx, int ny, int nu, int dim)
double maxnorm = 0; double maxnorm = 0;
for (int ydim = 0; ydim <= dim; ydim++) for (int ydim = 0; ydim <= dim; ydim++)
{ {
Symmetry s(ydim, dim-ydim); Symmetry s{ydim, dim-ydim};
printf("\tSymmetry: "); s.print(); printf("\tSymmetry: "); s.print();
UGSTensor res(ng, TensorDimens(s, nvs)); UGSTensor res(ng, TensorDimens(s, nvs));
res.getData().zeros(); res.getData().zeros();
@ -602,7 +602,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3); Symmetry s{2, 3};
IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2; IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2;
return index_forward<FGSTensor>(s, nvs); return index_forward<FGSTensor>(s, nvs);
} }
@ -618,7 +618,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3); Symmetry s{2, 3};
IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2; IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2;
return index_forward<UGSTensor>(s, nvs); return index_forward<UGSTensor>(s, nvs);
} }
@ -634,7 +634,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2;
return index_forward<FGSTensor>(s, nvs); return index_forward<FGSTensor>(s, nvs);
} }
@ -650,7 +650,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2;
return index_forward<UGSTensor>(s, nvs); return index_forward<UGSTensor>(s, nvs);
} }
@ -666,7 +666,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(1, 1, 3); Symmetry s{1, 1, 3};
IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2;
return index_backward<FGSTensor>(s, nvs); return index_backward<FGSTensor>(s, nvs);
} }
@ -682,7 +682,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 4; nvs[1] = 2; nvs[2] = 4; IntSequence nvs(3); nvs[0] = 4; nvs[1] = 2; nvs[2] = 4;
return index_backward<FGSTensor>(s, nvs); return index_backward<FGSTensor>(s, nvs);
} }
@ -698,7 +698,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(1, 1, 3); Symmetry s{1, 1, 3};
IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2;
return index_backward<UGSTensor>(s, nvs); return index_backward<UGSTensor>(s, nvs);
} }
@ -714,7 +714,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 4; nvs[1] = 2; nvs[2] = 4; IntSequence nvs(3); nvs[0] = 4; nvs[1] = 2; nvs[2] = 4;
return index_backward<UGSTensor>(s, nvs); return index_backward<UGSTensor>(s, nvs);
} }
@ -730,7 +730,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3); Symmetry s{2, 3};
IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2; IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2;
return index_offset<FGSTensor>(s, nvs); return index_offset<FGSTensor>(s, nvs);
} }
@ -746,7 +746,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3); Symmetry s{2, 3};
IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2; IntSequence nvs(2); nvs[0] = 4; nvs[1] = 2;
return index_offset<UGSTensor>(s, nvs); return index_offset<UGSTensor>(s, nvs);
} }
@ -762,7 +762,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2;
return index_offset<FGSTensor>(s, nvs); return index_offset<FGSTensor>(s, nvs);
} }
@ -778,7 +778,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 3, 2); Symmetry s{2, 3, 2};
IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 5; nvs[1] = 2; nvs[2] = 2;
return index_offset<UGSTensor>(s, nvs); return index_offset<UGSTensor>(s, nvs);
} }
@ -808,7 +808,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(1, 2, 2); Symmetry s{1, 2, 2};
IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2; IntSequence nvs(3); nvs[0] = 3; nvs[1] = 3; nvs[2] = 2;
return gs_fold_unfold(5, s, nvs); return gs_fold_unfold(5, s, nvs);
} }
@ -838,7 +838,7 @@ public:
bool bool
run() const override run() const override
{ {
Symmetry s(2, 1, 2); Symmetry s{2, 1, 2};
IntSequence nvs(3); nvs[0] = 6; nvs[1] = 2; nvs[2] = 6; IntSequence nvs(3); nvs[0] = 6; nvs[1] = 2; nvs[2] = 6;
return gs_fold_unfold(5, s, nvs); return gs_fold_unfold(5, s, nvs);
} }
@ -883,7 +883,7 @@ public:
run() const override run() const override
{ {
IntSequence bnvs(2); bnvs[0] = 3; bnvs[1] = 2; IntSequence bnvs(2); bnvs[0] = 3; bnvs[1] = 2;
return dense_prod(Symmetry(1, 2), bnvs, 2, 3, 2); return dense_prod(Symmetry{1, 2}, bnvs, 2, 3, 2);
} }
}; };
@ -898,7 +898,7 @@ public:
run() const override run() const override
{ {
IntSequence bnvs(2); bnvs[0] = 10; bnvs[1] = 7; IntSequence bnvs(2); bnvs[0] = 10; bnvs[1] = 7;
return dense_prod(Symmetry(2, 3), bnvs, 3, 15, 10); return dense_prod(Symmetry{2, 3}, bnvs, 3, 15, 10);
} }
}; };
@ -913,7 +913,7 @@ public:
run() const override run() const override
{ {
IntSequence bnvs(2); bnvs[0] = 13; bnvs[1] = 11; IntSequence bnvs(2); bnvs[0] = 13; bnvs[1] = 11;
return dense_prod(Symmetry(3, 2), bnvs, 3, 20, 20); return dense_prod(Symmetry{3, 2}, bnvs, 3, 20, 20);
} }
}; };

2
mex/sources/k_order_perturbation/k_ord_dynare.cc
View File

@ -248,7 +248,7 @@ KordpDynare::populateDerivativesContainer(const TwoDMatrix &g, int ord, const st
} }
// md container // md container
md.remove(Symmetry(ord)); md.remove(Symmetry{ord});
md.insert(mdTi); md.insert(mdTi);
// No need to delete mdTi, it will be deleted by TensorContainer destructor // No need to delete mdTi, it will be deleted by TensorContainer destructor
} }

24
mex/sources/k_order_perturbation/k_order_perturbation.cc
View File

@ -294,18 +294,18 @@ extern "C" {
for (int i = 0; i < 12; ++i) for (int i = 0; i < 12; ++i)
c_fieldnames[i] = fieldnames[i].c_str(); c_fieldnames[i] = fieldnames[i].c_str();
plhs[ii] = mxCreateStructMatrix(1, 1, 12, c_fieldnames); plhs[ii] = mxCreateStructMatrix(1, 1, 12, c_fieldnames);
copy_derivatives(plhs[ii], Symmetry(1, 0, 0, 0), derivs, "gy"); copy_derivatives(plhs[ii], Symmetry{1, 0, 0, 0}, derivs, "gy");
copy_derivatives(plhs[ii], Symmetry(0, 1, 0, 0), derivs, "gu"); copy_derivatives(plhs[ii], Symmetry{0, 1, 0, 0}, derivs, "gu");
copy_derivatives(plhs[ii], Symmetry(2, 0, 0, 0), derivs, "gyy"); copy_derivatives(plhs[ii], Symmetry{2, 0, 0, 0}, derivs, "gyy");
copy_derivatives(plhs[ii], Symmetry(0, 2, 0, 0), derivs, "guu"); copy_derivatives(plhs[ii], Symmetry{0, 2, 0, 0}, derivs, "guu");
copy_derivatives(plhs[ii], Symmetry(1, 1, 0, 0), derivs, "gyu"); copy_derivatives(plhs[ii], Symmetry{1, 1, 0, 0}, derivs, "gyu");
copy_derivatives(plhs[ii], Symmetry(0, 0, 0, 2), derivs, "gss"); copy_derivatives(plhs[ii], Symmetry{0, 0, 0, 2}, derivs, "gss");
copy_derivatives(plhs[ii], Symmetry(3, 0, 0, 0), derivs, "gyyy"); copy_derivatives(plhs[ii], Symmetry{3, 0, 0, 0}, derivs, "gyyy");
copy_derivatives(plhs[ii], Symmetry(0, 3, 0, 0), derivs, "guuu"); copy_derivatives(plhs[ii], Symmetry{0, 3, 0, 0}, derivs, "guuu");
copy_derivatives(plhs[ii], Symmetry(2, 1, 0, 0), derivs, "gyyu"); copy_derivatives(plhs[ii], Symmetry{2, 1, 0, 0}, derivs, "gyyu");
copy_derivatives(plhs[ii], Symmetry(1, 2, 0, 0), derivs, "gyuu"); copy_derivatives(plhs[ii], Symmetry{1, 2, 0, 0}, derivs, "gyuu");
copy_derivatives(plhs[ii], Symmetry(1, 0, 0, 2), derivs, "gyss"); copy_derivatives(plhs[ii], Symmetry{1, 0, 0, 2}, derivs, "gyss");
copy_derivatives(plhs[ii], Symmetry(0, 1, 0, 2), derivs, "guss"); copy_derivatives(plhs[ii], Symmetry{0, 1, 0, 2}, derivs, "guss");
} }
} }
} }