CreateProcessData.cpp

Go to the documentation of this file.

// SPDX-FileCopyrightText: Copyright (c) OpenGeoSys Community (opengeosys.org)
// SPDX-License-Identifier: BSD-3-Clause
 
#include <range/v3/algorithm/contains.hpp>
 
#include "BaseLib/Algorithm.h"
#include "BaseLib/ConfigTree.h"
#include "NumLib/ODESolver/TimeDiscretizationBuilder.h"
#ifdef USE_PETSC
#include "NumLib/ODESolver/PETScNonlinearSolver.h"
#endif  // USE_PETSC
#include "CreateProcessData.h"
#include "NumLib/TimeStepping/CreateTimeStepper.h"
 
namespace ProcessLib
{
static std::unique_ptr<ProcessData> makeProcessData(
    std::unique_ptr<NumLib::TimeStepAlgorithm>&& timestepper,
    NumLib::NonlinearSolverBase& nonlinear_solver, int const process_id,
    std::string process_name, Process& process,
    std::unique_ptr<NumLib::TimeDiscretization>&& time_disc,
    std::unique_ptr<NumLib::ConvergenceCriterion>&& conv_crit,
    bool const compensate_non_equilibrium_initial_residuum)
{
    using Tag = NumLib::NonlinearSolverTag;
 
    if (auto* nonlinear_solver_picard =
            dynamic_cast<NumLib::NonlinearSolver<Tag::Picard>*>(
                &nonlinear_solver))
    {
        nonlinear_solver_picard->compensateNonEquilibriumInitialResiduum(
            compensate_non_equilibrium_initial_residuum);
        return std::make_unique<ProcessData>(
            std::move(timestepper), Tag::Picard, *nonlinear_solver_picard,
            std::move(conv_crit), std::move(time_disc), process_id,
            std::move(process_name), process);
    }
    if (auto* nonlinear_solver_newton =
            dynamic_cast<NumLib::NonlinearSolver<Tag::Newton>*>(
                &nonlinear_solver))
    {
        nonlinear_solver_newton->compensateNonEquilibriumInitialResiduum(
            compensate_non_equilibrium_initial_residuum);
        return std::make_unique<ProcessData>(
            std::move(timestepper), Tag::Newton, *nonlinear_solver_newton,
            std::move(conv_crit), std::move(time_disc), process_id,
            std::move(process_name), process);
    }
#ifdef USE_PETSC
    if (auto* nonlinear_solver_petsc =
            dynamic_cast<NumLib::PETScNonlinearSolver*>(&nonlinear_solver))
    {
        for (auto const& pv : process.getProcessVariables(process_id))
        {
            if (pv.get().compensateNonEquilibriumInitialResiduum())
            {
                OGS_FATAL(
                    "The PETSc SNES nonlinear solver does not support "
                    "'compensate_non_equilibrium_initial_residuum', but it "
                    "is set for process variable '{}'. Please use the "
                    "Picard or Newton-Raphson solver instead.",
                    pv.get().getName());
            }
        }
        return std::make_unique<ProcessData>(
            std::move(timestepper), Tag::Newton, *nonlinear_solver_petsc,
            std::move(conv_crit), std::move(time_disc), process_id,
            std::move(process_name), process);
    }
#endif  // USE_PETSC
 
    OGS_FATAL("Encountered unknown nonlinear solver type. Aborting");
}
 
std::vector<std::unique_ptr<ProcessData>> createPerProcessData(
    BaseLib::ConfigTree const& config,
    std::vector<std::unique_ptr<Process>> const& processes,
    std::map<std::string, std::unique_ptr<NumLib::NonlinearSolverBase>> const&
        nonlinear_solvers,
    bool const compensate_non_equilibrium_initial_residuum,
    std::vector<double> const& fixed_times_for_output)
{
    std::vector<std::unique_ptr<ProcessData>> per_process_data;
    std::vector<std::string> process_names;
    int process_id = 0;

    for (auto pcs_config : config.getConfigSubtreeList("process"))
    {
        auto const pcs_name = pcs_config.getConfigAttribute<std::string>("ref");
        auto& pcs = *BaseLib::getIfOrError(
            processes,
            [&pcs_name](std::unique_ptr<Process> const& p)
            { return p->name == pcs_name; },
            "A process with the given name has not been defined.");
 
        auto const process_name =
            pcs_config.getConfigParameter<std::string>("process_name", "");
        if (process_name != "")
        {
            if (ranges::contains(process_names, process_name))
            {
                OGS_FATAL(
                    "The given process name is not unique! Please check the "
                    "element "
                    "'time_loop/process/name' in the project file. Found "
                    "duplicate "
                    "process name '{:s}'.",
                    process_name);
            }
            process_names.emplace_back(process_name);
        }
 
        auto const nl_slv_name =
            pcs_config.getConfigParameter<std::string>("nonlinear_solver");
        auto& nl_slv = *BaseLib::getOrError(
            nonlinear_solvers, nl_slv_name,
            "A nonlinear solver with the given name has not been defined.");
 
        auto time_disc = NumLib::createTimeDiscretization(
            pcs_config.getConfigSubtree("time_discretization"));
 
        auto timestepper = NumLib::createTimeStepper(
            pcs_config.getConfigSubtree("time_stepping"),
            fixed_times_for_output);
 
        auto conv_crit = NumLib::createConvergenceCriterion(
            pcs_config.getConfigSubtree("convergence_criterion"));

        auto output = pcs_config.getConfigSubtreeOptional("output");
        if (output)
        {
            OGS_FATAL(
                "In order to make the specification of output in the project "
                "file consistent, the variables output tags were moved from "
                "xpath "
                "'//OpenGeoSysProject/time_loop/processes/process/output' to "
                "the global output section, i.e., to the xpath "
                "'//OpenGeoSysProject/time_loop/output'. This has to be done "
                "in the current project file!");
        }
 
        per_process_data.emplace_back(makeProcessData(
            std::move(timestepper), nl_slv, process_id, std::move(process_name),
            pcs, std::move(time_disc), std::move(conv_crit),
            compensate_non_equilibrium_initial_residuum));
        ++process_id;
    }
 
    if (per_process_data.size() != processes.size())
    {
        if (processes.size() > 1)
        {
            OGS_FATAL(
                "Some processes have not been configured to be solved by this  "
                "time loop.");
        }
        else
        {
            INFO(
                "The equations of the coupled processes will be solved by the "
                "staggered scheme.");
        }
    }
 
    return per_process_data;
}
}  // namespace ProcessLib