From b32c212e5615f472f606bc54948a0cd40cd64b67 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Thu, 4 Oct 2012 11:25:42 +0200 Subject: [PATCH] Fix mass balance report. Now reports surface volumes for water and oil, instead of useless (for mass balance) reservoir volumes. Also makes the case with multiple transport steps more accurate. --- opm/polymer/SimulatorCompressiblePolymer.cpp | 89 +++++++++++--------- 1 file changed, 51 insertions(+), 38 deletions(-) diff --git a/opm/polymer/SimulatorCompressiblePolymer.cpp b/opm/polymer/SimulatorCompressiblePolymer.cpp index efb09ea61..8e07c55a3 100644 --- a/opm/polymer/SimulatorCompressiblePolymer.cpp +++ b/opm/polymer/SimulatorCompressiblePolymer.cpp @@ -265,22 +265,16 @@ namespace Opm double ttime = 0.0; Opm::time::StopWatch total_timer; total_timer.start(); - double init_satvol[2] = { 0.0 }; + double init_surfvol[2] = { 0.0 }; double init_polymass = 0.0; - double satvol[2] = { 0.0 }; + double inplace_surfvol[2] = { 0.0 }; double polymass = 0.0; double polymass_adsorbed = 0.0; - double injected[2] = { 0.0 }; - double produced[2] = { 0.0 }; - double polyinj = 0.0; - double polyprod = 0.0; double tot_injected[2] = { 0.0 }; double tot_produced[2] = { 0.0 }; double tot_polyinj = 0.0; double tot_polyprod = 0.0; - Opm::computeSaturatedVol(porevol, state.saturation(), init_satvol); - std::cout << "\nInitial saturations are " << init_satvol[0]/tot_porevol_init - << " " << init_satvol[1]/tot_porevol_init << std::endl; + Opm::computeSaturatedVol(porevol, state.surfacevol(), init_surfvol); Opm::Watercut watercut; watercut.push(0.0, 0.0, 0.0); Opm::WellReport wellreport; @@ -339,17 +333,32 @@ namespace Opm stepsize /= double(num_transport_substeps_); std::cout << "Making " << num_transport_substeps_ << " transport substeps." << std::endl; } + double injected[2] = { 0.0 }; + double produced[2] = { 0.0 }; + double polyinj = 0.0; + double polyprod = 0.0; for (int tr_substep = 0; tr_substep < num_transport_substeps_; ++tr_substep) { tsolver_.solve(&state.faceflux()[0], initial_pressure, state.pressure(), &initial_porevol[0], &porevol[0], &transport_src[0], stepsize, inflow_c, state.saturation(), state.surfacevol(), state.concentration(), state.maxconcentration()); + double substep_injected[2] = { 0.0 }; + double substep_produced[2] = { 0.0 }; + double substep_polyinj = 0.0; + double substep_polyprod = 0.0; Opm::computeInjectedProduced(props_, poly_props_, state.pressure(), state.surfacevol(), state.saturation(), state.concentration(), state.maxconcentration(), - transport_src, stepsize, inflow_c, injected, produced, - polyinj, polyprod); + transport_src, stepsize, inflow_c, + substep_injected, substep_produced, + substep_polyinj, substep_polyprod); + injected[0] += substep_injected[0]; + injected[1] += substep_injected[1]; + produced[0] += substep_produced[0]; + produced[1] += substep_produced[1]; + polyinj += substep_polyinj; + polyprod += substep_polyprod; if (gravity_ != 0 && use_segregation_split_) { tsolver_.solveGravity(columns_, stepsize, state.saturation(), state.surfacevol(), @@ -362,7 +371,7 @@ namespace Opm ttime += tt; // Report volume balances. - Opm::computeSaturatedVol(porevol, state.saturation(), satvol); + Opm::computeSaturatedVol(porevol, state.surfacevol(), inplace_surfvol); polymass = Opm::computePolymerMass(porevol, state.saturation(), state.concentration(), poly_props_.deadPoreVol()); polymass_adsorbed = Opm::computePolymerAdsorbed(grid_, props_, poly_props_, state, rock_comp_props_); @@ -374,40 +383,44 @@ namespace Opm tot_polyprod += polyprod; std::cout.precision(5); const int width = 18; - std::cout << "\nVolume and polymer mass balance: " - " water(pv) oil(pv) polymer(kg)\n"; - std::cout << " Saturated volumes: " - << std::setw(width) << satvol[0]/tot_porevol_init - << std::setw(width) << satvol[1]/tot_porevol_init + std::cout << "\nMass balance: " + " water(surfvol) oil(surfvol) polymer(kg)\n"; + std::cout << " In-place: " + << std::setw(width) << inplace_surfvol[0] + << std::setw(width) << inplace_surfvol[1] << std::setw(width) << polymass << std::endl; - std::cout << " Adsorbed volumes: " + std::cout << " Adsorbed: " << std::setw(width) << 0.0 << std::setw(width) << 0.0 << std::setw(width) << polymass_adsorbed << std::endl; - std::cout << " Injected volumes: " - << std::setw(width) << injected[0]/tot_porevol_init - << std::setw(width) << injected[1]/tot_porevol_init + std::cout << " Injected: " + << std::setw(width) << injected[0] + << std::setw(width) << injected[1] << std::setw(width) << polyinj << std::endl; - std::cout << " Produced volumes: " - << std::setw(width) << produced[0]/tot_porevol_init - << std::setw(width) << produced[1]/tot_porevol_init + std::cout << " Produced: " + << std::setw(width) << produced[0] + << std::setw(width) << produced[1] << std::setw(width) << polyprod << std::endl; - std::cout << " Total inj volumes: " - << std::setw(width) << tot_injected[0]/tot_porevol_init - << std::setw(width) << tot_injected[1]/tot_porevol_init + std::cout << " Total inj: " + << std::setw(width) << tot_injected[0] + << std::setw(width) << tot_injected[1] << std::setw(width) << tot_polyinj << std::endl; - std::cout << " Total prod volumes: " - << std::setw(width) << tot_produced[0]/tot_porevol_init - << std::setw(width) << tot_produced[1]/tot_porevol_init + std::cout << " Total prod: " + << std::setw(width) << tot_produced[0] + << std::setw(width) << tot_produced[1] << std::setw(width) << tot_polyprod << std::endl; - std::cout << " In-place + prod - inj: " - << std::setw(width) << (satvol[0] + tot_produced[0] - tot_injected[0])/tot_porevol_init - << std::setw(width) << (satvol[1] + tot_produced[1] - tot_injected[1])/tot_porevol_init - << std::setw(width) << (polymass + tot_polyprod - tot_polyinj + polymass_adsorbed) << std::endl; - std::cout << " Init - now - pr + inj: " - << std::setw(width) << (init_satvol[0] - satvol[0] - tot_produced[0] + tot_injected[0])/tot_porevol_init - << std::setw(width) << (init_satvol[1] - satvol[1] - tot_produced[1] + tot_injected[1])/tot_porevol_init - << std::setw(width) << (init_polymass - polymass - tot_polyprod + tot_polyinj - polymass_adsorbed) + const double balance[3] = { init_surfvol[0] - inplace_surfvol[0] - tot_produced[0] + tot_injected[0], + init_surfvol[1] - inplace_surfvol[1] - tot_produced[1] + tot_injected[1], + init_polymass - polymass - tot_polyprod + tot_polyinj - polymass_adsorbed }; + std::cout << " Initial - inplace + inj - prod: " + << std::setw(width) << balance[0] + << std::setw(width) << balance[1] + << std::setw(width) << balance[2] + << std::endl; + std::cout << " Relative mass error: " + << std::setw(width) << balance[0]/(init_surfvol[0] + tot_injected[0]) + << std::setw(width) << balance[1]/(init_surfvol[1] + tot_injected[1]) + << std::setw(width) << balance[2]/(init_polymass + tot_polyinj) << std::endl; std::cout.precision(8);