From 564b13d0c11dd48b7792ada3b5ea97fbe4605ca1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Thu, 3 Jan 2013 14:07:51 +0100 Subject: [PATCH] Added solveTofTracer() method. --- .../reorder/TransportModelTracerTof.cpp | 76 ++++++++++++++++++- .../reorder/TransportModelTracerTof.hpp | 19 +++++ 2 files changed, 94 insertions(+), 1 deletion(-) diff --git a/opm/core/transport/reorder/TransportModelTracerTof.cpp b/opm/core/transport/reorder/TransportModelTracerTof.cpp index 4e0f6c0c..27562844 100644 --- a/opm/core/transport/reorder/TransportModelTracerTof.cpp +++ b/opm/core/transport/reorder/TransportModelTracerTof.cpp @@ -32,7 +32,14 @@ namespace Opm /// \param[in] grid A 2d or 3d grid. TransportModelTracerTof::TransportModelTracerTof(const UnstructuredGrid& grid, const bool use_multidim_upwind) - : grid_(grid), use_multidim_upwind_(use_multidim_upwind) + : grid_(grid), + darcyflux_(0), + porevolume_(0), + source_(0), + tof_(0), + tracer_(0), + num_tracers_(0), + use_multidim_upwind_(use_multidim_upwind) { } @@ -68,6 +75,62 @@ namespace Opm face_tof_.resize(grid_.number_of_faces); std::fill(face_tof_.begin(), face_tof_.end(), 0.0); } + num_tracers_ = 0; + reorderAndTransport(grid_, darcyflux); + } + + + + + /// Solve for time-of-flight and a number of tracers. + /// One tracer will be used for each inflow flux specified in + /// the source parameter. + /// \param[in] darcyflux Array of signed face fluxes. + /// \param[in] porevolume Array of pore volumes. + /// \param[in] source Source term. Sign convention is: + /// (+) inflow flux, + /// (-) outflow flux. + /// \param[out] tof Array of time-of-flight values (1 per cell). + /// \param[out] tracer Array of tracer values (N per cell, where N is + /// the number of cells c for which source[c] > 0.0). + void TransportModelTracerTof::solveTofTracer(const double* darcyflux, + const double* porevolume, + const double* source, + std::vector& tof, + std::vector& tracer) + { + darcyflux_ = darcyflux; + porevolume_ = porevolume; + source_ = source; +#ifndef NDEBUG + // Sanity check for sources. + const double cum_src = std::accumulate(source, source + grid_.number_of_cells, 0.0); + if (std::fabs(cum_src) > *std::max_element(source, source + grid_.number_of_cells)*1e-2) { + THROW("Sources do not sum to zero: " << cum_src); + } +#endif + tof.resize(grid_.number_of_cells); + std::fill(tof.begin(), tof.end(), 0.0); + tof_ = &tof[0]; + // Find the tracer heads (injectors). + std::vector tracerheads; + for (int c = 0; c < grid_.number_of_cells; ++c) { + if (source[c] > 0.0) { + tracerheads.push_back(c); + } + } + num_tracers_ = tracerheads.size(); + tracer.resize(grid_.number_of_cells*num_tracers_); + std::fill(tracer.begin(), tracer.end(), 0.0); + for (int tr = 0; tr < num_tracers_; ++tr) { + tracer[tracerheads[tr]*num_tracers_ + tr] = 1.0; + } + tracer_ = &tracer[0]; + if (use_multidim_upwind_) { + face_tof_.resize(grid_.number_of_faces); + std::fill(face_tof_.begin(), face_tof_.end(), 0.0); + THROW("Multidimensional upwind not yet implemented for tracer."); + } reorderAndTransport(grid_, darcyflux); } @@ -107,6 +170,9 @@ namespace Opm // face. if (other != -1) { upwind_term += flux*tof_[other]; + for (int tr = 0; tr < num_tracers_; ++tr) { + tracer_[num_tracers_*cell + tr] += flux*tracer_[num_tracers_*other + tr]; + } } } else { downwind_flux += flux; @@ -115,6 +181,14 @@ namespace Opm // Compute tof. tof_[cell] = (porevolume_[cell] - upwind_term)/downwind_flux; + + // Compute tracers (if any). + // Do not change tracer solution in source cells. + if (source_[cell] <= 0.0) { + for (int tr = 0; tr < num_tracers_; ++tr) { + tracer_[num_tracers_*cell + tr] *= -1.0/downwind_flux; + } + } } diff --git a/opm/core/transport/reorder/TransportModelTracerTof.hpp b/opm/core/transport/reorder/TransportModelTracerTof.hpp index d9540abd..04983981 100644 --- a/opm/core/transport/reorder/TransportModelTracerTof.hpp +++ b/opm/core/transport/reorder/TransportModelTracerTof.hpp @@ -59,6 +59,23 @@ namespace Opm const double* source, std::vector& tof); + /// Solve for time-of-flight and a number of tracers. + /// One tracer will be used for each inflow flux specified in + /// the source parameter. + /// \param[in] darcyflux Array of signed face fluxes. + /// \param[in] porevolume Array of pore volumes. + /// \param[in] source Source term. Sign convention is: + /// (+) inflow flux, + /// (-) outflow flux. + /// \param[out] tof Array of time-of-flight values (1 per cell). + /// \param[out] tracer Array of tracer values (N per cell, where N is + /// the number of cells c for which source[c] > 0.0). + void solveTofTracer(const double* darcyflux, + const double* porevolume, + const double* source, + std::vector& tof, + std::vector& tracer); + private: virtual void solveSingleCell(const int cell); void solveSingleCellMultidimUpwind(const int cell); @@ -73,6 +90,8 @@ namespace Opm const double* porevolume_; // one volume per cell const double* source_; // one volumetric source term per cell double* tof_; + double* tracer_; + int num_tracers_; bool use_multidim_upwind_; std::vector face_tof_; // For multidim upwind face tofs. mutable std::vector adj_faces_; // For multidim upwind logic.