Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
16 changes: 8 additions & 8 deletions integration/VODE/vode_dvjust.H
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@ void dvjust (int IORD, BurnT& state, DvodeT& vstate)
const int NQM1 = vstate.NQ - 1;
const int NQM2 = vstate.NQ - 2;

amrex::Real HSUM{}, XI{}, ALPH0{}, ALPH1{}, PROD{}, XIOLD{}, T1{};
amrex::Real XI{};

// Check to see if the order is being increased or decreased.

Expand All @@ -35,7 +35,7 @@ void dvjust (int IORD, BurnT& state, DvodeT& vstate)
vstate.el(j) = 0.0_rt;
}
vstate.el(3) = 1.0_rt;
HSUM = 0.0_rt;
amrex::Real HSUM = 0.0_rt;
for (int j = 1; j <= NQM2; ++j) {
// Construct coefficients of x*x*(x+xi(1))*...*(x+xi(j)).
HSUM += vstate.tau(j);
Expand All @@ -62,11 +62,11 @@ void dvjust (int IORD, BurnT& state, DvodeT& vstate)
}

vstate.el(3) = 1.0_rt;
ALPH0 = -1.0_rt;
ALPH1 = 1.0_rt;
PROD = 1.0_rt;
XIOLD = 1.0_rt;
HSUM = vstate.HSCAL;
amrex::Real ALPH0 = -1.0_rt;
amrex::Real ALPH1 = 1.0_rt;
amrex::Real PROD = 1.0_rt;
amrex::Real XIOLD = 1.0_rt;
amrex::Real HSUM = vstate.HSCAL;

if (vstate.NQ != 1) {
for (int j = 1; j <= NQM1; ++j) {
Expand All @@ -84,7 +84,7 @@ void dvjust (int IORD, BurnT& state, DvodeT& vstate)
}
}

T1 = (-ALPH0 - ALPH1) / PROD;
amrex::Real T1 = (-ALPH0 - ALPH1) / PROD;
// Load column L + 1 in YH array.
for (int i = 1; i <= int_neqs; ++i) {
vstate.yh(i,vstate.L+1) = T1 * vstate.yh(i,VODE_LMAX);
Expand Down
29 changes: 14 additions & 15 deletions integration/VODE/vode_dvstep.H
Original file line number Diff line number Diff line change
Expand Up @@ -82,14 +82,13 @@ int dvstep (BurnT& state, DvodeT& vstate)
constexpr amrex::Real ONEPSM = 1.00001e0_rt;
constexpr amrex::Real THRESH = 1.5e0_rt;

amrex::Real CNQUOT{}, DDN{}, DSM{}, DUP{}, TOLD{};
amrex::Real FLOTL{}, R{};
amrex::Real DSM{};
int NCF{}, NFLAG{};

constexpr int int_neqs = integrator_neqs<BurnT>();

int kflag = 0;
TOLD = vstate.tn;
amrex::Real TOLD = vstate.tn;
NCF = 0;
vstate.JCUR = 0;
NFLAG = 0;
Expand Down Expand Up @@ -117,7 +116,7 @@ int dvstep (BurnT& state, DvodeT& vstate)

// Rescale the history array for a change in H by a factor of ETA.

R = 1.0_rt;
amrex::Real R = 1.0_rt;

for (int j = 2; j <= vstate.L; ++j) {
R *= vstate.ETA;
Expand Down Expand Up @@ -201,7 +200,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
vstate.ETA = amrex::max(vstate.ETA, HMIN / std::abs(vstate.H));

// Rescale the history array for a change in H by a factor of ETA.
R = 1.0_rt;
amrex::Real R = 1.0_rt;

for (int j = 2; j <= vstate.L; ++j) {
R *= vstate.ETA;
Expand Down Expand Up @@ -357,7 +356,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
vstate.ETAMAX = ETAMX2;
}

R = 1.0_rt / vstate.tq(2);
amrex::Real R = 1.0_rt / vstate.tq(2);
for (int i = 1; i <= int_neqs; ++i) {
vstate.acor(i) *= R;
}
Expand Down Expand Up @@ -389,15 +388,15 @@ int dvstep (BurnT& state, DvodeT& vstate)
if (kflag > KFC) {

// Compute ratio of new H to current H at the current order.
FLOTL = vstate.L;
amrex::Real FLOTL = vstate.L;
vstate.ETA = 1.0_rt / (std::pow(BIAS2 * DSM, 1.0_rt / FLOTL) + ADDON);
vstate.ETA = amrex::max(vstate.ETA, HMIN / std::abs(vstate.H), ETAMIN);
if ((kflag <= -2) && (vstate.ETA > ETAMXF)) {
vstate.ETA = ETAMXF;
}

// Rescale the history array for a change in H by a factor of ETA.
R = 1.0_rt;
amrex::Real R = 1.0_rt;

for (int j = 2; j <= vstate.L; ++j) {
R *= vstate.ETA;
Expand Down Expand Up @@ -434,7 +433,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
vstate.NQWAIT = vstate.L;

// Rescale the history array for a change in H by a factor of ETA.
R = 1.0_rt;
amrex::Real R = 1.0_rt;

for (int j = 2; j <= vstate.L; ++j) {
R *= vstate.ETA;
Expand Down Expand Up @@ -477,7 +476,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
bool already_set_eta = false;

// Compute ratio of new H to current H at the current order.
FLOTL = vstate.L;
amrex::Real FLOTL = vstate.L;
const amrex::Real ETAQ = 1.0_rt / (std::pow(BIAS2 * DSM, 1.0_rt / FLOTL) + ADDON);
if (vstate.NQWAIT == 0) {
vstate.NQWAIT = 2;
Expand All @@ -486,7 +485,7 @@ int dvstep (BurnT& state, DvodeT& vstate)

if (vstate.NQ != 1) {
// Compute ratio of new H to current H at the current order less one.
DDN = 0.0_rt;
amrex::Real DDN = 0.0_rt;
for (int i = 1; i <= int_neqs; ++i) {
DDN += (vstate.yh(i,vstate.L) * vstate.ewt(i)) * (vstate.yh(i,vstate.L) * vstate.ewt(i));
}
Expand All @@ -498,11 +497,11 @@ int dvstep (BurnT& state, DvodeT& vstate)

if (vstate.L != VODE_LMAX) {
// Compute ratio of new H to current H at current order plus one.
CNQUOT = (vstate.tq(5) / vstate.CONP) * std::pow(vstate.H / vstate.tau(2), vstate.L);
amrex::Real CNQUOT = (vstate.tq(5) / vstate.CONP) * std::pow(vstate.H / vstate.tau(2), vstate.L);
for (int i = 1; i <= int_neqs; ++i) {
vstate.savf(i) = vstate.acor(i) - CNQUOT * vstate.yh(i,VODE_LMAX);
}
DUP = 0.0_rt;
amrex::Real DUP = 0.0_rt;
for (int i = 1; i <= int_neqs; ++i) {
DUP += (vstate.savf(i) * vstate.ewt(i)) * (vstate.savf(i) * vstate.ewt(i));
}
Expand Down Expand Up @@ -547,7 +546,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
if (vstate.n_step <= 10) {
vstate.ETAMAX = ETAMX2;
}
R = 1.0_rt / vstate.tq(2);
amrex::Real R = 1.0_rt / vstate.tq(2);
for (int i = 1; i <= int_neqs; ++i) {
vstate.acor(i) *= R;
}
Expand All @@ -561,7 +560,7 @@ int dvstep (BurnT& state, DvodeT& vstate)
if (vstate.n_step <= 10) {
vstate.ETAMAX = ETAMX2;
}
R = 1.0_rt / vstate.tq(2);
amrex::Real R = 1.0_rt / vstate.tq(2);
for (int i = 1; i <= int_neqs; ++i) {
vstate.acor(i) *= R;
}
Expand Down
Loading