From 1fefb44867eca11d2ad7911ae71a71e84ad0dbd4 Mon Sep 17 00:00:00 2001 From: Thomas Fetzer Date: Wed, 8 Jun 2016 14:34:19 +0200 Subject: [PATCH 1/4] [plotmateriallaw] Remove superfluous checking warning Now the values are checked directly and only finite values are passed to gnuplot --- dumux/io/plotmateriallaw.hh | 143 +++++++++++++++++++++++----------- dumux/io/plotmateriallaw3p.hh | 112 ++++++++++++++++++-------- 2 files changed, 174 insertions(+), 81 deletions(-) diff --git a/dumux/io/plotmateriallaw.hh b/dumux/io/plotmateriallaw.hh index 7271fe152f..9be34713eb 100644 --- a/dumux/io/plotmateriallaw.hh +++ b/dumux/io/plotmateriallaw.hh @@ -24,6 +24,8 @@ #ifndef DUMUX_PLOT_FLUID_MATRIX_LAW_HH #define DUMUX_PLOT_FLUID_MATRIX_LAW_HH +#include + #include #include @@ -73,19 +75,24 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sw(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector sw; + std::vector pc; Scalar satInterval = upperSat - lowerSat; Scalar pcMin = 0.0; Scalar pcMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar swTemp, pcTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - pc[i] = MaterialLaw::pc(params, sw[i]); - pcMin = std::min(pcMin, pc[i]); - pcMax = std::max(pcMax, pc[i]); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + pcTemp = MaterialLaw::pc(params, swTemp); + if (checkValues(swTemp, pcTemp)) + { + sw.push_back(swTemp); + pc.push_back(pcTemp); + pcMin = std::min(pcMin, pcTemp); + pcMax = std::max(pcMax, pcTemp); + } } gnuplotpcsw_.setXRange(lowerSat, upperSat); @@ -109,25 +116,31 @@ public: Scalar upperpc = 5000.0, std::string plotName = "") { - std::vector sat(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector sw; + std::vector pc; Scalar pcInterval = upperpc - lowerpc; Scalar swMin = 1e100; Scalar swMax = -1e100; + Scalar pcTemp, swTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - pc[i] = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); - sat[i] = MaterialLaw::sw(params, pc[i]); - swMin = std::min(swMin, sat[i]); - swMax = std::max(swMax, sat[i]); + pcTemp = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); + swTemp = MaterialLaw::sw(params, pcTemp); + if (checkValues(pcTemp, swTemp)) + { + pc.push_back(pcTemp); + sw.push_back(swTemp); + swMin = std::min(swMin, swTemp); + swMax = std::max(swMax, swTemp); + } } gnuplotswpc_.setXRange(lowerpc, upperpc); gnuplotswpc_.setYRange(swMin, swMax); gnuplotswpc_.setXlabel("capillary pressure [Pa]"); gnuplotswpc_.setYlabel("wetting phase saturation [-]"); - gnuplotswpc_.addDataSetToPlot(pc, sat, plotName + "_Sw-pc"); + gnuplotswpc_.addDataSetToPlot(pc, sw, plotName + "_Sw-pc"); gnuplotswpc_.plot("sw-pc"); } @@ -144,26 +157,31 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sat(numIntervals_ + 1); - std::vector dpcdsw(numIntervals_ + 1); + std::vector sw; + std::vector dpcdsw; Scalar satInterval = upperSat - lowerSat; Scalar dpcdswMin = 1e100; Scalar dpcdswMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar swTemp, dpcdswTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sat[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - dpcdsw[i] = MaterialLaw::dpc_dsw(params, sat[i]); - dpcdswMin = std::min(dpcdswMin, dpcdsw[i]); - dpcdswMax = std::max(dpcdswMax, dpcdsw[i]); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + dpcdswTemp = MaterialLaw::dpc_dsw(params, swTemp); + if (checkValues(swTemp, dpcdsw)) + { + sw.push_back(swTemp); + dpcdsw.push_back(dpcdswTemp); + dpcdswMin = std::min(dpcdswMin, dpcdswTemp); + dpcdswMax = std::max(dpcdswMax, dpcdswTemp); + } } gnuplotdpcdsw_.setXRange(lowerSat, upperSat); gnuplotdpcdsw_.setYRange(dpcdswMin, dpcdswMax); gnuplotdpcdsw_.setXlabel("wetting phase saturation [-]"); gnuplotdpcdsw_.setYlabel("gradient of the pc-Sw curve [Pa]"); - gnuplotdpcdsw_.addDataSetToPlot(sat, dpcdsw, plotName + "_dpcdSw-Sw"); + gnuplotdpcdsw_.addDataSetToPlot(sw, dpcdsw, plotName + "_dpcdSw-Sw"); gnuplotdpcdsw_.plot("dpcdsw"); } @@ -180,18 +198,24 @@ public: Scalar upperpc = 5000.0, std::string plotName = "") { - std::vector dswdpc(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector pc; + std::vector dswdpc; Scalar pcInterval = upperpc - lowerpc; Scalar dswdpcMin = 1e100; Scalar dswdpcMax = -1e100; + Scalar dswdpcTemp, pcTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - pc[i] = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); - dswdpc[i] = MaterialLaw::dsw_dpc(params, pc[i]); - dswdpcMin = std::min(dswdpcMin, dswdpc[i]); - dswdpcMax = std::max(dswdpcMax, dswdpc[i]); + pcTemp = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); + dswdpcTemp = MaterialLaw::dsw_dpc(params, pcTemp); + if (checkValues(pcTemp, dswdpcTemp)) + { + pc.push_back(pcTemp); + dswdpc.push_back(dswdpcTemp); + dswdpcMin = std::min(dswdpcMin, dswdpcTemp); + dswdpcMax = std::max(dswdpcMax, dswdpcTemp); + } } gnuplotdswdpc_.setXRange(lowerpc, upperpc); @@ -215,21 +239,27 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sw(numIntervals_ + 1); - std::vector krw(numIntervals_ + 1); - std::vector krn(numIntervals_ + 1); + std::vector sw; + std::vector krw; + std::vector krn; Scalar satInterval = upperSat - lowerSat; Scalar krMin = 1e100; Scalar krMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName + "_kr"); + Scalar swTemp, krwTemp, krnTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - krw[i] = MaterialLaw::krw(params, sw[i]); - krn[i] = MaterialLaw::krn(params, sw[i]); - krMin = std::min(krMin, std::min(krw[i], krn[i])); - krMax = std::max(krMax, std::max(krw[i], krn[i])); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + krwTemp = MaterialLaw::krw(params, swTemp); + krnTemp = MaterialLaw::krn(params, swTemp); + if (checkValues(swTemp, krwTemp) && checkValues(swTemp, krnTemp)) + { + sw.push_back(swTemp); + krw.push_back(krwTemp); + krn.push_back(krnTemp); + krMin = std::min({krMin, krwTemp, krnTemp}); + krMax = std::max({krMax, krwTemp, krnTemp}); + } } gnuplotkr_.setXRange(lowerSat, upperSat); @@ -254,21 +284,27 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sw(numIntervals_+1); - std::vector dkrw_dsw(numIntervals_+1); - std::vector dkrn_dsw(numIntervals_+1); + std::vector sw; + std::vector dkrw_dsw; + std::vector dkrn_dsw; Scalar satInterval = upperSat - lowerSat; Scalar dkrdswMin = 1e100; Scalar dkrdswMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName + "_dkr_dsw"); + Scalar swTemp, dkrwdswTemp, dkrndswTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - dkrw_dsw[i] = MaterialLaw::dkrw_dsw(params, sw[i]); - dkrn_dsw[i] = MaterialLaw::dkrn_dsw(params, sw[i]); - dkrdswMin = std::min(dkrdswMin, std::min(dkrw_dsw[i], dkrn_dsw[i])); - dkrdswMax = std::max(dkrdswMax, std::max(dkrw_dsw[i], dkrn_dsw[i])); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + dkrwdswTemp = MaterialLaw::dkrw_dsw(params, swTemp); + dkrndswTemp = MaterialLaw::dkrn_dsw(params, swTemp); + if (checkValues(swTemp, dkrwdswTemp) && checkValues(swTemp, dkrndswTemp)) + { + sw.push_back(swTemp); + dkrw_dsw.push_back(dkrwdswTemp); + dkrn_dsw.push_back(dkrndswTemp); + dkrdswMin = std::min({dkrdswMin, dkrwdswTemp, dkrndswTemp}); + dkrdswMax = std::max({dkrdswMax, dkrwdswTemp, dkrndswTemp}); + } } gnuplotkrdsw_.setXRange(lowerSat, upperSat); @@ -280,6 +316,7 @@ public: gnuplotkrdsw_.plot("dkrndsw"); } +private: /*! * \brief Check the validity range for wetting saturation, to avoid an * assert of the used material laws @@ -289,6 +326,7 @@ public: * \param upperSat Maximum x-value * \param plotName Name of the plotted curve */ + DUNE_DEPRECATED_MSG("checkEffectiveSaturation is deprecated.") void checkEffectiveSaturation(const MaterialLawParams ¶ms, Scalar lowerSat, Scalar upperSat, @@ -300,7 +338,18 @@ public: Dune::dwarn << "warning: fluid-matrix law " << plotName << " can only be plotted for sw < 1.0 - snr" << std::endl; } -private: + /*! + * \brief Check the values for occurrences of nan and inf + * + * \param value1 A data point value + * \param value2 An other data point value + */ + bool checkValues(Scalar value1, Scalar value2) + { + return !std::isnan(value1) && !std::isinf(value1) + && !std::isnan(value2) && !std::isinf(value2); + } + int numIntervals_; GnuplotInterface gnuplotpcsw_; GnuplotInterface gnuplotswpc_; diff --git a/dumux/io/plotmateriallaw3p.hh b/dumux/io/plotmateriallaw3p.hh index 1cb520bb75..3c9f025287 100644 --- a/dumux/io/plotmateriallaw3p.hh +++ b/dumux/io/plotmateriallaw3p.hh @@ -24,6 +24,8 @@ #ifndef DUMUX_PLOT_FLUID_MATRIX_LAW_HH #define DUMUX_PLOT_FLUID_MATRIX_LAW_HH +#include + #include #include @@ -91,19 +93,24 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sw(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector sw; + std::vector pc; Scalar satInterval = upperSat - lowerSat; Scalar pcMin = 0.0; Scalar pcMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar swTemp, pcTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - pc[i] = MaterialLaw::pcgw(params, sw[i]); - pcMin = std::min(pcMin, pc[i]); - pcMax = std::max(pcMax, pc[i]); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + pcTemp = MaterialLaw::pcgw(params, swTemp); + if (checkValues(swTemp, pcTemp)) + { + sw.push_back(swTemp); + pc.push_back(pcTemp); + pcMin = std::min(pcMin, pcTemp); + pcMax = std::max(pcMax, pcTemp); + } } gnuplotpc_.setXRange(lowerSat, upperSat); @@ -127,19 +134,24 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector sw(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector sw; + std::vector pc; Scalar satInterval = upperSat - lowerSat; Scalar pcMin = 0.0; Scalar pcMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar swTemp, pcTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - pc[i] = MaterialLaw::pcnw(params, sw[i]); - pcMin = std::min(pcMin, pc[i]); - pcMax = std::max(pcMax, pc[i]); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + pcTemp = MaterialLaw::pcnw(params, swTemp); + if (checkValues(swTemp, pcTemp)) + { + sw.push_back(swTemp); + pc.push_back(pcTemp); + pcMin = std::min(pcMin, pcTemp); + pcMax = std::max(pcMax, pcTemp); + } } gnuplotpc_.setXRange(lowerSat, upperSat); @@ -163,19 +175,24 @@ public: Scalar upperSat = 1.0, std::string plotName = "") { - std::vector st(numIntervals_+1); - std::vector pc(numIntervals_+1); + std::vector st; + std::vector pc; Scalar satInterval = upperSat - lowerSat; Scalar pcMin = 0.0; Scalar pcMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar stTemp, pcTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - st[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - pc[i] = MaterialLaw::pcgn(params, st[i]); - pcMin = std::min(pcMin, pc[i]); - pcMax = std::max(pcMax, pc[i]); + stTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + pcTemp = MaterialLaw::pcgn(params, stTemp); + if (checkValues(stTemp, pcTemp)) + { + st.push_back(stTemp); + pc.push_back(pcTemp); + pcMin = std::min(pcMin, pcTemp); + pcMax = std::max(pcMax, pcTemp); + } } gnuplotpc_.setXRange(lowerSat, upperSat); @@ -207,16 +224,25 @@ public: Scalar satInterval = upperSat - lowerSat; Scalar krMin = 1e100; Scalar krMax = -1e100; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName + "_kr"); + Scalar swTemp, krwTemp, krnTemp, krgTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sw[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - krw[i] = MaterialLaw::krw(params, sw[i], 0.0); - krn[i] = MaterialLaw::krn(params, sw[i], 1-sw[i]); - krg[i] = MaterialLaw::krg(params, sw[i], 0.0); - krMin = std::min(krMin, std::min({krw[i], krn[i], krg[i]})); - krMax = std::max(krMax, std::max({krw[i], krn[i], krg[i]})); + swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + krwTemp = MaterialLaw::krw(params, swTemp); + krnTemp = MaterialLaw::krn(params, swTemp); + krgTemp = MaterialLaw::krg(params, swTemp); + if (checkValues(swTemp, krwTemp) + && checkValues(swTemp, krnTemp) + && checkValues(swTemp, krgTemp)) + { + sw.push_back(swTemp); + krw.push_back(krwTemp); + krn.push_back(krnTemp); + krn.push_back(krgTemp); + krMin = std::min({krMin, krwTemp, krnTemp, krgTemp}); + krMax = std::max({krMax, krwTemp, krnTemp, krgTemp}); + } } gnuplotkr_.setXRange(lowerSat, upperSat); @@ -247,14 +273,19 @@ public: Scalar satInterval = upperSat - lowerSat; Scalar alphaMin = -2; Scalar alphaMax = 2; - checkEffectiveSaturation(params, lowerSat, upperSat, plotName); + Scalar snTemp, alphaTemp = 0.0; for (int i = 0; i <= numIntervals_; i++) { - sn[i] = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - alpha[i] = MaterialLaw::pcAlpha(params, sn[i]); - alphaMin = std::min(alphaMin, alpha[i]); - alphaMax = std::max(alphaMax, alpha[i]); + snTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); + alphaTemp = MaterialLaw::pcAlpha(params, snTemp); + if (checkValues(snTemp, alphaTemp)) + { + sn.push_back(snTemp); + alpha.push_back(alphaTemp); + alphaMin = std::min(alphaMin, alphaTemp); + alphaMax = std::max(alphaMax, alphaTemp); + } } gnuplotpcAlpha_.setXRange(lowerSat, upperSat); @@ -265,6 +296,7 @@ public: gnuplotpcAlpha_.plot("alpha"); } +private: /*! * \brief Check the validity range for wetting saturation, to avoid an * assert of the used material laws @@ -274,6 +306,7 @@ public: * \param upperSat Maximum x-value * \param plotName Name of the plotted curve */ + DUNE_DEPRECATED_MSG("checkEffectiveSaturation is deprecated.") void checkEffectiveSaturation(const MaterialLawParams ¶ms, Scalar lowerSat, Scalar upperSat, @@ -285,7 +318,18 @@ public: Dune::dwarn << "warning: fluid-matrix law " << plotName << " can only be plotted for sw < 1.0 - snr" << std::endl; } -private: + /*! + * \brief Check the values for occurrences of nan and inf + * + * \param value1 A data point value + * \param value2 An other data point value + */ + bool checkValues(Scalar value1, Scalar value2) + { + return !std::isnan(value1) && !std::isinf(value1) + && !std::isnan(value2) && !std::isinf(value2); + } + int numIntervals_; GnuplotInterface gnuplotpc_; GnuplotInterface gnuplotpcAlpha_; -- GitLab From 9391d1b312bc2adba8cc5c911f735facc8b4a571 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Wed, 8 Jun 2016 17:40:39 +0200 Subject: [PATCH 2/4] [plotmateriallaw] Fix function arguments and typo * The 3p kr() requires different arguments as the 2p kr() * Wrong vector was pushed back --- dumux/io/plotmateriallaw3p.hh | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/dumux/io/plotmateriallaw3p.hh b/dumux/io/plotmateriallaw3p.hh index 3c9f025287..ef882d3950 100644 --- a/dumux/io/plotmateriallaw3p.hh +++ b/dumux/io/plotmateriallaw3p.hh @@ -229,9 +229,9 @@ public: for (int i = 0; i <= numIntervals_; i++) { swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); - krwTemp = MaterialLaw::krw(params, swTemp); - krnTemp = MaterialLaw::krn(params, swTemp); - krgTemp = MaterialLaw::krg(params, swTemp); + krwTemp = MaterialLaw::krw(params, swTemp, 0.0); + krnTemp = MaterialLaw::krn(params, swTemp, 1.0 - swTemp); + krgTemp = MaterialLaw::krg(params, swTemp, 0.0); if (checkValues(swTemp, krwTemp) && checkValues(swTemp, krnTemp) && checkValues(swTemp, krgTemp)) @@ -239,7 +239,7 @@ public: sw.push_back(swTemp); krw.push_back(krwTemp); krn.push_back(krnTemp); - krn.push_back(krgTemp); + krg.push_back(krgTemp); krMin = std::min({krMin, krwTemp, krnTemp, krgTemp}); krMax = std::max({krMax, krwTemp, krnTemp, krgTemp}); } -- GitLab From f717aba3aa433c8bdabe4b7265aa891f16442114 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Thu, 9 Jun 2016 08:13:15 +0200 Subject: [PATCH 3/4] [plotmateriallaw] Correct private / public stuff --- dumux/io/plotmateriallaw.hh | 16 ++++++++-------- dumux/io/plotmateriallaw3p.hh | 18 +++++++++--------- 2 files changed, 17 insertions(+), 17 deletions(-) diff --git a/dumux/io/plotmateriallaw.hh b/dumux/io/plotmateriallaw.hh index 9be34713eb..ef87efa516 100644 --- a/dumux/io/plotmateriallaw.hh +++ b/dumux/io/plotmateriallaw.hh @@ -86,7 +86,7 @@ public: { swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); pcTemp = MaterialLaw::pc(params, swTemp); - if (checkValues(swTemp, pcTemp)) + if (checkValues_(swTemp, pcTemp)) { sw.push_back(swTemp); pc.push_back(pcTemp); @@ -127,7 +127,7 @@ public: { pcTemp = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); swTemp = MaterialLaw::sw(params, pcTemp); - if (checkValues(pcTemp, swTemp)) + if (checkValues_(pcTemp, swTemp)) { pc.push_back(pcTemp); sw.push_back(swTemp); @@ -168,7 +168,7 @@ public: { swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); dpcdswTemp = MaterialLaw::dpc_dsw(params, swTemp); - if (checkValues(swTemp, dpcdsw)) + if (checkValues_(swTemp, dpcdsw)) { sw.push_back(swTemp); dpcdsw.push_back(dpcdswTemp); @@ -209,7 +209,7 @@ public: { pcTemp = lowerpc + pcInterval * Scalar(i) / Scalar(numIntervals_); dswdpcTemp = MaterialLaw::dsw_dpc(params, pcTemp); - if (checkValues(pcTemp, dswdpcTemp)) + if (checkValues_(pcTemp, dswdpcTemp)) { pc.push_back(pcTemp); dswdpc.push_back(dswdpcTemp); @@ -252,7 +252,7 @@ public: swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); krwTemp = MaterialLaw::krw(params, swTemp); krnTemp = MaterialLaw::krn(params, swTemp); - if (checkValues(swTemp, krwTemp) && checkValues(swTemp, krnTemp)) + if (checkValues_(swTemp, krwTemp) && checkValues_(swTemp, krnTemp)) { sw.push_back(swTemp); krw.push_back(krwTemp); @@ -297,7 +297,7 @@ public: swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); dkrwdswTemp = MaterialLaw::dkrw_dsw(params, swTemp); dkrndswTemp = MaterialLaw::dkrn_dsw(params, swTemp); - if (checkValues(swTemp, dkrwdswTemp) && checkValues(swTemp, dkrndswTemp)) + if (checkValues_(swTemp, dkrwdswTemp) && checkValues_(swTemp, dkrndswTemp)) { sw.push_back(swTemp); dkrw_dsw.push_back(dkrwdswTemp); @@ -316,7 +316,6 @@ public: gnuplotkrdsw_.plot("dkrndsw"); } -private: /*! * \brief Check the validity range for wetting saturation, to avoid an * assert of the used material laws @@ -338,13 +337,14 @@ private: Dune::dwarn << "warning: fluid-matrix law " << plotName << " can only be plotted for sw < 1.0 - snr" << std::endl; } +private: /*! * \brief Check the values for occurrences of nan and inf * * \param value1 A data point value * \param value2 An other data point value */ - bool checkValues(Scalar value1, Scalar value2) + bool checkValues_(Scalar value1, Scalar value2) { return !std::isnan(value1) && !std::isinf(value1) && !std::isnan(value2) && !std::isinf(value2); diff --git a/dumux/io/plotmateriallaw3p.hh b/dumux/io/plotmateriallaw3p.hh index ef882d3950..7da7ab693d 100644 --- a/dumux/io/plotmateriallaw3p.hh +++ b/dumux/io/plotmateriallaw3p.hh @@ -104,7 +104,7 @@ public: { swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); pcTemp = MaterialLaw::pcgw(params, swTemp); - if (checkValues(swTemp, pcTemp)) + if (checkValues_(swTemp, pcTemp)) { sw.push_back(swTemp); pc.push_back(pcTemp); @@ -145,7 +145,7 @@ public: { swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); pcTemp = MaterialLaw::pcnw(params, swTemp); - if (checkValues(swTemp, pcTemp)) + if (checkValues_(swTemp, pcTemp)) { sw.push_back(swTemp); pc.push_back(pcTemp); @@ -186,7 +186,7 @@ public: { stTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); pcTemp = MaterialLaw::pcgn(params, stTemp); - if (checkValues(stTemp, pcTemp)) + if (checkValues_(stTemp, pcTemp)) { st.push_back(stTemp); pc.push_back(pcTemp); @@ -232,9 +232,9 @@ public: krwTemp = MaterialLaw::krw(params, swTemp, 0.0); krnTemp = MaterialLaw::krn(params, swTemp, 1.0 - swTemp); krgTemp = MaterialLaw::krg(params, swTemp, 0.0); - if (checkValues(swTemp, krwTemp) - && checkValues(swTemp, krnTemp) - && checkValues(swTemp, krgTemp)) + if (checkValues_(swTemp, krwTemp) + && checkValues_(swTemp, krnTemp) + && checkValues_(swTemp, krgTemp)) { sw.push_back(swTemp); krw.push_back(krwTemp); @@ -279,7 +279,7 @@ public: { snTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_); alphaTemp = MaterialLaw::pcAlpha(params, snTemp); - if (checkValues(snTemp, alphaTemp)) + if (checkValues_(snTemp, alphaTemp)) { sn.push_back(snTemp); alpha.push_back(alphaTemp); @@ -296,7 +296,6 @@ public: gnuplotpcAlpha_.plot("alpha"); } -private: /*! * \brief Check the validity range for wetting saturation, to avoid an * assert of the used material laws @@ -318,13 +317,14 @@ private: Dune::dwarn << "warning: fluid-matrix law " << plotName << " can only be plotted for sw < 1.0 - snr" << std::endl; } +private: /*! * \brief Check the values for occurrences of nan and inf * * \param value1 A data point value * \param value2 An other data point value */ - bool checkValues(Scalar value1, Scalar value2) + bool checkValues_(Scalar value1, Scalar value2) { return !std::isnan(value1) && !std::isinf(value1) && !std::isnan(value2) && !std::isinf(value2); -- GitLab From 5caa39b34c7dc0db95294d206174b96c84bfd816 Mon Sep 17 00:00:00 2001 From: Kilian Weishaupt Date: Thu, 9 Jun 2016 08:14:23 +0200 Subject: [PATCH 4/4] [test] Include plotMaterialLaw in 3p test --- .../3p/implicit/infiltration3pproblem.hh | 2 ++ .../3p/implicit/infiltration3pspatialparams.hh | 18 ++++++++++++++++++ .../3p/implicit/test_box3p.input | 3 +++ .../3p/implicit/test_cc3p.input | 3 +++ 4 files changed, 26 insertions(+) diff --git a/test/porousmediumflow/3p/implicit/infiltration3pproblem.hh b/test/porousmediumflow/3p/implicit/infiltration3pproblem.hh index 7a73d40545..4dc9f78780 100644 --- a/test/porousmediumflow/3p/implicit/infiltration3pproblem.hh +++ b/test/porousmediumflow/3p/implicit/infiltration3pproblem.hh @@ -143,6 +143,8 @@ public: FluidSystem::init(282.15, 284.15, 3, 8e4, 3e5, 200); name_ = GET_RUNTIME_PARAM(TypeTag, std::string, Problem.Name); + + this->spatialParams().plotMaterialLaw(); } /*! diff --git a/test/porousmediumflow/3p/implicit/infiltration3pspatialparams.hh b/test/porousmediumflow/3p/implicit/infiltration3pspatialparams.hh index ea8b9fcc01..c7bc3e18a8 100644 --- a/test/porousmediumflow/3p/implicit/infiltration3pspatialparams.hh +++ b/test/porousmediumflow/3p/implicit/infiltration3pspatialparams.hh @@ -30,6 +30,7 @@ #include #include #include +#include namespace Dumux { @@ -117,11 +118,26 @@ public: // parameters for adsorption materialParams_.setKdNAPL(0.); materialParams_.setRhoBulk(1500.); + + plotFluidMatrixInteractions_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, bool, Output, + PlotFluidMatrixInteractions); } ~InfiltrationThreePSpatialParams() {} + /*! + * \brief This is called from the problem and creates a gnuplot output + * of e.g the pc-Sw curve + */ + void plotMaterialLaw() + { + PlotMaterialLaw plotMaterialLaw(plotFluidMatrixInteractions_); + + plotMaterialLaw.plotpc(materialParams_); + plotMaterialLaw.plotkr(materialParams_); + } + /*! * \brief Intrinsic permability * @@ -188,6 +204,8 @@ private: Scalar porosity_; MaterialLawParams materialParams_; + + bool plotFluidMatrixInteractions_; }; } diff --git a/test/porousmediumflow/3p/implicit/test_box3p.input b/test/porousmediumflow/3p/implicit/test_box3p.input index edd8a86f23..6445dff4ae 100644 --- a/test/porousmediumflow/3p/implicit/test_box3p.input +++ b/test/porousmediumflow/3p/implicit/test_box3p.input @@ -15,3 +15,6 @@ porosity = 0.40 vanGenuchtenAlpha = 0.0005 vanGenuchtenN = 4.0 +[Output] +PlotFluidMatrixInteractions = false + diff --git a/test/porousmediumflow/3p/implicit/test_cc3p.input b/test/porousmediumflow/3p/implicit/test_cc3p.input index 08652caa0f..764e1a3f0c 100644 --- a/test/porousmediumflow/3p/implicit/test_cc3p.input +++ b/test/porousmediumflow/3p/implicit/test_cc3p.input @@ -15,3 +15,6 @@ porosity = 0.40 vanGenuchtenAlpha = 0.0005 vanGenuchtenN = 4.0 +[Output] +PlotFluidMatrixInteractions = false + -- GitLab