diff --git a/.clang-format b/.clang-format
index ff37c38a4c657078c863bd6afcb5c0702b32e1b9..b2ce7c0e6df0152c3ab0697053d14af385407722 100644
--- a/.clang-format
+++ b/.clang-format
@@ -19,6 +19,7 @@ BreakInheritanceList: AfterColon
 Cpp11BracedListStyle: true
 DeriveLineEnding: false
 IndentPPDirectives: AfterHash
+PackConstructorInitializers: NextLine
 PPIndentWidth: 1
 PointerAlignment: Left
 SpaceAfterTemplateKeyword: false
diff --git a/mex/sources/bytecode/Evaluate.cc b/mex/sources/bytecode/Evaluate.cc
index 1adf1cf5faacce227a11da0daa29a8ed6543306f..a09ace9d14503cb5bb4b1f8f8df01ab1884da5d3 100644
--- a/mex/sources/bytecode/Evaluate.cc
+++ b/mex/sources/bytecode/Evaluate.cc
@@ -32,8 +32,7 @@
 
 Evaluate::Evaluate(const filesystem::path& codfile, bool steady_state_arg,
                    const BasicSymbolTable& symbol_table_arg) :
-    symbol_table {symbol_table_arg},
-    steady_state {steady_state_arg}
+    symbol_table {symbol_table_arg}, steady_state {steady_state_arg}
 {
   ifstream CompiledCode {codfile, ios::in | ios::binary | ios::ate};
   if (!CompiledCode.is_open())
diff --git a/mex/sources/bytecode/Interpreter.cc b/mex/sources/bytecode/Interpreter.cc
index a9edd5930c90dde33a4ec71e204624b0448b0f73..151cb6d0fb388658674da35cc145c25d7c86300b 100644
--- a/mex/sources/bytecode/Interpreter.cc
+++ b/mex/sources/bytecode/Interpreter.cc
@@ -42,9 +42,14 @@ Interpreter::Interpreter(Evaluate& evaluator_arg, double* params_arg, double* y_
                          const BasicSymbolTable& symbol_table_arg, int verbosity_arg) :
     symbol_table {symbol_table_arg},
     steady_state {steady_state_arg},
-    block_decomposed {block_decomposed_arg}, evaluator {evaluator_arg},
-    minimal_solving_periods {minimal_solving_periods_arg}, y_size {y_size_arg}, y_kmin {y_kmin_arg},
-    y_kmax {y_kmax_arg}, periods {periods_arg}, verbosity {verbosity_arg}
+    block_decomposed {block_decomposed_arg},
+    evaluator {evaluator_arg},
+    minimal_solving_periods {minimal_solving_periods_arg},
+    y_size {y_size_arg},
+    y_kmin {y_kmin_arg},
+    y_kmax {y_kmax_arg},
+    periods {periods_arg},
+    verbosity {verbosity_arg}
 {
   pivotva = nullptr;
   mem_mngr.init_Mem();
diff --git a/mex/sources/k_order_welfare/approximation_welfare.cc b/mex/sources/k_order_welfare/approximation_welfare.cc
index c33855ff8ee0230b89d2f28fcc9985dd3cfa13e6..bb56d6e49b3705750f82bf3736acac8d7f9b9309 100644
--- a/mex/sources/k_order_welfare/approximation_welfare.cc
+++ b/mex/sources/k_order_welfare/approximation_welfare.cc
@@ -27,8 +27,8 @@ ApproximationWelfare::ApproximationWelfare(KordwDynare& w, double discount_facto
                                            const FGSContainer& rule_ders_arg,
                                            const FGSContainer& rule_ders_s_arg, Journal& j) :
     welfare {w},
-    discount_factor(discount_factor_arg), nvs {welfare.getModel().nys(), welfare.getModel().nexog(),
-                                               welfare.getModel().nexog(), 1},
+    discount_factor(discount_factor_arg),
+    nvs {welfare.getModel().nys(), welfare.getModel().nexog(), welfare.getModel().nexog(), 1},
     journal {j}
 {
   rule_ders = std::make_unique<FGSContainer>(rule_ders_arg);
diff --git a/mex/sources/k_order_welfare/k_ord_objective.cc b/mex/sources/k_order_welfare/k_ord_objective.cc
index e432767b0d0726daf4096ac7cbb1b8fc54d8db47..d3b2ffcf5299ac0a43bce3d51064260f811690ee 100644
--- a/mex/sources/k_order_welfare/k_ord_objective.cc
+++ b/mex/sources/k_order_welfare/k_ord_objective.cc
@@ -27,8 +27,12 @@ KordwDynare::KordwDynare(KordpDynare& m, ConstVector& NNZD_arg, Journal& jr, Vec
                          std::unique_ptr<ObjectiveAC> objectiveFile_arg,
                          const std::vector<int>& dr_order) :
     model {m},
-    NNZD {NNZD_arg}, journal {jr}, params {inParams},
-    resid(1), ud {1}, objectiveFile {std::move(objectiveFile_arg)}
+    NNZD {NNZD_arg},
+    journal {jr},
+    params {inParams},
+    resid(1),
+    ud {1},
+    objectiveFile {std::move(objectiveFile_arg)}
 {
   dynppToDyn = dr_order;
   dynToDynpp.resize(model.ny());
@@ -309,11 +313,30 @@ KOrderWelfare::KOrderWelfare(int num_stat, int num_pred, int num_both, int num_f
                              const FGSContainer& g_arg, const FGSContainer& gs_arg,
                              const TwoDMatrix& v, Journal& jr) :
     ypart(num_stat, num_pred, num_both, num_forw),
-    ny(ypart.ny()), nu(nu), maxk(ucont.getMaxDim()),
-    order(ord), discount_factor {discount_factor}, nvs {ypart.nys(), nu, nu, 1}, _uU(4), _fU(4),
-    _uW(4), _fW(4), _uWrond(4), _fWrond(4), _ug(4), _fg(g_arg), _ugs(4), _fgs(gs_arg),
-    _uXstack(&_ug, ny), _fXstack(&_fg, ny), _uGstack(&_ugs, ypart.nys(), nu),
-    _fGstack(&_fgs, ypart.nys(), nu), _um(maxk, v), _fm(_um), u(ucont), journal(jr)
+    ny(ypart.ny()),
+    nu(nu),
+    maxk(ucont.getMaxDim()),
+    order(ord),
+    discount_factor {discount_factor},
+    nvs {ypart.nys(), nu, nu, 1},
+    _uU(4),
+    _fU(4),
+    _uW(4),
+    _fW(4),
+    _uWrond(4),
+    _fWrond(4),
+    _ug(4),
+    _fg(g_arg),
+    _ugs(4),
+    _fgs(gs_arg),
+    _uXstack(&_ug, ny),
+    _fXstack(&_fg, ny),
+    _uGstack(&_ugs, ypart.nys(), nu),
+    _fGstack(&_fgs, ypart.nys(), nu),
+    _um(maxk, v),
+    _fm(_um),
+    u(ucont),
+    journal(jr)
 {
   KORD_RAISE_IF(v.ncols() != nu, "Wrong number of columns of Vcov in KOrderWelfare constructor");
   KORD_RAISE_IF(nu != v.nrows(), "Wrong number of rows of Vcov in KOrderWelfare constructor");
diff --git a/mex/sources/libkorder/k_ord_dynare.cc b/mex/sources/libkorder/k_ord_dynare.cc
index 6d231442d503c6cfff31fe97f604aa8619acbde2..a2cae057244f1906b9294e099d2c281910c95db8 100644
--- a/mex/sources/libkorder/k_ord_dynare.cc
+++ b/mex/sources/libkorder/k_ord_dynare.cc
@@ -33,11 +33,28 @@ KordpDynare::KordpDynare(const std::vector<std::string>& endo, const std::vector
                          std::unique_ptr<DynamicModelAC> dynamicModelFile_arg,
                          const std::vector<int>& dr_order, const ConstTwoDMatrix& llincidence) :
     nStat {nstat},
-    nBoth {nboth}, nPred {npred}, nForw {nforw}, nExog {nexog}, nPar {npar}, nYs {npred + nboth},
-    nYss {nboth + nforw}, nY {nstat + npred + nboth + nforw}, nJcols {nExog + nY + nYs + nYss},
-    NNZD {nnzd}, nSteps {nsteps}, nOrder {norder}, journal {jr}, ySteady {ysteady},
-    params {inParams}, vCov {vcov}, md {1}, dnl {endo}, denl {exo}, dsnl {*this, dnl, denl},
-    ll_Incidence {llincidence}, dynamicModelFile {std::move(dynamicModelFile_arg)}
+    nBoth {nboth},
+    nPred {npred},
+    nForw {nforw},
+    nExog {nexog},
+    nPar {npar},
+    nYs {npred + nboth},
+    nYss {nboth + nforw},
+    nY {nstat + npred + nboth + nforw},
+    nJcols {nExog + nY + nYs + nYss},
+    NNZD {nnzd},
+    nSteps {nsteps},
+    nOrder {norder},
+    journal {jr},
+    ySteady {ysteady},
+    params {inParams},
+    vCov {vcov},
+    md {1},
+    dnl {endo},
+    denl {exo},
+    dsnl {*this, dnl, denl},
+    ll_Incidence {llincidence},
+    dynamicModelFile {std::move(dynamicModelFile_arg)}
 {
   computeJacobianPermutation(dr_order);
 }
diff --git a/mex/sources/libkorder/kord/approximation.cc b/mex/sources/libkorder/kord/approximation.cc
index e375c513211ecf49d16811a7e9d3596350b21cb1..49dbf4b9979c0fab4d724af4236749f74a16f6e7 100644
--- a/mex/sources/libkorder/kord/approximation.cc
+++ b/mex/sources/libkorder/kord/approximation.cc
@@ -52,10 +52,14 @@ ZAuxContainer::getType(int i, const Symmetry& s) const
 Approximation::Approximation(DynamicModel& m, Journal& j, int ns, bool dr_centr, bool pruned_dr,
                              double qz_crit) :
     model(m),
-    journal(j), ypart(model.nstat(), model.npred(), model.nboth(), model.nforw()),
-    mom(UNormalMoments(model.order(), model.getVcov())), nvs {ypart.nys(), model.nexog(),
-                                                              model.nexog(), 1},
-    steps(ns), dr_centralize(dr_centr), pruning(pruned_dr), qz_criterium(qz_crit),
+    journal(j),
+    ypart(model.nstat(), model.npred(), model.nboth(), model.nforw()),
+    mom(UNormalMoments(model.order(), model.getVcov())),
+    nvs {ypart.nys(), model.nexog(), model.nexog(), 1},
+    steps(ns),
+    dr_centralize(dr_centr),
+    pruning(pruned_dr),
+    qz_criterium(qz_crit),
     ss(ypart.ny(), steps + 1)
 {
   ss.nans();
diff --git a/mex/sources/libkorder/kord/decision_rule.hh b/mex/sources/libkorder/kord/decision_rule.hh
index 80ce73016ac1a35a85afefeed6d73bef9ca0edc3..084a4c01f6ad79014070e773b0c5ab1daa158568 100644
--- a/mex/sources/libkorder/kord/decision_rule.hh
+++ b/mex/sources/libkorder/kord/decision_rule.hh
@@ -127,15 +127,13 @@ public:
   }
   DecisionRuleImpl(const _Tg& g, const PartitionY& yp, int nuu, const ConstVector& ys,
                    double sigma) :
-      ctraits<t>::Tpol(yp.ny(), yp.nys() + nuu),
-      ysteady(ys), ypart(yp), nu(nuu)
+      ctraits<t>::Tpol(yp.ny(), yp.nys() + nuu), ysteady(ys), ypart(yp), nu(nuu)
   {
     fillTensors(g, sigma);
   }
   DecisionRuleImpl(const _Tg& g, const PartitionY& yp, int nuu, const ConstVector& ys, double sigma,
                    bool pruning) :
-      ctraits<t>::Tpol(yp.ny(), yp.nys() + nuu),
-      ysteady(ys), ypart(yp), nu(nuu)
+      ctraits<t>::Tpol(yp.ny(), yp.nys() + nuu), ysteady(ys), ypart(yp), nu(nuu)
   {
     if (pruning)
       fillTensorsPruning(g);
@@ -149,7 +147,9 @@ public:
     fillTensors(W, nys);
   }
   DecisionRuleImpl(const DecisionRuleImpl<t>& dr, const ConstVector& fixpoint) :
-      ctraits<t>::Tpol(dr.ypart.ny(), dr.ypart.nys() + dr.nu), ysteady(fixpoint), ypart(dr.ypart),
+      ctraits<t>::Tpol(dr.ypart.ny(), dr.ypart.nys() + dr.nu),
+      ysteady(fixpoint),
+      ypart(dr.ypart),
       nu(dr.nu)
   {
     centralize(dr);
diff --git a/mex/sources/libkorder/kord/first_order.hh b/mex/sources/libkorder/kord/first_order.hh
index 37ec26ef9e74188802da484610ed87c030e84eb7..d22d1a09895edf5fc19863b73379b19c77a311ef 100644
--- a/mex/sources/libkorder/kord/first_order.hh
+++ b/mex/sources/libkorder/kord/first_order.hh
@@ -61,8 +61,13 @@ public:
   FirstOrder(int num_stat, int num_pred, int num_both, int num_forw, int num_u,
              const FSSparseTensor& f, Journal& jr, double qz_crit) :
       ypart(num_stat, num_pred, num_both, num_forw),
-      nu(num_u), gy(ypart.ny(), ypart.nys()), gu(ypart.ny(), nu), alphar(ypart.ny() + ypart.nboth),
-      alphai(ypart.ny() + ypart.nboth), beta(ypart.ny() + ypart.nboth), qz_criterium(qz_crit),
+      nu(num_u),
+      gy(ypart.ny(), ypart.nys()),
+      gu(ypart.ny(), nu),
+      alphar(ypart.ny() + ypart.nboth),
+      alphai(ypart.ny() + ypart.nboth),
+      beta(ypart.ny() + ypart.nboth),
+      qz_criterium(qz_crit),
       journal(jr)
   {
     solve(FFSTensor(f));
diff --git a/mex/sources/libkorder/kord/korder.cc b/mex/sources/libkorder/kord/korder.cc
index 040e30050f98ed2d1c5c2ce00718477094a63ffc..adbcaeaba7f964c946759aba43b287ea15818925 100644
--- a/mex/sources/libkorder/kord/korder.cc
+++ b/mex/sources/libkorder/kord/korder.cc
@@ -293,14 +293,29 @@ KOrder::KOrder(int num_stat, int num_pred, int num_both, int num_forw,
                const TensorContainer<FSSparseTensor>& fcont, const TwoDMatrix& gy,
                const TwoDMatrix& gu, const TwoDMatrix& v, Journal& jr) :
     ypart(num_stat, num_pred, num_both, num_forw),
-    ny(ypart.ny()), nu(gu.ncols()), maxk(fcont.getMaxDim()), nvs {ypart.nys(), nu, nu, 1}, _ug(4),
-    _fg(4), _ugs(4), _fgs(4), _ugss(4), _fgss(4), _uG(4), _fG(4),
+    ny(ypart.ny()),
+    nu(gu.ncols()),
+    maxk(fcont.getMaxDim()),
+    nvs {ypart.nys(), nu, nu, 1},
+    _ug(4),
+    _fg(4),
+    _ugs(4),
+    _fgs(4),
+    _ugss(4),
+    _fgss(4),
+    _uG(4),
+    _fG(4),
     _uZstack(&_uG, ypart.nyss(), &_ug, ny, ypart.nys(), nu),
-    _fZstack(&_fG, ypart.nyss(), &_fg, ny, ypart.nys(), nu), _uGstack(&_ugs, ypart.nys(), nu),
-    _fGstack(&_fgs, ypart.nys(), nu), _um(maxk, v), _fm(_um), f(fcont),
+    _fZstack(&_fG, ypart.nyss(), &_fg, ny, ypart.nys(), nu),
+    _uGstack(&_ugs, ypart.nys(), nu),
+    _fGstack(&_fgs, ypart.nys(), nu),
+    _um(maxk, v),
+    _fm(_um),
+    f(fcont),
     matA(f.get(Symmetry {1}), _uZstack.getStackSizes(), gy, ypart),
     matS(f.get(Symmetry {1}), _uZstack.getStackSizes(), gy, ypart),
-    matB(f.get(Symmetry {1}), _uZstack.getStackSizes()), journal(jr)
+    matB(f.get(Symmetry {1}), _uZstack.getStackSizes()),
+    journal(jr)
 {
   KORD_RAISE_IF(gy.ncols() != ypart.nys(), "Wrong number of columns in gy in KOrder constructor");
   KORD_RAISE_IF(v.ncols() != nu, "Wrong number of columns of Vcov in KOrder constructor");
diff --git a/mex/sources/libkorder/kord/korder_stoch.cc b/mex/sources/libkorder/kord/korder_stoch.cc
index 418765326b1f730ce236aa62719565a18c80cbc0..b18f3c4d30233844defd3c9a1054b47da53b3a13 100644
--- a/mex/sources/libkorder/kord/korder_stoch.cc
+++ b/mex/sources/libkorder/kord/korder_stoch.cc
@@ -45,10 +45,21 @@ MatrixAA::MatrixAA(const FSSparseTensor& f, const IntSequence& ss, const TwoDMat
 KOrderStoch::KOrderStoch(const PartitionY& yp, int nu, const TensorContainer<FSSparseTensor>& fcont,
                          const FGSContainer& hh, Journal& jr) :
     nvs {yp.nys(), nu, nu, 1},
-    ypart(yp), journal(jr), _ug(4), _fg(4), _ugs(4), _fgs(4), _uG(4), _fG(4), _uh(nullptr),
-    _fh(&hh), _uZstack(&_uG, ypart.nyss(), &_ug, ypart.ny(), ypart.nys(), nu),
+    ypart(yp),
+    journal(jr),
+    _ug(4),
+    _fg(4),
+    _ugs(4),
+    _fgs(4),
+    _uG(4),
+    _fG(4),
+    _uh(nullptr),
+    _fh(&hh),
+    _uZstack(&_uG, ypart.nyss(), &_ug, ypart.ny(), ypart.nys(), nu),
     _fZstack(&_fG, ypart.nyss(), &_fg, ypart.ny(), ypart.nys(), nu),
-    _uGstack(&_ugs, ypart.nys(), nu), _fGstack(&_fgs, ypart.nys(), nu), f(fcont),
+    _uGstack(&_ugs, ypart.nys(), nu),
+    _fGstack(&_fgs, ypart.nys(), nu),
+    f(fcont),
     matA(fcont.get(Symmetry {1}), _uZstack.getStackSizes(), hh.get(Symmetry {1, 0, 0, 0}), ypart)
 {
 }
@@ -57,10 +68,21 @@ KOrderStoch::KOrderStoch(const PartitionY& yp, int nu, const TensorContainer<FSS
 KOrderStoch::KOrderStoch(const PartitionY& yp, int nu, const TensorContainer<FSSparseTensor>& fcont,
                          const UGSContainer& hh, Journal& jr) :
     nvs {yp.nys(), nu, nu, 1},
-    ypart(yp), journal(jr), _ug(4), _fg(4), _ugs(4), _fgs(4), _uG(4), _fG(4), _uh(&hh),
-    _fh(nullptr), _uZstack(&_uG, ypart.nyss(), &_ug, ypart.ny(), ypart.nys(), nu),
+    ypart(yp),
+    journal(jr),
+    _ug(4),
+    _fg(4),
+    _ugs(4),
+    _fgs(4),
+    _uG(4),
+    _fG(4),
+    _uh(&hh),
+    _fh(nullptr),
+    _uZstack(&_uG, ypart.nyss(), &_ug, ypart.ny(), ypart.nys(), nu),
     _fZstack(&_fG, ypart.nyss(), &_fg, ypart.ny(), ypart.nys(), nu),
-    _uGstack(&_ugs, ypart.nys(), nu), _fGstack(&_fgs, ypart.nys(), nu), f(fcont),
+    _uGstack(&_ugs, ypart.nys(), nu),
+    _fGstack(&_fgs, ypart.nys(), nu),
+    f(fcont),
     matA(fcont.get(Symmetry {1}), _uZstack.getStackSizes(), hh.get(Symmetry {1, 0, 0, 0}), ypart)
 {
 }
diff --git a/mex/sources/libkorder/sylv/GeneralMatrix.cc b/mex/sources/libkorder/sylv/GeneralMatrix.cc
index 99d1e6be196125518cd8cb7a7087c2d3c6a29451..812661522552854845f60aecda7cc427adbe1fb2 100644
--- a/mex/sources/libkorder/sylv/GeneralMatrix.cc
+++ b/mex/sources/libkorder/sylv/GeneralMatrix.cc
@@ -339,7 +339,9 @@ GeneralMatrix::gemm_partial_right(const std::string& trans, const ConstGeneralMa
 }
 
 ConstGeneralMatrix::ConstGeneralMatrix(const GeneralMatrix& m, int i, int j, int nrows, int ncols) :
-    data(m.getData(), j * m.getLD() + i, (ncols - 1) * m.getLD() + nrows), rows(nrows), cols(ncols),
+    data(m.getData(), j * m.getLD() + i, (ncols - 1) * m.getLD() + nrows),
+    rows(nrows),
+    cols(ncols),
     ld(m.getLD())
 {
   // FIXME: check that the submatrix is fully in the matrix
@@ -348,7 +350,9 @@ ConstGeneralMatrix::ConstGeneralMatrix(const GeneralMatrix& m, int i, int j, int
 ConstGeneralMatrix::ConstGeneralMatrix(const ConstGeneralMatrix& m, int i, int j, int nrows,
                                        int ncols) :
     data(m.getData(), j * m.getLD() + i, (ncols - 1) * m.getLD() + nrows),
-    rows(nrows), cols(ncols), ld(m.getLD())
+    rows(nrows),
+    cols(ncols),
+    ld(m.getLD())
 {
   // FIXME: check that the submatrix is fully in the matrix
 }
diff --git a/mex/sources/libkorder/sylv/GeneralMatrix.hh b/mex/sources/libkorder/sylv/GeneralMatrix.hh
index a9168c51a6aa7f2cbc0195694076f8efff9953fc..498f2bd5c1020466bc764ac1bcaf838c85f345a1 100644
--- a/mex/sources/libkorder/sylv/GeneralMatrix.hh
+++ b/mex/sources/libkorder/sylv/GeneralMatrix.hh
@@ -586,8 +586,11 @@ protected:
 
 public:
   SVDDecomp(const GeneralMatrix& A) :
-      minmn(std::min<int>(A.nrows(), A.ncols())), sigma(minmn), U(A.nrows(), A.nrows()),
-      VT(A.ncols(), A.ncols()), conv(false)
+      minmn(std::min<int>(A.nrows(), A.ncols())),
+      sigma(minmn),
+      U(A.nrows(), A.nrows()),
+      VT(A.ncols(), A.ncols()),
+      conv(false)
   {
     construct(A);
   }
diff --git a/mex/sources/libkorder/sylv/GeneralSylvester.cc b/mex/sources/libkorder/sylv/GeneralSylvester.cc
index 785de3e274fcbb7e82ab4b75aada4b7d5ff0390b..43548103e8a0ba0da18c348639b4e3f019eea670 100644
--- a/mex/sources/libkorder/sylv/GeneralSylvester.cc
+++ b/mex/sources/libkorder/sylv/GeneralSylvester.cc
@@ -31,8 +31,12 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols, const C
                                    const ConstVector& db, const ConstVector& dc,
                                    const ConstVector& dd, const SylvParams& ps) :
     pars(ps),
-    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)
+    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)
 {
   init();
 }
@@ -41,8 +45,12 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols, const C
                                    const ConstVector& db, const ConstVector& dc, Vector& dd,
                                    const SylvParams& ps) :
     pars(ps),
-    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)
+    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)
 {
   init();
 }
@@ -51,8 +59,12 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols, const C
                                    const ConstVector& db, const ConstVector& dc,
                                    const ConstVector& dd, bool alloc_for_check) :
     pars(alloc_for_check),
-    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)
+    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)
 {
   init();
 }
@@ -61,8 +73,12 @@ GeneralSylvester::GeneralSylvester(int ord, int n, int m, int zero_cols, const C
                                    const ConstVector& db, const ConstVector& dc, Vector& dd,
                                    bool alloc_for_check) :
     pars(alloc_for_check),
-    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)
+    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)
 {
   init();
 }
diff --git a/mex/sources/libkorder/sylv/KronVector.cc b/mex/sources/libkorder/sylv/KronVector.cc
index 939bea7e674ddd1171182d9f9d68d11a2f4d3f70..bf1a696ff2346127078fd3925feb4fa65f024f25 100644
--- a/mex/sources/libkorder/sylv/KronVector.cc
+++ b/mex/sources/libkorder/sylv/KronVector.cc
@@ -35,7 +35,9 @@ KronVector::KronVector(Vector& v, int mm, int nn, int dp) : Vector(v), m(mm), n(
 }
 
 KronVector::KronVector(KronVector& v, int i) :
-    Vector(v, i * power(v.m, v.depth - 1) * v.n, power(v.m, v.depth - 1) * v.n), m(v.m), n(v.n),
+    Vector(v, i * power(v.m, v.depth - 1) * v.n, power(v.m, v.depth - 1) * v.n),
+    m(v.m),
+    n(v.n),
     depth(v.depth - 1)
 {
   if (depth < 0)
@@ -87,15 +89,19 @@ ConstKronVector::ConstKronVector(ConstVector v, int mm, int nn, int dp) :
 ConstKronVector::ConstKronVector(const KronVector& v, int i) :
     ConstVector(v, i * power(v.getM(), v.getDepth() - 1) * v.getN(),
                 power(v.getM(), v.getDepth() - 1) * v.getN()),
-    m(v.getM()), n(v.getN()), depth(v.getDepth() - 1)
+    m(v.getM()),
+    n(v.getN()),
+    depth(v.getDepth() - 1)
 {
   if (depth < 0)
     throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
 }
 
 ConstKronVector::ConstKronVector(const ConstKronVector& v, int i) :
-    ConstVector(v, i * power(v.m, v.depth - 1) * v.n, power(v.m, v.depth - 1) * v.n), m(v.getM()),
-    n(v.getN()), depth(v.getDepth() - 1)
+    ConstVector(v, i * power(v.m, v.depth - 1) * v.n, power(v.m, v.depth - 1) * v.n),
+    m(v.getM()),
+    n(v.getN()),
+    depth(v.getDepth() - 1)
 {
   if (depth < 0)
     throw SYLV_MES_EXCEPTION("Bad KronVector pick, depth < 0.");
diff --git a/mex/sources/libkorder/sylv/QuasiTriangular.cc b/mex/sources/libkorder/sylv/QuasiTriangular.cc
index cc68318c791f87ab8124f716b58cf2456c489f49..80f045b489e021e41eb339e883e630de880c00ad 100644
--- a/mex/sources/libkorder/sylv/QuasiTriangular.cc
+++ b/mex/sources/libkorder/sylv/QuasiTriangular.cc
@@ -395,8 +395,7 @@ QuasiTriangular::QuasiTriangular(double r, const QuasiTriangular& t) :
 
 QuasiTriangular::QuasiTriangular(double r, const QuasiTriangular& t, double r2,
                                  const QuasiTriangular& t2) :
-    SqSylvMatrix(t.nrows()),
-    diagonal(getData().base(), t.diagonal)
+    SqSylvMatrix(t.nrows()), diagonal(getData().base(), t.diagonal)
 {
   setMatrix(r, t);
   addMatrix(r2, t2);
@@ -414,8 +413,7 @@ QuasiTriangular::QuasiTriangular(const ConstVector& d, int d_size) :
 
 QuasiTriangular::QuasiTriangular([[maybe_unused]] const std::string& dummy,
                                  const QuasiTriangular& t) :
-    SqSylvMatrix(t.nrows()),
-    diagonal(getData().base(), t.diagonal)
+    SqSylvMatrix(t.nrows()), diagonal(getData().base(), t.diagonal)
 {
   Vector aux(t.getData());
   blas_int d_size = diagonal.getSize();
diff --git a/mex/sources/libkorder/sylv/QuasiTriangularZero.cc b/mex/sources/libkorder/sylv/QuasiTriangularZero.cc
index 1e43924a84d1fcdcce7419ee7a1b68b0c193a86b..06eb2794c18f081c5688dfbc49c50b167c2e69bc 100644
--- a/mex/sources/libkorder/sylv/QuasiTriangularZero.cc
+++ b/mex/sources/libkorder/sylv/QuasiTriangularZero.cc
@@ -43,8 +43,7 @@ QuasiTriangularZero::QuasiTriangularZero(double r, const QuasiTriangularZero& t)
 
 QuasiTriangularZero::QuasiTriangularZero(double r, const QuasiTriangularZero& t, double r2,
                                          const QuasiTriangularZero& t2) :
-    QuasiTriangular(r, t, r2, t2),
-    nz(t.nz), ru(t.ru)
+    QuasiTriangular(r, t, r2, t2), nz(t.nz), ru(t.ru)
 {
   ru.mult(r);
   ru.add(r2, t2.ru);
@@ -57,7 +56,8 @@ QuasiTriangularZero::QuasiTriangularZero(const std::string& dummy, const QuasiTr
 }
 
 QuasiTriangularZero::QuasiTriangularZero(const SchurDecompZero& decomp) :
-    QuasiTriangular(decomp.getT().getData(), decomp.getT().nrows()), nz(decomp.getZeroCols()),
+    QuasiTriangular(decomp.getT().getData(), decomp.getT().nrows()),
+    nz(decomp.getZeroCols()),
     ru(decomp.getRU())
 {
 }
diff --git a/mex/sources/libkorder/sylv/SylvParams.hh b/mex/sources/libkorder/sylv/SylvParams.hh
index 8289cfc7ff52ce4b92282f1b23b0a1673a090440..f9c111be1b6886390c641b3238ab771dd324bdc9 100644
--- a/mex/sources/libkorder/sylv/SylvParams.hh
+++ b/mex/sources/libkorder/sylv/SylvParams.hh
@@ -210,7 +210,10 @@ public:
   DoubleParamItem cpu_time;       // time of the job in CPU seconds
 
   SylvParams(bool wc = false) :
-      method(solve_method::recurse), convergence_tol(1.e-30), max_num_iter(15), bs_norm(1.3),
+      method(solve_method::recurse),
+      convergence_tol(1.e-30),
+      max_num_iter(15),
+      bs_norm(1.3),
       want_check(wc)
   {
   }
diff --git a/mex/sources/libkorder/sylv/TriangularSylvester.cc b/mex/sources/libkorder/sylv/TriangularSylvester.cc
index 95aec463514b8913d6eb8c56313f36b22ab5d256..6b12fae22a6d1511f0330f29f945f2f9bf374043 100644
--- a/mex/sources/libkorder/sylv/TriangularSylvester.cc
+++ b/mex/sources/libkorder/sylv/TriangularSylvester.cc
@@ -33,15 +33,13 @@ TriangularSylvester::TriangularSylvester(const QuasiTriangular& k, const QuasiTr
 
 TriangularSylvester::TriangularSylvester(const SchurDecompZero& kdecomp,
                                          const SchurDecomp& fdecomp) :
-    SylvesterSolver(kdecomp, fdecomp),
-    matrixKK {matrixK->square()}, matrixFF {matrixF->square()}
+    SylvesterSolver(kdecomp, fdecomp), matrixKK {matrixK->square()}, matrixFF {matrixF->square()}
 {
 }
 
 TriangularSylvester::TriangularSylvester(const SchurDecompZero& kdecomp,
                                          const SimilarityDecomp& fdecomp) :
-    SylvesterSolver(kdecomp, fdecomp),
-    matrixKK {matrixK->square()}, matrixFF {matrixF->square()}
+    SylvesterSolver(kdecomp, fdecomp), matrixKK {matrixK->square()}, matrixFF {matrixF->square()}
 {
 }
 
diff --git a/mex/sources/libkorder/sylv/Vector.hh b/mex/sources/libkorder/sylv/Vector.hh
index ecf596c6cfd6c83fedfb7994653e8f233ee3fc6b..6c148300edfaca2540af335ab4019e2a081ff1d5 100644
--- a/mex/sources/libkorder/sylv/Vector.hh
+++ b/mex/sources/libkorder/sylv/Vector.hh
@@ -57,7 +57,9 @@ public:
   }
   Vector(const Vector& v);
   Vector(Vector&& v) :
-      len {std::exchange(v.len, 0)}, s {v.s}, data {std::exchange(v.data, nullptr)},
+      len {std::exchange(v.len, 0)},
+      s {v.s},
+      data {std::exchange(v.data, nullptr)},
       destroy {std::exchange(v.destroy, false)}
   {
   }
diff --git a/mex/sources/libkorder/tl/fine_container.hh b/mex/sources/libkorder/tl/fine_container.hh
index 6e57f5c02c28382f61bd4325719d0481d4df6ae9..8569002ca2271cb17ddbd9b373f036120564bfa7 100644
--- a/mex/sources/libkorder/tl/fine_container.hh
+++ b/mex/sources/libkorder/tl/fine_container.hh
@@ -117,9 +117,10 @@ public:
      conts. */
 
   FineContainer(const _Stype& sc, int max) :
-      SizeRefinement(sc.getStackSizes(), sc.numConts(), max), StackContainer<_Ttype>(
-                                                                  numRefinements(), getNC()),
-      ref_conts(getNC()), stack_cont(sc)
+      SizeRefinement(sc.getStackSizes(), sc.numConts(), max),
+      StackContainer<_Ttype>(numRefinements(), getNC()),
+      ref_conts(getNC()),
+      stack_cont(sc)
   {
     for (int i = 0; i < numRefinements(); i++)
       _Stype::stack_sizes[i] = getRefSize(i);
diff --git a/mex/sources/libkorder/tl/int_sequence.hh b/mex/sources/libkorder/tl/int_sequence.hh
index f37815ac3cf69079cf6fb04e223e82fa50416534..fb95e5b3d8c49d5ae85730920eefe32335be6e71 100644
--- a/mex/sources/libkorder/tl/int_sequence.hh
+++ b/mex/sources/libkorder/tl/int_sequence.hh
@@ -87,7 +87,8 @@ public:
   }
   // Move constructor
   IntSequence(IntSequence&& s) noexcept :
-      data {std::exchange(s.data, nullptr)}, length {std::exchange(s.length, 0)},
+      data {std::exchange(s.data, nullptr)},
+      length {std::exchange(s.length, 0)},
       destroy {std::exchange(s.destroy, false)}
   {
   }
diff --git a/mex/sources/libkorder/tl/pyramid_prod2.cc b/mex/sources/libkorder/tl/pyramid_prod2.cc
index d4f6bb32a7398979cef313a3ed76b72816bdc495..1ce556f8ffef7c97415c02b8108bd9bb4cc527a4 100644
--- a/mex/sources/libkorder/tl/pyramid_prod2.cc
+++ b/mex/sources/libkorder/tl/pyramid_prod2.cc
@@ -26,7 +26,9 @@
    ‘end_seq’ according to ‘unit_flag’ and columns lengths. */
 
 IrregTensorHeader::IrregTensorHeader(const StackProduct<FGSTensor>& sp, const IntSequence& c) :
-    nv(sp.getAllSize()), unit_flag(sp.dimen()), cols(sp.createPackedColumns(c, unit_flag)),
+    nv(sp.getAllSize()),
+    unit_flag(sp.dimen()),
+    cols(sp.createPackedColumns(c, unit_flag)),
     end_seq(sp.dimen())
 {
   for (int i = 0; i < sp.dimen(); i++)
diff --git a/mex/sources/libkorder/tl/sparse_tensor.cc b/mex/sources/libkorder/tl/sparse_tensor.cc
index f8930f12bef0ef7696e842e36108ead4d848ce30..2730a33c8f0f009e7ddff3bba119dafdfd449281 100644
--- a/mex/sources/libkorder/tl/sparse_tensor.cc
+++ b/mex/sources/libkorder/tl/sparse_tensor.cc
@@ -204,8 +204,7 @@ FSSparseTensor::print() const
 /* This is the same as FGSTensor slicing constructor from FSSparseTensor. */
 GSSparseTensor::GSSparseTensor(const FSSparseTensor& t, const IntSequence& ss,
                                const IntSequence& coor, TensorDimens td) :
-    SparseTensor(td.dimen(), t.nrows(), td.calcFoldMaxOffset()),
-    tdims(std::move(td))
+    SparseTensor(td.dimen(), t.nrows(), td.calcFoldMaxOffset()), tdims(std::move(td))
 {
   // set ‘lb’ and ‘ub’ to lower and upper bounds of slice indices
   /* The same code is present in FGSTensor slicing constructor, see it for
diff --git a/mex/sources/libkorder/tl/stack_container.cc b/mex/sources/libkorder/tl/stack_container.cc
index 220ed8db087bb794f4b3a47a0e0d0fbc4f4a9e44..747f1d34a31e05fcfad845ea2dca6bb6becb179f 100644
--- a/mex/sources/libkorder/tl/stack_container.cc
+++ b/mex/sources/libkorder/tl/stack_container.cc
@@ -76,8 +76,7 @@ WorkerFoldMAADense::operator()(std::mutex& mut)
 
 WorkerFoldMAADense::WorkerFoldMAADense(const FoldedStackContainer& container, Symmetry s,
                                        const FGSContainer& dcontainer, FGSTensor& outten) :
-    cont(container),
-    sym(std::move(s)), dense_cont(dcontainer), out(outten)
+    cont(container), sym(std::move(s)), dense_cont(dcontainer), out(outten)
 {
 }
 
@@ -140,8 +139,7 @@ WorkerFoldMAASparse1::operator()(std::mutex& mut)
 WorkerFoldMAASparse1::WorkerFoldMAASparse1(const FoldedStackContainer& container,
                                            const FSSparseTensor& ten, FGSTensor& outten,
                                            IntSequence c) :
-    cont(container),
-    t(ten), out(outten), coor(std::move(c))
+    cont(container), t(ten), out(outten), coor(std::move(c))
 {
 }
 
@@ -196,8 +194,7 @@ WorkerFoldMAASparse2::operator()(std::mutex& mut)
 WorkerFoldMAASparse2::WorkerFoldMAASparse2(const FoldedStackContainer& container,
                                            const FSSparseTensor& ten, FGSTensor& outten,
                                            IntSequence c) :
-    cont(container),
-    t(ten), out(outten), coor(std::move(c))
+    cont(container), t(ten), out(outten), coor(std::move(c))
 {
 }
 
@@ -270,8 +267,7 @@ WorkerFoldMAASparse4::operator()(std::mutex& mut)
 WorkerFoldMAASparse4::WorkerFoldMAASparse4(const FoldedStackContainer& container,
                                            const FSSparseTensor& ten, FGSTensor& outten,
                                            IntSequence c) :
-    cont(container),
-    t(ten), out(outten), coor(std::move(c))
+    cont(container), t(ten), out(outten), coor(std::move(c))
 {
 }
 
@@ -405,8 +401,7 @@ WorkerUnfoldMAADense::operator()(std::mutex& mut)
 
 WorkerUnfoldMAADense::WorkerUnfoldMAADense(const UnfoldedStackContainer& container, Symmetry s,
                                            const UGSContainer& dcontainer, UGSTensor& outten) :
-    cont(container),
-    sym(std::move(s)), dense_cont(dcontainer), out(outten)
+    cont(container), sym(std::move(s)), dense_cont(dcontainer), out(outten)
 {
 }
 
@@ -496,8 +491,7 @@ WorkerUnfoldMAASparse1::operator()(std::mutex& mut)
 WorkerUnfoldMAASparse1::WorkerUnfoldMAASparse1(const UnfoldedStackContainer& container,
                                                const FSSparseTensor& ten, UGSTensor& outten,
                                                IntSequence c) :
-    cont(container),
-    t(ten), out(outten), coor(std::move(c))
+    cont(container), t(ten), out(outten), coor(std::move(c))
 {
 }
 
@@ -560,8 +554,7 @@ WorkerUnfoldMAASparse2::operator()(std::mutex& mut)
 WorkerUnfoldMAASparse2::WorkerUnfoldMAASparse2(const UnfoldedStackContainer& container,
                                                const FSSparseTensor& ten, UGSTensor& outten,
                                                IntSequence c) :
-    cont(container),
-    t(ten), out(outten), coor(std::move(c))
+    cont(container), t(ten), out(outten), coor(std::move(c))
 {
 }
 
diff --git a/mex/sources/libkorder/tl/tensor.hh b/mex/sources/libkorder/tl/tensor.hh
index cc4506636130a1960769a5fc174bad46986fe7c7..1a02940b767a7561e91de2cf8483b222fb55465e 100644
--- a/mex/sources/libkorder/tl/tensor.hh
+++ b/mex/sources/libkorder/tl/tensor.hh
@@ -180,13 +180,17 @@ protected:
 
 public:
   Tensor(indor io, IntSequence last, int r, int c, int d) :
-      TwoDMatrix(r, c), in_beg(*this, d),
-      in_end(*this, std::move(last), (io == indor::along_row) ? r : c), dim(d)
+      TwoDMatrix(r, c),
+      in_beg(*this, d),
+      in_end(*this, std::move(last), (io == indor::along_row) ? r : c),
+      dim(d)
   {
   }
   Tensor(indor io, IntSequence first, IntSequence last, int r, int c, int d) :
-      TwoDMatrix(r, c), in_beg(*this, std::move(first), 0),
-      in_end(*this, std::move(last), (io == indor::along_row) ? r : c), dim(d)
+      TwoDMatrix(r, c),
+      in_beg(*this, std::move(first), 0),
+      in_end(*this, std::move(last), (io == indor::along_row) ? r : c),
+      dim(d)
   {
   }
   Tensor(int first_row, int num, Tensor& t) :
@@ -194,8 +198,10 @@ public:
   {
   }
   Tensor(const Tensor& t) :
-      TwoDMatrix(t), in_beg(*this, t.in_beg.getCoor(), *(t.in_beg)),
-      in_end(*this, t.in_end.getCoor(), *(t.in_end)), dim(t.dim)
+      TwoDMatrix(t),
+      in_beg(*this, t.in_beg.getCoor(), *(t.in_beg)),
+      in_end(*this, t.in_end.getCoor(), *(t.in_end)),
+      dim(t.dim)
   {
   }
   Tensor(Tensor&&) = default;
diff --git a/mex/sources/libkorder/tl/tests/monoms.cc b/mex/sources/libkorder/tl/tests/monoms.cc
index 547ec250fcc8e8d98a24e11f447a99594c4bbbe0..a8262f37aa0bbf77577b2e7fe63c6280a0fd6a88 100644
--- a/mex/sources/libkorder/tl/tests/monoms.cc
+++ b/mex/sources/libkorder/tl/tests/monoms.cc
@@ -381,8 +381,7 @@ Monom4Vector::print() const
 
 SparseDerivGenerator::SparseDerivGenerator(int nf, int ny, int nu, int nup, int nbigg, int ng,
                                            int mx, double prob, int maxdim) :
-    maxdimen(maxdim),
-    bigg(4), g(4), rcont(4)
+    maxdimen(maxdim), bigg(4), g(4), rcont(4)
 {
   intgen.init(nf, ny, nu, nup, nbigg, mx, prob);
 
@@ -407,8 +406,7 @@ SparseDerivGenerator::SparseDerivGenerator(int nf, int ny, int nu, int nup, int
 
 DenseDerivGenerator::DenseDerivGenerator(int ng, int nx, int ny, int nu, int mx, double prob,
                                          int maxdim) :
-    maxdimen(maxdim),
-    xcont(0), rcont(0), ts(maxdimen), uxcont(0), uts(maxdimen)
+    maxdimen(maxdim), xcont(0), rcont(0), ts(maxdimen), uxcont(0), uts(maxdimen)
 {
   intgen.init(ng, nx, ny, nu, nu, mx, prob);
   Monom1Vector g(nx, ng);
diff --git a/mex/sources/perfect_foresight_problem/DynamicModelCaller.cc b/mex/sources/perfect_foresight_problem/DynamicModelCaller.cc
index 283c93321079070fa0c3fdd4ca57025b65fac1d3..fd8558e758b540d5f0d55f671b927a499ebdc4e2 100644
--- a/mex/sources/perfect_foresight_problem/DynamicModelCaller.cc
+++ b/mex/sources/perfect_foresight_problem/DynamicModelCaller.cc
@@ -95,7 +95,9 @@ DynamicModelDllCaller::DynamicModelDllCaller(size_t ntt, mwIndex ny, mwIndex nx,
                                              const int32_T* g1_sparse_colptr_arg, bool linear_arg,
                                              bool compute_jacobian_arg) :
     DynamicModelCaller {linear_arg, compute_jacobian_arg},
-    params {params_arg}, steady_state {steady_state_arg}, g1_sparse_colptr {g1_sparse_colptr_arg}
+    params {params_arg},
+    steady_state {steady_state_arg},
+    g1_sparse_colptr {g1_sparse_colptr_arg}
 {
   tt = std::make_unique<double[]>(ntt);
   y_p = std::make_unique<double[]>(3 * ny);
@@ -134,11 +136,13 @@ DynamicModelMatlabCaller::DynamicModelMatlabCaller(std::string basename_arg, mwI
                                                    const mxArray* g1_sparse_colptr_mx_arg,
                                                    bool linear_arg, bool compute_jacobian_arg) :
     DynamicModelCaller {linear_arg, compute_jacobian_arg},
-    basename {std::move(basename_arg)}, y_mx {mxCreateDoubleMatrix(3 * ny, 1, mxREAL)},
-    x_mx {mxCreateDoubleMatrix(nx, 1, mxREAL)}, jacobian_mx {nullptr}, params_mx {mxDuplicateArray(
-                                                                           params_mx_arg)},
-    steady_state_mx {mxDuplicateArray(steady_state_mx_arg)}, g1_sparse_rowval_mx {mxDuplicateArray(
-                                                                 g1_sparse_rowval_mx_arg)},
+    basename {std::move(basename_arg)},
+    y_mx {mxCreateDoubleMatrix(3 * ny, 1, mxREAL)},
+    x_mx {mxCreateDoubleMatrix(nx, 1, mxREAL)},
+    jacobian_mx {nullptr},
+    params_mx {mxDuplicateArray(params_mx_arg)},
+    steady_state_mx {mxDuplicateArray(steady_state_mx_arg)},
+    g1_sparse_rowval_mx {mxDuplicateArray(g1_sparse_rowval_mx_arg)},
     g1_sparse_colval_mx {mxDuplicateArray(g1_sparse_colval_mx_arg)},
     g1_sparse_colptr_mx {mxDuplicateArray(g1_sparse_colptr_mx_arg)}
 {