225doInit(
bool rst, std::size_t numGridDof,
226 std::size_t gasPhaseIdx, std::size_t oilPhaseIdx, std::size_t waterPhaseIdx)
228 const auto& tracers = eclState_.tracer();
230 if (tracers.size() == 0) {
235 const std::size_t numTracers = tracers.size();
236 enableSolTracers_.resize(numTracers);
237 tracerConcentration_.resize(numTracers);
238 freeTracerConcentration_.resize(numTracers);
239 solTracerConcentration_.resize(numTracers);
242 tracerPhaseIdx_.resize(numTracers);
243 for (std::size_t tracerIdx = 0; tracerIdx < numTracers; tracerIdx++) {
244 const auto& tracer = tracers[tracerIdx];
246 if (tracer.phase == Phase::WATER)
247 tracerPhaseIdx_[tracerIdx] = waterPhaseIdx;
248 else if (tracer.phase == Phase::OIL)
249 tracerPhaseIdx_[tracerIdx] = oilPhaseIdx;
250 else if (tracer.phase == Phase::GAS)
251 tracerPhaseIdx_[tracerIdx] = gasPhaseIdx;
253 tracerConcentration_[tracerIdx].resize(numGridDof);
254 freeTracerConcentration_[tracerIdx].resize(numGridDof);
255 solTracerConcentration_[tracerIdx].resize(numGridDof);
262 if (tracer.free_concentration.has_value()){
263 const auto& free_concentration = tracer.free_concentration.value();
264 int tblkDatasize = free_concentration.size();
265 if (tblkDatasize < cartMapper_.cartesianSize()){
266 throw std::runtime_error(
"Wrong size of TBLKF for" + tracer.name);
268 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
269 int cartDofIdx = cartMapper_.cartesianIndex(globalDofIdx);
270 tracerConcentration_[tracerIdx][globalDofIdx][0] = free_concentration[cartDofIdx];
271 freeTracerConcentration_[tracerIdx][globalDofIdx] = free_concentration[cartDofIdx];
275 else if (tracer.free_tvdp.has_value()) {
276 const auto& free_tvdp = tracer.free_tvdp.value();
277 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
278 tracerConcentration_[tracerIdx][globalDofIdx][0] =
279 free_tvdp.evaluate(
"TRACER_CONCENTRATION",
280 centroids_(globalDofIdx)[2]);
281 freeTracerConcentration_[tracerIdx][globalDofIdx] =
282 free_tvdp.evaluate(
"TRACER_CONCENTRATION",
283 centroids_(globalDofIdx)[2]);
287 OpmLog::warning(fmt::format(
"No TBLKF or TVDPF given for free tracer {}. "
288 "Initial values set to zero. ", tracer.name));
289 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
290 tracerConcentration_[tracerIdx][globalDofIdx][0] = 0.0;
291 freeTracerConcentration_[tracerIdx][globalDofIdx] = 0.0;
296 if (tracer.phase != Phase::WATER &&
297 ((tracer.phase == Phase::GAS && FluidSystem::enableDissolvedGas()) ||
298 (tracer.phase == Phase::OIL && FluidSystem::enableVaporizedOil()))) {
300 if (tracer.solution_concentration.has_value()){
301 enableSolTracers_[tracerIdx] =
true;
302 const auto& solution_concentration = tracer.solution_concentration.value();
303 int tblkDatasize = solution_concentration.size();
304 if (tblkDatasize < cartMapper_.cartesianSize()){
305 throw std::runtime_error(
"Wrong size of TBLKS for" + tracer.name);
307 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
308 int cartDofIdx = cartMapper_.cartesianIndex(globalDofIdx);
309 tracerConcentration_[tracerIdx][globalDofIdx][1] = solution_concentration[cartDofIdx];
310 solTracerConcentration_[tracerIdx][globalDofIdx] = solution_concentration[cartDofIdx];
314 else if (tracer.solution_tvdp.has_value()) {
315 enableSolTracers_[tracerIdx] =
true;
316 const auto& solution_tvdp = tracer.solution_tvdp.value();
317 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
318 tracerConcentration_[tracerIdx][globalDofIdx][1] =
319 solution_tvdp.evaluate(
"TRACER_CONCENTRATION",
320 centroids_(globalDofIdx)[2]);
321 solTracerConcentration_[tracerIdx][globalDofIdx] =
322 solution_tvdp.evaluate(
"TRACER_CONCENTRATION",
323 centroids_(globalDofIdx)[2]);
328 enableSolTracers_[tracerIdx] =
false;
329 OpmLog::warning(fmt::format(
"No TBLKS or TVDPS given for solution tracer {}. "
330 "Initial values set to zero. ", tracer.name));
331 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
332 tracerConcentration_[tracerIdx][globalDofIdx][1] = 0.0;
333 solTracerConcentration_[tracerIdx][globalDofIdx] = 0.0;
339 enableSolTracers_[tracerIdx] =
false;
340 for (std::size_t globalDofIdx = 0; globalDofIdx < numGridDof; ++globalDofIdx) {
341 tracerConcentration_[tracerIdx][globalDofIdx][1] = 0.0;
342 solTracerConcentration_[tracerIdx][globalDofIdx] = 0.0;
348 tracerMatrix_ = std::make_unique<TracerMatrix>(numGridDof, numGridDof, TracerMatrix::random);
351 using NeighborSet = std::set<unsigned>;
352 std::vector<NeighborSet> neighbors(numGridDof);
354 Stencil stencil(gridView_, dofMapper_);
355 for (
const auto& elem : elements(gridView_)) {
356 stencil.update(elem);
358 for (
unsigned primaryDofIdx = 0; primaryDofIdx < stencil.numPrimaryDof(); ++primaryDofIdx) {
359 unsigned myIdx = stencil.globalSpaceIndex(primaryDofIdx);
361 for (
unsigned dofIdx = 0; dofIdx < stencil.numDof(); ++dofIdx) {
362 unsigned neighborIdx = stencil.globalSpaceIndex(dofIdx);
363 neighbors[myIdx].insert(neighborIdx);
369 for (
unsigned dofIdx = 0; dofIdx < numGridDof; ++ dofIdx) {
370 tracerMatrix_->setrowsize(dofIdx, neighbors[dofIdx].size());
372 tracerMatrix_->endrowsizes();
377 for (
unsigned dofIdx = 0; dofIdx < numGridDof; ++ dofIdx) {
378 for (
const auto& index : neighbors[dofIdx]) {
379 tracerMatrix_->addindex(dofIdx, index);
382 tracerMatrix_->endindices();