diff --git a/examples/sim_simple.cpp b/examples/sim_simple.cpp index b8403b7ef..caa824ee3 100644 --- a/examples/sim_simple.cpp +++ b/examples/sim_simple.cpp @@ -90,8 +90,8 @@ phaseMobility(const Opm::IncompPropertiesInterface& props, std::vector dmw = { krwjac/mu[0] }; std::vector dmo = { krojac/mu[1] }; - std::vector pmobc = { ADB::function(krw / mu[0], dmw) , - ADB::function(kro / mu[1], dmo) }; + std::vector pmobc = { ADB::function(krw / mu[0], std::move(dmw)) , + ADB::function(kro / mu[1], std::move(dmo)) }; return pmobc; } diff --git a/opm/autodiff/AutoDiffBlock.hpp b/opm/autodiff/AutoDiffBlock.hpp index aaf9789d6..5ba57b32f 100644 --- a/opm/autodiff/AutoDiffBlock.hpp +++ b/opm/autodiff/AutoDiffBlock.hpp @@ -28,6 +28,7 @@ #include +#include #include #include #include @@ -98,9 +99,14 @@ namespace Opm /// Construct an empty AutoDiffBlock. static AutoDiffBlock null() { - V val; - std::vector jac; - return AutoDiffBlock(val, jac); + return AutoDiffBlock(V(), {}); + } + + /// Create an AutoDiffBlock representing a constant. + /// \param[in] val values + static AutoDiffBlock constant(V&& val) + { + return AutoDiffBlock(std::move(val)); } /// Create an AutoDiffBlock representing a constant. @@ -126,7 +132,8 @@ namespace Opm for (int i = 0; i < num_blocks; ++i) { jac[i] = M(num_elem, blocksizes[i]); } - return AutoDiffBlock(val, jac); + V val_copy(val); + return AutoDiffBlock(std::move(val_copy), std::move(jac)); } /// Create an AutoDiffBlock representing a single variable block. @@ -136,7 +143,7 @@ namespace Opm /// The resulting object will have size() equal to block_pattern[index]. /// Its jacobians will all be zero, except for derivative()[index], which /// will be an identity matrix. - static AutoDiffBlock variable(const int index, const V& val, const std::vector& blocksizes) + static AutoDiffBlock variable(const int index, V&& val, const std::vector& blocksizes) { std::vector jac; const int num_elem = val.size(); @@ -164,15 +171,42 @@ namespace Opm assert(jac[i].size()==0); } } - return AutoDiffBlock(val, jac); + return AutoDiffBlock(std::move(val), std::move(jac)); + } + + /// Create an AutoDiffBlock representing a single variable block. + /// \param[in] index index of the variable you are constructing + /// \param[in] val values + /// \param[in] blocksizes block pattern + /// The resulting object will have size() equal to block_pattern[index]. + /// Its jacobians will all be zero, except for derivative()[index], which + /// will be an identity matrix. + static AutoDiffBlock variable(const int index, const V& val, const std::vector& blocksizes) + { + V value = val; + return variable(index, std::move(value), blocksizes); } /// Create an AutoDiffBlock by directly specifying values and jacobians. + /// This version of function() moves its arguments and is therefore + /// quite efficient, but leaves the argument variables empty (but valid). + /// \param[in] val values + /// \param[in] jac vector of jacobians + static AutoDiffBlock function(V&& val, std::vector&& jac) + { + return AutoDiffBlock(std::move(val), std::move(jac)); + } + + /// Create an AutoDiffBlock by directly specifying values and jacobians. + /// This version of function() copies its arguments and is therefore + /// less efficient than the other (moving) overload. /// \param[in] val values /// \param[in] jac vector of jacobians static AutoDiffBlock function(const V& val, const std::vector& jac) { - return AutoDiffBlock(val, jac); + V val_copy(val); + std::vector jac_copy(jac); + return AutoDiffBlock(std::move(val_copy), std::move(jac_copy)); } /// Construct a set of primary variables, each initialized to @@ -259,7 +293,7 @@ namespace Opm assert(jac[block].cols() == rhs.jac_[block].cols()); jac[block] += rhs.jac_[block]; } - return function(val_ + rhs.val_, jac); + return function(val_ + rhs.val_, std::move(jac)); } /// Elementwise operator - @@ -282,7 +316,7 @@ namespace Opm assert(jac[block].cols() == rhs.jac_[block].cols()); jac[block] -= rhs.jac_[block]; } - return function(val_ - rhs.val_, jac); + return function(val_ - rhs.val_, std::move(jac)); } /// Elementwise operator * @@ -318,7 +352,7 @@ namespace Opm jac[block] = D2*jac_[block] + D1*rhs.jac_[block]; } } - return function(val_ * rhs.val_, jac); + return function(val_ * rhs.val_, std::move(jac)); } /// Elementwise operator / @@ -358,7 +392,7 @@ namespace Opm jac[block] = D3 * (D2*jac_[block] - D1*rhs.jac_[block]); } } - return function(val_ / rhs.val_, jac); + return function(val_ / rhs.val_, std::move(jac)); } /// I/O. @@ -374,6 +408,13 @@ namespace Opm return os; } + /// Efficient swap function. + void swap(AutoDiffBlock& other) + { + val_.swap(other.val_); + jac_.swap(other.jac_); + } + /// Number of elements int size() const { @@ -411,13 +452,17 @@ namespace Opm private: AutoDiffBlock(const V& val) - : val_(val) + : val_(val) { } - AutoDiffBlock(const V& val, - const std::vector& jac) - : val_(val), jac_(jac) + AutoDiffBlock(V&& val) + : val_(std::move(val)) + { + } + + AutoDiffBlock(V&& val, std::vector&& jac) + : val_(std::move(val)), jac_(std::move(jac)) { #ifndef NDEBUG const int num_elem = val_.size(); @@ -457,7 +502,7 @@ namespace Opm fastSparseProduct(lhs, rhs.derivative()[block], jac[block]); } typename AutoDiffBlock::V val = lhs*rhs.value().matrix(); - return AutoDiffBlock::function(val, jac); + return AutoDiffBlock::function(std::move(val), std::move(jac)); } @@ -549,7 +594,7 @@ namespace Opm for (int block=0; block::function( lhs.value() * rhs, jac ); + return AutoDiffBlock::function( lhs.value() * rhs, std::move(jac) ); } diff --git a/opm/autodiff/AutoDiffHelpers.hpp b/opm/autodiff/AutoDiffHelpers.hpp index 8036582d9..c49c6c65a 100644 --- a/opm/autodiff/AutoDiffHelpers.hpp +++ b/opm/autodiff/AutoDiffHelpers.hpp @@ -403,7 +403,8 @@ collapseJacs(const AutoDiffBlock& x) // Build final jacobian. std::vector jacs(1); collapseJacs( x, jacs[ 0 ] ); - return ADB::function(x.value(), jacs); + ADB::V val = x.value(); + return ADB::function(std::move(val), std::move(jacs)); } diff --git a/opm/autodiff/BlackoilPropsAdFromDeck.cpp b/opm/autodiff/BlackoilPropsAdFromDeck.cpp index 9d328ce88..2f1075aab 100644 --- a/opm/autodiff/BlackoilPropsAdFromDeck.cpp +++ b/opm/autodiff/BlackoilPropsAdFromDeck.cpp @@ -335,7 +335,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& props_[phase_usage_.phase_pos[Water]]->mu(n, pvt_region_.data(), pw.value().data(), T.value().data(), rs, mu.data(), dmudp.data(), dmudr.data()); if (pw.derivative().empty()) { - return ADB::constant(mu); + return ADB::constant(std::move(mu)); } else { ADB::M dmudp_diag = spdiag(dmudp); const int num_blocks = pw.numBlocks(); @@ -343,7 +343,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(dmudp_diag, pw.derivative()[block], jacs[block]); } - return ADB::function(mu, jacs); + return ADB::function(std::move(mu), std::move(jacs)); } } @@ -383,7 +383,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& fastSparseProduct(dmudr_diag, rs.derivative()[block], temp); jacs[block] += temp; } - return ADB::function(mu, jacs); + return ADB::function(std::move(mu), std::move(jacs)); } /// Gas viscosity. @@ -422,7 +422,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& fastSparseProduct(dmudr_diag, rv.derivative()[block], temp); jacs[block] += temp; } - return ADB::function(mu, jacs); + return ADB::function(std::move(mu), std::move(jacs)); } @@ -459,7 +459,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(dbdp_diag, pw.derivative()[block], jacs[block]); } - return ADB::function(b, jacs); + return ADB::function(std::move(b), std::move(jacs)); } /// Oil formation volume factor. @@ -499,7 +499,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& fastSparseProduct(dbdr_diag, rs.derivative()[block], temp); jacs[block] += temp; } - return ADB::function(b, jacs); + return ADB::function(std::move(b), std::move(jacs)); } /// Gas formation volume factor. @@ -539,7 +539,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& fastSparseProduct(dbdr_diag, rv.derivative()[block], temp); jacs[block] += temp; } - return ADB::function(b, jacs); + return ADB::function(std::move(b), std::move(jacs)); } @@ -568,7 +568,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(drbubdp_diag, po.derivative()[block], jacs[block]); } - return ADB::function(rbub, jacs); + return ADB::function(std::move(rbub), std::move(jacs)); } /// Bubble point curve for Rs as function of oil pressure. @@ -609,7 +609,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(drvdp_diag, po.derivative()[block], jacs[block]); } - return ADB::function(rv, jacs); + return ADB::function(std::move(rv), std::move(jacs)); } /// Condensation curve for Rv as function of oil pressure. @@ -686,7 +686,8 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& jacs[block] += temp; } } - relperms.emplace_back(ADB::function(kr.col(phase1_pos), jacs)); + ADB::V val = kr.col(phase1_pos); + relperms.emplace_back(ADB::function(std::move(val), std::move(jacs))); } else { relperms.emplace_back(ADB::null()); } @@ -746,7 +747,8 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& jacs[block] += temp; } } - adbCapPressures.emplace_back(ADB::function(pc.col(phase1_pos), jacs)); + ADB::V val = pc.col(phase1_pos); + adbCapPressures.emplace_back(ADB::function(std::move(val), std::move(jacs))); } else { adbCapPressures.emplace_back(ADB::null()); } @@ -849,7 +851,7 @@ BlackoilPropsAdFromDeck::BlackoilPropsAdFromDeck(const BlackoilPropsAdFromDeck& for (int block = 0; block < num_blocks; ++block) { jacs[block] = dfactor_dso_diag * so.derivative()[block]; } - r = ADB::function(factor, jacs)*r; + r = ADB::function(std::move(factor), std::move(jacs))*r; } } diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.hpp b/opm/autodiff/FullyImplicitBlackoilSolver.hpp index 71d6d90a8..6ed224c6c 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.hpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.hpp @@ -281,11 +281,6 @@ namespace Opm { std::vector computeRelPerm(const SolutionState& state) const; - std::vector - computeRelPermWells(const SolutionState& state, - const DataBlock& well_s, - const std::vector& well_cells) const; - void computeMassFlux(const int actph , const V& transi, diff --git a/opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp b/opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp index d0625ee32..5d7e63e09 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp @@ -542,21 +542,20 @@ namespace detail { // Pressure. int nextvar = 0; - state.pressure = vars[ nextvar++ ]; + state.pressure = std::move(vars[ nextvar++ ]); // temperature const V temp = Eigen::Map(& x.temperature()[0], x.temperature().size()); state.temperature = ADB::constant(temp); // Saturations - const std::vector& bpat = vars[0].blockPattern(); { - ADB so = ADB::constant(V::Ones(nc, 1), bpat); + ADB so = ADB::constant(V::Ones(nc, 1)); if (active_[ Water ]) { - ADB& sw = vars[ nextvar++ ]; - state.saturation[pu.phase_pos[ Water ]] = sw; + state.saturation[pu.phase_pos[ Water ]] = std::move(vars[ nextvar++ ]); + const ADB& sw = state.saturation[pu.phase_pos[ Water ]]; so -= sw; } @@ -572,7 +571,7 @@ namespace detail { // RS and RV is only defined if both oil and gas phase are active. const ADB& sw = (active_[ Water ] ? state.saturation[ pu.phase_pos[ Water ] ] - : ADB::constant(V::Zero(nc, 1), bpat)); + : ADB::constant(V::Zero(nc, 1))); state.canonical_phase_pressures = computePressures(state.pressure, sw, so, sg); const ADB rsSat = fluidRsSat(state.canonical_phase_pressures[ Oil ], so , cells_); if (has_disgas_) { @@ -591,15 +590,15 @@ namespace detail { if (active_[ Oil ]) { // Note that so is never a primary variable. - state.saturation[pu.phase_pos[ Oil ]] = so; + state.saturation[pu.phase_pos[ Oil ]] = std::move(so); } } // Qs. - state.qs = vars[ nextvar++ ]; + state.qs = std::move(vars[ nextvar++ ]); // Bhp. - state.bhp = vars[ nextvar++ ]; + state.bhp = std::move(vars[ nextvar++ ]); assert(nextvar == int(vars.size())); @@ -946,7 +945,7 @@ namespace detail { } // total rates at std - ADB qt_s = ADB::constant(V::Zero(nw), state.bhp.blockPattern()); + ADB qt_s = ADB::constant(V::Zero(nw)); for (int phase = 0; phase < np; ++phase) { qt_s += subset(state.qs, Span(nw, 1, phase*nw)); } @@ -954,7 +953,7 @@ namespace detail { // compute avg. and total wellbore phase volumetric rates at std. conds const DataBlock compi = Eigen::Map(wells().comp_frac, nw, np); std::vector wbq(np, ADB::null()); - ADB wbqt = ADB::constant(V::Zero(nw), state.pressure.blockPattern()); + ADB wbqt = ADB::constant(V::Zero(nw)); for (int phase = 0; phase < np; ++phase) { const int pos = pu.phase_pos[phase]; wbq[phase] = (isInj * compi.col(pos)) * qt_s - q_ps[phase]; @@ -995,7 +994,7 @@ namespace detail { ADB cqt_i = -(isInjInx * Tw) * (mt * drawdown); // compute volume ratio between connection at standard conditions - ADB volRat = ADB::constant(V::Zero(nperf), state.pressure.blockPattern()); + ADB volRat = ADB::constant(V::Zero(nperf)); std::vector cmix_s(np, ADB::null()); for (int phase = 0; phase < np; ++phase) { cmix_s[phase] = wops_.w2p * mix_s[phase]; @@ -1223,16 +1222,22 @@ namespace detail { // Update primary variables, if necessary. if (bhp_changed) { + // We will set the bhp primary variable to the new ones, + // but we do not change the derivatives here. ADB::V new_bhp = Eigen::Map(xw.bhp().data(), nw); - bhp = ADB::function(new_bhp, bhp.derivative()); + // Avoiding the copy below would require a value setter method + // in AutoDiffBlock. + std::vector old_derivs = bhp.derivative(); + bhp = ADB::function(std::move(new_bhp), std::move(old_derivs)); } if (rates_changed) { // Need to reshuffle well rates, from phase running fastest // to wells running fastest. // The transpose() below switches the ordering. const DataBlock wrates = Eigen::Map(xw.wellRates().data(), nw, np).transpose(); - const ADB::V new_qs = Eigen::Map(wrates.data(), nw*np); - well_phase_flow_rate = ADB::function(new_qs, well_phase_flow_rate.derivative()); + ADB::V new_qs = Eigen::Map(wrates.data(), nw*np); + std::vector old_derivs = well_phase_flow_rate.derivative(); + well_phase_flow_rate = ADB::function(std::move(new_qs), std::move(old_derivs)); } } @@ -1607,22 +1612,21 @@ namespace detail { { using namespace Opm::AutoDiffGrid; const int nc = numCells(grid_); - const std::vector& bpat = state.pressure.blockPattern(); - const ADB null = ADB::constant(V::Zero(nc, 1), bpat); + const ADB zero = ADB::constant(V::Zero(nc)); const Opm::PhaseUsage& pu = fluid_.phaseUsage(); - const ADB sw = (active_[ Water ] - ? state.saturation[ pu.phase_pos[ Water ] ] - : null); + const ADB& sw = (active_[ Water ] + ? state.saturation[ pu.phase_pos[ Water ] ] + : zero); - const ADB so = (active_[ Oil ] - ? state.saturation[ pu.phase_pos[ Oil ] ] - : null); + const ADB& so = (active_[ Oil ] + ? state.saturation[ pu.phase_pos[ Oil ] ] + : zero); - const ADB sg = (active_[ Gas ] - ? state.saturation[ pu.phase_pos[ Gas ] ] - : null); + const ADB& sg = (active_[ Gas ] + ? state.saturation[ pu.phase_pos[ Gas ] ] + : zero); return fluid_.relperm(sw, so, sg, cells_); } @@ -1637,9 +1641,8 @@ namespace detail { { using namespace Opm::AutoDiffGrid; const int nc = numCells(grid_); - const std::vector& bpat = state.pressure.blockPattern(); - const ADB null = ADB::constant(V::Zero(nc, 1), bpat); + const ADB null = ADB::constant(V::Zero(nc)); const Opm::PhaseUsage& pu = fluid_.phaseUsage(); const ADB& sw = (active_[ Water ] @@ -1718,38 +1721,6 @@ namespace detail { - template - std::vector - FullyImplicitBlackoilSolver::computeRelPermWells(const SolutionState& state, - const DataBlock& well_s, - const std::vector& well_cells) const - { - const int nw = wells().number_of_wells; - const int nperf = wells().well_connpos[nw]; - const std::vector& bpat = state.pressure.blockPattern(); - - const ADB null = ADB::constant(V::Zero(nperf), bpat); - - const Opm::PhaseUsage& pu = fluid_.phaseUsage(); - const ADB sw = (active_[ Water ] - ? ADB::constant(well_s.col(pu.phase_pos[ Water ]), bpat) - : null); - - const ADB so = (active_[ Oil ] - ? ADB::constant(well_s.col(pu.phase_pos[ Oil ]), bpat) - : null); - - const ADB sg = (active_[ Gas ] - ? ADB::constant(well_s.col(pu.phase_pos[ Gas ]), bpat) - : null); - - return fluid_.relperm(sw, so, sg, well_cells); - } - - - - - template void FullyImplicitBlackoilSolver::computeMassFlux(const int actph , @@ -2263,9 +2234,9 @@ namespace detail { for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(dpm_diag, p.derivative()[block], jacs[block]); } - return ADB::function(pm, jacs); + return ADB::function(std::move(pm), std::move(jacs)); } else { - return ADB::constant(V::Constant(n, 1.0), p.blockPattern()); + return ADB::constant(V::Constant(n, 1.0)); } } @@ -2291,9 +2262,9 @@ namespace detail { for (int block = 0; block < num_blocks; ++block) { fastSparseProduct(dtm_diag, p.derivative()[block], jacs[block]); } - return ADB::function(tm, jacs); + return ADB::function(std::move(tm), std::move(jacs)); } else { - return ADB::constant(V::Constant(n, 1.0), p.blockPattern()); + return ADB::constant(V::Constant(n, 1.0)); } } diff --git a/opm/autodiff/ImpesTPFAAD.cpp b/opm/autodiff/ImpesTPFAAD.cpp index 90fcaa3d3..e1ffaeae8 100644 --- a/opm/autodiff/ImpesTPFAAD.cpp +++ b/opm/autodiff/ImpesTPFAAD.cpp @@ -326,7 +326,6 @@ namespace { const ADB& p = vars[0]; const ADB T = ADB::constant(T0); const ADB& bhp = vars[1]; - std::vector bpat = p.blockPattern(); // Compute T_ij * (p_i - p_j). const ADB nkgradp = transi * (ops_.ngrad * p); @@ -351,10 +350,10 @@ namespace { const ADB nkgradp_well = transw * (p_perfcell - p_perfwell); const Selector cell_to_well_selector(nkgradp_well.value()); - cell_residual_ = ADB::constant(pv, bpat); - well_residual_ = ADB::constant(V::Zero(nw,1), bpat); - ADB divcontrib_sum = ADB::constant(V::Zero(nc,1), bpat); - qs_ = ADB::constant(V::Zero(nw*np, 1), bpat); + cell_residual_ = ADB::constant(pv); + well_residual_ = ADB::constant(V::Zero(nw,1)); + ADB divcontrib_sum = ADB::constant(V::Zero(nc,1)); + qs_ = ADB::constant(V::Zero(nw*np, 1)); for (int phase = 0; phase < np; ++phase) { const ADB cell_b = fluidFvf(phase, p, T, cells); const ADB cell_rho = fluidRho(phase, p, T, cells); @@ -381,7 +380,7 @@ namespace { const ADB well_contrib = superset(perf_flux*perf_b, well_cells, nc); const ADB divcontrib = delta_t * (ops_.div * (flux * face_b) + well_contrib); const V qcontrib = delta_t * q; - const ADB pvcontrib = ADB::constant(pv*z0, bpat); + const ADB pvcontrib = ADB::constant(pv*z0); const ADB component_contrib = pvcontrib + qcontrib; divcontrib_sum = divcontrib_sum - divcontrib/cell_b; cell_residual_ = cell_residual_ - (component_contrib/cell_b); diff --git a/opm/autodiff/NewtonIterationBlackoilCPR.cpp b/opm/autodiff/NewtonIterationBlackoilCPR.cpp index 1df51a984..0a014c5f6 100644 --- a/opm/autodiff/NewtonIterationBlackoilCPR.cpp +++ b/opm/autodiff/NewtonIterationBlackoilCPR.cpp @@ -314,7 +314,7 @@ namespace Opm if (eq == n) { continue; } - retval.push_back(ADB::function(vals[eq], jacs[eq])); + retval.push_back(ADB::function(std::move(vals[eq]), std::move(jacs[eq]))); } return retval; } @@ -338,7 +338,8 @@ namespace Opm // Build C. std::vector C_jacs = equation.derivative(); C_jacs.erase(C_jacs.begin() + n); - ADB eq_coll = collapseJacs(ADB::function(equation.value(), C_jacs)); + V equation_value = equation.value(); + ADB eq_coll = collapseJacs(ADB::function(std::move(equation_value), std::move(C_jacs))); const M& C = eq_coll.derivative()[0]; // Use sparse LU to solve the block submatrices @@ -395,7 +396,7 @@ namespace Opm // A concession to MRST, to obtain more similar behaviour: // swap the first two equations, so that oil is first, then water. auto eqs = eqs_in; - std::swap(eqs[0], eqs[1]); + eqs[0].swap(eqs[1]); // Characterize the material balance equations. const int n = eqs[0].size(); @@ -462,7 +463,7 @@ namespace Opm L.setFromTriplets(t.begin(), t.end()); // Combine in single block. - ADB total_residual = eqs[0]; + ADB total_residual = std::move(eqs[0]); for (int phase = 1; phase < num_phases; ++phase) { total_residual = vertcat(total_residual, eqs[phase]); } diff --git a/opm/autodiff/TransportSolverTwophaseAd.cpp b/opm/autodiff/TransportSolverTwophaseAd.cpp index 9ac7402ca..117ca51a7 100644 --- a/opm/autodiff/TransportSolverTwophaseAd.cpp +++ b/opm/autodiff/TransportSolverTwophaseAd.cpp @@ -129,8 +129,8 @@ namespace Opm std::vector dmw = { krwjac/mu[0] }; std::vector dmo = { krojac/mu[1] }; - std::vector pmobc = { ADB::function(krw / mu[0], dmw) , - ADB::function(kro / mu[1], dmo) }; + std::vector pmobc = { ADB::function(krw / mu[0], std::move(dmw)) , + ADB::function(kro / mu[1], std::move(dmo)) }; return pmobc; } diff --git a/tests/test_block.cpp b/tests/test_block.cpp index 102237924..51da5f0dd 100644 --- a/tests/test_block.cpp +++ b/tests/test_block.cpp @@ -77,12 +77,10 @@ BOOST_AUTO_TEST_CASE(ConstantInitialisation) { typedef AutoDiffBlock ADB; - std::vector blocksizes = { 3, 1, 2 }; - ADB::V v(3); v << 0.2, 1.2, 13.4; - ADB a = ADB::constant(v, blocksizes); + ADB a = ADB::constant(v); BOOST_REQUIRE(a.value().matrix() == v.matrix()); const std::vector& J = a.derivative(); @@ -137,7 +135,9 @@ BOOST_AUTO_TEST_CASE(FunctionInitialisation) jacs[j].insert(0,0) = -1.0; } - ADB f = ADB::function(v, jacs); + ADB::V v_copy(v); + std::vector jacs_copy(jacs); + ADB f = ADB::function(std::move(v_copy), std::move(jacs_copy)); BOOST_REQUIRE(f.value().matrix() == v.matrix());