#3042 OPM flowdiag apps update from repo. Commit 24ff768dc509b6c6bbd0121ef46a5932fae92961

This commit is contained in:
Bjørn Erik Jensen
2018-06-12 14:36:08 +02:00
parent 4491b1de63
commit 6ee84b598c
5 changed files with 1384 additions and 538 deletions

View File

@@ -31,9 +31,17 @@
#include <opm/utility/ECLEndPointScaling.hpp>
#include <opm/utility/ECLPropTable.hpp>
#include <opm/utility/ECLPvtCommon.hpp>
#include <opm/utility/ECLUnitHandling.hpp>
#include <opm/utility/imported/Units.hpp>
#include <exception>
#include <initializer_list>
#include <stdexcept>
#include <type_traits>
#include <utility>
#include <vector>
namespace {
@@ -70,6 +78,23 @@ namespace {
return sp;
}
std::vector<double>
makeTable(const std::size_t ncol,
std::initializer_list<double> data)
{
auto result = std::vector<double>(data.size(), 0.0);
const auto nrows = data.size() / ncol;
auto di = std::begin(data);
for (auto i = 0*nrows; i < nrows; ++i) {
for (auto j = 0*ncol; j < ncol; ++j, ++di) {
result[i + j*nrows] = *di;
}
}
return result;
}
std::vector<double> sw_core1d_example()
{
return {
@@ -179,6 +204,192 @@ namespace {
1.000000000000000e+00,
};
}
std::vector<double> sgfn_raw()
{
return makeTable(5, {
0, 0, 0, 0, 0,
5.000e-02, 1.655000000000000e-03, 0, 3.310000000000000e-02, 0,
1.000e-01, 6.913000000000000e-03, 0, 1.051600000000000e-01, 0,
1.500e-01, 1.621300000000000e-02, 0, 1.860000000000001e-01, 0,
2.000e-01, 2.999000000000000e-02, 0, 2.755399999999998e-01, 0,
2.500e-01, 4.865500000000000e-02, 0, 3.733000000000000e-01, 0,
3.000e-01, 7.257300000000000e-02, 0, 4.783600000000001e-01, 0,
3.500e-01, 1.020460000000000e-01, 0, 5.894600000000001e-01, 0,
4.000e-01, 1.372870000000000e-01, 0, 7.048199999999992e-01, 0,
4.500e-01, 1.784020000000000e-01, 0, 8.223000000000005e-01, 0,
5.000e-01, 2.253680000000000e-01, 0, 9.393200000000004e-01, 0,
5.500e-01, 2.780300000000000e-01, 0, 1.053239999999999e+00, 0,
6.000e-01, 3.360930000000000e-01, 0, 1.161260000000001e+00, 0,
6.500e-01, 3.991350000000000e-01, 0, 1.260840000000000e+00, 0,
7.000e-01, 4.666310000000000e-01, 0, 1.349920000000002e+00, 0,
7.500e-01, 5.380000000000000e-01, 0, 1.427379999999999e+00, 0,
8.000e-01, 6.126650000000000e-01, 0, 1.493299999999998e+00, 0,
8.500e-01, 6.901690000000000e-01, 0, 1.550080000000002e+00, 0,
9.000e-01, 7.703950000000001e-01, 0, 1.604519999999999e+00, 0,
9.500e-01, 8.542180000000000e-01, 0, 1.676460000000002e+00, 0,
9.999e-01, 9.499000000000000e-01, 0, 1.917474949899796e+00, 0,
1.000e+00, 9.500000000000000e-01, 0, 1.000000000000000e+00, 0,
});
}
std::vector<double> sofn_raw()
{
return makeTable(5, {
0, 0, 0, 0, 0,
9.999999999998899e-05, 0, 0, 0, 0,
5.000000000000004e-02, 1.000000000000000e-05, 7.000000000000000e-06, 1.999999999999998e-04, 1.402805611222443e-04,
9.999999999999998e-02, 1.200000000000000e-04, 8.800000000000000e-05, 2.200000000000003e-03, 1.620000000000002e-03,
1.500000000000000e-01, 5.100000000000000e-04, 3.840000000000000e-04, 7.799999999999994e-03, 5.919999999999996e-03,
2.000000000000000e-01, 1.490000000000000e-03, 1.117000000000000e-03, 1.960000000000003e-02, 1.466000000000002e-02,
2.500000000000000e-01, 3.460000000000000e-03, 2.597000000000000e-03, 3.939999999999996e-02, 2.959999999999997e-02,
3.000000000000000e-01, 6.990000000000000e-03, 5.254000000000000e-03, 7.059999999999993e-02, 5.313999999999995e-02,
3.500000000000000e-01, 1.284000000000000e-02, 9.662000000000000e-03, 1.170000000000002e-01, 8.816000000000013e-02,
4.000000000000000e-01, 2.199000000000000e-02, 1.658600000000000e-02, 1.829999999999998e-01, 1.384799999999999e-01,
4.500000000000000e-01, 3.572000000000000e-02, 2.703500000000000e-02, 2.746000000000004e-01, 2.089800000000003e-01,
5.000000000000000e-01, 5.565000000000000e-02, 4.232400000000000e-02, 3.985999999999996e-01, 3.057799999999997e-01,
5.500000000000000e-01, 8.373999999999999e-02, 6.415100000000000e-02, 5.617999999999994e-01, 4.365399999999996e-01,
6.000000000000000e-01, 1.223700000000000e-01, 9.467100000000001e-02, 7.726000000000013e-01, 6.104000000000009e-01,
6.500000000000000e-01, 1.741500000000000e-01, 1.365540000000000e-01, 1.035599999999999e+00, 8.376599999999993e-01,
7.000000000000000e-01, 2.417700000000000e-01, 1.929920000000000e-01, 1.352400000000002e+00, 1.128760000000001e+00,
7.500000000000000e-01, 3.275700000000000e-01, 2.675890000000000e-01, 1.715999999999998e+00, 1.491939999999999e+00,
8.000000000000000e-01, 4.328600000000000e-01, 3.640430000000000e-01, 2.105799999999999e+00, 1.929079999999998e+00,
8.500000000000000e-01, 5.571700000000001e-01, 4.855060000000000e-01, 2.486200000000004e+00, 2.429260000000003e+00,
9.000000000000000e-01, 6.974600000000000e-01, 6.335620000000000e-01, 2.805799999999996e+00, 2.961119999999997e+00,
9.500000000000000e-01, 8.478200000000000e-01, 8.068880000000000e-01, 3.007200000000005e+00, 3.466520000000006e+00,
9.999000000000000e-01, 9.990000000000000e-01, 9.996137760000001e-01, 3.029659318637271e+00, 3.862239999999997e+00,
1.000000000000000e+00, 1.000000000000000e+00, 1.000000000000000e+00, 1.000000000000111e+01, 3.862239999999221e+00,
});
}
std::vector<double> swfn_raw()
{
return makeTable(5, {
0, 0, 3.756330000000000e+00, 0, 0,
1.00e-04, 0, 3.756320000000000e+00, 0, -1.000000000006551e-01,
5.00e-02, 8.6000e-04, 1.869810000000000e+00, 1.723446893787575e-02, -3.780581162324650e+01,
1.00e-01, 2.6300e-03, 1.237310000000000e+00, 3.539999999999999e-02, -1.265000000000000e+01,
1.50e-01, 5.2400e-03, 9.182100000000000e-01, 5.220000000000001e-02, -6.382000000000001e+00,
2.00e-01, 8.7700e-03, 7.245100000000000e-01, 7.059999999999998e-02, -3.873999999999998e+00,
2.50e-01, 1.3380e-02, 5.934100000000000e-01, 9.220000000000000e-02, -2.622000000000000e+00,
3.00e-01, 1.9270e-02, 4.981100000000000e-01, 1.178000000000000e-01, -1.906000000000000e+00,
3.50e-01, 2.6720e-02, 4.251100000000000e-01, 1.490000000000001e-01, -1.460000000000000e+00,
4.00e-01, 3.6080e-02, 3.669100000000000e-01, 1.871999999999998e-01, -1.163999999999998e+00,
4.50e-01, 4.7810e-02, 3.191100000000000e-01, 2.346000000000000e-01, -9.560000000000004e-01,
5.00e-01, 6.2500e-02, 2.788100000000000e-01, 2.938000000000001e-01, -8.060000000000003e-01,
5.50e-01, 8.0900e-02, 2.440100000000000e-01, 3.679999999999997e-01, -6.959999999999993e-01,
6.00e-01, 1.0394e-01, 2.135100000000000e-01, 4.608000000000007e-01, -6.100000000000008e-01,
6.50e-01, 1.3277e-01, 1.863100000000000e-01, 5.765999999999993e-01, -5.439999999999996e-01,
7.00e-01, 1.6869e-01, 1.616100000000000e-01, 7.184000000000011e-01, -4.940000000000007e-01,
7.50e-01, 2.1302e-01, 1.390100000000000e-01, 8.865999999999988e-01, -4.519999999999998e-01,
8.00e-01, 2.6667e-01, 1.180100000000000e-01, 1.073000000000000e+00, -4.199999999999994e-01,
8.50e-01, 3.2918e-01, 9.830999999999999e-02, 1.250200000000001e+00, -3.940000000000007e-01,
9.00e-01, 3.9706e-01, 7.961000000000000e-02, 1.357600000000000e+00, -3.739999999999996e-01,
9.50e-01, 4.6103e-01, 6.161000000000000e-02, 1.279400000000001e+00, -3.600000000000005e-01,
1.00e+00, 5.0000e-01, 4.408000000000000e-02, 7.793999999999994e-01, -3.505999999999996e-01,
});
}
Opm::ECLPropTableRawData makeSatFunc(std::vector<double> data)
{
auto table = Opm::ECLPropTableRawData{};
table.data = std::move(data);
table.numTables = 1;
table.numPrimary = 1;
table.numCols = 5;
table.numRows = table.data.size() / table.numCols;
return table;
}
Opm::SatFuncInterpolant sgfn()
{
// Input tables follow METRIC unit conventions.
const auto usys = Opm::ECLUnits::createUnitSystem(1);
auto convert = Opm::SatFuncInterpolant::ConvertUnits{};
convert.indep = [](const double s) { return s; };
// Krg(Sg)
convert.column.emplace_back(
[](const double kr) { return kr; });
// Pcgo(Sg)
convert.column.push_back(
Opm::ECLPVT::CreateUnitConverter::ToSI::pressure(*usys));
// dKrg/dSg
convert.column.emplace_back(
[](const double dkr) { return dkr; });
// dPcgo/dSg
convert.column.push_back(
Opm::ECLPVT::CreateUnitConverter::ToSI::pressure(*usys));
return {
makeSatFunc(sgfn_raw()), convert
};
}
Opm::SatFuncInterpolant sofn()
{
auto convert = Opm::SatFuncInterpolant::ConvertUnits{};
convert.indep = [](const double s) { return s; };
// Krow(So)
convert.column.emplace_back(
[](const double kr) { return kr; });
// Krog(So)
convert.column.emplace_back(
[](const double kr) { return kr; });
// dKrow/dSo
convert.column.emplace_back(
[](const double dkr) { return dkr; });
// dKrow/dSo
convert.column.emplace_back(
[](const double dkr) { return dkr; });
return {
makeSatFunc(sofn_raw()), convert
};
}
Opm::SatFuncInterpolant swfn()
{
// Input tables follow METRIC unit conventions.
const auto usys = Opm::ECLUnits::createUnitSystem(1);
auto convert = Opm::SatFuncInterpolant::ConvertUnits{};
convert.indep = [](const double s) { return s; };
// Krw(Sw)
convert.column.emplace_back(
[](const double kr) { return kr; });
// Pcow(Sw)
convert.column.push_back(
Opm::ECLPVT::CreateUnitConverter::ToSI::pressure(*usys));
// dKrw/dSw
convert.column.emplace_back(
[](const double dkr) { return dkr; });
// dPcow/dSw
convert.column.push_back(
Opm::ECLPVT::CreateUnitConverter::ToSI::pressure(*usys));
return {
makeSatFunc(swfn_raw()), convert
};
}
} // Namespace Anonymous
// =====================================================================
@@ -625,6 +836,178 @@ BOOST_AUTO_TEST_CASE (ScaledBoth)
BOOST_AUTO_TEST_SUITE_END ()
// =====================================================================
BOOST_AUTO_TEST_SUITE (HorizontalScaling_SatFuncCurves)
BOOST_AUTO_TEST_CASE (KrGas_2Pt)
{
const auto interp = sgfn();
using Interp = std::remove_cv<
std::remove_reference<decltype(interp)>::type
>::type;
const auto SG = std::vector<double>{
0.0, 0.025, 0.05, 0.075, 0.1 , 0.125, 0.15 , 0.175,
0.2, 0.225, 0.25, 0.275, 0.3 , 0.325, 0.35 , 0.375,
0.4, 0.425, 0.45, 0.475, 0.49995, 0.5,
};
// Live Sg in [0, 0.5]. Map to [0, 1] for table lookup
const auto eps = ::Opm::SatFunc::
TwoPointScaling{ { 0.0 }, { 0.5 } };
const auto tep = ::Opm::SatFunc::EPSEvalInterface::TableEndPoints {
interp.connateSat()[0] , // low
interp.connateSat()[0] , // disp (ignored)
interp.maximumSat()[0] // max
};
// Forward: Scaled SG -> Lookup/Table sg
{
const auto expect = interp.saturationPoints(Interp::InTable{0});
const auto sg = eps.eval(tep, associate(SG));
check_is_close(sg, expect);
}
// Reverse: Lookup/Table sg -> Scaled SG
{
const auto sg = interp.saturationPoints(Interp::InTable{0});
const auto scaled = eps.reverse(tep, associate(sg));
check_is_close(scaled, SG);
}
}
BOOST_AUTO_TEST_CASE (KrOW_3pt)
{
const auto interp = sofn();
using Interp = std::remove_cv<
std::remove_reference<decltype(interp)>::type
>::type;
const auto swl = 0.1;
const auto swcr = 0.1; // Swcr == Swco in table.
const auto sgl = 0.0;
// Scaled end-points:
// SOWCR = 0.2, SWCR = 0.25, SWL = 0.15, SGL = 0
// => Mobile So in [ 0.2, 0.85 ], Sr = 0.75.
const auto eps = ::Opm::SatFunc::
ThreePointScaling{ { 0.2 }, { 0.75 }, { 0.85 } };
const auto tep = ::Opm::SatFunc::EPSEvalInterface::TableEndPoints {
interp.criticalSat(Interp::ResultColumn{0})[0], // low
1.0 - (swcr + sgl), // disp (Sr)
1.0 - (swl + sgl), // max
};
// Forward: Scaled SO -> Lookup/Table so.
{
const auto SO = std::vector<double> {
0.0, 0.05, 0.1, 0.15,
// ------------ Sowcr -----------------
0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
0.5, 0.55, 0.6, 0.65, 0.7, 0.725, 0.749,
// ------------ Sr --------------------
0.75, 0.7501, 0.8, 0.825, 0.85,
};
const auto so = eps.eval(tep, associate(SO));
const auto expect = std::vector<double>{
// so = Sowcr. (SO < SLO)
9.9999999999988987e-05, // SO = 0.0
9.9999999999988987e-05, // SO = 0.05
9.9999999999988987e-05, // SO = 0.1
9.9999999999988987e-05, // SO = 0.15
// ------------ SLO ------------------
// so = tep.low + t*(tep.disp - tep.low)
// = 1.0e-4 + t*(0.9 - 1.0e-4)
//
// t = (SO - SLO) / (SR - SLO)
// = (SO - 0.2) / (0.75 - 0.2)
9.9999999999988987e-05, // SO = 0.2
0.081909090909090876, // SO = 0.25
0.16371818181818176, // SO = 0.3
0.24552727272727265, // SO = 0.35
0.32733636363636365, // SO = 0.4
0.40914545454545453, // SO = 0.45
0.49095454545454542, // SO = 0.5
0.57276363636363636, // SO = 0.55
0.6545727272727272, // SO = 0.6
0.73638181818181814, // SO = 0.65
0.81819090909090897, // SO = 0.7
0.85909545454545433, // SO = 0.725
0.89836381818181799, // SO = 0.749
// ------------ Sr --------------------
// Everything below here maps to Somax (=1-Swco) because the
// 'tep' data says that Sdisp == Smax in the (purported) input
// table. The actual SOFN data (sofn_raw()) does not account
// for Swco = 0.1 (> 0) and therefore provides satfunc values
// for So > 0.9.
0.90000000000000002, // SO = 0.75
0.90000000000000002, // SO = 0.751
0.90000000000000002, // SO = 0.8
0.90000000000000002, // SO = 0.825
0.90000000000000002, // SO = 0.85
};
check_is_close(so, expect);
}
// Reverse: Lookup/Table so -> Scaled SO.
{
const auto SO = std::vector<double> {
0.20000000000000001, // so = 0
0.20000000000000001, // so = 9.9999999999988987e-05
0.23049783309256588, // so = 0.050000000000000037
0.2610567840871208, // so = 0.099999999999999978
0.29161573508167576, // so = 0.14999999999999999
0.32217468607623073, // so = 0.20000000000000001
0.35273363707078564, // so = 0.25
0.38329258806534061, // so = 0.29999999999999999
0.41385153905989558, // so = 0.34999999999999998
0.44441049005445055, // so = 0.40000000000000002
0.47496944104900546, // so = 0.45000000000000001
0.50552839204356048, // so = 0.5
0.53608734303811545, // so = 0.55000000000000004
0.56664629403267042, // so = 0.59999999999999998
0.59720524502722527, // so = 0.65000000000000002
0.62776419602178024, // so = 0.69999999999999996
0.6583231470163351, // so = 0.75
0.68888209801089018, // so = 0.80000000000000004
0.71944104900544503, // so = 0.84999999999999998
0.84999999999999998, // so = 0.90000000000000002
0.84999999999999998, // so = 0.94999999999999996
0.84999999999999998, // so = 0.99990000000000001
0.84999999999999998, // so = 1
};
const auto so = interp.saturationPoints(Interp::InTable{0});
const auto scaled = eps.reverse(tep, associate(so));
check_is_close(scaled, SO);
}
}
BOOST_AUTO_TEST_SUITE_END ()
// =====================================================================
// Three-point (alternative) scaling, applicable to relperm only.
// ---------------------------------------------------------------------
@@ -955,11 +1338,12 @@ BOOST_AUTO_TEST_CASE (Sw2_Regular)
SFPt{ 0.9, 0.81 }, // Maximum
};
// Scaled residual displacement sat: 0.7
// Scaled residual displacement sat: 0.7 (unchanged)
// Scaled Kr at Scaled Sr (KRxR): 0.6
// Scaled maximum saturation: 0.9 (unchanged)
// Scaled Kr at Smax: (KRx) 0.7
const auto scaler = Opm::SatFunc::
CritSatVerticalScaling({ 0.71 }, { 0.6 }, { 0.7 });
CritSatVerticalScaling({ 0.71 }, { 0.6 }, { 0.9 }, { 0.7 });
const auto y = scaler.vertScale(f, sw, kr);
@@ -1005,7 +1389,7 @@ BOOST_AUTO_TEST_CASE (Core1D_Coincident_FVal)
};
const auto scaler = Opm::SatFunc::
CritSatVerticalScaling({ 0.93 }, { 0.9 }, { 1.0 });
CritSatVerticalScaling({ 0.93 }, { 0.9 }, { 1.0 }, { 1.0 });
const auto y = scaler.vertScale(f, s, kr);
@@ -1065,3 +1449,149 @@ BOOST_AUTO_TEST_CASE (Core1D_Coincident_FVal)
}
BOOST_AUTO_TEST_SUITE_END ()
// =====================================================================
BOOST_AUTO_TEST_SUITE (VerticalScaling_SatFuncCurves)
BOOST_AUTO_TEST_CASE (Pcow_2Pt)
{
using SFPt = ::Opm::SatFunc::VerticalScalingInterface::FunctionValues::Point;
using SFVal = ::Opm::SatFunc::VerticalScalingInterface::FunctionValues;
const auto interp = swfn();
using Interp = std::remove_cv<
std::remove_reference<decltype(interp)>::type
>::type;
using ImportedOpm::unit::barsa;
// Maximum Pcow reset to 10 barsa.
const auto vscale = ::Opm::SatFunc::PureVerticalScaling {
{ 10.0*barsa }
};
const auto sw = interp.saturationPoints(Interp::InTable{0});
const auto pc = interp // Pc_ow = ResultColumn{1}
.interpolate(Interp::InTable{0}, Interp::ResultColumn{1}, sw);
const auto f = SFVal{
SFPt{ 1.0, 0.04408*barsa }, // Displacement
SFPt{ 0.0, 3.75633*barsa }, // Maximum
};
const auto PC = vscale.vertScale(f, associate(sw), pc);
const auto expect = std::vector<double> { // pc * (10.0 / 3.75633)
1000000,
999997.33782708121,
497775.75452635949,
329393.31741353922,
244443.37957527695,
192877.09013851287,
157976.00317331014,
132605.49525733895,
113171.63295024665,
97677.786562948415,
84952.60001118113,
74224.043148498662,
64959.681391145081,
56840.053988866792,
49598.943649785826,
43023.376540399804,
37006.865743957533,
31416.302614520024,
26171.821964523886,
21193.558606405721,
16401.64735260214,
11734.85822598121,
};
}
BOOST_AUTO_TEST_CASE (KrOG_3Pt)
{
using SFPt = ::Opm::SatFunc::VerticalScalingInterface::FunctionValues::Point;
using SFVal = ::Opm::SatFunc::VerticalScalingInterface::FunctionValues;
const auto interp = sofn();
using Interp = std::remove_cv<
std::remove_reference<decltype(interp)>::type
>::type;
const auto so = interp.saturationPoints(Interp::InTable{0});
const auto krog = interp // krog = ResultColumn{1}
.interpolate(Interp::InTable{0}, Interp::ResultColumn{1}, so);
const auto f = SFVal{
SFPt{ 0.6, 0.094671 }, // Displacement
SFPt{ 1.0, 1.0 }, // Maximum
};
const auto vscale = ::Opm::SatFunc::CritSatVerticalScaling {
// s , Kr
{ 0.6 }, { 0.25 },
{ 1.0 }, { 0.6 }
};
const auto KROG = vscale.vertScale(f, associate(so), krog);
const auto expect = std::vector<double> {
// Left interval (So <= Sr = f.disp.sat = 0.6)
//
// Pure vertical scaling (multiply by a constant factor).
//
// Kr = Kr(So) * (KRORG / Kr(Sr; table))
// = Kr(So) * (KRORG / f.disp.val)
// = Kr(So) * (0.25 / 0.094671)
//
0, // So = 0, Kr(So) = 0,
0, // So = 1.0e-4 Kr(So) = 0,
1.8485069345417287e-05, // So = 0.05 Kr(So) = 7.0000e-06
0.00023238372891381731, // So = 0.1 Kr(So) = 8.8000e-05
0.0010140380898057484, // So = 0.15 Kr(So) = 3.8400e-04
0.0029496889226901584, // So = 0.2 Kr(So) = 1.1170e-03
0.0068579607271498132, // So = 0.25 Kr(So) = 2.5970e-03
0.013874364905831774, // So = 0.3 Kr(So) = 5.2540e-03
0.025514677145060262, // So = 0.35 Kr(So) = 9.6620e-03
0.043799051451870158, // So = 0.4 Kr(So) = 1.6586e-02
0.071391978536193765, // So = 0.45 Kr(So) = 2.7035e-02
0.11176601071077732, // So = 0.5 Kr(So) = 4.2324e-02
0.16940509765398062, // So = 0.55 Kr(So) = 6.4151e-02
0.25, // So = 0.6 Kr(So) = 9.4671e-02
// ------------------------------------------
//
// Right interval (So >= Sr)
//
// Scaled Kr values defined by linear function between relperm
// points (Kr(Sr),KRORG) and (Krmax,KRO)
//
// Kr(So) - Kr(Sr)
// Kr = KRORG + ----------------- (KRO - KRORG)
// Krmax - Kr(Sr)
//
// Kr(So) - f.disp.val
// = KRORG + ------------------------ (KRO - KRORG)
// f.max.val - f.disp.val
//
// Kr(So) - 0.094671
// = 0.25 + ------------------- (0.6 - 0.25)
// 1 - 0.094671
//
0.26619195894531161, // So = 0.65 Kr(So) = 0.136554
0.28801087781347995, // So = 0.7 Kr(So) = 0.192992
0.31685006224256596, // So = 0.75 Kr(So) = 0.267589
0.35413915825075742, // So = 0.8 Kr(So) = 0.364043
0.40109672837167476, // So = 0.85 Kr(So) = 0.485506
0.45833514667043695, // So = 0.9 Kr(So) = 0.633562
0.52534294162674566, // So = 0.95 Kr(So) = 0.806888
0.59985068588325352, // So = 0.9999 Kr(So) = 0.999613776
0.59999999999999998, // So = 1.0 Kr(So) = 1.0
};
check_is_close(KROG, expect);
}
BOOST_AUTO_TEST_SUITE_END ()