OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #240 from andlaus/temperature_dependent_PVT

Browse Source 
PVT properties: allow them to be temperature dependent
This commit is contained in:
Atgeirr Flø Rasmussen 2014-12-02 08:49:00 +01:00
commit 710297ecd6
12 changed files with 324 additions and 93 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

66
opm/autodiff/BlackoilPropsAd.cpp
View File

@ -103,9 +103,11 @@ namespace Opm
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAd::muWat(const V& pw, V BlackoilPropsAd::muWat(const V& pw,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Water]) { if (!pu_.phase_used[Water]) {
@ -116,17 +118,19 @@ namespace Opm
const int np = props_.numPhases(); const int np = props_.numPhases();
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block mu(n, np); Block mu(n, np);
props_.viscosity(n, pw.data(), z.data(), cells.data(), mu.data(), 0); props_.viscosity(n, pw.data(), T.data(), z.data(), cells.data(), mu.data(), 0);
return mu.col(pu_.phase_pos[Water]); return mu.col(pu_.phase_pos[Water]);
} }
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAd::muOil(const V& po, V BlackoilPropsAd::muOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -145,15 +149,17 @@ namespace Opm
z.col(pu_.phase_pos[Gas]) = rs; z.col(pu_.phase_pos[Gas]) = rs;
} }
Block mu(n, np); Block mu(n, np);
props_.viscosity(n, po.data(), z.data(), cells.data(), mu.data(), 0); props_.viscosity(n, po.data(), T.data(), z.data(), cells.data(), mu.data(), 0);
return mu.col(pu_.phase_pos[Oil]); return mu.col(pu_.phase_pos[Oil]);
} }
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAd::muGas(const V& pg, V BlackoilPropsAd::muGas(const V& pg,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Gas]) { if (!pu_.phase_used[Gas]) {
@ -164,17 +170,19 @@ namespace Opm
const int np = props_.numPhases(); const int np = props_.numPhases();
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block mu(n, np); Block mu(n, np);
props_.viscosity(n, pg.data(), z.data(), cells.data(), mu.data(), 0); props_.viscosity(n, pg.data(), T.data(), z.data(), cells.data(), mu.data(), 0);
return mu.col(pu_.phase_pos[Gas]); return mu.col(pu_.phase_pos[Gas]);
} }
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAd::muGas(const V& pg, V BlackoilPropsAd::muGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -193,19 +201,21 @@ namespace Opm
z.col(pu_.phase_pos[Gas]) = V::Ones(n, 1); z.col(pu_.phase_pos[Gas]) = V::Ones(n, 1);
} }
Block mu(n, np); Block mu(n, np);
props_.viscosity(n, pg.data(), z.data(), cells.data(), mu.data(), 0); props_.viscosity(n, pg.data(), T.data(), z.data(), cells.data(), mu.data(), 0);
return mu.col(pu_.phase_pos[Gas]); return mu.col(pu_.phase_pos[Gas]);
} }
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAd::muWat(const ADB& pw, ADB BlackoilPropsAd::muWat(const ADB& pw,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
#if 1 #if 1
return ADB::constant(muWat(pw.value(), cells), pw.blockPattern()); return ADB::constant(muWat(pw.value(), T.value(), cells), pw.blockPattern());
#else #else
if (!pu_.phase_used[Water]) { if (!pu_.phase_used[Water]) {
OPM_THROW(std::runtime_error, "Cannot call muWat(): water phase not present."); OPM_THROW(std::runtime_error, "Cannot call muWat(): water phase not present.");
@ -216,7 +226,7 @@ namespace Opm
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block mu(n, np); Block mu(n, np);
Block dmu(n, np); Block dmu(n, np);
props_.viscosity(n, pw.value().data(), z.data(), cells.data(), mu.data(), dmu.data()); props_.viscosity(n, pw.value().data(), T.data(), z.data(), cells.data(), mu.data(), dmu.data());
ADB::M dmu_diag = spdiag(dmu.col(pu_.phase_pos[Water])); ADB::M dmu_diag = spdiag(dmu.col(pu_.phase_pos[Water]));
const int num_blocks = pw.numBlocks(); const int num_blocks = pw.numBlocks();
std::vector<ADB::M> jacs(num_blocks); std::vector<ADB::M> jacs(num_blocks);
@ -229,17 +239,19 @@ namespace Opm
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAd::muOil(const ADB& po, ADB BlackoilPropsAd::muOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
{ {
#if 1 #if 1
return ADB::constant(muOil(po.value(), rs.value(), cond, cells), po.blockPattern()); return ADB::constant(muOil(po.value(), T.value(), rs.value(), cond, cells), po.blockPattern());
#else #else
if (!pu_.phase_used[Oil]) { if (!pu_.phase_used[Oil]) {
OPM_THROW(std::runtime_error, "Cannot call muOil(): oil phase not present."); OPM_THROW(std::runtime_error, "Cannot call muOil(): oil phase not present.");
@ -272,13 +284,15 @@ namespace Opm
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAd::muGas(const ADB& pg, ADB BlackoilPropsAd::muGas(const ADB& pg,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
#if 1 #if 1
return ADB::constant(muGas(pg.value(), cells), pg.blockPattern()); return ADB::constant(muGas(pg.value(), T.value(), cells), pg.blockPattern());
#else #else
if (!pu_.phase_used[Gas]) { if (!pu_.phase_used[Gas]) {
OPM_THROW(std::runtime_error, "Cannot call muGas(): gas phase not present."); OPM_THROW(std::runtime_error, "Cannot call muGas(): gas phase not present.");
@ -301,17 +315,19 @@ namespace Opm
} }
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAd::muGas(const ADB& pg, ADB BlackoilPropsAd::muGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
{ {
#if 1 #if 1
return ADB::constant(muGas(pg.value(), rv.value(),cond,cells), pg.blockPattern()); return ADB::constant(muGas(pg.value(), T.value(), rv.value(),cond,cells), pg.blockPattern());
#else #else
if (!pu_.phase_used[Gas]) { if (!pu_.phase_used[Gas]) {
OPM_THROW(std::runtime_error, "Cannot call muGas(): gas phase not present."); OPM_THROW(std::runtime_error, "Cannot call muGas(): gas phase not present.");
@ -359,9 +375,11 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAd::bWat(const V& pw, V BlackoilPropsAd::bWat(const V& pw,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Water]) { if (!pu_.phase_used[Water]) {
@ -372,18 +390,20 @@ namespace Opm
const int np = props_.numPhases(); const int np = props_.numPhases();
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block matrix(n, np*np); Block matrix(n, np*np);
props_.matrix(n, pw.data(), z.data(), cells.data(), matrix.data(), 0); props_.matrix(n, pw.data(), T.data(), z.data(), cells.data(), matrix.data(), 0);
const int wi = pu_.phase_pos[Water]; const int wi = pu_.phase_pos[Water];
return matrix.col(wi*np + wi); return matrix.col(wi*np + wi);
} }
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAd::bOil(const V& po, V BlackoilPropsAd::bOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -402,16 +422,18 @@ namespace Opm
z.col(pu_.phase_pos[Gas]) = rs; z.col(pu_.phase_pos[Gas]) = rs;
} }
Block matrix(n, np*np); Block matrix(n, np*np);
props_.matrix(n, po.data(), z.data(), cells.data(), matrix.data(), 0); props_.matrix(n, po.data(), T.data(), z.data(), cells.data(), matrix.data(), 0);
const int oi = pu_.phase_pos[Oil]; const int oi = pu_.phase_pos[Oil];
return matrix.col(oi*np + oi); return matrix.col(oi*np + oi);
} }
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAd::bGas(const V& pg, V BlackoilPropsAd::bGas(const V& pg,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Gas]) { if (!pu_.phase_used[Gas]) {
@ -422,18 +444,20 @@ namespace Opm
const int np = props_.numPhases(); const int np = props_.numPhases();
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block matrix(n, np*np); Block matrix(n, np*np);
props_.matrix(n, pg.data(), z.data(), cells.data(), matrix.data(), 0); props_.matrix(n, pg.data(), pg.data(), z.data(), cells.data(), matrix.data(), 0);
const int gi = pu_.phase_pos[Gas]; const int gi = pu_.phase_pos[Gas];
return matrix.col(gi*np + gi); return matrix.col(gi*np + gi);
} }
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAd::bGas(const V& pg, V BlackoilPropsAd::bGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -452,16 +476,18 @@ namespace Opm
z.col(pu_.phase_pos[Gas]) = V::Ones(n, 1); z.col(pu_.phase_pos[Gas]) = V::Ones(n, 1);
} }
Block matrix(n, np*np); Block matrix(n, np*np);
props_.matrix(n, pg.data(), z.data(), cells.data(), matrix.data(), 0); props_.matrix(n, pg.data(), T.data(), z.data(), cells.data(), matrix.data(), 0);
const int gi = pu_.phase_pos[Gas]; const int gi = pu_.phase_pos[Gas];
return matrix.col(gi*np + gi); return matrix.col(gi*np + gi);
} }
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAd::bWat(const ADB& pw, ADB BlackoilPropsAd::bWat(const ADB& pw,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Water]) { if (!pu_.phase_used[Water]) {
@ -473,7 +499,7 @@ namespace Opm
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block matrix(n, np*np); Block matrix(n, np*np);
Block dmatrix(n, np*np); Block dmatrix(n, np*np);
props_.matrix(n, pw.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data()); props_.matrix(n, pw.value().data(), T.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data());
const int phase_ind = pu_.phase_pos[Water]; const int phase_ind = pu_.phase_pos[Water];
const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column. const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column.
ADB::M db_diag = spdiag(dmatrix.col(column)); ADB::M db_diag = spdiag(dmatrix.col(column));
@ -487,11 +513,13 @@ namespace Opm
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAd::bOil(const ADB& po, ADB BlackoilPropsAd::bOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -511,7 +539,7 @@ namespace Opm
} }
Block matrix(n, np*np); Block matrix(n, np*np);
Block dmatrix(n, np*np); Block dmatrix(n, np*np);
props_.matrix(n, po.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data()); props_.matrix(n, po.value().data(), T.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data());
const int phase_ind = pu_.phase_pos[Oil]; const int phase_ind = pu_.phase_pos[Oil];
const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column. const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column.
ADB::M db_diag = spdiag(dmatrix.col(column)); ADB::M db_diag = spdiag(dmatrix.col(column));
@ -528,9 +556,11 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAd::bGas(const ADB& pg, ADB BlackoilPropsAd::bGas(const ADB& pg,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!pu_.phase_used[Gas]) { if (!pu_.phase_used[Gas]) {
@ -542,7 +572,7 @@ namespace Opm
Block z = Block::Zero(n, np); Block z = Block::Zero(n, np);
Block matrix(n, np*np); Block matrix(n, np*np);
Block dmatrix(n, np*np); Block dmatrix(n, np*np);
props_.matrix(n, pg.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data()); props_.matrix(n, pg.value().data(), T.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data());
const int phase_ind = pu_.phase_pos[Gas]; const int phase_ind = pu_.phase_pos[Gas];
const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column. const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column.
ADB::M db_diag = spdiag(dmatrix.col(column)); ADB::M db_diag = spdiag(dmatrix.col(column));
@ -556,11 +586,13 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAd::bGas(const ADB& pg, ADB BlackoilPropsAd::bGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& /*cond*/, const std::vector<PhasePresence>& /*cond*/,
const Cells& cells) const const Cells& cells) const
@ -580,7 +612,7 @@ namespace Opm
} }
Block matrix(n, np*np); Block matrix(n, np*np);
Block dmatrix(n, np*np); Block dmatrix(n, np*np);
props_.matrix(n, pg.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data()); props_.matrix(n, pg.value().data(), T.value().data(), z.data(), cells.data(), matrix.data(), dmatrix.data());
const int phase_ind = pu_.phase_pos[Gas]; const int phase_ind = pu_.phase_pos[Gas];
const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column. const int column = phase_ind*np + phase_ind; // Index of our sought diagonal column.
ADB::M db_diag = spdiag(dmatrix.col(column)); ADB::M db_diag = spdiag(dmatrix.col(column));

32
opm/autodiff/BlackoilPropsAd.hpp
View File

@ -101,72 +101,88 @@ namespace Opm
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muWat(const V& pw, V muWat(const V& pw,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muOil(const V& po, V muOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muGas(const V& pg, V muGas(const V& pg,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n gas solution factor values. /// \param[in] rv Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muGas(const V& pg, V muGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muWat(const ADB& pw, ADB muWat(const ADB& pw,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muOil(const ADB& po, ADB muOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n gas solution factor values. /// \param[in] rv Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
@ -175,73 +191,89 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bWat(const V& pw, V bWat(const V& pw,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bOil(const V& po, V bOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bWat(const ADB& pw, ADB bWat(const ADB& pw,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bOil(const ADB& po, ADB bOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;

64
opm/autodiff/BlackoilPropsAdFromDeck.cpp
View File

@ -269,9 +269,11 @@ namespace Opm
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAdFromDeck::muWat(const V& pw, V BlackoilPropsAdFromDeck::muWat(const V& pw,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Water]) { if (!phase_usage_.phase_used[Water]) {
@ -284,18 +286,20 @@ namespace Opm
V dmudr(n); V dmudr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Water]]->mu(n, &pvtTableIdx_[0], pw.data(), rs, props_[phase_usage_.phase_pos[Water]]->mu(n, &pvtTableIdx_[0], pw.data(), T.data(), rs,
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
return mu; return mu;
} }
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAdFromDeck::muOil(const V& po, V BlackoilPropsAdFromDeck::muOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -309,16 +313,18 @@ namespace Opm
V dmudp(n); V dmudp(n);
V dmudr(n); V dmudr(n);
props_[phase_usage_.phase_pos[Oil]]->mu(n, &pvtTableIdx_[0], po.data(), rs.data(), &cond[0], props_[phase_usage_.phase_pos[Oil]]->mu(n, &pvtTableIdx_[0], po.data(), T.data(), rs.data(), &cond[0],
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
return mu; return mu;
} }
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAdFromDeck::muGas(const V& pg, V BlackoilPropsAdFromDeck::muGas(const V& pg,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Gas]) { if (!phase_usage_.phase_used[Gas]) {
@ -331,16 +337,18 @@ namespace Opm
V dmudr(n); V dmudr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.data(), rs, props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.data(), T.data(), rs,
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
return mu; return mu;
} }
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V BlackoilPropsAdFromDeck::muGas(const V& pg, V BlackoilPropsAdFromDeck::muGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -354,16 +362,18 @@ namespace Opm
V dmudp(n); V dmudp(n);
V dmudr(n); V dmudr(n);
props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.data(), rv.data(),&cond[0], props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.data(), T.data(), rv.data(),&cond[0],
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
return mu; return mu;
} }
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAdFromDeck::muWat(const ADB& pw, ADB BlackoilPropsAdFromDeck::muWat(const ADB& pw,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Water]) { if (!phase_usage_.phase_used[Water]) {
@ -376,7 +386,7 @@ namespace Opm
V dmudr(n); V dmudr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Water]]->mu(n, &pvtTableIdx_[0], pw.value().data(), rs, props_[phase_usage_.phase_pos[Water]]->mu(n, &pvtTableIdx_[0], pw.value().data(), T.value().data(), rs,
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
ADB::M dmudp_diag = spdiag(dmudp); ADB::M dmudp_diag = spdiag(dmudp);
const int num_blocks = pw.numBlocks(); const int num_blocks = pw.numBlocks();
@ -389,11 +399,13 @@ namespace Opm
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAdFromDeck::muOil(const ADB& po, ADB BlackoilPropsAdFromDeck::muOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -407,7 +419,7 @@ namespace Opm
V dmudp(n); V dmudp(n);
V dmudr(n); V dmudr(n);
props_[phase_usage_.phase_pos[Oil]]->mu(n, &pvtTableIdx_[0], po.value().data(), rs.value().data(), props_[phase_usage_.phase_pos[Oil]]->mu(n, &pvtTableIdx_[0], po.value().data(), T.value().data(), rs.value().data(),
&cond[0], mu.data(), dmudp.data(), dmudr.data()); &cond[0], mu.data(), dmudp.data(), dmudr.data());
ADB::M dmudp_diag = spdiag(dmudp); ADB::M dmudp_diag = spdiag(dmudp);
@ -422,9 +434,11 @@ namespace Opm
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAdFromDeck::muGas(const ADB& pg, ADB BlackoilPropsAdFromDeck::muGas(const ADB& pg,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Gas]) { if (!phase_usage_.phase_used[Gas]) {
@ -437,7 +451,7 @@ namespace Opm
V dmudr(n); V dmudr(n);
const double* rv = 0; const double* rv = 0;
props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.value().data(), rv, props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.value().data(), T.value().data(), rv,
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
ADB::M dmudp_diag = spdiag(dmudp); ADB::M dmudp_diag = spdiag(dmudp);
@ -451,11 +465,13 @@ namespace Opm
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB BlackoilPropsAdFromDeck::muGas(const ADB& pg, ADB BlackoilPropsAdFromDeck::muGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -469,7 +485,7 @@ namespace Opm
V dmudp(n); V dmudp(n);
V dmudr(n); V dmudr(n);
props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.value().data(), rv.value().data(),&cond[0], props_[phase_usage_.phase_pos[Gas]]->mu(n, &pvtTableIdx_[0], pg.value().data(), T.value().data(), rv.value().data(),&cond[0],
mu.data(), dmudp.data(), dmudr.data()); mu.data(), dmudp.data(), dmudr.data());
ADB::M dmudp_diag = spdiag(dmudp); ADB::M dmudp_diag = spdiag(dmudp);
@ -502,9 +518,11 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAdFromDeck::bWat(const V& pw, V BlackoilPropsAdFromDeck::bWat(const V& pw,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Water]) { if (!phase_usage_.phase_used[Water]) {
@ -518,7 +536,7 @@ namespace Opm
V dbdr(n); V dbdr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Water]]->b(n, &pvtTableIdx_[0], pw.data(), rs, props_[phase_usage_.phase_pos[Water]]->b(n, &pvtTableIdx_[0], pw.data(), T.data(), rs,
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
return b; return b;
@ -526,11 +544,13 @@ namespace Opm
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAdFromDeck::bOil(const V& po, V BlackoilPropsAdFromDeck::bOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -545,7 +565,7 @@ namespace Opm
V dbdp(n); V dbdp(n);
V dbdr(n); V dbdr(n);
props_[phase_usage_.phase_pos[Oil]]->b(n, &pvtTableIdx_[0], po.data(), rs.data(), &cond[0], props_[phase_usage_.phase_pos[Oil]]->b(n, &pvtTableIdx_[0], po.data(), T.data(), rs.data(), &cond[0],
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
return b; return b;
@ -553,9 +573,11 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAdFromDeck::bGas(const V& pg, V BlackoilPropsAdFromDeck::bGas(const V& pg,
const V& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Gas]) { if (!phase_usage_.phase_used[Gas]) {
@ -569,7 +591,7 @@ namespace Opm
V dbdr(n); V dbdr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.data(), rs, props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.data(), T.data(), rs,
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
return b; return b;
@ -577,11 +599,13 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V BlackoilPropsAdFromDeck::bGas(const V& pg, V BlackoilPropsAdFromDeck::bGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -596,7 +620,7 @@ namespace Opm
V dbdp(n); V dbdp(n);
V dbdr(n); V dbdr(n);
props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.data(), rv.data(), &cond[0], props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.data(), T.data(), rv.data(), &cond[0],
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
return b; return b;
@ -604,9 +628,11 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAdFromDeck::bWat(const ADB& pw, ADB BlackoilPropsAdFromDeck::bWat(const ADB& pw,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Water]) { if (!phase_usage_.phase_used[Water]) {
@ -620,7 +646,7 @@ namespace Opm
V dbdr(n); V dbdr(n);
const double* rs = 0; const double* rs = 0;
props_[phase_usage_.phase_pos[Water]]->b(n, &pvtTableIdx_[0], pw.value().data(), rs, props_[phase_usage_.phase_pos[Water]]->b(n, &pvtTableIdx_[0], pw.value().data(), T.value().data(), rs,
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
ADB::M dbdp_diag = spdiag(dbdp); ADB::M dbdp_diag = spdiag(dbdp);
@ -634,11 +660,13 @@ namespace Opm
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n taxonomies classifying fluid condition. /// \param[in] cond Array of n taxonomies classifying fluid condition.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAdFromDeck::bOil(const ADB& po, ADB BlackoilPropsAdFromDeck::bOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -653,7 +681,7 @@ namespace Opm
V dbdp(n); V dbdp(n);
V dbdr(n); V dbdr(n);
props_[phase_usage_.phase_pos[Oil]]->b(n, &pvtTableIdx_[0], po.value().data(), rs.value().data(), props_[phase_usage_.phase_pos[Oil]]->b(n, &pvtTableIdx_[0], po.value().data(), T.value().data(), rs.value().data(),
&cond[0], b.data(), dbdp.data(), dbdr.data()); &cond[0], b.data(), dbdp.data(), dbdr.data());
ADB::M dbdp_diag = spdiag(dbdp); ADB::M dbdp_diag = spdiag(dbdp);
@ -668,9 +696,11 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAdFromDeck::bGas(const ADB& pg, ADB BlackoilPropsAdFromDeck::bGas(const ADB& pg,
const ADB& T,
const Cells& cells) const const Cells& cells) const
{ {
if (!phase_usage_.phase_used[Gas]) { if (!phase_usage_.phase_used[Gas]) {
@ -684,7 +714,7 @@ namespace Opm
V dbdr(n); V dbdr(n);
const double* rv = 0; const double* rv = 0;
props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.value().data(), rv, props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.value().data(), T.value().data(), rv,
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
ADB::M dbdp_diag = spdiag(dbdp); ADB::M dbdp_diag = spdiag(dbdp);
@ -698,11 +728,13 @@ namespace Opm
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB BlackoilPropsAdFromDeck::bGas(const ADB& pg, ADB BlackoilPropsAdFromDeck::bGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const const Cells& cells) const
@ -717,7 +749,7 @@ namespace Opm
V dbdp(n); V dbdp(n);
V dbdr(n); V dbdr(n);
props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.value().data(), rv.value().data(), &cond[0], props_[phase_usage_.phase_pos[Gas]]->b(n, &pvtTableIdx_[0], pg.value().data(), T.value().data(), rv.value().data(), &cond[0],
b.data(), dbdp.data(), dbdr.data()); b.data(), dbdp.data(), dbdr.data());
ADB::M dbdp_diag = spdiag(dbdp); ADB::M dbdp_diag = spdiag(dbdp);

32
opm/autodiff/BlackoilPropsAdFromDeck.hpp
View File

@ -126,70 +126,86 @@ namespace Opm
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muWat(const V& pw, V muWat(const V& pw,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muOil(const V& po, V muOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muGas(const V& pg, V muGas(const V& pg,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
V muGas(const V& po, V muGas(const V& po,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muWat(const ADB& pw, ADB muWat(const ADB& pw,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muOil(const ADB& po, ADB muOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
@ -198,72 +214,88 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bWat(const V& pw, V bWat(const V& pw,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bOil(const V& po, V bOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bWat(const ADB& pw, ADB bWat(const ADB& pw,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bOil(const ADB& po, ADB bOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const Cells& cells) const; const Cells& cells) const;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const; const Cells& cells) const;

30
opm/autodiff/BlackoilPropsAdInterface.hpp
View File

@ -91,66 +91,80 @@ namespace Opm
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
V muWat(const V& pw, V muWat(const V& pw,
const V& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
V muOil(const V& po, V muOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
V muGas(const V& pg, V muGas(const V& pg,
const V& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Water viscosity. /// Water viscosity.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
ADB muWat(const ADB& pw, ADB muWat(const ADB& pw,
const ADB& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Oil viscosity. /// Oil viscosity.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
ADB muOil(const ADB& po, ADB muOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas viscosity. /// Gas viscosity.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n viscosity values. /// \return Array of n viscosity values.
virtual virtual
ADB muGas(const ADB& pg, ADB muGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
@ -159,80 +173,96 @@ namespace Opm
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
V bWat(const V& pw, V bWat(const V& pw,
const V& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
V bOil(const V& po, V bOil(const V& po,
const V& T,
const V& rs, const V& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
V bGas(const V& pg, V bGas(const V& pg,
const V& T,
const V& rv, const V& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Water formation volume factor. /// Water formation volume factor.
/// \param[in] pw Array of n water pressure values. /// \param[in] pw Array of n water pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
ADB bWat(const ADB& pw, ADB bWat(const ADB& pw,
const ADB& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Oil formation volume factor. /// Oil formation volume factor.
/// \param[in] po Array of n oil pressure values. /// \param[in] po Array of n oil pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rs Array of n gas solution factor values. /// \param[in] rs Array of n gas solution factor values.
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
ADB bOil(const ADB& po, ADB bOil(const ADB& po,
const ADB& T,
const ADB& rs, const ADB& rs,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const Cells& cells) const = 0; const Cells& cells) const = 0;
/// Gas formation volume factor. /// Gas formation volume factor.
/// \param[in] pg Array of n gas pressure values. /// \param[in] pg Array of n gas pressure values.
/// \param[in] T Array of n temperature values.
/// \param[in] rv Array of n vapor oil/gas ratio /// \param[in] rv Array of n vapor oil/gas ratio
/// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell. /// \param[in] cond Array of n objects, each specifying which phases are present with non-zero saturation in a cell.
/// \param[in] cells Array of n cell indices to be associated with the pressure values. /// \param[in] cells Array of n cell indices to be associated with the pressure values.
/// \return Array of n formation volume factor values. /// \return Array of n formation volume factor values.
virtual virtual
ADB bGas(const ADB& pg, ADB bGas(const ADB& pg,
const ADB& T,
const ADB& rv, const ADB& rv,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const Cells& cells) const = 0; const Cells& cells) const = 0;

4
opm/autodiff/FullyImplicitBlackoilSolver.hpp
View File

@ -142,6 +142,7 @@ namespace Opm {
struct SolutionState { struct SolutionState {
SolutionState(const int np); SolutionState(const int np);
ADB pressure; ADB pressure;
ADB temperature;
std::vector<ADB> saturation; std::vector<ADB> saturation;
ADB rs; ADB rs;
ADB rv; ADB rv;
@ -259,6 +260,7 @@ namespace Opm {
ADB ADB
fluidViscosity(const int phase, fluidViscosity(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
@ -267,6 +269,7 @@ namespace Opm {
ADB ADB
fluidReciprocFVF(const int phase, fluidReciprocFVF(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
@ -275,6 +278,7 @@ namespace Opm {
ADB ADB
fluidDensity(const int phase, fluidDensity(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,

41
opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp
View File

@ -357,6 +357,7 @@ namespace {
template<class T> template<class T>
FullyImplicitBlackoilSolver<T>::SolutionState::SolutionState(const int np) FullyImplicitBlackoilSolver<T>::SolutionState::SolutionState(const int np)
: pressure ( ADB::null()) : pressure ( ADB::null())
, temperature( ADB::null())
, saturation(np, ADB::null()) , saturation(np, ADB::null())
, rs ( ADB::null()) , rs ( ADB::null())
, rv ( ADB::null()) , rv ( ADB::null())
@ -413,6 +414,7 @@ namespace {
// automatically consistent with variableState() (and doing // automatically consistent with variableState() (and doing
// things automatically is all the rage in this module ;) // things automatically is all the rage in this module ;)
state.pressure = ADB::constant(state.pressure.value()); state.pressure = ADB::constant(state.pressure.value());
state.temperature = ADB::constant(state.temperature.value());
state.rs = ADB::constant(state.rs.value()); state.rs = ADB::constant(state.rs.value());
state.rv = ADB::constant(state.rv.value()); state.rv = ADB::constant(state.rv.value());
for (int phaseIdx= 0; phaseIdx < x.numPhases(); ++ phaseIdx) for (int phaseIdx= 0; phaseIdx < x.numPhases(); ++ phaseIdx)
@ -512,6 +514,10 @@ namespace {
int nextvar = 0; int nextvar = 0;
state.pressure = vars[ nextvar++ ]; state.pressure = vars[ nextvar++ ];
// temperature
const V temp = Eigen::Map<const V>(& x.temperature()[0], x.temperature().size());
state.temperature = ADB::constant(temp);
// Saturations // Saturations
const std::vector<int>& bpat = vars[0].blockPattern(); const std::vector<int>& bpat = vars[0].blockPattern();
{ {
@ -573,6 +579,7 @@ namespace {
const Opm::PhaseUsage& pu = fluid_.phaseUsage(); const Opm::PhaseUsage& pu = fluid_.phaseUsage();
const ADB& press = state.pressure; const ADB& press = state.pressure;
const ADB& temp = state.temperature;
const std::vector<ADB>& sat = state.saturation; const std::vector<ADB>& sat = state.saturation;
const ADB& rs = state.rs; const ADB& rs = state.rs;
const ADB& rv = state.rv; const ADB& rv = state.rv;
@ -586,7 +593,7 @@ namespace {
for (int phase = 0; phase < maxnp; ++phase) { for (int phase = 0; phase < maxnp; ++phase) {
if (active_[ phase ]) { if (active_[ phase ]) {
const int pos = pu.phase_pos[ phase ]; const int pos = pu.phase_pos[ phase ];
rq_[pos].b = fluidReciprocFVF(phase, pressures[phase], rs, rv, cond, cells_); rq_[pos].b = fluidReciprocFVF(phase, pressures[phase], temp, rs, rv, cond, cells_);
rq_[pos].accum[aix] = pv_mult * rq_[pos].b * sat[pos]; rq_[pos].accum[aix] = pv_mult * rq_[pos].b * sat[pos];
// DUMP(rq_[pos].b); // DUMP(rq_[pos].b);
// DUMP(rq_[pos].accum[aix]); // DUMP(rq_[pos].accum[aix]);
@ -620,6 +627,7 @@ namespace {
const std::vector<int> well_cells(wells_.well_cells, wells_.well_cells + nperf); const std::vector<int> well_cells(wells_.well_cells, wells_.well_cells + nperf);
// Compute b, rsmax, rvmax values for perforations. // Compute b, rsmax, rvmax values for perforations.
const ADB perf_press = subset(state.pressure, well_cells); const ADB perf_press = subset(state.pressure, well_cells);
const ADB perf_temp = subset(state.temperature, well_cells);
std::vector<PhasePresence> perf_cond(nperf); std::vector<PhasePresence> perf_cond(nperf);
const std::vector<PhasePresence>& pc = phaseCondition(); const std::vector<PhasePresence>& pc = phaseCondition();
for (int perf = 0; perf < nperf; ++perf) { for (int perf = 0; perf < nperf; ++perf) {
@ -630,21 +638,21 @@ namespace {
std::vector<double> rssat_perf(nperf, 0.0); std::vector<double> rssat_perf(nperf, 0.0);
std::vector<double> rvsat_perf(nperf, 0.0); std::vector<double> rvsat_perf(nperf, 0.0);
if (pu.phase_used[BlackoilPhases::Aqua]) { if (pu.phase_used[BlackoilPhases::Aqua]) {
const ADB bw = fluid_.bWat(perf_press, well_cells); const ADB bw = fluid_.bWat(perf_press, perf_temp, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Aqua]) = bw.value(); b.col(pu.phase_pos[BlackoilPhases::Aqua]) = bw.value();
} }
assert(active_[Oil]); assert(active_[Oil]);
const ADB perf_so = subset(state.saturation[pu.phase_pos[Oil]], well_cells); const ADB perf_so = subset(state.saturation[pu.phase_pos[Oil]], well_cells);
if (pu.phase_used[BlackoilPhases::Liquid]) { if (pu.phase_used[BlackoilPhases::Liquid]) {
const ADB perf_rs = subset(state.rs, well_cells); const ADB perf_rs = subset(state.rs, well_cells);
const ADB bo = fluid_.bOil(perf_press, perf_rs, perf_cond, well_cells); const ADB bo = fluid_.bOil(perf_press, perf_temp, perf_rs, perf_cond, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Liquid]) = bo.value(); b.col(pu.phase_pos[BlackoilPhases::Liquid]) = bo.value();
const V rssat = fluidRsSat(perf_press.value(), perf_so.value(), well_cells); const V rssat = fluidRsSat(perf_press.value(), perf_so.value(), well_cells);
rssat_perf.assign(rssat.data(), rssat.data() + nperf); rssat_perf.assign(rssat.data(), rssat.data() + nperf);
} }
if (pu.phase_used[BlackoilPhases::Vapour]) { if (pu.phase_used[BlackoilPhases::Vapour]) {
const ADB perf_rv = subset(state.rv, well_cells); const ADB perf_rv = subset(state.rv, well_cells);
const ADB bg = fluid_.bGas(perf_press, perf_rv, perf_cond, well_cells); const ADB bg = fluid_.bGas(perf_press, perf_temp, perf_rv, perf_cond, well_cells);
b.col(pu.phase_pos[BlackoilPhases::Vapour]) = bg.value(); b.col(pu.phase_pos[BlackoilPhases::Vapour]) = bg.value();
const V rvsat = fluidRvSat(perf_press.value(), perf_so.value(), well_cells); const V rvsat = fluidRvSat(perf_press.value(), perf_so.value(), well_cells);
rvsat_perf.assign(rvsat.data(), rvsat.data() + nperf); rvsat_perf.assign(rvsat.data(), rvsat.data() + nperf);
@ -1611,11 +1619,11 @@ namespace {
const std::vector<PhasePresence> cond = phaseCondition(); const std::vector<PhasePresence> cond = phaseCondition();
const ADB tr_mult = transMult(state.pressure); const ADB tr_mult = transMult(state.pressure);
const ADB mu = fluidViscosity(canonicalPhaseIdx, phasePressure, state.rs, state.rv,cond, cells_); const ADB mu = fluidViscosity(canonicalPhaseIdx, phasePressure, state.temperature, state.rs, state.rv,cond, cells_);
rq_[ actph ].mob = tr_mult * kr / mu; rq_[ actph ].mob = tr_mult * kr / mu;
const ADB rho = fluidDensity(canonicalPhaseIdx, phasePressure, state.rs, state.rv,cond, cells_); const ADB rho = fluidDensity(canonicalPhaseIdx, phasePressure, state.temperature, state.rs, state.rv,cond, cells_);
ADB& head = rq_[ actph ].head; ADB& head = rq_[ actph ].head;
@ -1910,7 +1918,8 @@ namespace {
template<class T> template<class T>
ADB ADB
FullyImplicitBlackoilSolver<T>::fluidViscosity(const int phase, FullyImplicitBlackoilSolver<T>::fluidViscosity(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
@ -1918,12 +1927,12 @@ namespace {
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.muWat(p, cells); return fluid_.muWat(p, temp, cells);
case Oil: { case Oil: {
return fluid_.muOil(p, rs, cond, cells); return fluid_.muOil(p, temp, rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.muGas(p, rv, cond, cells); return fluid_.muGas(p, temp, rv, cond, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -1937,6 +1946,7 @@ namespace {
ADB ADB
FullyImplicitBlackoilSolver<T>::fluidReciprocFVF(const int phase, FullyImplicitBlackoilSolver<T>::fluidReciprocFVF(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
@ -1944,12 +1954,12 @@ namespace {
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.bWat(p, cells); return fluid_.bWat(p, temp, cells);
case Oil: { case Oil: {
return fluid_.bOil(p, rs, cond, cells); return fluid_.bOil(p, temp, rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.bGas(p, rv, cond, cells); return fluid_.bGas(p, temp, rv, cond, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -1962,14 +1972,15 @@ namespace {
template<class T> template<class T>
ADB ADB
FullyImplicitBlackoilSolver<T>::fluidDensity(const int phase, FullyImplicitBlackoilSolver<T>::fluidDensity(const int phase,
const ADB& p , const ADB& p ,
const ADB& temp ,
const ADB& rs , const ADB& rs ,
const ADB& rv , const ADB& rv ,
const std::vector<PhasePresence>& cond, const std::vector<PhasePresence>& cond,
const std::vector<int>& cells) const const std::vector<int>& cells) const
{ {
const double* rhos = fluid_.surfaceDensity(); const double* rhos = fluid_.surfaceDensity();
ADB b = fluidReciprocFVF(phase, p, rs, rv, cond, cells); ADB b = fluidReciprocFVF(phase, p, temp, rs, rv, cond, cells);
ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b; ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b;
if (phase == Oil && active_[Gas]) { if (phase == Oil && active_[Gas]) {
// It is correct to index into rhos with canonical phase indices. // It is correct to index into rhos with canonical phase indices.

64
opm/autodiff/ImpesTPFAAD.cpp
View File

@ -280,8 +280,9 @@ namespace {
} }
V cell_rho_total = V::Zero(nc,1); V cell_rho_total = V::Zero(nc,1);
const Eigen::Map<const V> p(state.pressure().data(), nc, 1); const Eigen::Map<const V> p(state.pressure().data(), nc, 1);
const Eigen::Map<const V> T(state.temperature().data(), nc, 1);
for (int phase = 0; phase < np; ++phase) { for (int phase = 0; phase < np; ++phase) {
const V cell_rho = fluidRho(phase, p, cells); const V cell_rho = fluidRho(phase, p, T, cells);
const V cell_s = s.col(phase); const V cell_s = s.col(phase);
cell_rho_total += cell_s * cell_rho; cell_rho_total += cell_s * cell_rho;
} }
@ -318,10 +319,12 @@ namespace {
// Initialize AD variables: p (cell pressures) and bhp (well bhp). // Initialize AD variables: p (cell pressures) and bhp (well bhp).
const V p0 = Eigen::Map<const V>(&state.pressure()[0], nc, 1); const V p0 = Eigen::Map<const V>(&state.pressure()[0], nc, 1);
const V T0 = Eigen::Map<const V>(&state.temperature()[0], nc, 1);
const V bhp0 = Eigen::Map<const V>(&well_state.bhp()[0], nw, 1); const V bhp0 = Eigen::Map<const V>(&well_state.bhp()[0], nw, 1);
std::vector<V> vars0 = { p0, bhp0 }; std::vector<V> vars0 = { p0, bhp0 };
std::vector<ADB> vars = ADB::variables(vars0); std::vector<ADB> vars = ADB::variables(vars0);
const ADB& p = vars[0]; const ADB& p = vars[0];
const ADB T = ADB::constant(T0);
const ADB& bhp = vars[1]; const ADB& bhp = vars[1];
std::vector<int> bpat = p.blockPattern(); std::vector<int> bpat = p.blockPattern();
@ -331,6 +334,7 @@ namespace {
// Extract variables for perforation cell pressures // Extract variables for perforation cell pressures
// and corresponding perforation well pressures. // and corresponding perforation well pressures.
const ADB p_perfcell = subset(p, well_cells); const ADB p_perfcell = subset(p, well_cells);
const ADB T_perfcell = subset(T, well_cells);
// Construct matrix to map wells->perforations. // Construct matrix to map wells->perforations.
M well_to_perf(well_cells.size(), nw); M well_to_perf(well_cells.size(), nw);
typedef Eigen::Triplet<double> Tri; typedef Eigen::Triplet<double> Tri;
@ -352,20 +356,20 @@ namespace {
ADB divcontrib_sum = ADB::constant(V::Zero(nc,1), bpat); ADB divcontrib_sum = ADB::constant(V::Zero(nc,1), bpat);
qs_ = ADB::constant(V::Zero(nw*np, 1), bpat); qs_ = ADB::constant(V::Zero(nw*np, 1), bpat);
for (int phase = 0; phase < np; ++phase) { for (int phase = 0; phase < np; ++phase) {
const ADB cell_b = fluidFvf(phase, p, cells); const ADB cell_b = fluidFvf(phase, p, T, cells);
const ADB cell_rho = fluidRho(phase, p, cells); const ADB cell_rho = fluidRho(phase, p, T, cells);
const ADB well_b = fluidFvf(phase, p_perfwell, well_cells); const ADB well_b = fluidFvf(phase, p_perfwell, T_perfcell, well_cells);
const V kr = fluidKr(phase); const V kr = fluidKr(phase);
// Explicitly not asking for derivatives of viscosity, // Explicitly not asking for derivatives of viscosity,
// since they are not available yet. // since they are not available yet.
const V mu = fluidMu(phase, p.value(), cells); const V mu = fluidMu(phase, p.value(), T.value(), cells);
const V cell_mob = kr / mu; const V cell_mob = kr / mu;
const ADB head_diff_grav = (grav_ * cell_rho); const ADB head_diff_grav = (grav_ * cell_rho);
const ADB head = nkgradp + (grav_ * cell_rho); const ADB head = nkgradp + (grav_ * cell_rho);
const UpwindSelector<double> upwind(grid_, ops_, head.value()); const UpwindSelector<double> upwind(grid_, ops_, head.value());
const V face_mob = upwind.select(cell_mob); const V face_mob = upwind.select(cell_mob);
const V well_kr = fluidKrWell(phase); const V well_kr = fluidKrWell(phase);
const V well_mu = fluidMu(phase, p_perfwell.value(), well_cells); const V well_mu = fluidMu(phase, p_perfwell.value(), T_perfcell.value(), well_cells);
const V well_mob = well_kr / well_mu; const V well_mob = well_kr / well_mu;
const V perf_mob = cell_to_well_selector.select(subset(cell_mob, well_cells), well_mob); const V perf_mob = cell_to_well_selector.select(subset(cell_mob, well_cells), well_mob);
const ADB flux = face_mob * head; const ADB flux = face_mob * head;
@ -499,9 +503,11 @@ namespace {
const V transw = Eigen::Map<const V>(wells_.WI, nperf, 1); const V transw = Eigen::Map<const V>(wells_.WI, nperf, 1);
const V p = Eigen::Map<const V>(&state.pressure()[0], nc, 1); const V p = Eigen::Map<const V>(&state.pressure()[0], nc, 1);
const V T = Eigen::Map<const V>(&state.temperature()[0], nc, 1);
const V bhp = Eigen::Map<const V>(&well_state.bhp()[0], nw, 1); const V bhp = Eigen::Map<const V>(&well_state.bhp()[0], nw, 1);
const V p_perfcell = subset(p, well_cells); const V p_perfcell = subset(p, well_cells);
const V T_perfcell = subset(T, well_cells);
const V transi = subset(geo_.transmissibility(), const V transi = subset(geo_.transmissibility(),
ops_.internal_faces); ops_.internal_faces);
@ -515,15 +521,15 @@ namespace {
V perf_flux = V::Zero(nperf, 1); V perf_flux = V::Zero(nperf, 1);
for (int phase = 0; phase < np; ++phase) { for (int phase = 0; phase < np; ++phase) {
const V cell_rho = fluidRho(phase, p, cells); const V cell_rho = fluidRho(phase, p, T, cells);
const V head = nkgradp + (grav_ * cell_rho.matrix()).array(); const V head = nkgradp + (grav_ * cell_rho.matrix()).array();
const UpwindSelector<double> upwind(grid_, ops_, head); const UpwindSelector<double> upwind(grid_, ops_, head);
const V kr = fluidKr(phase); const V kr = fluidKr(phase);
const V mu = fluidMu(phase, p, cells); const V mu = fluidMu(phase, p, T, cells);
const V cell_mob = kr / mu; const V cell_mob = kr / mu;
const V face_mob = upwind.select(cell_mob); const V face_mob = upwind.select(cell_mob);
const V well_kr = fluidKrWell(phase); const V well_kr = fluidKrWell(phase);
const V well_mu = fluidMu(phase, p_perfwell, well_cells); const V well_mu = fluidMu(phase, p_perfwell, T_perfcell, well_cells);
const V well_mob = well_kr / well_mu; const V well_mob = well_kr / well_mu;
const V perf_mob = cell_to_well_selector.select(subset(cell_mob, well_cells), well_mob); const V perf_mob = cell_to_well_selector.select(subset(cell_mob, well_cells), well_mob);
@ -545,19 +551,19 @@ namespace {
V ImpesTPFAAD::fluidMu(const int phase, const V& p, const std::vector<int>& cells) const V ImpesTPFAAD::fluidMu(const int phase, const V& p, const V& T, const std::vector<int>& cells) const
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.muWat(p, cells); return fluid_.muWat(p, T, cells);
case Oil: { case Oil: {
V dummy_rs = V::Zero(p.size(), 1) * p; V dummy_rs = V::Zero(p.size(), 1) * p;
std::vector<PhasePresence> cond(dummy_rs.size()); std::vector<PhasePresence> cond(dummy_rs.size());
return fluid_.muOil(p, dummy_rs, cond, cells); return fluid_.muOil(p, T, dummy_rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.muGas(p, cells); return fluid_.muGas(p, T, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -567,19 +573,19 @@ namespace {
ADB ImpesTPFAAD::fluidMu(const int phase, const ADB& p, const std::vector<int>& cells) const ADB ImpesTPFAAD::fluidMu(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.muWat(p, cells); return fluid_.muWat(p, T, cells);
case Oil: { case Oil: {
ADB dummy_rs = V::Zero(p.size(), 1) * p; ADB dummy_rs = V::Zero(p.size(), 1) * p;
std::vector<PhasePresence> cond(dummy_rs.size()); std::vector<PhasePresence> cond(dummy_rs.size());
return fluid_.muOil(p, dummy_rs, cond, cells); return fluid_.muOil(p, T, dummy_rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.muGas(p, cells); return fluid_.muGas(p, T, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -589,19 +595,19 @@ namespace {
V ImpesTPFAAD::fluidFvf(const int phase, const V& p, const std::vector<int>& cells) const V ImpesTPFAAD::fluidFvf(const int phase, const V& p, const V& T, const std::vector<int>& cells) const
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.bWat(p, cells); return fluid_.bWat(p, T, cells);
case Oil: { case Oil: {
V dummy_rs = V::Zero(p.size(), 1) * p; V dummy_rs = V::Zero(p.size(), 1) * p;
std::vector<PhasePresence> cond(dummy_rs.size()); std::vector<PhasePresence> cond(dummy_rs.size());
return fluid_.bOil(p, dummy_rs, cond, cells); return fluid_.bOil(p, T, dummy_rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.bGas(p, cells); return fluid_.bGas(p, T, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -611,19 +617,19 @@ namespace {
ADB ImpesTPFAAD::fluidFvf(const int phase, const ADB& p, const std::vector<int>& cells) const ADB ImpesTPFAAD::fluidFvf(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const
{ {
switch (phase) { switch (phase) {
case Water: case Water:
return fluid_.bWat(p, cells); return fluid_.bWat(p, T, cells);
case Oil: { case Oil: {
ADB dummy_rs = V::Zero(p.size(), 1) * p; ADB dummy_rs = V::Zero(p.size(), 1) * p;
std::vector<PhasePresence> cond(dummy_rs.size()); std::vector<PhasePresence> cond(dummy_rs.size());
return fluid_.bOil(p, dummy_rs, cond, cells); return fluid_.bOil(p, T, dummy_rs, cond, cells);
} }
case Gas: case Gas:
return fluid_.bGas(p, cells); return fluid_.bGas(p, T, cells);
default: default:
OPM_THROW(std::runtime_error, "Unknown phase index " << phase); OPM_THROW(std::runtime_error, "Unknown phase index " << phase);
} }
@ -633,10 +639,10 @@ namespace {
V ImpesTPFAAD::fluidRho(const int phase, const V& p, const std::vector<int>& cells) const V ImpesTPFAAD::fluidRho(const int phase, const V& p, const V& T, const std::vector<int>& cells) const
{ {
const double* rhos = fluid_.surfaceDensity(); const double* rhos = fluid_.surfaceDensity();
V b = fluidFvf(phase, p, cells); V b = fluidFvf(phase, p, T, cells);
V rho = V::Constant(p.size(), 1, rhos[phase]) * b; V rho = V::Constant(p.size(), 1, rhos[phase]) * b;
return rho; return rho;
} }
@ -645,10 +651,10 @@ namespace {
ADB ImpesTPFAAD::fluidRho(const int phase, const ADB& p, const std::vector<int>& cells) const ADB ImpesTPFAAD::fluidRho(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const
{ {
const double* rhos = fluid_.surfaceDensity(); const double* rhos = fluid_.surfaceDensity();
ADB b = fluidFvf(phase, p, cells); ADB b = fluidFvf(phase, p, T, cells);
ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b; ADB rho = V::Constant(p.size(), 1, rhos[phase]) * b;
return rho; return rho;
} }

12
opm/autodiff/ImpesTPFAAD.hpp
View File

@ -104,12 +104,12 @@ namespace Opm {
void computeFluxes(BlackoilState& state, WellState& well_state) const; void computeFluxes(BlackoilState& state, WellState& well_state) const;
// Fluid interface forwarding calls to correct methods of fluid_. // Fluid interface forwarding calls to correct methods of fluid_.
V fluidMu(const int phase, const V& p, const std::vector<int>& cells) const; V fluidMu(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
ADB fluidMu(const int phase, const ADB& p, const std::vector<int>& cells) const; ADB fluidMu(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
V fluidFvf(const int phase, const V& p, const std::vector<int>& cells) const; V fluidFvf(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
ADB fluidFvf(const int phase, const ADB& p, const std::vector<int>& cells) const; ADB fluidFvf(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
V fluidRho(const int phase, const V& p, const std::vector<int>& cells) const; V fluidRho(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
ADB fluidRho(const int phase, const ADB& p, const std::vector<int>& cells) const; ADB fluidRho(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
V fluidKr(const int phase) const; V fluidKr(const int phase) const;
V fluidKrWell(const int phase) const; V fluidKrWell(const int phase) const;
}; };

53
opm/autodiff/RateConverter.hpp
View File

@ -267,6 +267,7 @@ namespace Opm {
, repcells_(Details::representative<typename Property::Cells>(rmap_)) , repcells_(Details::representative<typename Property::Cells>(rmap_))
, ncells_ (Details::countCells(rmap_)) , ncells_ (Details::countCells(rmap_))
, p_avg_ (rmap_.numRegions()) , p_avg_ (rmap_.numRegions())
, T_avg_ (rmap_.numRegions())
, Rmax_ (rmap_.numRegions(), props.numPhases()) , Rmax_ (rmap_.numRegions(), props.numPhases())
{} {}
@ -327,6 +328,7 @@ namespace Opm {
const PhaseUsage& pu = props_.phaseUsage(); const PhaseUsage& pu = props_.phaseUsage();
const V& p = getRegPress(r); const V& p = getRegPress(r);
const V& T = getRegTemp(r);
const typename Property::Cells& c = getRegCell (r); const typename Property::Cells& c = getRegCell (r);
const int iw = Details::PhasePos::water(pu); const int iw = Details::PhasePos::water(pu);
@ -338,7 +340,7 @@ namespace Opm {
if (Details::PhaseUsed::water(pu)) { if (Details::PhaseUsed::water(pu)) {
// q[w]_r = q[w]_s / bw // q[w]_r = q[w]_s / bw
const V& bw = props_.bWat(p, c); const V& bw = props_.bWat(p, T, c);
coeff[iw] = 1.0 / bw(0); coeff[iw] = 1.0 / bw(0);
} }
@ -351,7 +353,7 @@ namespace Opm {
if (Details::PhaseUsed::oil(pu)) { if (Details::PhaseUsed::oil(pu)) {
// q[o]_r = 1/(bo * (1 - rs*rv)) * (q[o]_s - rv*q[g]_s) // q[o]_r = 1/(bo * (1 - rs*rv)) * (q[o]_s - rv*q[g]_s)
const V& bo = props_.bOil(p, m.rs, m.cond, c); const V& bo = props_.bOil(p, T, m.rs, m.cond, c);
const double den = bo(0) * detR; const double den = bo(0) * detR;
coeff[io] += 1.0 / den; coeff[io] += 1.0 / den;
@ -364,7 +366,7 @@ namespace Opm {
if (Details::PhaseUsed::gas(pu)) { if (Details::PhaseUsed::gas(pu)) {
// q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s) // q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s)
const V& bg = props_.bGas(p, m.rv, m.cond, c); const V& bg = props_.bGas(p, T, m.rv, m.cond, c);
const double den = bg(0) * detR; const double den = bg(0) * detR;
coeff[ig] += 1.0 / den; coeff[ig] += 1.0 / den;
@ -404,6 +406,11 @@ namespace Opm {
*/ */
Eigen::ArrayXd p_avg_; Eigen::ArrayXd p_avg_;
/**
* Average temperature in each FIP region.
*/
Eigen::ArrayXd T_avg_;
/** /**
* Maximum dissolution and evaporation ratios at average * Maximum dissolution and evaporation ratios at average
* hydrocarbon pressure. * hydrocarbon pressure.
@ -474,6 +481,26 @@ namespace Opm {
p_avg_ /= ncells_; p_avg_ /= ncells_;
} }
/**
* Compute average temperature in all regions.
*
* \param[in] state Dynamic reservoir state.
*/
void
averageTemperature(const BlackoilState& state)
{
T_avg_.setZero();
const std::vector<double>& T = state.temperature();
for (std::vector<double>::size_type
i = 0, n = T.size(); i < n; ++i)
{
T_avg_(rmap_.region(i)) += T[i];
}
T_avg_ /= ncells_;
}
/** /**
* Compute maximum dissolution and evaporation ratios at * Compute maximum dissolution and evaporation ratios at
* average hydrocarbon pressure. * average hydrocarbon pressure.
@ -499,8 +526,8 @@ namespace Opm {
// pressure into account. This facility uses the // pressure into account. This facility uses the
// average *hydrocarbon* pressure rather than // average *hydrocarbon* pressure rather than
// average phase pressure. // average phase pressure.
Rmax_.col(io) = props_.rsSat(p_avg_, repcells_); Rmax_.col(io) = props_.rsSat(p_avg_, T_avg_, repcells_);
Rmax_.col(ig) = props_.rvSat(p_avg_, repcells_); Rmax_.col(ig) = props_.rvSat(p_avg_, T_avg_, repcells_);
} }
} }
@ -591,6 +618,22 @@ namespace Opm {
return p; return p;
} }
/**
* Retrieve average temperature in region.
*
* \param[in] r Particular region.
*
* \return Average temperature in region \c r.
*/
typename Property::V
getRegTemp(const RegionId r) const
{
typename Property::V T(1);
T << T_avg_(r);
return T;
}
/** /**
* Retrieve representative cell of region * Retrieve representative cell of region
* *

4
opm/autodiff/SimulatorCompressibleAd.cpp
View File

@ -357,7 +357,7 @@ namespace Opm
// Solve pressure equation. // Solve pressure equation.
if (check_well_controls_) { if (check_well_controls_) {
computeFractionalFlow(props_, allcells_, computeFractionalFlow(props_, allcells_,
state.pressure(), state.surfacevol(), state.saturation(), state.pressure(), state.temperature(), state.surfacevol(), state.saturation(),
fractional_flows); fractional_flows);
wells_manager_.applyExplicitReinjectionControls(well_resflows_phase, well_resflows_phase); wells_manager_.applyExplicitReinjectionControls(well_resflows_phase, well_resflows_phase);
} }
@ -446,7 +446,7 @@ namespace Opm
double injected[2] = { 0.0 }; double injected[2] = { 0.0 };
double produced[2] = { 0.0 }; double produced[2] = { 0.0 };
for (int tr_substep = 0; tr_substep < num_transport_substeps_; ++tr_substep) { for (int tr_substep = 0; tr_substep < num_transport_substeps_; ++tr_substep) {
tsolver_.solve(&state.faceflux()[0], &state.pressure()[0], tsolver_.solve(&state.faceflux()[0], &state.pressure()[0], &state.temperature()[0],
&initial_porevol[0], &porevol[0], &transport_src[0], stepsize, &initial_porevol[0], &porevol[0], &transport_src[0], stepsize,
state.saturation(), state.surfacevol()); state.saturation(), state.surfacevol());
double substep_injected[2] = { 0.0 }; double substep_injected[2] = { 0.0 };

15
tests/test_boprops_ad.cpp
View File

@ -124,7 +124,11 @@ BOOST_FIXTURE_TEST_CASE(ViscosityValue, TestFixture<SetupSimple>)
Vpw[3] = 8*Opm::unit::barsa; Vpw[3] = 8*Opm::unit::barsa;
Vpw[4] = 16*Opm::unit::barsa; Vpw[4] = 16*Opm::unit::barsa;
const Opm::BlackoilPropsAd::V VmuWat = boprops_ad.muWat(Vpw, cells); // standard temperature
V T;
T.resize(cells.size(), 273.15+20);
const Opm::BlackoilPropsAd::V VmuWat = boprops_ad.muWat(Vpw, T, cells);
// Zero pressure dependence in water viscosity // Zero pressure dependence in water viscosity
for (V::Index i = 0, n = VmuWat.size(); i < n; ++i) { for (V::Index i = 0, n = VmuWat.size(); i < n; ++i) {
@ -149,16 +153,21 @@ BOOST_FIXTURE_TEST_CASE(ViscosityAD, TestFixture<SetupSimple>)
Vpw[3] = 8*Opm::unit::barsa; Vpw[3] = 8*Opm::unit::barsa;
Vpw[4] = 16*Opm::unit::barsa; Vpw[4] = 16*Opm::unit::barsa;
// standard temperature
V T;
T.resize(cells.size(), 273.15+20);
typedef Opm::BlackoilPropsAd::ADB ADB; typedef Opm::BlackoilPropsAd::ADB ADB;
const V VmuWat = boprops_ad.muWat(Vpw, cells); const V VmuWat = boprops_ad.muWat(Vpw, T, cells);
for (V::Index i = 0, n = Vpw.size(); i < n; ++i) { for (V::Index i = 0, n = Vpw.size(); i < n; ++i) {
const std::vector<int> bp(1, grid.c_grid()->number_of_cells); const std::vector<int> bp(1, grid.c_grid()->number_of_cells);
const Opm::BlackoilPropsAd::Cells c(1, 0); const Opm::BlackoilPropsAd::Cells c(1, 0);
const V pw = V(1, 1) * Vpw[i]; const V pw = V(1, 1) * Vpw[i];
const ADB Apw = ADB::variable(0, pw, bp); const ADB Apw = ADB::variable(0, pw, bp);
const ADB AmuWat = boprops_ad.muWat(Apw, c); const ADB AT = ADB::constant(T);
const ADB AmuWat = boprops_ad.muWat(Apw, AT, c);
BOOST_CHECK_EQUAL(AmuWat.value()[0], VmuWat[i]); BOOST_CHECK_EQUAL(AmuWat.value()[0], VmuWat[i]);
} }