OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-07 23:13:01 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'master' into gravity-in-wells

Browse Source
This commit is contained in:
Bård Skaflestad 2012-10-17 13:30:37 +02:00
commit 8e38cd5f5c
2 changed files with 29 additions and 31 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

46
opm/core/transport/reorder/TransportModelTracerTofDiscGal.cpp
View File

@ -20,10 +20,11 @@
#include <opm/core/transport/reorder/TransportModelTracerTofDiscGal.hpp>
#include <opm/core/grid.h>
#include <opm/core/utility/ErrorMacros.hpp>
#include <opm/core/utility/VelocityInterpolation.hpp>
#include <opm/core/linalg/blas_lapack.h>
#include <algorithm>
#include <numeric>
#include <cmath>
#include <numeric>
namespace Opm
{
@ -403,32 +404,6 @@ namespace Opm
// Initial version: only a constant interpolation.
static void interpolateVelocity(const UnstructuredGrid& grid,
const int cell,
const double* darcyflux,
const double* /*x*/,
double* v)
{
const int dim = grid.dimensions;
std::fill(v, v + dim, 0.0);
const double* cc = grid.cell_centroids + cell*dim;
for (int hface = grid.cell_facepos[cell]; hface < grid.cell_facepos[cell+1]; ++hface) {
const int face = grid.cell_faces[hface];
const double* fc = grid.face_centroids + face*dim;
double flux = 0.0;
if (cell == grid.face_cells[2*face]) {
flux = darcyflux[face];
} else {
flux = -darcyflux[face];
}
for (int dd = 0; dd < dim; ++dd) {
v[dd] += flux * (fc[dd] - cc[dd]) / grid.cell_volumes[cell];
}
}
}
// --------------- Methods of TransportModelTracerTofDiscGal ---------------
@ -436,11 +411,19 @@ namespace Opm
/// Construct solver.
/// \param[in] grid A 2d or 3d grid.
TransportModelTracerTofDiscGal::TransportModelTracerTofDiscGal(const UnstructuredGrid& grid)
/// \param[in] use_cvi If true, use corner point velocity interpolation.
/// Otherwise, use the basic constant interpolation.
TransportModelTracerTofDiscGal::TransportModelTracerTofDiscGal(const UnstructuredGrid& grid,
const bool use_cvi)
: grid_(grid),
coord_(grid.dimensions),
velocity_(grid.dimensions)
{
if (use_cvi) {
velocity_interpolation_.reset(new VelocityInterpolationECVI(grid));
} else {
velocity_interpolation_.reset(new VelocityInterpolationConstant(grid));
}
}
@ -482,9 +465,11 @@ namespace Opm
tof_coeff_ = &tof_coeff[0];
rhs_.resize(num_basis);
jac_.resize(num_basis*num_basis);
orig_jac_.resize(num_basis*num_basis);
basis_.resize(num_basis);
basis_nb_.resize(num_basis);
grad_basis_.resize(num_basis*grid_.dimensions);
velocity_interpolation_->setupFluxes(darcyflux);
reorderAndTransport(grid_, darcyflux);
}
@ -560,7 +545,7 @@ namespace Opm
quad.quadPtCoord(quad_pt, &coord_[0]);
DGBasis::eval(grid_, cell, degree_, &coord_[0], &basis_[0]);
DGBasis::evalGrad(grid_, cell, degree_, &coord_[0], &grad_basis_[0]);
interpolateVelocity(grid_, cell, darcyflux_, &coord_[0], &velocity_[0]);
velocity_interpolation_->interpolate(cell, &coord_[0], &velocity_[0]);
const double w = quad.quadPtWeight(quad_pt);
for (int j = 0; j < num_basis; ++j) {
for (int i = 0; i < num_basis; ++i) {
@ -633,6 +618,7 @@ namespace Opm
std::vector<MAT_SIZE_T> piv(num_basis);
MAT_SIZE_T ldb = num_basis;
MAT_SIZE_T info = 0;
orig_jac_ = jac_;
dgesv_(&n, &nrhs, &jac_[0], &lda, &piv[0], &rhs_[0], &ldb, &info);
if (info != 0) {
// Print the local matrix and rhs.
@ -640,7 +626,7 @@ namespace Opm
<< " with A = \n";
for (int row = 0; row < n; ++row) {
for (int col = 0; col < n; ++col) {
std::cerr << " " << jac_[row + n*col];
std::cerr << " " << orig_jac_[row + n*col];
}
std::cerr << '\n';
}

14
opm/core/transport/reorder/TransportModelTracerTofDiscGal.hpp
View File

@ -21,15 +21,18 @@
#define OPM_TRANSPORTMODELTRACERTOFDISCGAL_HEADER_INCLUDED
#include <opm/core/transport/reorder/TransportModelInterface.hpp>
#include <boost/shared_ptr.hpp>
#include <vector>
#include <map>
#include <ostream>
struct UnstructuredGrid;
namespace Opm
{
class IncompPropertiesInterface;
class VelocityInterpolationInterface;
/// Implements a discontinuous Galerkin solver for
/// (single-phase) time-of-flight using reordering.
@ -45,7 +48,10 @@ namespace Opm
public:
/// Construct solver.
/// \param[in] grid A 2d or 3d grid.
TransportModelTracerTofDiscGal(const UnstructuredGrid& grid);
/// \param[in] use_cvi If true, use corner point velocity interpolation.
/// Otherwise, use the basic constant interpolation.
TransportModelTracerTofDiscGal(const UnstructuredGrid& grid,
const bool use_cvi);
/// Solve for time-of-flight.
@ -72,7 +78,12 @@ namespace Opm
virtual void solveMultiCell(const int num_cells, const int* cells);
private:
// Disable copying and assignment.
TransportModelTracerTofDiscGal(const TransportModelTracerTofDiscGal&);
TransportModelTracerTofDiscGal& operator=(const TransportModelTracerTofDiscGal&);
const UnstructuredGrid& grid_;
boost::shared_ptr<VelocityInterpolationInterface> velocity_interpolation_;
const double* darcyflux_; // one flux per grid face
const double* porevolume_; // one volume per cell
const double* source_; // one volumetric source term per cell
@ -80,6 +91,7 @@ namespace Opm
double* tof_coeff_;
std::vector<double> rhs_; // single-cell right-hand-side
std::vector<double> jac_; // single-cell jacobian
std::vector<double> orig_jac_; // single-cell jacobian (copy)
// Below: storage for quantities needed by solveSingleCell().
std::vector<double> coord_;
std::vector<double> basis_;