From 1fbe8e3b6d855dd4cd23f8d1f249f28e6d0d0269 Mon Sep 17 00:00:00 2001 From: Markus Blatt Date: Thu, 7 May 2015 11:21:09 +0200 Subject: [PATCH] Distribute the geology information. As it turns out initializing the Geology on a distributed grid result in wrong values for e.g. saturation. Therefore with this commit we resort to initializing the global geology and distribute it using communication. --- examples/flow_cp.cpp | 22 +++++-- opm/autodiff/GeoProps.hpp | 2 + opm/autodiff/RedistributeDataHandles.hpp | 76 +++++++++++++++++++++++- 3 files changed, 93 insertions(+), 7 deletions(-) diff --git a/examples/flow_cp.cpp b/examples/flow_cp.cpp index 6e8c1714b..c3937f59c 100644 --- a/examples/flow_cp.cpp +++ b/examples/flow_cp.cpp @@ -263,6 +263,13 @@ try new_props->setSwatInitScaling(state.saturation(),pc); } + bool use_gravity = (gravity[0] != 0.0 || gravity[1] != 0.0 || gravity[2] != 0.0); + const double *grav = use_gravity ? &gravity[0] : 0; + + std::shared_ptr geoprops; + std::shared_ptr distributed_geology; + geoprops.reset(new Opm::DerivedGeology(*grid, *new_props, eclipseState, false, grav)); + distributed_geology=geoprops; // copy for serial case. BlackoilState distributed_state; std::shared_ptr distributed_props = new_props; Dune::CpGrid distributed_grid = *grid; @@ -295,11 +302,16 @@ try *distributed_props); grid->scatterData(state_handle); grid->scatterData(props_handle); + // Create a distributed Geology. Some values will be updated using communication + // below + distributed_geology.reset(new Opm::DerivedGeology(distributed_grid, + *distributed_props, eclipseState, + false, grav)); + Opm::GeologyDataHandle geo_handle(global_grid, distributed_grid, + *geoprops, *distributed_geology); + grid->scatterData(geo_handle); } - bool use_gravity = (gravity[0] != 0.0 || gravity[1] != 0.0 || gravity[2] != 0.0); - const double *grav = use_gravity ? &gravity[0] : 0; - // Solver for Newton iterations. std::unique_ptr fis_solver; @@ -318,13 +330,11 @@ try // initialize variables simtimer.init(timeMap); - Opm::DerivedGeology geology(distributed_grid, *distributed_props, eclipseState, false, grav); - std::vector threshold_pressures = thresholdPressures(eclipseState, distributed_grid); SimulatorFullyImplicitBlackoil simulator(param, distributed_grid, - geology, + *distributed_geology, *distributed_props, rock_comp->isActive() ? rock_comp.get() : 0, *fis_solver, diff --git a/opm/autodiff/GeoProps.hpp b/opm/autodiff/GeoProps.hpp index 358befab5..54792137e 100644 --- a/opm/autodiff/GeoProps.hpp +++ b/opm/autodiff/GeoProps.hpp @@ -178,6 +178,8 @@ namespace Opm const Vector& gravityPotential() const { return gpot_ ;} const Vector& z() const { return z_ ;} const double* gravity() const { return gravity_;} + Vector& poreVolume() { return pvol_ ;} + Vector& transmissibility() { return trans_ ;} private: template diff --git a/opm/autodiff/RedistributeDataHandles.hpp b/opm/autodiff/RedistributeDataHandles.hpp index fcba01bc9..29617014e 100644 --- a/opm/autodiff/RedistributeDataHandles.hpp +++ b/opm/autodiff/RedistributeDataHandles.hpp @@ -24,6 +24,80 @@ namespace Opm { + +/// \brief a data handle to distribute Derived Geology +class GeologyDataHandle +{ +public: + /// \brief type of the data we send + typedef double DataType; + /// \brief Constructor + /// \param sendGrid The grid that the data is attached to when sending. + /// \param recvGrid The grid that the data is attached to when receiving. + /// \param sendGeology The state where we will retieve the values to be sent. + /// \param recvGeology The state where we will store the received values. + GeologyDataHandle(const Dune::CpGrid& sendGrid, + const Dune::CpGrid& recvGrid, + const DerivedGeology& sendGeology, + DerivedGeology& recvGeology) + : sendGrid_(sendGrid), recvGrid_(recvGrid), sendGeology_(sendGeology), + recvGeology_(recvGeology) + {} + + bool fixedsize(int /*dim*/, int /*codim*/) + { + return false; + } + template + std::size_t size(const T& e) + { + if ( T::codimension == 0) + { + return 1 + sendGrid_.numCellFaces(e.index()); + } + else + { + OPM_THROW(std::logic_error, "Data handle can only be used for elements"); + } + } + template + void gather(B& buffer, const T& e) + { + assert( T::codimension == 0); + buffer.write(sendGeology_.poreVolume()[e.index()]); + for ( int i=0; i< sendGrid_.numCellFaces(e.index()); ++i ) + { + buffer.write(sendGeology_.transmissibility()[sendGrid_.cellFace(e.index(), i)]); + } + } + template + void scatter(B& buffer, const T& e, std::size_t size) + { + assert( T::codimension == 0); + double val; + buffer.read(val); + recvGeology_.poreVolume()[e.index()]=val; + for ( int i=0; i< recvGrid_.numCellFaces(e.index()); ++i ) + { + buffer.read(val); + recvGeology_.transmissibility()[recvGrid_.cellFace(e.index(), i)]=val; + } + } + bool contains(int dim, int codim) + { + return dim==3 && codim==0; + } +private: + /// \brief The grid that the data we send is associated with. + const Dune::CpGrid& sendGrid_; + /// \brief The grid that the data we receive is associated with. + const Dune::CpGrid& recvGrid_; + /// \brief The data to send. + const DerivedGeology& sendGeology_; + /// \brief The data to receive. + DerivedGeology& recvGeology_; +}; + /// \brief a data handle to distribute the BlackoilState class BlackoilStateDataHandle { @@ -34,7 +108,7 @@ public: /// \param sendGrid The grid that the data is attached to when sending. /// \param recvGrid The grid that the data is attached to when receiving. /// \param sendState The state where we will retieve the values to be sent. - /// \parame recvState The state where we will store the received values. + /// \param recvState The state where we will store the received values. BlackoilStateDataHandle(const Dune::CpGrid& sendGrid, const Dune::CpGrid& recvGrid, const BlackoilState& sendState,