Merge remote-tracking branch 'upstream/opm-parser-integrate' into wellsmanager-rates

Conflicts:
	opm/core/wells/WellsManager.cpp
This commit is contained in:
Joakim Hove 2014-01-30 09:50:09 +01:00
commit 3e2df33a65
9 changed files with 265 additions and 163 deletions

View File

@ -37,7 +37,7 @@ enum WellControlType {
struct WellControls;
bool
well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2);
well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2 , bool verbose);
struct WellControls *
well_controls_create(void);

View File

@ -263,7 +263,7 @@ struct Wells *
clone_wells(const struct Wells *W);
bool
wells_equal(const struct Wells *W1, const struct Wells *W2);
wells_equal(const struct Wells *W1, const struct Wells *W2 , bool verbose);

View File

@ -274,6 +274,24 @@ namespace Opm
return;
}
// global_cell is a map from compressed cells to Cartesian grid cells.
// We must make the inverse lookup.
const int* global_cell = grid.global_cell;
const int* cpgdim = grid.cartdims;
std::map<int,int> cartesian_to_compressed;
if (global_cell) {
for (int i = 0; i < grid.number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(global_cell[i], i));
}
}
else {
for (int i = 0; i < grid.number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(i, i));
}
}
// Obtain phase usage data.
PhaseUsage pu = phaseUsageFromDeck(eclipseState);
@ -294,129 +312,68 @@ namespace Opm
std::vector<WellConstPtr> wells = schedule->getWells(timeStep);
well_names.reserve(wells.size());
well_data.reserve(wells.size());
wellperf_data.resize(wells.size());
int well_index = 0;
for (auto wellIter= wells.begin(); wellIter != wells.end(); ++wellIter) {
WellConstPtr well = (*wellIter);
well_names_to_index[well->name()] = well_index;
well_names.push_back(well->name());
{
WellData wd;
// If negative (defaulted), set refdepth to a marker
// value, will be changed after getting perforation
// data to the centroid of the cell of the top well
// perforation.
wd.reference_bhp_depth = (well->getRefDepth() < 0.0) ? -1e100 : well->getRefDepth();
wd.welspecsline = -1;
{ // WELSPECS handling
well_names_to_index[well->name()] = well_index;
well_names.push_back(well->name());
{
WellData wd;
// If negative (defaulted), set refdepth to a marker
// value, will be changed after getting perforation
// data to the centroid of the cell of the top well
// perforation.
wd.reference_bhp_depth = (well->getRefDepth() < 0.0) ? -1e100 : well->getRefDepth();
wd.welspecsline = -1;
if (well->isInjector( timeStep ))
wd.type = INJECTOR;
else
wd.type = PRODUCER;
well_data.push_back(wd);
}
}
if (well->isProducer(timeStep))
wd.type = PRODUCER;
else
wd.type = INJECTOR;
well_data.push_back(wd);
{ // COMPDAT handling
CompletionSetConstPtr completionSet = well->getCompletions(timeStep);
for (size_t c=0; c<completionSet->size(); c++) {
CompletionConstPtr completion = completionSet->get(c);
int i = completion->getI();
int j = completion->getJ();
int k = completion->getK();
int cart_grid_indx = i + cpgdim[0]*(j + cpgdim[1]*k);
std::map<int, int>::const_iterator cgit = cartesian_to_compressed.find(cart_grid_indx);
if (cgit == cartesian_to_compressed.end()) {
OPM_THROW(std::runtime_error, "Cell with i,j,k indices " << i << ' ' << j << ' '
<< k << " not found in grid (well = " << well->name() << ')');
}
int cell = cgit->second;
PerfData pd;
pd.cell = cell;
if (completion->getCF() > 0.0) {
pd.well_index = completion->getCF();
} else {
double radius = 0.5*completion->getDiameter();
if (radius <= 0.0) {
radius = 0.5*unit::feet;
OPM_MESSAGE("**** Warning: Well bore internal radius set to " << radius);
}
std::array<double, 3> cubical = getCubeDim(grid, cell);
const double* cell_perm = &permeability[grid.dimensions*grid.dimensions*cell];
pd.well_index = computeWellIndex(radius, cubical, cell_perm, completion->getDiameter());
}
wellperf_data[well_index].push_back(pd);
}
}
well_index++;
}
const int num_wells = well_data.size();
wellperf_data.resize(num_wells);
// global_cell is a map from compressed cells to Cartesian grid cells.
// We must make the inverse lookup.
const int* global_cell = grid.global_cell;
const int* cpgdim = grid.cartdims;
std::map<int,int> cartesian_to_compressed;
if (global_cell) {
for (int i = 0; i < grid.number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(global_cell[i], i));
}
}
else {
for (int i = 0; i < grid.number_of_cells; ++i) {
cartesian_to_compressed.insert(std::make_pair(i, i));
}
}
// Get COMPDAT data
// It is *not* allowed to have multiple lines corresponding to
// the same perforation!
const COMPDAT& compdat = deck.getCOMPDAT();
const int num_compdat = compdat.compdat.size();
for (int kw = 0; kw < num_compdat; ++kw) {
// Extract well name, or the part of the well name that
// comes before the '*'.
std::string name = compdat.compdat[kw].well_;
std::string::size_type len = name.find('*');
if (len != std::string::npos) {
name = name.substr(0, len);
}
// Look for well with matching name.
bool found = false;
for (int wix = 0; wix < num_wells; ++wix) {
if (well_names[wix].compare(0,len, name) == 0) { // equal
// We have a matching name.
int ix = compdat.compdat[kw].grid_ind_[0] - 1;
int jy = compdat.compdat[kw].grid_ind_[1] - 1;
int kz1 = compdat.compdat[kw].grid_ind_[2] - 1;
int kz2 = compdat.compdat[kw].grid_ind_[3] - 1;
WellConstPtr well = schedule->getWell(well_names[wix]);
if (ix < 0) {
// Defaulted I location. Extract from WELSPECS.
ix = well->getHeadI() - 1;
}
if (jy < 0) {
// Defaulted J location. Extract from WELSPECS.
jy = well->getHeadJ() - 1;
}
if (kz1 < 0) {
// Defaulted KZ1. Use top layer.
kz1 = 0;
}
if (kz2 < 0) {
// Defaulted KZ2. Use bottom layer.
kz2 = cpgdim[2] - 1;
}
for (int kz = kz1; kz <= kz2; ++kz) {
int cart_grid_indx = ix + cpgdim[0]*(jy + cpgdim[1]*kz);
std::map<int, int>::const_iterator cgit =
cartesian_to_compressed.find(cart_grid_indx);
if (cgit == cartesian_to_compressed.end()) {
OPM_THROW(std::runtime_error, "Cell with i,j,k indices " << ix << ' ' << jy << ' '
<< kz << " not found in grid (well = " << name << ')');
}
int cell = cgit->second;
PerfData pd;
pd.cell = cell;
if (compdat.compdat[kw].connect_trans_fac_ > 0.0) {
pd.well_index = compdat.compdat[kw].connect_trans_fac_;
} else {
double radius = 0.5*compdat.compdat[kw].diameter_;
if (radius <= 0.0) {
radius = 0.5*unit::feet;
OPM_MESSAGE("**** Warning: Well bore internal radius set to " << radius);
}
std::array<double, 3> cubical = getCubeDim(grid, cell);
const double* cell_perm = &permeability[grid.dimensions*grid.dimensions*cell];
pd.well_index = computeWellIndex(radius, cubical, cell_perm,
compdat.compdat[kw].skin_factor_);
}
wellperf_data[wix].push_back(pd);
}
found = true;
break;
}
}
if (!found) {
OPM_THROW(std::runtime_error, "Undefined well name: " << compdat.compdat[kw].well_
<< " in COMPDAT");
}
}
// Set up reference depths that were defaulted. Count perfs.
const int num_wells = well_data.size();
int num_perfs = 0;
assert(grid.dimensions == 3);
for (int w = 0; w < num_wells; ++w) {
@ -563,8 +520,8 @@ namespace Opm
}
std::copy(cf, cf + pu.num_phases, w_->comp_frac + well_index*pu.num_phases);
well_index++;
}
well_index++;
}
// Get WCONPROD data
@ -587,6 +544,7 @@ namespace Opm
well_found = true;
assert(well_data[wix].type == w_->type[wix]);
if (well_data[wix].type != PRODUCER) {
std::cout << "Looking at well: " << well_names[wix] << std::endl;
OPM_THROW(std::runtime_error, "Found WCONPROD entry for a non-producer well: " << well_names[wix]);
}
// Add all controls that are present in well.

View File

@ -124,7 +124,7 @@ namespace Opm
// Data
Wells* w_;
WellCollection well_collection_;
WellCollection well_collection_;

View File

@ -293,19 +293,38 @@ well_controls_add_new(enum WellControlType type , double target , const double *
bool
well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2)
well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2 , bool verbose)
/* ---------------------------------------------------------------------- */
{
bool are_equal = true;
are_equal = (ctrls1->num == ctrls2->num);
are_equal = are_equal && (ctrls1->number_of_phases == ctrls2->number_of_phases);
if (ctrls1->num != ctrls2->num) {
are_equal = false;
if (verbose)
printf("ctrls1->num:%d ctrls2->num:%d \n",ctrls1->num , ctrls2->num);
}
if (ctrls1->number_of_phases != ctrls2->number_of_phases) {
are_equal = false;
if (verbose)
printf("ctrls1->number_of_phases:%d ctrls2->number_of_phases:%d \n",ctrls1->number_of_phases , ctrls2->number_of_phases);
}
if (!are_equal) {
return are_equal;
}
are_equal = are_equal && (memcmp(ctrls1->type, ctrls2->type, ctrls1->num * sizeof *ctrls1->type ) == 0);
are_equal = are_equal && (memcmp(ctrls1->target, ctrls2->target, ctrls1->num * sizeof *ctrls1->target ) == 0);
are_equal = are_equal && (memcmp(ctrls1->distr, ctrls2->distr, ctrls1->num * ctrls1->number_of_phases * sizeof *ctrls1->distr ) == 0);
if (memcmp(ctrls1->type, ctrls2->type, ctrls1->num * sizeof *ctrls1->type ) != 0) {
are_equal = false;
if (verbose)
printf("The ->type vectors are different \n");
}
if (memcmp(ctrls1->target, ctrls2->target, ctrls1->num * sizeof *ctrls1->target ) != 0) {
are_equal = false;
if (verbose)
printf("The ->target vectors are different \n");
}
are_equal = are_equal && (ctrls1->cpty == ctrls2->cpty);
return are_equal;

View File

@ -541,7 +541,7 @@ clone_wells(const struct Wells *W)
/* ---------------------------------------------------------------------- */
bool
wells_equal(const struct Wells *W1, const struct Wells *W2)
wells_equal(const struct Wells *W1, const struct Wells *W2 , bool verbose)
/* ---------------------------------------------------------------------- */
{
bool are_equal = true;
@ -567,10 +567,25 @@ wells_equal(const struct Wells *W1, const struct Wells *W2)
are_equal = are_equal && (strcmp(W1->name[i], W2->name[i]) == 0);
else
are_equal = are_equal && (W1->name[i] == W2->name[i]);
if (verbose && !are_equal)
printf("Well name[%d] %s and %s are different \n", i , W1->name[i] , W2->name[i]);
}
if (W1->type[i] != W2->type[i]) {
are_equal = false;
if (verbose)
printf("Well->type[%d] different %d %d \n",i , W1->type[i] , W2->type[i] );
}
if (W1->depth_ref[i] != W2->depth_ref[i]) {
are_equal = false;
if (verbose)
printf("Well->depth_ref[%d] different %g %g \n",i , W1->depth_ref[i] , W2->depth_ref[i] );
}
if (!well_controls_equal(W1->ctrls[i], W2->ctrls[i],verbose)) {
are_equal = false;
if (verbose)
printf("Well controls are different for well[%d]:%s \n",i,W1->name[i]);
}
are_equal = are_equal && (W1->type[i] == W2->type[i]);
are_equal = are_equal && (W1->depth_ref[i] == W2->depth_ref[i]);
are_equal = are_equal && (well_controls_equal(W1->ctrls[i], W2->ctrls[i]));
}
@ -581,10 +596,16 @@ wells_equal(const struct Wells *W1, const struct Wells *W2)
are_equal = are_equal && (mgmt1->well_cpty == mgmt2->well_cpty);
}
are_equal = are_equal && (memcmp(W1->comp_frac, W2->comp_frac, W1->number_of_wells * W1->number_of_phases * sizeof *W1->comp_frac ) == 0);
are_equal = are_equal && (memcmp(W1->well_connpos, W2->well_connpos, (1 + W1->number_of_wells) * sizeof *W1->well_connpos ) == 0);
if (!are_equal) {
return are_equal;
if (memcmp(W1->comp_frac, W2->comp_frac, W1->number_of_wells * W1->number_of_phases * sizeof *W1->comp_frac ) != 0) {
are_equal = false;
if (verbose)
printf("Component fractions different \n");
}
if (memcmp(W1->well_connpos, W2->well_connpos, (1 + W1->number_of_wells) * sizeof *W1->well_connpos ) != 0) {
are_equal = false;
if (verbose)
printf("perforation position map difference \n");
}
{

View File

@ -204,7 +204,7 @@ BOOST_AUTO_TEST_CASE(Copy)
BOOST_CHECK_EQUAL( dist1[p] , dist2[p]);
}
}
BOOST_CHECK( well_controls_equal( c1 , c2 ));
BOOST_CHECK( well_controls_equal( c1 , c2 , false) );
}
}
}
@ -219,7 +219,7 @@ BOOST_AUTO_TEST_CASE(Equals_WellsEqual_ReturnsTrue) {
std::shared_ptr<Wells> W2(create_wells(nphases, nwells, nperfs),
destroy_wells);
BOOST_CHECK(wells_equal(W1.get(), W2.get()));
BOOST_CHECK(wells_equal(W1.get(), W2.get() , false));
}
@ -232,5 +232,5 @@ BOOST_AUTO_TEST_CASE(Equals_WellsDiffer_ReturnsFalse) {
std::shared_ptr<Wells> W2(create_wells(nphases, 3, nperfs),
destroy_wells);
BOOST_CHECK(!wells_equal(W1.get(), W2.get()));
BOOST_CHECK(!wells_equal(W1.get(), W2.get() , false ));
}

View File

@ -174,8 +174,6 @@ void check_controls_epoch1( struct WellControls ** ctrls) {
}
BOOST_AUTO_TEST_CASE(Constructor_Works) {
Opm::EclipseGridParser Deck("wells_manager_data.data");
Opm::GridManager gridManager(Deck);
@ -212,11 +210,15 @@ BOOST_AUTO_TEST_CASE(New_Constructor_Works) {
Opm::WellsManager wellsManager(eclipseState, 0, Deck, *gridManager.c_grid(), NULL);
Opm::WellsManager oldWellsManager(Deck, *gridManager.c_grid(), NULL);
const Wells* wells = wellsManager.c_wells();
wells_static_check( wells );
check_controls_epoch0( wells->ctrls );
std::cout << "Checking new well structure, epoch 0" << std::endl;
wells_static_check( wellsManager.c_wells() );
BOOST_CHECK(wells_equal(wells, oldWellsManager.c_wells()));
std::cout << "Checking old well structure, epoch 0" << std::endl;
wells_static_check( oldWellsManager.c_wells() );
check_controls_epoch0( wellsManager.c_wells()->ctrls );
BOOST_CHECK(wells_equal(wellsManager.c_wells(), oldWellsManager.c_wells() , false));
}
Deck.setCurrentEpoch(1);
@ -224,11 +226,41 @@ BOOST_AUTO_TEST_CASE(New_Constructor_Works) {
Opm::WellsManager wellsManager(eclipseState, 1,Deck, *gridManager.c_grid(), NULL);
Opm::WellsManager oldWellsManager(Deck, *gridManager.c_grid(), NULL);
const Wells* wells = wellsManager.c_wells();
wells_static_check( wells );
check_controls_epoch1( wells->ctrls );
std::cout << "Checking new well structure, epoch 1" << std::endl;
wells_static_check( wellsManager.c_wells() );
BOOST_CHECK(wells_equal(wells, oldWellsManager.c_wells()));
std::cout << "Checking old well structure, epoch 1" << std::endl;
wells_static_check( oldWellsManager.c_wells() );
check_controls_epoch1( wellsManager.c_wells()->ctrls );
BOOST_CHECK(wells_equal( wellsManager.c_wells(), oldWellsManager.c_wells(),false));
}
}
BOOST_AUTO_TEST_CASE(New_Constructor_Works_ExpandedData) {
Opm::ParserPtr parser(new Opm::Parser());
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(parser->parseFile("wells_manager_data_expanded.data")));
Opm::EclipseGridParser Deck("wells_manager_data_expanded.data");
Opm::GridManager gridManager(Deck);
Deck.setCurrentEpoch(0);
{
Opm::WellsManager wellsManager(eclipseState, 0, Deck, *gridManager.c_grid(), NULL);
Opm::WellsManager oldWellsManager(Deck, *gridManager.c_grid(), NULL);
BOOST_CHECK(wells_equal(wellsManager.c_wells(), oldWellsManager.c_wells(),false));
}
Deck.setCurrentEpoch(1);
{
Opm::WellsManager wellsManager(eclipseState, 1,Deck, *gridManager.c_grid(), NULL);
Opm::WellsManager oldWellsManager(Deck, *gridManager.c_grid(), NULL);
BOOST_CHECK(wells_equal( wellsManager.c_wells(), oldWellsManager.c_wells(), true));
}
}
@ -244,8 +276,8 @@ BOOST_AUTO_TEST_CASE(WellsEqual) {
Deck.setCurrentEpoch(1);
Opm::WellsManager wellsManager1(Deck, *gridManager.c_grid(), NULL);
BOOST_CHECK( wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells()) );
BOOST_CHECK( !wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells()) );
BOOST_CHECK( wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells(),false) );
BOOST_CHECK( !wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells(),false) );
}
@ -259,15 +291,15 @@ BOOST_AUTO_TEST_CASE(ControlsEqual) {
Deck.setCurrentEpoch(1);
Opm::WellsManager wellsManager1(Deck, *gridManager.c_grid(), NULL);
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0]));
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1]));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager1.c_wells()->ctrls[0]));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager1.c_wells()->ctrls[1]));
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager1.c_wells()->ctrls[0] , false));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager1.c_wells()->ctrls[1] , false));
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[1]));
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[0]));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0]));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1]));
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
}

View File

@ -0,0 +1,72 @@
OIL
GAS
WATER
DIMENS
10 10 5 /
GRID
DXV
10*1000.0 /
DYV
10*1000.0 /
DZV
10.0 20.0 30.0 10.0 5.0 /
SCHEDULE
WELSPECS
'INJ1' 'G' 1 1 8335 'GAS' /
'PROD1' 'G' 10 10 8400 'OIL' /
/
COMPDAT
'INJ1' 1 1 1 2 'OPEN' 1 10.6092 0.5 /
'INJ1' 1 1 3 5 'OPEN' 1 12.6092 0.5 /
'INJ1' 2 2 1 1 'OPEN' 1 14.6092 0.5 /
'PROD1' 10 3 3 3 'OPEN' 0 10.6092 0.5 /
/
WCONPROD
'PROD1' 'OPEN' 'ORAT' 20000 4* 1000 /
/
WCONINJE
'INJ1' 'GAS' 'OPEN' 'RATE' 100 200 400 /
/
DATES
1 'FEB' 2000 /
/
WELSPECS
'INJ1' 'G3' 1 1 8335 'GAS' /
'QNJ2' 'G3' 1 1 8335 'GAS' /
/
COMPDAT
'QNJ2' 3 4 1 2 'OPEN' 1 10.6092 0.5 /
'QNJ2' 4 4 3 5 'OPEN' 1 12.6092 0.5 /
/
WCONPROD
'PROD1' 'OPEN' 'RESV' 999 3* 123 100 /
/
WCONINJE
'INJ1' 'WATER' 'OPEN' 'RESV' 10 20 40 /
'QNJ2' 'WATER' 'OPEN' 'RESV' 7 33 39 /
/
TSTEP
14.0
/
END