C++ SimulatorTimer::simulationTimeElapsed方法代码示例

//.........这里部分代码省略.........
            ptime += pt;

            // Optionally, check if well controls are satisfied.
            if (check_well_controls_) {
                Opm::computePhaseFlowRatesPerWell(*wells_,
                                                  well_state.perfRates(),
                                                  fractional_flows,
                                                  well_resflows_phase);
                std::cout << "Checking well conditions." << std::endl;
                // For testing we set surface := reservoir
                well_control_passed = wells_manager_.conditionsMet(well_state.bhp(), well_resflows_phase, well_resflows_phase);
                ++well_control_iteration;
                if (!well_control_passed && well_control_iteration > max_well_control_iterations_) {
                    OPM_THROW(std::runtime_error, "Could not satisfy well conditions in " << max_well_control_iterations_ << " tries.");
                }
                if (!well_control_passed) {
                    std::cout << "Well controls not passed, solving again." << std::endl;
                } else {
                    std::cout << "Well conditions met." << std::endl;
                }
            }
        } while (!well_control_passed);

        // Update pore volumes if rock is compressible.
        if (rock_comp_props_ && rock_comp_props_->isActive()) {
            initial_porevol = porevol;
            computePorevolume(grid_, props_.porosity(), *rock_comp_props_, state.pressure(), porevol);
        }

        // Process transport sources (to include bdy terms and well flows).
        Opm::computeTransportSource(props_, wells_, well_state, transport_src);

        // Find inflow rate.
        const double current_time = timer.simulationTimeElapsed();
        double stepsize = timer.currentStepLength();
        polymer_inflow_.getInflowValues(current_time, current_time + stepsize, polymer_inflow_c);


        // Solve transport.
        transport_timer.start();
        if (num_transport_substeps_ != 1) {
            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(), state.temperature(), &initial_porevol[0], &porevol[0],
                           &transport_src[0], &polymer_inflow_c[0], stepsize,
                           state.saturation(), state.surfacevol(),
                           state.getCellData( state.CONCENTRATION ), state.getCellData( state.CMAX ));
            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,
                                         transport_src, polymer_inflow_c, stepsize,
                                         substep_injected, substep_produced,
                                         substep_polyinj, substep_polyprod);
            injected[0] += substep_injected[0];
            injected[1] += substep_injected[1];
            produced[0] += substep_produced[0];

    SimulatorReport SimulatorFullyImplicitCompressiblePolymer::Impl::run(SimulatorTimer& timer,
                                                                         PolymerBlackoilState& state)
    {
        WellStateFullyImplicitBlackoil prev_well_state;
        // Initialisation.
        std::vector<double> porevol;
        if (rock_comp_props_ && rock_comp_props_->isActive()) {
            computePorevolume(grid_, props_.porosity(), *rock_comp_props_, state.pressure(), porevol);
        } else {
            computePorevolume(grid_, props_.porosity(), porevol);
        }
        std::vector<double> initial_porevol = porevol;

        std::vector<double> polymer_inflow_c(grid_.number_of_cells);
        // Main simulation loop.
        Opm::time::StopWatch solver_timer;
        double stime = 0.0;
        Opm::time::StopWatch step_timer;
        Opm::time::StopWatch total_timer;
        total_timer.start();
        std::string tstep_filename = output_dir_ + "/step_timing.txt";
        std::ofstream tstep_os(tstep_filename.c_str());

        //Main simulation loop.
        while (!timer.done()) {
#if 0
            double tot_injected[2] = { 0.0 };
            double tot_produced[2] = { 0.0 };
            Opm::Watercut watercut;
            watercut.push(0.0, 0.0, 0.0);
            std::vector<double> fractional_flows;
            std::vector<double> well_resflows_phase;
            if (wells_) {
                well_resflows_phase.resize((wells_->number_of_phases)*(wells_->number_of_wells), 0.0);
            }
            std::fstream tstep_os;
            if (output_) {
                std::string filename = output_dir_ + "/step_timing.param";
                tstep_os.open(filename.c_str(), std::fstream::out | std::fstream::app);
            }
#endif
            // Report timestep and (optionally) write state to disk.

            step_timer.start();
            timer.report(std::cout);

            WellsManager wells_manager(eclipse_state_,
                                       timer.currentStepNum(),
                                       Opm::UgGridHelpers::numCells(grid_),
                                       Opm::UgGridHelpers::globalCell(grid_),
                                       Opm::UgGridHelpers::cartDims(grid_),
                                       Opm::UgGridHelpers::dimensions(grid_),
                                       Opm::UgGridHelpers::cell2Faces(grid_),
                                       Opm::UgGridHelpers::beginFaceCentroids(grid_),
                                       props_.permeability());
            const Wells* wells = wells_manager.c_wells();
            WellStateFullyImplicitBlackoil well_state;
            well_state.init(wells, state.blackoilState(), prev_well_state);
            //Compute polymer inflow.
            std::unique_ptr<PolymerInflowInterface> polymer_inflow_ptr;
            if (deck_->hasKeyword("WPOLYMER")) {
                if (wells_manager.c_wells() == 0) {
                    OPM_THROW(std::runtime_error, "Cannot control polymer injection via WPOLYMER without wells.");
                }
                polymer_inflow_ptr.reset(new PolymerInflowFromDeck(deck_, eclipse_state_, *wells, Opm::UgGridHelpers::numCells(grid_), timer.currentStepNum()));
            } else {
                polymer_inflow_ptr.reset(new PolymerInflowBasic(0.0*Opm::unit::day,
                                                                1.0*Opm::unit::day,
                                                                0.0));
            }
            std::vector<double> polymer_inflow_c(Opm::UgGridHelpers::numCells(grid_));
            polymer_inflow_ptr->getInflowValues(timer.simulationTimeElapsed(),
                                                timer.simulationTimeElapsed() + timer.currentStepLength(),
                                                polymer_inflow_c);

            if (output_ && (timer.currentStepNum() % output_interval_ == 0)) {
                if (output_vtk_) {
                    outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_);
                }
                outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_);
            }
            if (output_) {
                if (timer.currentStepNum() == 0) {
                    output_writer_.writeInit(timer);
                }
                output_writer_.writeTimeStep(timer, state.blackoilState(), well_state);
            }
            // Run solver.
            solver_timer.start();
            FullyImplicitCompressiblePolymerSolver solver(grid_, props_, geo_, rock_comp_props_, polymer_props_, *wells_manager.c_wells(), linsolver_);
            solver.step(timer.currentStepLength(), state, well_state, polymer_inflow_c);
            // Stop timer and report.
            solver_timer.stop();
            const double st = solver_timer.secsSinceStart();
            std::cout << "Fully implicit solver took:  " << st << " seconds." << std::endl;

            stime += st;
            // Update pore volumes if rock is compressible.
            if (rock_comp_props_ && rock_comp_props_->isActive()) {
                initial_porevol = porevol;
//.........这里部分代码省略.........