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

Added (fluid) compressible pressure solver for polymer.

Browse Source
This commit is contained in:
Xavier Raynaud
2012-08-27 13:16:10 +02:00
parent 72658ac91c
commit 9a271a6038
4 changed files with 374 additions and 1 deletions
Download Patch File Download Diff File Expand all files Collapse all files

176
opm/polymer/CompressibleTpfaPolymer.cpp Normal file
View File

@@ -0,0 +1,176 @@
/*
Copyright 2012 SINTEF ICT, Applied Mathematics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <opm/polymer/PolymerBlackoilState.hpp>
#include <opm/polymer/CompressibleTpfaPolymer.hpp>
#include <opm/core/fluid/BlackoilPropertiesInterface.hpp>
#include <opm/core/fluid/RockCompressibility.hpp>
#include <opm/core/pressure/tpfa/ifs_tpfa.h>
#include <opm/core/pressure/tpfa/trans_tpfa.h>
#include <opm/core/pressure/mimetic/mimetic.h>
#include <opm/core/pressure/flow_bc.h>
#include <opm/core/linalg/LinearSolverInterface.hpp>
#include <opm/core/linalg/sparse_sys.h>
#include <opm/polymer/polymerUtilities.hpp>
#include <opm/core/simulator/WellState.hpp>
#include <opm/core/utility/ErrorMacros.hpp>
#include <opm/core/utility/miscUtilities.hpp>
#include <opm/core/newwells.h>
#include <iomanip>
#include <cmath>
#include <algorithm>
namespace Opm
{
/// Construct solver
/// \param[in] grid A 2d or 3d grid.
/// \param[in] props Rock and fluid properties.
/// \param[in] linsolver Linear solver to use.
/// \param[in] residual_tol Solution accepted if inf-norm of residual is smaller.
/// \param[in] change_tol Solution accepted if inf-norm of change in pressure is smaller.
/// \param[in] maxiter Maximum acceptable number of iterations.
/// \param[in] gravity Gravity vector. If non-null, the array should
/// have D elements.
/// \param[in] wells The wells argument. Will be used in solution,
/// is ignored if NULL.
/// Note: this class observes the well object, and
/// makes the assumption that the well topology
/// and completions does not change during the
/// run. However, controls (only) are allowed
/// to change.
CompressibleTpfaPolymer::CompressibleTpfaPolymer(const UnstructuredGrid& grid,
const BlackoilPropertiesInterface& props,
const RockCompressibility* rock_comp_props,
const PolymerProperties& poly_props,
const LinearSolverInterface& linsolver,
const double residual_tol,
const double change_tol,
const int maxiter,
const double* gravity,
const Wells* wells)
: CompressibleTpfa(grid, props, rock_comp_props, linsolver,
residual_tol, change_tol, maxiter,
gravity, wells),
poly_props_(poly_props),
c_(0),
cmax_(0)
{
}
/// Solve the pressure equation. The nonlinear equations ares
/// solved by a Newton-Raphson scheme. May throw an exception if
/// the number of iterations exceed maxiter (set in constructor).
void CompressibleTpfaPolymer::solve(const double dt,
PolymerBlackoilState& state,
WellState& well_state)
{
c_ = &state.concentration();
cmax_ = &state.maxconcentration();
CompressibleTpfa::solve(dt, state.blackoilState(), well_state);
}
/// Compute per-solve dynamic properties.
void CompressibleTpfaPolymer::computePerSolveDynamicData(const double /* dt */,
const BlackoilState& state,
const WellState& /* well_state */)
{
// std::vector<double> cell_relperm__;
// std::vector<double> cell_eff_relperm_;
const int nc = grid_.number_of_cells;
const int np = props_.numPhases();
cell_relperm_.resize(nc*np);
cell_eff_relperm_.resize(nc*np);
const double* cell_s = &state.saturation()[0];
props_.relperm(nc, cell_s, &allcells_[0], &cell_relperm_[0], 0);
std::copy(cell_relperm_.begin(), cell_relperm_.end(), cell_eff_relperm_.begin());
for (int cell; cell < grid_.number_of_cells; ++cell) {
// only the water phase is modified by the presence og polymer.
poly_props_.effectiveRelperm((*c_)[cell], (*cmax_)[cell], &cell_relperm_[nc + 0], cell_eff_relperm_[nc + 0]);
}
computeWellPotentials(state);
if (rock_comp_props_ && rock_comp_props_->isActive()) {
computePorevolume(grid_, props_.porosity(), *rock_comp_props_, state.pressure(), initial_porevol_);
}
}
/// Compute per-iteration dynamic properties for cells.
void CompressibleTpfaPolymer::computeCellDynamicData(const double /*dt*/,
const BlackoilState& state,
const WellState& /*well_state*/)
{
// These are the variables that get computed by this function:
//
// std::vector<double> cell_A_;
// std::vector<double> cell_dA_;
// std::vector<double> cell_viscosity_;
// std::vector<double> cell_eff_viscosity_;
// std::vector<double> cell_phasemob_;
// std::vector<double> cell_voldisc_;
// std::vector<double> porevol_; // Only modified if rock_comp_props_ is non-null.
// std::vector<double> rock_comp_; // Empty unless rock_comp_props_ is non-null.
const int nc = grid_.number_of_cells;
const int np = props_.numPhases();
const double* cell_p = &state.pressure()[0];
const double* cell_z = &state.surfacevol()[0];
const double* cell_s = &state.saturation()[0];
cell_A_.resize(nc*np*np);
cell_dA_.resize(nc*np*np);
props_.matrix(nc, cell_p, cell_z, &allcells_[0], &cell_A_[0], &cell_dA_[0]);
cell_viscosity_.resize(nc*np);
props_.viscosity(nc, cell_p, cell_z, &allcells_[0], &cell_viscosity_[0], 0);
cell_phasemob_.resize(nc*np);
for (int cell; cell < nc; ++cell) {
// Only the water viscosity is modified by polymer.
poly_props_.effectiveVisc((*c_)[cell], &cell_viscosity_[nc + 0], cell_eff_viscosity_[nc + 0]);
}
std::transform(cell_eff_relperm_.begin(), cell_eff_relperm_.end(),
cell_eff_viscosity_.begin(),
cell_phasemob_.begin(),
std::divides<double>());
// Volume discrepancy: we have that
// z = Au, voldiscr = sum(u) - 1,
// but I am not sure it is actually needed.
// Use zero for now.
// TODO: Check this!
cell_voldisc_.clear();
cell_voldisc_.resize(nc, 0.0);
if (rock_comp_props_ && rock_comp_props_->isActive()) {
computePorevolume(grid_, props_.porosity(), *rock_comp_props_, state.pressure(), porevol_);
rock_comp_.resize(nc);
for (int cell = 0; cell < nc; ++cell) {
rock_comp_[cell] = rock_comp_props_->rockComp(state.pressure()[cell]);
}
}
}
} // namespace Opm

114
opm/polymer/CompressibleTpfaPolymer.hpp Normal file
View File

@@ -0,0 +1,114 @@
/*
Copyright 2012 SINTEF ICT, Applied Mathematics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_INCOMPTPFAPOLYMER_HEADER_INCLUDED
#define OPM_INCOMPTPFAPOLYMER_HEADER_INCLUDED
#include <opm/core/pressure/CompressibleTpfa.hpp>
#include <vector>
struct UnstructuredGrid;
struct Wells;
struct FlowBoundaryConditions;
namespace Opm
{
class BlackoilState;
class PolymerBlackoilState;
class RockCompressibility;
class PolymerProperties;
class LinearSolverInterface;
class PolymerState;
class WellState;
/// Encapsulating a tpfa pressure solver for the compressible-fluid case with polymer.
/// Supports gravity, wells controlled by bhp or reservoir rates,
/// boundary conditions and simple sources as driving forces.
/// Rock compressibility can be included, and necessary nonlinear
/// iterations are handled.
/// Below we use the shortcuts D for the number of dimensions, N
/// for the number of cells and F for the number of faces.
class CompressibleTpfaPolymer : public CompressibleTpfa
{
public:
/// Construct solver, possibly with rock compressibility.
/// \param[in] grid A 2d or 3d grid.
/// \param[in] props Rock and fluid properties.
/// \param[in] rock_comp_props Rock compressibility properties. May be null.
/// \param[in] poly_props Polymer properties.
/// \param[in] linsolver Linear solver to use.
/// \param[in] residual_tol Solution accepted if inf-norm of residual is smaller.
/// \param[in] change_tol Solution accepted if inf-norm of change in pressure is smaller.
/// \param[in] maxiter Maximum acceptable number of iterations.
/// \param[in] gravity Gravity vector. If non-null, the array should
/// have D elements.
/// \param[in] wells The wells argument. Will be used in solution,
/// is ignored if NULL.
/// Note: this class observes the well object, and
/// makes the assumption that the well topology
/// and completions does not change during the
/// run. However, controls (only) are allowed
/// to change.
CompressibleTpfaPolymer(const UnstructuredGrid& grid,
const BlackoilPropertiesInterface& props,
const RockCompressibility* rock_comp_props,
const PolymerProperties& poly_props,
const LinearSolverInterface& linsolver,
const double residual_tol,
const double change_tol,
const int maxiter,
const double* gravity,
const Wells* wells);
/// Solve the pressure equation. The nonlinear equations ares solved by a
/// Newton-Raphson scheme.
/// May throw an exception if the number of iterations
/// exceed maxiter (set in constructor).
void solve(const double dt,
PolymerBlackoilState& state,
WellState& well_state);
private:
virtual void computeCellDynamicData(const double dt,
const BlackoilState& state,
const WellState& well_state);
virtual void computePerSolveDynamicData(const double dt,
const BlackoilState& state,
const WellState& well_state);
private:
// ------ Data that will remain unmodified after construction. ------
const PolymerProperties& poly_props_;
// ------ Data that will be updated every solve() call. ------
const std::vector<double>* c_;
const std::vector<double>* cmax_;
std::vector<double> cell_eff_viscosity_;
std::vector<double> cell_relperm_;
std::vector<double> cell_eff_relperm_;
};
} // namespace Opm
#endif // OPM_INCOMPTPFAPOLYMER_HEADER_INCLUDED

82
opm/polymer/PolymerBlackoilState.hpp Normal file
View File

@@ -0,0 +1,82 @@
/*
Copyright 2012 SINTEF ICT, Applied Mathematics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_POLYMERSTATE_HEADER_INCLUDED
#define OPM_POLYMERSTATE_HEADER_INCLUDED
#include <opm/core/simulator/BlackoilState.hpp>
#include <opm/core/grid.h>
#include <vector>
namespace Opm
{
/// Simulator state for a compressible two-phase simulator with polymer.
/// We use the Blackoil state parameters.
class PolymerBlackoilState
{
public:
void init(const UnstructuredGrid& g, int num_phases)
{
state_blackoil_.init(g, num_phases);
concentration_.resize(g.number_of_cells, 0.0);
cmax_.resize(g.number_of_cells, 0.0);
}
enum ExtremalSat { MinSat = BlackoilState::MinSat, MaxSat = BlackoilState::MaxSat };
void setFirstSat(const std::vector<int>& cells,
const Opm::BlackoilPropertiesInterface& props,
ExtremalSat es)
{
// A better solution for embedding BlackoilState::ExtremalSat could perhaps
// be found, to avoid the cast.
state_blackoil_.setFirstSat(cells, props, static_cast<BlackoilState::ExtremalSat>(es));
}
std::vector<double>& pressure () { return state_blackoil_.pressure(); }
std::vector<double>& facepressure() { return state_blackoil_.facepressure(); }
std::vector<double>& faceflux () { return state_blackoil_.faceflux(); }
std::vector<double>& saturation () { return state_blackoil_.saturation(); }
std::vector<double>& concentration() { return concentration_; }
std::vector<double>& maxconcentration() { return cmax_; }
const std::vector<double>& pressure () const { return state_blackoil_.pressure(); }
const std::vector<double>& facepressure() const { return state_blackoil_.facepressure(); }
const std::vector<double>& faceflux () const { return state_blackoil_.faceflux(); }
const std::vector<double>& saturation () const { return state_blackoil_.saturation(); }
const std::vector<double>& concentration() const { return concentration_; }
const std::vector<double>& maxconcentration() const { return cmax_; }
BlackoilState& blackoilState() { return state_blackoil_; }
const BlackoilState& blackoilState() const { return state_blackoil_; }
private:
BlackoilState state_blackoil_;
std::vector<double> concentration_;
std::vector<double> cmax_;
};
} // namespace Opm
#endif // OPM_POLYMERSTATE_HEADER_INCLUDED