non_embedded_formula_wrapper.h

Go to the documentation of this file.

/*
 * Copyright 2021 MusicScience37 (Kenta Kabashima)
 *
 * Licensed under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License.
 * You may obtain a copy of the License at
 *
 *     http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
#pragma once
 
#include "num_collect/base/index_type.h"
#include "num_collect/constants/half.h"  // IWYU pragma: keep
#include "num_collect/logging/log_tag_view.h"
#include "num_collect/ode/concepts/formula.h"
#include "num_collect/ode/embedded_solver.h"
 
namespace num_collect::ode {
 
template <concepts::formula BaseFormula>
class non_embedded_formula_wrapper {
public:
    using base_formula_type = BaseFormula;
 
    using problem_type = typename base_formula_type::problem_type;
 
    using variable_type = typename problem_type::variable_type;
 
    using scalar_type = typename problem_type::scalar_type;
 
    static constexpr index_type stages = base_formula_type::stages;
 
    static constexpr index_type order = base_formula_type::order;
 
    static constexpr index_type lesser_order = base_formula_type::order;
 
    static constexpr logging::log_tag_view log_tag = base_formula_type::log_tag;
 
    explicit non_embedded_formula_wrapper(
        const problem_type& problem = problem_type())
        : formula_(problem) {}
 
    void step(scalar_type time, scalar_type step_size,
        const variable_type& current, variable_type& estimate) {
        formula_.step(time, step_size, current, estimate);
    }
 
    void step_embedded(scalar_type time, scalar_type step_size,
        const variable_type& current, variable_type& estimate,
        variable_type& error) {
        variable_type rough_estimate;
        step(time, step_size, current, rough_estimate);
        variable_type half_estimate;
        const scalar_type half_step_size =
            constants::half<scalar_type> * step_size;
        step(time, half_step_size, current, half_estimate);
        step(time, half_step_size, half_estimate, estimate);
        error = rough_estimate - estimate;
    }
 
    [[nodiscard]] auto problem() -> problem_type& { return formula_.problem(); }
 
    [[nodiscard]] auto problem() const -> const problem_type& {
        return formula_.problem();
    }
 
private:
    base_formula_type formula_;
};
 
template <typename Formula>
using non_embedded_auto_solver =
    embedded_solver<non_embedded_formula_wrapper<Formula>>;
 
}  // namespace num_collect::ode