mirror of
https://github.com/OPM/opm-upscaling.git
synced 2025-02-25 18:45:23 -06:00
add a cpchop-variant for studying depth trends
This commit is contained in:
Håvard Berland
committed by
Havard Berland
Havard Berland
parent
8182c85b49
commit
7fbbec6b47
@@ -9,7 +9,8 @@ noinst_PROGRAMS = upscale_perm \
|
||||
upscale_cap \
|
||||
cpchop \
|
||||
exp_variogram \
|
||||
cpregularize
|
||||
cpregularize \
|
||||
cpchop_depthtrend
|
||||
|
||||
upscale_perm_SOURCES = upscale_perm.cpp
|
||||
|
||||
@@ -27,6 +28,8 @@ exp_variogram_SOURCES = exp_variogram.cpp
|
||||
|
||||
cpregularize_SOURCES = cpregularize.cpp
|
||||
|
||||
cpchop_depthtrend_SOURCES = cpchop_depthtrend.cpp
|
||||
|
||||
AM_CPPFLAGS += $(DUNEMPICPPFLAGS) $(BOOST_CPPFLAGS) $(SUPERLU_CPPFLAGS)
|
||||
AM_LDFLAGS += $(DUNEMPILDFLAGS) $(BOOST_LDFLAGS) $(SUPERLU_LDFLAGS)
|
||||
LDADD = $(DUNEMPILIBS) $(BOOST_UNIT_TEST_FRAMEWORK_LIB) \
|
||||
|
||||
211
examples/cpchop_depthtrend.cpp
Normal file
211
examples/cpchop_depthtrend.cpp
Normal file
@@ -0,0 +1,211 @@
|
||||
/*
|
||||
Copyright 2010 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 <dune/common/CornerpointChopper.hpp>
|
||||
#include <dune/upscaling/SinglePhaseUpscaler.hpp>
|
||||
#include <dune/porsol/common/setupBoundaryConditions.hpp>
|
||||
#include <dune/common/Units.hpp>
|
||||
|
||||
#include <boost/random/mersenne_twister.hpp>
|
||||
#include <boost/random/uniform_int.hpp>
|
||||
#include <boost/random/uniform_real.hpp>
|
||||
#include <boost/random/variate_generator.hpp>
|
||||
|
||||
#include <ios>
|
||||
#include <iomanip>
|
||||
#include <sys/utsname.h>
|
||||
#include <cmath> // for min()
|
||||
#include <ctime>
|
||||
#include <sstream>
|
||||
#include <fstream>
|
||||
#include <iostream>
|
||||
|
||||
/**
|
||||
This program is a variant of cpchop. Instead of subsampling randomly,
|
||||
it picks subsamples downwards in a model. It is specifically designed
|
||||
for extracting laterally extensive subsamples in specified intervals
|
||||
downwards through a model, in order to analyze depth trends in
|
||||
porosity in particular
|
||||
|
||||
Default is to pick a new subsample every meter (zresolution=1)
|
||||
The subsample height might be smaller or larger than zresolution.
|
||||
*/
|
||||
|
||||
|
||||
int main(int argc, char** argv)
|
||||
{
|
||||
Dune::parameter::ParameterGroup param(argc, argv);
|
||||
std::string gridfilename = param.get<std::string>("gridfilename");
|
||||
Dune::CornerPointChopper ch(gridfilename);
|
||||
|
||||
// The cells with i coordinate in [imin, imax) are included, similar for j.
|
||||
// The z limits may be changed inside the chopper to match actual min/max z.
|
||||
const int* dims = ch.dimensions();
|
||||
int imin = param.getDefault("imin", 0);
|
||||
int imax = param.getDefault("imax", dims[0]);
|
||||
int jmin = param.getDefault("jmin", 0);
|
||||
int jmax = param.getDefault("jmax", dims[1]);
|
||||
double zmin = param.getDefault("zmin", ch.zLimits().first);
|
||||
double zmax = param.getDefault("zmax", ch.zLimits().second);
|
||||
int ilen = param.getDefault("ilen", imax - imin);
|
||||
int jlen = param.getDefault("jlen", jmax - jmin);
|
||||
double zlen = param.getDefault("zlen", zmax - zmin);
|
||||
double zresolution = param.getDefault("zresolution", 1.0);
|
||||
bool upscale = param.getDefault("upscale", true);
|
||||
bool resettoorigin = param.getDefault("resettoorigin", true);
|
||||
boost::mt19937::result_type userseed = param.getDefault("seed", 0);
|
||||
|
||||
int outputprecision = param.getDefault("outputprecision", 8);
|
||||
std::string filebase = param.getDefault<std::string>("filebase", "");
|
||||
std::string resultfile = param.getDefault<std::string>("resultfile", "");
|
||||
|
||||
double minperm = param.getDefault("minperm", 1e-9);
|
||||
double minpermSI = Dune::unit::convert::from(minperm, Dune::prefix::milli*Dune::unit::darcy);
|
||||
|
||||
double z_tolerance = param.getDefault("z_tolerance", 0.0);
|
||||
double residual_tolerance = param.getDefault("residual_tolerance", 1e-8);
|
||||
double linsolver_verbosity = param.getDefault("linsolver_verbosity", 0);
|
||||
double linsolver_type = param.getDefault("linsolver_type", 1);
|
||||
|
||||
// Check that we do not have any user input
|
||||
// that goes outside the coordinates described in
|
||||
// the cornerpoint file (runtime-exception will be thrown in case of error)
|
||||
ch.verifyInscribedShoebox(imin, ilen, imax,
|
||||
jmin, jlen, jmax,
|
||||
zmin, zlen, zmax);
|
||||
|
||||
// Random number generator from boost.
|
||||
boost::mt19937 gen;
|
||||
|
||||
// Seed the random number generators with the current time, unless specified on command line
|
||||
// Warning: Current code does not allow 0 for the seed!!
|
||||
if (userseed == 0) {
|
||||
gen.seed(time(NULL));
|
||||
}
|
||||
else {
|
||||
gen.seed(userseed);
|
||||
}
|
||||
|
||||
|
||||
// Note that end is included in interval for uniform_int.
|
||||
boost::uniform_int<> disti(imin, imax - ilen);
|
||||
boost::uniform_int<> distj(jmin, jmax - jlen);
|
||||
boost::variate_generator<boost::mt19937&, boost::uniform_int<> > ri(gen, disti);
|
||||
boost::variate_generator<boost::mt19937&, boost::uniform_int<> > rj(gen, distj);
|
||||
|
||||
// Storage for results
|
||||
std::vector<double> zstarts;
|
||||
std::vector<double> porosities;
|
||||
std::vector<double> permxs;
|
||||
std::vector<double> permys;
|
||||
std::vector<double> permzs;
|
||||
|
||||
|
||||
/* z_start is the topmost point of the subsample to extract */
|
||||
for (double zstart = 0.0; zstart <= zmax-zlen; zstart += zresolution) {
|
||||
/* Horizontally, we pick by random, even though default behaviour is
|
||||
to have ilen=imax-min so that there is no randomness */
|
||||
int istart = ri();
|
||||
int jstart = rj();
|
||||
ch.chop(istart, istart + ilen, jstart, jstart + jlen, zstart, std::min(zstart + zlen,zmax), resettoorigin);
|
||||
std::string subsampledgrdecl = filebase;
|
||||
|
||||
// Output grdecl-data to file if a filebase is supplied.
|
||||
if (filebase != "") {
|
||||
std::ostringstream oss;
|
||||
oss << 'Z' << std::setw(4) << std::setfill('0') << zstart;
|
||||
subsampledgrdecl += oss.str();
|
||||
subsampledgrdecl += ".grdecl";
|
||||
ch.writeGrdecl(subsampledgrdecl);
|
||||
}
|
||||
|
||||
try { /* The upscaling may fail to converge on icky grids, lets just pass by those */
|
||||
if (upscale) {
|
||||
Dune::EclipseGridParser subparser = ch.subparser();
|
||||
subparser.convertToSI();
|
||||
Dune::SinglePhaseUpscaler upscaler;
|
||||
upscaler.init(subparser, Dune::SinglePhaseUpscaler::Fixed, minpermSI, z_tolerance,
|
||||
residual_tolerance, linsolver_verbosity, linsolver_type, false);
|
||||
|
||||
Dune::SinglePhaseUpscaler::permtensor_t upscaled_K = upscaler.upscaleSinglePhase();
|
||||
upscaled_K *= (1.0/(Dune::prefix::milli*Dune::unit::darcy));
|
||||
|
||||
|
||||
zstarts.push_back(zstart);
|
||||
porosities.push_back(upscaler.upscalePorosity());
|
||||
permxs.push_back(upscaled_K(0,0));
|
||||
permys.push_back(upscaled_K(1,1));
|
||||
permzs.push_back(upscaled_K(2,2));
|
||||
|
||||
}
|
||||
}
|
||||
catch (...) {
|
||||
std::cerr << "Warning: Upscaling chopped subsample at z=" << zstart << "failed, proceeding to next subsample\n";
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
if (upscale) {
|
||||
|
||||
// Make stream of output data, to be outputted to screen and optionally to file
|
||||
std::stringstream outputtmp;
|
||||
|
||||
outputtmp << "################################################################################################" << std::endl;
|
||||
outputtmp << "# Results from depth thrend analysis on subsamples" << std::endl;
|
||||
outputtmp << "#" << std::endl;
|
||||
time_t now = time(NULL);
|
||||
outputtmp << "# Finished: " << asctime(localtime(&now));
|
||||
|
||||
utsname hostname; uname(&hostname);
|
||||
outputtmp << "# Hostname: " << hostname.nodename << std::endl;
|
||||
outputtmp << "#" << std::endl;
|
||||
outputtmp << "# Options used:" << std::endl;
|
||||
outputtmp << "# gridfilename: " << gridfilename << std::endl;
|
||||
outputtmp << "# i; min,len,max: " << imin << " " << ilen << " " << imax << std::endl;
|
||||
outputtmp << "# j; min,len,max: " << jmin << " " << jlen << " " << jmax << std::endl;
|
||||
outputtmp << "# z; min,len,max: " << zmin << " " << zlen << " " << zmax << std::endl;
|
||||
outputtmp << "# zresolution: " << zresolution << std::endl;
|
||||
outputtmp << "################################################################################################" << std::endl;
|
||||
outputtmp << "# zstart porosity permx permy permz" << std::endl;
|
||||
|
||||
const int fieldwidth = outputprecision + 8;
|
||||
for (size_t sample = 1; sample <= porosities.size(); ++sample) {
|
||||
outputtmp <<
|
||||
std::showpoint << std::setw(fieldwidth) << std::setprecision(outputprecision) << zstarts[sample-1] << '\t' <<
|
||||
std::showpoint << std::setw(fieldwidth) << std::setprecision(outputprecision) << porosities[sample-1] << '\t' <<
|
||||
std::showpoint << std::setw(fieldwidth) << std::setprecision(outputprecision) << permxs[sample-1] << '\t' <<
|
||||
std::showpoint << std::setw(fieldwidth) << std::setprecision(outputprecision) << permys[sample-1] << '\t' <<
|
||||
std::showpoint << std::setw(fieldwidth) << std::setprecision(outputprecision) << permzs[sample-1] << '\t' <<
|
||||
std::endl;
|
||||
}
|
||||
|
||||
if (resultfile != "") {
|
||||
std::cout << "Writing results to " << resultfile << std::endl;
|
||||
std::ofstream outfile;
|
||||
outfile.open(resultfile.c_str(), std::ios::out | std::ios::trunc);
|
||||
outfile << outputtmp.str();
|
||||
outfile.close();
|
||||
}
|
||||
|
||||
|
||||
|
||||
std::cout << outputtmp.str();
|
||||
}
|
||||
}
|
||||
Reference in New Issue
Block a user