From 8ccf1b3dbe89fea04b61fe2e2c0da453d38f5caf Mon Sep 17 00:00:00 2001 From: Andreas Lauser Date: Mon, 29 Jun 2015 15:28:30 +0200 Subject: [PATCH] multiplexer saturation function: get rid of the pesky union and the reinterpret_cast basically, we now store a pointer to the base class and convert this to the concrete incarnation of the saturation function. the only disadvantage of this is that it requires to allocate the concrete objects dynamically, i.e., there's another layer of indirection. On the plus side it saves a few bytes per object and it makes things quite a bit simpler, so it is a win IMO. thanks to [at]atgeirr for suggesting this. --- opm/core/props/satfunc/SatFuncBase.hpp | 2 + opm/core/props/satfunc/SatFuncMultiplexer.hpp | 147 +++++++++--------- 2 files changed, 72 insertions(+), 77 deletions(-) diff --git a/opm/core/props/satfunc/SatFuncBase.hpp b/opm/core/props/satfunc/SatFuncBase.hpp index 95e13ab0..88b42b9b 100644 --- a/opm/core/props/satfunc/SatFuncBase.hpp +++ b/opm/core/props/satfunc/SatFuncBase.hpp @@ -77,6 +77,8 @@ namespace Opm class SatFuncBase : public BlackoilPhases { public: + virtual ~SatFuncBase() {}; + void init(Opm::EclipseStateConstPtr eclipseState, const int table_num, const PhaseUsage phase_usg, diff --git a/opm/core/props/satfunc/SatFuncMultiplexer.hpp b/opm/core/props/satfunc/SatFuncMultiplexer.hpp index 4538e8a2..c882dcab 100644 --- a/opm/core/props/satfunc/SatFuncMultiplexer.hpp +++ b/opm/core/props/satfunc/SatFuncMultiplexer.hpp @@ -31,6 +31,8 @@ #include +#include + namespace Opm { @@ -50,43 +52,41 @@ public: }; // this is a helper macro which helps to save us from RSI in the following -#define OPM_MULTIPLEXER_SATFUNC_CALL(codeToCall) \ - switch (satFuncType_) { \ - case Gwseg: { \ - typedef SatFuncGwseg SatFunc; \ - if (false) SatFunc dummyForWarning; \ - auto& satFunc = getGwseg(); \ - codeToCall; \ - break; \ - } \ - case Stone2: { \ - typedef SatFuncStone2 SatFunc; \ - if (false) SatFunc dummyForWarning; \ - auto& satFunc = getStone2(); \ - codeToCall; \ - break; \ - } \ - case Simple: { \ - typedef SatFuncSimple SatFunc; \ - if (false) SatFunc dummyForWarning; \ - auto& satFunc = getSimple(); \ - codeToCall; \ - break; \ - } \ +#define OPM_MULTIPLEXER_CONST const +#define OPM_MULTIPLEXER_NON_CONST +#define OPM_MULTIPLEXER_SATFUNC_CALL(codeToCall, CONST) \ + switch (satFuncType_) { \ + case Gwseg: { \ + typedef SatFuncGwseg SatFunc; \ + __attribute__((unused)) CONST SatFunc& satFunc = \ + *static_cast(satFunc_.get()); \ + codeToCall; \ + break; \ + } \ + case Stone2: { \ + typedef SatFuncStone2 SatFunc; \ + __attribute__((unused)) CONST SatFunc& satFunc = \ + *static_cast(satFunc_.get()); \ + codeToCall; \ + break; \ + } \ + case Simple: { \ + typedef SatFuncSimple SatFunc; \ + __attribute__((unused)) CONST SatFunc& satFunc = \ + *static_cast(satFunc_.get()); \ + codeToCall; \ + break; \ + } \ }; SatFuncMultiplexer() - : wasInitialized_(false) + {} + + SatFuncMultiplexer(const SatFuncMultiplexer&) {} ~SatFuncMultiplexer() - { - if (!wasInitialized_) - return; - - // call the destructor of the selected saturation function - OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.~SatFunc()); - } + {} /*! * \brief Pick the correct saturation function type and initialize the object using @@ -103,7 +103,7 @@ public: OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.init(eclState, tableIdx, phaseUsage, - /*samples=*/0)); + /*samples=*/0), OPM_MULTIPLEXER_NON_CONST); } /*! @@ -111,14 +111,27 @@ public: */ void setSatFuncType(SatFuncType newSatFuncType) { - assert(!wasInitialized_); - // set the type of the saturation function satFuncType_ = newSatFuncType; // call the default constructor for the selected saturation function - OPM_MULTIPLEXER_SATFUNC_CALL(new(&satFunc) SatFunc()); - wasInitialized_ = true; + switch (satFuncType_) { + case Gwseg: { + typedef SatFuncGwseg SatFunc; + satFunc_.reset(new SatFunc); + break; + } + case Stone2: { + typedef SatFuncStone2 SatFunc; + satFunc_.reset(new SatFunc); + break; + } + case Simple: { + typedef SatFuncSimple SatFunc; + satFunc_.reset(new SatFunc); + break; + } + } } /*! @@ -135,27 +148,11 @@ public: * parameters). */ const SatFuncBase& getSatFuncBase() const - { - switch (satFuncType()) { - case Gwseg: - return getGwseg(); - case Stone2: - return getStone2(); - case Simple: - return getSimple(); - } - OPM_THROW(std::logic_error, "Unhandled saturation function type"); - } + { return *satFunc_; } SatFuncBase& getSatFuncBase() - { - return const_cast&>(static_cast(this)->getSatFuncBase()); - } + { return *satFunc_; } -// the strict aliasing compiler warnings need to be ignored here. This is safe because we -// deal with the alignment of the data ourselfs. -#pragma GCC diagnostic push -#pragma GCC diagnostic ignored "-Wstrict-aliasing" /*! * \brief Return the raw "Gwseg" saturation function object. * @@ -165,13 +162,13 @@ public: SatFuncGwseg& getGwseg() { assert(satFuncType_ == Gwseg); - return *reinterpret_cast*>(satFuncsData_.gwseg); + return *static_cast*>(satFunc_.get()); } const SatFuncGwseg& getGwseg() const { assert(satFuncType_ == Gwseg); - return *reinterpret_cast*>(satFuncsData_.gwseg); + return *static_cast*>(satFunc_.get()); } /*! @@ -183,13 +180,13 @@ public: SatFuncSimple& getSimple() { assert(satFuncType_ == Simple); - return *reinterpret_cast*>(satFuncsData_.simple); + return *static_cast*>(satFunc_.get()); } const SatFuncSimple& getSimple() const { assert(satFuncType_ == Simple); - return *reinterpret_cast*>(satFuncsData_.simple); + return *static_cast*>(satFunc_.get()); } /*! @@ -201,67 +198,63 @@ public: SatFuncStone2& getStone2() { assert(satFuncType_ == Stone2); - return *reinterpret_cast*>(satFuncsData_.stone2); + return *static_cast*>(satFunc_.get()); } const SatFuncStone2& getStone2() const { assert(satFuncType_ == Stone2); - return *reinterpret_cast*>(satFuncsData_.stone2); + return *static_cast*>(satFunc_.get()); } -#pragma GCC diagnostic pop template void evalKr(const FluidState& fluidState, double* kr) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr), OPM_MULTIPLEXER_CONST); } template void evalKrDeriv(const FluidState& fluidState, double* kr, double* dkrds) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds), OPM_MULTIPLEXER_CONST); } template void evalPc(const FluidState& fluidState, double* pc) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPc(fluidState, pc)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPc(fluidState, pc), OPM_MULTIPLEXER_CONST); } template void evalPcDeriv(const FluidState& fluidState, double* pc, double* dpcds) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPcDeriv(fluidState, pc, dpcds)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPcDeriv(fluidState, pc, dpcds), OPM_MULTIPLEXER_CONST); } template void evalKr(const FluidState& fluidState, double* kr, const EPSTransforms* epst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr, epst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr, epst), OPM_MULTIPLEXER_CONST); } template void evalKr(const FluidState& fluidState, double* kr, const EPSTransforms* epst, const EPSTransforms* epst_hyst, const SatHyst* sat_hyst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr, epst, epst_hyst, sat_hyst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKr(fluidState, kr, epst, epst_hyst, sat_hyst), OPM_MULTIPLEXER_CONST); } template void evalKrDeriv(const FluidState& fluidState, double* kr, double* dkrds, const EPSTransforms* epst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds, epst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds, epst), OPM_MULTIPLEXER_CONST); } template void evalKrDeriv(const FluidState& fluidState, double* kr, double* dkrds, const EPSTransforms* epst, const EPSTransforms* epst_hyst, const SatHyst* sat_hyst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds, epst, epst_hyst, sat_hyst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalKrDeriv(fluidState, kr, dkrds, epst, epst_hyst, sat_hyst), OPM_MULTIPLEXER_CONST); } template void evalPc(const FluidState& fluidState, double* pc, const EPSTransforms* epst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPc(fluidState, pc, epst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPc(fluidState, pc, epst), OPM_MULTIPLEXER_CONST); } template void evalPcDeriv(const FluidState& fluidState, double* pc, double* dpcds, const EPSTransforms* epst) const - { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPcDeriv(fluidState, pc, dpcds, epst)); } + { OPM_MULTIPLEXER_SATFUNC_CALL(satFunc.evalPcDeriv(fluidState, pc, dpcds, epst), OPM_MULTIPLEXER_CONST); } // undefine the helper macro here because we don't need it anymore #undef OPM_MULTIPLEXER_SATFUNC_CALL +#undef OPM_MULTIPLEXER_NON_CONST +#undef OPM_MULTIPLEXER_CONST private: - bool wasInitialized_; SatFuncType satFuncType_; - union { - char gwseg[sizeof(SatFuncGwseg)]; - char simple[sizeof(SatFuncSimple)]; - char stone2[sizeof(SatFuncStone2)]; - } satFuncsData_ __attribute__((aligned(sizeof(double)))); + std::unique_ptr > satFunc_; }; } // namespace Opm