OilfieldToolsNet/opm-core
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge remote-tracking branch 'upstream/master'

Browse Source
This commit is contained in:
Atgeirr Flø Rasmussen
2012-07-19 15:29:14 +02:00
parent 2dc1744bce 8fefee7295
commit 342a2f2668
4 changed files with 324 additions and 86 deletions
Download Patch File Download Diff File Expand all files Collapse all files

127
opm/core/grid/cart_grid.h
View File

@@ -30,25 +30,136 @@
#ifndef OPM_CART_GRID_H_HEADER
#define OPM_CART_GRID_H_HEADER
/**
* \file
* Routines to construct fully formed grid structures from a simple Cartesian
* (i.e., tensor product) description.
*
* The cells are lexicographically ordered with the @c i index cycling the most
* rapidly, followed by the @c j index and then, in three space dimensions, the
* @c k (`layer') index as the least rapidly cycling index.
*/
#ifdef __cplusplus
extern "C" {
#endif
struct UnstructuredGrid;
struct UnstructuredGrid *create_grid_cart2d(int nx, int ny);
struct UnstructuredGrid *create_grid_cart3d(int nx, int ny, int nz);
struct UnstructuredGrid *create_grid_hexa3d(int nx, int ny, int nz,
double dx, double dy, double dz);
/**
* Form geometrically Cartesian grid in two space dimensions with unit-sized
* cells.
*
* @param[in] nx Number of cells in @c x direction.
* @param[in] ny Number of cells in @c y direction.
*
* @return Fully formed grid structure containing valid geometric primitives.
* Must be destroyed using function destroy_grid().
*/
struct UnstructuredGrid *
create_grid_cart2d(int nx, int ny);
struct UnstructuredGrid *create_grid_tensor2d(int nx, int ny,
double x[], double y[]);
/**
* Form geometrically Cartesian grid in three space dimensions with unit-sized
* cells.
*
* @param[in] nx Number of cells in @c x direction.
* @param[in] ny Number of cells in @c y direction.
* @param[in] nz Number of cells in @c z direction.
*
* @return Fully formed grid structure containing valid geometric primitives.
* Must be destroyed using function destroy_grid().
*/
struct UnstructuredGrid *
create_grid_cart3d(int nx, int ny, int nz);
/**
* Form geometrically Cartesian grid in three space dimensions with equally
* sized cells.
*
* Each cell has physical size (volume) \f$\mathit{dx}\times \mathit{dy}\times
* \mathit{dz}\f$.
*
* @param[in] nx Number of cells in @c x direction.
* @param[in] ny Number of cells in @c y direction.
* @param[in] nz Number of cells in @c z direction.
*
* @param[in] dx Length, in meters, of each cell's @c x extent.
* @param[in] dy Length, in meters, of each cell's @c y extent.
* @param[in] dz Length, in meters, of each cell's @c z extent.
*
* @return Fully formed grid structure containing valid geometric primitives.
* Must be destroyed using function destroy_grid().
*/
struct UnstructuredGrid *
create_grid_hexa3d(int nx, int ny, int nz,
double dx, double dy, double dz);
/**
* Form tensor product (Cartesian) grid in two space dimensions.
*
* The size (volume) of cell \f$(i,j)\f$ is
* \f[
* v_{ij} = (x_{i+1} - x_i)\cdot (y_{j+1} - y_j)
* \f]
* Similar relations hold for the cell and interface centroids as well as the
* interface areas and normal vectors. In other words, cell \f$(i,j)\f$ is the
* convex hull bounded by the tensor product of nodes \f$x_i\f$, \f$x_{i+1}\f$,
* \f$y_j\f$, and \f$y_{j+1}\f$.
*
* @param[in] nx Number of cells in @c x direction.
* @param[in] ny Number of cells in @c y direction.
*
* @param[in] x Position along @c x axis of each grid line with constant @c x
* coordinate. Array of size <CODE>nx + 1</CODE>.
* @param[in] y Position along @c y axis of each grid line with constant @c y
* coordinate. Array of size <CODE>ny + 1</CODE>.
*
* @return Fully formed grid structure containing valid geometric primitives.
* Must be destroyed using function destroy_grid().
*/
struct UnstructuredGrid *
create_grid_tensor2d(int nx, int ny,
double x[], double y[]);
/**
* Form tensor product (i.e., topologically Cartesian) grid in three space
* dimensions--possibly with a variable top-layer topography.
*
* If @c depthz is @c NULL, then geometric information such as volumes and
* centroids is calculated from analytic expressions. Otherwise, these values
* are computed using function compute_geometry().
*
* @param[in] nx Number of cells in @c x direction.
* @param[in] ny Number of cells in @c y direction.
* @param[in] nz Number of cells in @c z direction.
*
* @param[in] x Position along @c x axis of each grid line with constant @c x
* coordinate. Array of size <CODE>nx + 1</CODE>.
* @param[in] y Position along @c y axis of each grid line with constant @c y
* coordinate. Array of size <CODE>ny + 1</CODE>.
* @param[in] z Distance (depth) from top-layer measured along the @c z axis of
* each grid line with constant @c z coordinate. Array of size
* <CODE>nz + 1</CODE>.
*
* @param[in] depthz
* Top-layer topography specification. If @c NULL, interpreted as
* horizontal top-layer at <CODE>z=0</CODE>. Otherwise, must be
* an array of size <CODE>(nx + 1) * (ny + 1)</CODE>, ordered
* lexicographically, that defines the depth of each top-layer
* pillar vertex.
*
* @return Fully formed grid structure containing valid geometric primitives.
* Must be destroyed using function destroy_grid().
*/
struct UnstructuredGrid *
create_grid_tensor3d(int nx, int ny, int nz,
double x[], double y[], double z[],
const double depthz[]);
#ifdef __cplusplus
}
#endif

50
opm/core/grid/cornerpoint_grid.h
View File

@@ -34,6 +34,12 @@
#ifndef OPM_CORNERPOINT_GRID_HEADER_INCLUDED
#define OPM_CORNERPOINT_GRID_HEADER_INCLUDED
/**
* \file
* Routines to form a complete UnstructuredGrid from a corner-point
* specification.
*/
#include <opm/core/grid.h>
#include <opm/core/grid/cpgpreprocess/preprocess.h>
@@ -41,9 +47,53 @@
extern "C" {
#endif
/**
* Construct grid representation from corner-point specification of a
* particular geological model.
*
* Pinched cells will be removed irrespective of any explicit "active" map
* in the geological model input specification. Geometric primitives such
* as cell barycenters (i.e., centroids), volumes and interface areas are
* computed internally using function compute_geometry(). The caller does
* not need to compute this information separately.
*
* @param[in] in Corner-point specification. If "actnum" is NULL, then the
* specification is interpreted as if all cells are initially
* active.
*
* @param[in] tol Absolute tolerance of node-coincidence.
*
* @return Fully formed grid data structure that manages the grid defined by
* the input corner-point specification. Must be destroyed using function
* destroy_grid().
*/
struct UnstructuredGrid *
create_grid_cornerpoint(const struct grdecl *in, double tol);
/**
* Compute derived geometric primitives in a grid.
*
* This function computes values for each of the following quantities
* - Quantities pertaining to interfaces (connections, faces)
* -# Barycenters (centroids), <CODE>g->dimensions</CODE> scalars per face
* stored sequentially in <CODE>g->face_centroids</CODE>.
* -# Areas, one scalar per face stored sequentially in
* <CODE>g->face_areas</CODE>.
* -# Normals, <CODE>g->dimensions</CODE> scalars per face stored
* sequentially in <CODE>g->face_normals</CODE>. The Euclidian norm of
* each normal is equal to the corresponding face's area.
*
* - Quantities pertaining to cells (volumes)
* -# Barycenters (centroids), <CODE>g->dimensions</CODE> scalars per cell
* stored sequentially in <CODE>g->cell_centroids</CODE>.
* -# Volumes, one scalar per cell stored sequentially in
* <CODE>g->cell_volumes</CODE>.
*
* These fields must be allocated prior to calling compute_geometry().
*
* @param[in,out] g Grid structure.
*/
void compute_geometry(struct UnstructuredGrid *g);
#ifdef __cplusplus

124
opm/core/grid/cpgpreprocess/preprocess.h
View File

@@ -35,47 +35,111 @@
#ifndef OPENRS_PREPROCESS_HEADER
#define OPENRS_PREPROCESS_HEADER
/**
* \file
* Low-level corner-point processing routines and supporting data structures.
*
* User code should typically employ higher-level routines such as
* create_grid_cornerpoint() in order to construct fully formed UnstructuredGrid
* data structures from a corner-point specification. Incidentally, the routines
* provided by this module are used to implement function
* create_grid_cornerpoint().
*/
#ifdef __cplusplus
extern "C" {
#endif
/* Input structure holding raw cornerpoint spec. */
struct grdecl{
int dims[3];
const double *coord;
const double *zcorn;
const int *actnum;
/**
* Raw corner-point specification of a particular geological model.
*/
struct grdecl {
int dims[3]; /**< Cartesian box dimensions. */
const double *coord; /**< Pillar end-points. */
const double *zcorn; /**< Explicit "active" map. May be NULL.*/
const int *actnum; /**< Corner-point depths. */
};
/* Constant: I J K */
enum face_tag { LEFT, BACK, TOP };
/* Output structure holding grid topology */
struct processed_grid{
int m; /** Upper bound on "number_of_faces" */
int n; /** Upper bound on "number_of_nodes" */
int dimensions[3]; /* Cartesian dimension */
int number_of_faces;
int *face_nodes; /* Nodes numbers of each face sequentially. */
int *face_ptr; /* Start position for each face in face_nodes. */
int *face_neighbors; /* Global cell numbers. 2 ints per face sequentially */
enum face_tag *face_tag;
int number_of_nodes;
int number_of_nodes_on_pillars; /** Total number of unique cell vertices that lie on pillars. */
double *node_coordinates; /* 3 doubles per node, sequentially */
int number_of_cells; /* number of active cells */
int *local_cell_index; /* Global to local map */
/**
* Connection taxonomy.
*/
enum face_tag {
LEFT, /**< Connection topologically parallel to J-K plane. */
BACK, /**< Connection topologically parallel to I-K plane. */
TOP /**< Connection topologically parallel to I-J plane. */
};
void process_grdecl (const struct grdecl *g,
double tol,
struct processed_grid *out);
/**
* Result structure representing minimal derived topology and geometry of
* a geological model in corner-point format.
*/
struct processed_grid {
int m; /**< Upper bound on "number_of_faces". For internal use in
function process_grid()'s memory management. */
int n; /**< Upper bound on "number_of_nodes". For internal use in
function process_grid()'s memory management. */
int dimensions[3]; /**< Cartesian box dimensions. */
int number_of_faces; /**< Total number of unique grid faces
(i.e., connections). */
int *face_nodes; /**< Node (vertex) numbers of each face,
stored sequentially. */
int *face_ptr; /**< Start position for each face's
`face_nodes'. */
int *face_neighbors; /**< Global cell numbers. Two elements per
face, stored sequentially. */
enum face_tag *face_tag; /**< Classification of grid's individual
connections (faces). */
int number_of_nodes; /**< Number of unique grid vertices. */
int number_of_nodes_on_pillars; /**< Total number of unique cell
vertices that lie on pillars. */
double *node_coordinates; /**< Vertex coordinates. Three doubles
(\f$x\f$, \f$y\f$, \f$z\f$) per vertex,
stored sequentially. */
int number_of_cells; /**< Number of active grid cells. */
int *local_cell_index; /**< Deceptively named local-to-global cell
index mapping. */
};
/**
* Construct a prototypical grid representation from a corner-point
* specification.
*
* Pinched cells will be removed irrespective of any explicit "active" map
* in the geological model input specification. On input, the result
* structure "out" must point to a valid management structure. In other
* words, the result structure must point to a region of memory that is
* typically backed by automatic or allocated (dynamic) storage duration.
*
* @param[in] g Corner-point specification. If "actnum" is NULL, then
* the specification is interpreted as if all cells are
* initially active.
* @param[in] tol Absolute tolerance of node-coincidence.
* @param[in,out] out Minimal grid representation featuring face-to-cell
* neighbourship definition, vertex geometry, face's
* constituent vertices, and local-to-global cell
* mapping.
*/
void process_grdecl(const struct grdecl *g ,
double tol,
struct processed_grid *out);
/**
* Release memory resources acquired in previous grid processing using
* function process_grdecl().
*
* Note: This function releases the resources associated to the individual
* fields of the processed_grid, but does not free() the structure itself.
*
* @param[in,out] g Prototypical grid representation obtained in an earlier
* call to function process_grdecl().
*/
void free_processed_grid(struct processed_grid *g);
#ifdef __cplusplus
}
#endif

109
opm/core/utility/Units.hpp
View File

@@ -35,18 +35,24 @@
#ifndef OPENRS_UNITS_HEADER
#define OPENRS_UNITS_HEADER
/**
* \file
* Constants and routines to assist in handling units of measurement. These are
* geared towards handling common units in reservoir descriptions.
*/
namespace Opm
{
namespace prefix
/// Conversion prefix for units.
{
const double micro = 1.0e-6;
const double milli = 1.0e-3;
const double centi = 1.0e-2;
const double deci = 1.0e-1;
const double kilo = 1.0e3;
const double mega = 1.0e6;
const double giga = 1.0e9;
const double micro = 1.0e-6; /**< Unit prefix [\f$\mu\f$] */
const double milli = 1.0e-3; /**< Unit prefix [m] */
const double centi = 1.0e-2; /**< Non-standard unit prefix [c] */
const double deci = 1.0e-1; /**< Non-standard unit prefix [d] */
const double kilo = 1.0e3; /**< Unit prefix [k] */
const double mega = 1.0e6; /**< Unit prefix [M] */
const double giga = 1.0e9; /**< Unit prefix [G] */
} // namespace prefix
namespace unit
@@ -137,15 +143,6 @@ namespace Opm
const double Poise = prefix::deci*Pas;
/// @}
/// \name Permeability
/// @{
///
/// A porous medium with a permeability of 1 darcy permits a
/// flow (flux) of \f$1\,cm^3/s\f$ of a fluid with viscosity
/// \f$1\,cP\f$ (\f$1\,mPa\cdot s\f$) under a pressure
/// gradient of \f$1\,atm/cm\f$ acting across an area of
/// \f$1\,cm^2\f$.
///
namespace perm_details {
const double p_grad = atm / (prefix::centi*meter);
const double area = square(prefix::centi*meter);
@@ -156,47 +153,63 @@ namespace Opm
// == 1e-7 [m^2] / 101325
// == 9.869232667160130e-13 [m^2]
}
/// \name Permeability
/// @{
///
/// A porous medium with a permeability of 1 darcy permits a flow (flux)
/// of \f$1\,\mathit{cm}^3/s\f$ of a fluid with viscosity
/// \f$1\,\mathit{cP}\f$ (\f$1\,mPa\cdot s\f$) under a pressure gradient
/// of \f$1\,\mathit{atm}/\mathit{cm}\f$ acting across an area of
/// \f$1\,\mathit{cm}^2\f$.
///
const double darcy = perm_details::darcy;
/// @}
// Unit conversion support.
//
// Note: Under the penalty of treason will you be
//
// using namespace Opm::unit::convert;
//
// I mean it!
//
/**
* Unit conversion routines.
*/
namespace convert {
// Convert from external units of measurements to equivalent
// internal units of measurements. Note: The internal units
// of measurements are *ALWAYS*, and exclusively, SI.
//
// Example: Convert a double kx, containing a permeability
// value in units of milli-darcy (mD) to the equivalent
// value in SI units (m^2).
//
// using namespace Opm::unit;
// using namespace Opm::prefix;
// convert::from(kx, milli*darcy);
//
/**
* Convert from external units of measurements to equivalent
* internal units of measurements. Note: The internal units of
* measurements are *ALWAYS*, and exclusively, SI.
*
* Example: Convert a double @c kx, containing a permeability value
* in units of milli-darcy (mD) to the equivalent value in SI units
* (i.e., \f$m^2\f$).
* \code
* using namespace Opm::unit;
* using namespace Opm::prefix;
* convert::from(kx, milli*darcy);
* \endcode
*
* @param[in] q Physical quantity.
* @param[in] unit Physical unit of measurement.
* @return Value of @c q in equivalent SI units of measurements.
*/
inline double from(const double q, const double unit)
{
return q * unit;
}
// Convert from internal units of measurements to equivalent
// external units of measurements. Note: The internal units
// of measurements are *ALWAYS*, and exclusively, SI.
//
// Example: Convert a std::vector<double> p, containing
// pressure values in the SI unit Pascal (i.e., unit::Pascal)
// to the equivalent values in Psi (unit::psia).
//
// using namespace Opm::unit;
// std::transform(p.begin(), p.end(), p.begin(),
// boost::bind(convert::to, _1, psia));
//
/**
* Convert from internal units of measurements to equivalent
* external units of measurements. Note: The internal units of
* measurements are *ALWAYS*, and exclusively, SI.
*
* Example: Convert a <CODE>std::vector<double> p</CODE>, containing
* pressure values in the SI unit Pascal (i.e., unit::Pascal) to the
* equivalent values in Psi (unit::psia).
* \code
* using namespace Opm::unit;
* std::transform(p.begin(), p.end(), p.begin(),
* boost::bind(convert::to, _1, psia));
* \endcode
*
* @param[in] q Physical quantity, measured in SI units.
* @param[in] unit Physical unit of measurement.
* @return Value of @c q in unit <CODE>unit</CODE>.
*/
inline double to(const double q, const double unit)
{
return q / unit;