apfel/ome/rpd_distribution.h

0001 /******************************************************************************
0002  * This file is part of libome                                                *
0003  * Copyright (C) 2025 Arnd Behring, Kay Schoenwald                            *
0004  * SPDX-License-Identifier: GPL-3.0-or-later                                  *
0005  ******************************************************************************/
0006
0007 /**
0008  * \file
0009  * \brief Container for mixed distributions with regular, plus and delta part
0010  */
0011
0012 #ifndef LIBOME_RPD_DISTRIBUTION_H
0013 #define LIBOME_RPD_DISTRIBUTION_H
0014
0015 #include <optional>
0016
0017 namespace apfel
0018 {
0019   namespace ome
0020   {
0021     /**
0022      * \brief Mixed distribution with regular, plus and delta part
0023      *
0024      * \details
0025      * This is a container that models a distribution \f$f(\dots,x)\f$ which is a
0026      * sum of a regular, a plus and a delta part:
0027      * \f[
0028      *   f(\dots,x) = f_\mathrm{reg}(\dots,x)
0029      *     +\left[f_+(\dots,x)\right]_+
0030      *     +\delta(1-x) f_\delta(\dots)
0031      * \f]
0032      * Whether or not each of the three parts is present is optional. All three
0033      * parts should have the compatible evaluation signatures: the parameters
0034      * denoted by \f$\dots\f$ above should be the same for each of the three
0035      * parts.
0036      *
0037      * \tparam Tfuncreg Type for the wrapped regular part
0038      * \tparam Tfuncplus Type for the wrapped plus part
0039      * \tparam Tfuncdelta Type for the wrapped delta part
0040      */
0041     template<typename Tfuncreg,
0042              typename Tfuncplus,
0043              typename Tfuncdelta>
0044     class rpd_distribution
0045     {
0046     public:
0047       /// Type alias for the regular part template parameter
0048       using element_regular_type = Tfuncreg;
0049       /// Type alias for the plus part template parameter
0050       using element_plus_type = Tfuncplus;
0051       /// Type alias for the delta part template parameter
0052       using element_delta_type = Tfuncdelta;
0053
0054
0055       /**
0056        * \brief Default construtor
0057        *
0058        * \details
0059        * Constructs an "empty" container view with no regular, plus or delta
0060        * part.
0061        */
0062       rpd_distribution()
0063         : regular_(std::nullopt), plus_(std::nullopt), delta_(std::nullopt) {};
0064
0065       /**
0066        * \brief Construct distribution from regular, plus and delta parts
0067        *
0068        * \details
0069        * The individual parts are wrapped in std::optional so that their absence
0070        * can be signalled by passing std::nullopt.
0071        *
0072        * \param regular Regular part. Can be a numerical value or a callable
0073        * \param plus Plus part. Can be a numerical value or a callable
0074        * \param delta Delta part. Can be a numerical value or a callable
0075        */
0076       rpd_distribution(
0077         std::optional<element_regular_type> regular,
0078         std::optional<element_plus_type> plus,
0079         std::optional<element_delta_type> delta
0080       )
0081         : regular_(regular), plus_(plus), delta_(delta) {};
0082
0083       /**
0084        * \brief Check whether the regular part is present
0085        *
0086        * \return True if the regular part is present, false otherwise.
0087        */
0088       bool has_regular() const
0089       {
0090         return(regular_.has_value());
0091       }
0092
0093       /**
0094        * \brief Check whether the plus part is present
0095        *
0096        * \return True if the plus part is present, false otherwise.
0097        */
0098       bool has_plus() const
0099       {
0100         return(plus_.has_value());
0101       }
0102
0103       /**
0104        * \brief Check whether the delta part is present
0105        *
0106        * \return True if the delta part is present, false otherwise.
0107        */
0108       bool has_delta() const
0109       {
0110         return(delta_.has_value());
0111       }
0112
0113       /**
0114        * \brief Retrieve the regular part
0115        *
0116        * \return Constant reference to an std::optional wrapping the regular
0117        *         part. If the regular part is absent, the std::optional
0118        *         contains std::nullopt.
0119        */
0120       const std::optional<element_regular_type>& get_regular() const
0121       {
0122         return(regular_);
0123       };
0124
0125       /**
0126        * \brief Retrieve the plus part
0127        *
0128        * \return Constant reference to an std::optional wrapping the plus part.
0129        *         If the plus part is absent, the std::optional contains
0130        *         std::nullopt.
0131        */
0132       const std::optional<element_plus_type>& get_plus() const
0133       {
0134         return(plus_);
0135       };
0136
0137       /**
0138        * \brief Retrieve the delta part
0139        *
0140        * \return Constant reference to an std::optional wrapping the delta part.
0141        *         If the delta part is absent, the std::optional contains
0142        *         std::nullopt.
0143        */
0144       const std::optional<element_delta_type>& get_delta() const
0145       {
0146         return(delta_);
0147       };
0148
0149       /**
0150        * \brief Extract the n-th coefficient from the regular, plus and delta parts
0151        *        into a new rpd_distribution container
0152        *
0153        * \details
0154        * In order to be applicable, the regular, plus and delta parts must have a
0155        * coefficient_type, coefficient_has_view and coefficient_view_type member
0156        * types. Also the n-th coefficient must be extractable using operator[][n].
0157        * If coefficient_has_view is std::true_type, the function extracts the view
0158        * of the n-th coefficient instead of the coefficient itself.
0159        *
0160        * \param n Exponent of the coefficient to extract
0161        *
0162        * \return rpd_distribution container with the n-th coefficient (or a view
0163        *         on it, if available) of the regular, plus and delta parts
0164        */
0165       auto get_coefficient(int n) const
0166       {
0167         using reg_type
0168         = std::conditional_t<element_regular_type::coefficient_has_view::value,
0169         typename element_regular_type::coefficient_view_type,
0170         typename element_regular_type::coefficient_type>;
0171         using plus_type
0172         = std::conditional_t<element_plus_type::coefficient_has_view::value,
0173         typename element_plus_type::coefficient_view_type,
0174         typename element_plus_type::coefficient_type>;
0175         using delta_type
0176         = std::conditional_t<element_delta_type::coefficient_has_view::value,
0177         typename element_delta_type::coefficient_view_type,
0178         typename element_delta_type::coefficient_type>;
0179         std::optional<reg_type> reg;
0180         std::optional<plus_type> plus;
0181         std::optional<delta_type> delta;
0182
0183         if(has_regular())
0184           {
0185             if constexpr (element_regular_type::coefficient_has_view::value)
0186               reg = std::make_optional((*regular_)[n].get_view());
0187             else
0188               reg = std::make_optional((*regular_)[n]);
0189           }
0190         if(has_plus())
0191           {
0192             if constexpr (element_plus_type::coefficient_has_view::value)
0193               plus = std::make_optional((*plus_)[n].get_view());
0194             else
0195               plus = std::make_optional((*plus_)[n]);
0196           }
0197         if(has_delta())
0198           {
0199             if constexpr (element_delta_type::coefficient_has_view::value)
0200               delta = std::make_optional((*delta_)[n].get_view());
0201             else
0202               delta = std::make_optional((*delta_)[n]);
0203           }
0204
0205         return(
0206                 rpd_distribution<reg_type, plus_type, delta_type> {reg, plus, delta}
0207               );
0208       };
0209
0210       /**
0211        * \brief Apply the truncate operation on each of the regular, plus and
0212        *        delta parts (if supported)
0213        *
0214        * \param n Exponent of the last order to include
0215        *
0216        * \return An rpd_distribution with truncated regular, plus and delta
0217        *         parts
0218        */
0219       auto truncate(int n) const
0220       {
0221         using reg_type = decltype(std::declval<element_regular_type>().truncate(n));
0222         using plus_type = decltype(std::declval<element_plus_type>().truncate(n));
0223         using delta_type = decltype(std::declval<element_delta_type>().truncate(n));
0224
0225         std::optional<reg_type> reg;
0226         std::optional<plus_type> plus;
0227         std::optional<delta_type> delta;
0228
0229         if(has_regular())
0230           {
0231             reg = std::make_optional(regular_->truncate(n));
0232           }
0233         if(has_plus())
0234           {
0235             plus = std::make_optional(plus_->truncate(n));
0236           }
0237         if(has_delta())
0238           {
0239             delta = std::make_optional(delta_->truncate(n));
0240           }
0241
0242         return(
0243                 rpd_distribution<reg_type, plus_type, delta_type> {reg, plus, delta}
0244               );
0245       };
0246
0247     private:
0248       std::optional<element_regular_type> regular_;
0249       std::optional<element_plus_type> plus_;
0250       std::optional<element_delta_type> delta_;
0251     };
0252
0253   }
0254 }
0255
0256 #endif