12#ifndef DUMUX_PLOT_FLUID_MATRIX_LAW_HH 
   13#define DUMUX_PLOT_FLUID_MATRIX_LAW_HH 
   29template<
class Scalar, 
class MaterialLaw>
 
   32    using MaterialLawParams = 
typename MaterialLaw::Params;
 
   50                const MaterialLawParams ¶ms,
 
   51                Scalar lowerSat = 0.0,
 
   52                Scalar upperSat = 1.0,
 
   53                std::string curveTitle = 
"")
 
   55        addpcgw(gnuplot, params, lowerSat, upperSat, curveTitle);
 
   56        addpcnw(gnuplot, params, lowerSat, upperSat, curveTitle);
 
   57        addpcgn(gnuplot, params, lowerSat, upperSat, curveTitle);
 
 
   70                  const MaterialLawParams ¶ms,
 
   71                  Scalar lowerSat = 0.0,
 
   72                  Scalar upperSat = 1.0,
 
   73                  std::string curveTitle = 
"")
 
   75        std::vector<Scalar> sw;
 
   76        std::vector<Scalar> pc;
 
   77        Scalar satInterval = upperSat - lowerSat;
 
   79        Scalar pcMax = -1e100;
 
   81        Scalar swTemp, pcTemp = 0.0;
 
   82        for (
int i = 0; i <= numIntervals_; i++)
 
   84            swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_);
 
   85            pcTemp = MaterialLaw::pcgw(params, swTemp);
 
   88            if (checkValues_(swTemp, pcTemp))
 
   92                pcMin = min(pcMin, pcTemp);
 
   93                pcMax = max(pcMax, pcTemp);
 
   97        gnuplot.
setXlabel(
"wetting phase saturation [-]");
 
   98        gnuplot.
setYlabel(
"capillary pressure [Pa]");
 
 
  112                  const MaterialLawParams ¶ms,
 
  113                  Scalar lowerSat = 0.0,
 
  114                  Scalar upperSat = 1.0,
 
  115                  std::string curveTitle = 
"")
 
  117        std::vector<Scalar> sw;
 
  118        std::vector<Scalar> pc;
 
  119        Scalar satInterval = upperSat - lowerSat;
 
  121        Scalar pcMax = -1e100;
 
  123        Scalar swTemp, pcTemp = 0.0;
 
  124        for (
int i = 0; i <= numIntervals_; i++)
 
  126            swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_);
 
  127            pcTemp = MaterialLaw::pcnw(params, swTemp);
 
  130            if (checkValues_(swTemp, pcTemp))
 
  132                sw.push_back(swTemp);
 
  133                pc.push_back(pcTemp);
 
  134                pcMin = min(pcMin, pcTemp);
 
  135                pcMax = max(pcMax, pcTemp);
 
  139        gnuplot.
setXlabel(
"wetting phase saturation [-]");
 
  140        gnuplot.
setYlabel(
"capillary pressure [Pa]");
 
 
  154                  const MaterialLawParams ¶ms,
 
  155                  Scalar lowerSat = 0.0,
 
  156                  Scalar upperSat = 1.0,
 
  157                  std::string curveTitle = 
"")
 
  159        std::vector<Scalar> st;
 
  160        std::vector<Scalar> pc;
 
  161        Scalar satInterval = upperSat - lowerSat;
 
  163        Scalar pcMax = -1e100;
 
  165        Scalar stTemp, pcTemp = 0.0;
 
  166        for (
int i = 0; i <= numIntervals_; i++)
 
  168            stTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_);
 
  169            pcTemp = MaterialLaw::pcgn(params, stTemp);
 
  172            if (checkValues_(stTemp, pcTemp))
 
  174                st.push_back(stTemp);
 
  175                pc.push_back(pcTemp);
 
  176                pcMin = min(pcMin, pcTemp);
 
  177                pcMax = max(pcMax, pcTemp);
 
  181        gnuplot.
setXlabel(
"wetting phase saturation [-]");
 
  182        gnuplot.
setYlabel(
"capillary pressure [Pa]");
 
 
  197                const MaterialLawParams ¶ms,
 
  198                Scalar lowerSat = 0.0,
 
  199                Scalar upperSat = 1.0,
 
  200                std::string curveTitle = 
"")
 
  202        std::vector<Scalar> sw(numIntervals_ + 1);
 
  203        std::vector<Scalar> krw(numIntervals_ + 1);
 
  204        std::vector<Scalar> krn(numIntervals_ + 1);
 
  205        std::vector<Scalar> krg(numIntervals_ + 1);
 
  206        Scalar satInterval = upperSat - lowerSat;
 
  207        Scalar krMin = 1e100;
 
  208        Scalar krMax = -1e100;
 
  210        Scalar swTemp, krwTemp, krnTemp, krgTemp = 0.0;
 
  211        for (
int i = 0; i <= numIntervals_; i++)
 
  213            swTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_);
 
  214            krwTemp = MaterialLaw::krw(params, swTemp, 0.0);
 
  215            krnTemp = MaterialLaw::krn(params, swTemp, 1.0 - swTemp);
 
  216            krgTemp = MaterialLaw::krg(params, swTemp, 0.0);
 
  219            if (checkValues_(swTemp, krwTemp)
 
  220                && checkValues_(swTemp, krnTemp)
 
  221                && checkValues_(swTemp, krgTemp))
 
  223                sw.push_back(swTemp);
 
  224                krw.push_back(krwTemp);
 
  225                krn.push_back(krnTemp);
 
  226                krg.push_back(krgTemp);
 
  227                krMin = min({krMin, krwTemp, krnTemp, krgTemp});
 
  228                krMax = max({krMax, krwTemp, krnTemp, krgTemp});
 
  232        gnuplot.
setXlabel(
"wetting phase saturation [-]");
 
  233        gnuplot.
setYlabel(
"relative permeability [-]");
 
 
  249                     const MaterialLawParams ¶ms,
 
  250                     Scalar lowerSat = 0.0,
 
  251                     Scalar upperSat = 1.0,
 
  252                     std::string curveTitle = 
"")
 
  254        std::vector<Scalar> sn(numIntervals_ + 1);
 
  255        std::vector<Scalar> alpha(numIntervals_ + 1);
 
  256        Scalar satInterval = upperSat - lowerSat;
 
  257        Scalar alphaMin = -2;
 
  260        Scalar snTemp, alphaTemp = 0.0;
 
  261        for (
int i = 0; i <= numIntervals_; i++)
 
  263            snTemp = lowerSat + satInterval * Scalar(i) / Scalar(numIntervals_);
 
  264            alphaTemp = MaterialLaw::pcAlpha(params, snTemp);
 
  267            if (checkValues_(snTemp, alphaTemp))
 
  269                sn.push_back(snTemp);
 
  270                alpha.push_back(alphaTemp);
 
  271                alphaMin = min(alphaMin, alphaTemp);
 
  272                alphaMax = max(alphaMax, alphaTemp);
 
  276        gnuplot.
setXlabel(
"nonwetting phase saturation [-]");
 
  277        gnuplot.
setYlabel(
"transition function [-]");
 
 
  288    bool checkValues_(Scalar value1, Scalar value2)
 
  292        return !isnan(value1) && !isinf(value1)
 
  293               && !isnan(value2) && !isinf(value2);
 
 
Interface for passing data sets to a file and plotting them, if gnuplot is installed.
Definition gnuplotinterface.hh:45
void setXlabel(const std::string &label)
Sets the label for the x-axis.
Definition gnuplotinterface.hh:267
void addDataSetToPlot(const DataX &x, const DataY &y, const std::string &fileName, const std::string &options="with lines")
Adds a data set and writes a data file.
Definition gnuplotinterface.hh:231
void setYlabel(const std::string &label)
Sets the label for the y-axis.
Definition gnuplotinterface.hh:277
void addPcAlpha(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the transition (2P/3P) function.
Definition plotmateriallaw3p.hh:248
void addpc(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the capillary pressure-saturation curve for all phases.
Definition plotmateriallaw3p.hh:49
void addpcgw(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the capillary pressure-saturation curve for the water-gas interphase.
Definition plotmateriallaw3p.hh:69
void addpcgn(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the capillary pressure-saturation curve for the gas-NAPL interface.
Definition plotmateriallaw3p.hh:153
PlotMaterialLaw(bool interaction=true)
Constructor.
Definition plotmateriallaw3p.hh:36
void addkr(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the relative permeabilities.
Definition plotmateriallaw3p.hh:196
void addpcnw(GnuplotInterface< Scalar > &gnuplot, const MaterialLawParams ¶ms, Scalar lowerSat=0.0, Scalar upperSat=1.0, std::string curveTitle="")
Plot the capillary pressure-saturation curve for the water-NAPL interface.
Definition plotmateriallaw3p.hh:111