boostorg / odeint Goto Github PK

View Code? Open in Web Editor NEW

39.0 14.0 57.0 4.03 MB

Boost.odeint

License: Boost Software License 1.0

C++ 98.86% Python 0.34% Makefile 0.23% HTML 0.05% Shell 0.05% C 0.12% Fortran 0.14% CMake 0.07% Cuda 0.15%

odeint's Introduction

odeint is a highly flexible library for solving ordinary differential equations.

odeint's People

Contributors

Stargazers

Watchers

Forkers

sanderbog danieljames simongibbons marcelraad jzmaddock akumta limbolily nasailja crossbuild sehe autodesk-forks jason-marshall danielae xltools lwang-astro asomers j0ran corelli-diy myestro hmito mike-devel lubkoll thecount2a kojoley rcv-legado madcatx cryoris zoq jinhang2008 f3rm4rf3r mephphisto samd2 lhofmann eldiener opensourceeda vardhanyadav1714 fieryswampshire longgang zhaofeng-shu33 ekra-ltd angel6708 jingliang2005 financialengineerlab fanquake isabella232 lunarflu axylp gtton siddhantkadwe ka-de seanpm2001 mborland georgthegreat sam-caldwell grafikrobot infinityofficialnetwork magni-mar

odeint's Issues

`static_cast<time_type>` in `integrate_times`

Hi,

it seems to me that there are a number of static_cast<time_type> in the integrate_times implementation, e.g.

odeint/include/boost/numeric/odeint/integrate/detail/integrate_times.hpp

Line 64 in 8b4768e

    
           while( less_with_sign( current_time , static_cast<time_type>(*start_time) , current_dt ) )

odeint/include/boost/numeric/odeint/integrate/detail/integrate_times.hpp

Line 97 in 8b4768e

    
           while( less_with_sign( current_time , static_cast<time_type>(*start_time) , dt ) )

odeint/include/boost/numeric/odeint/integrate/detail/integrate_times.hpp

Line 139 in 8b4768e

Time last_time_point = static_cast<time_type>(*last_time_iterator);

(and potentially more; L139 is where I encountered a compiler error) that don't fit the the target variables' types. Say, for instance, Time in these cases is a Boost Units quantity, then time_type, if I understand correctly, is a plain double (or other floating point format). If I further understand the TimeIterator type correctly, that type should return a Time instance on valid dereference.

Now the static_cast tries to cast a quantity to, say, double (which fails), which is then reassigned to a quantity (which should
also fail).

It seems to me that simply removing the static casts might fix the issue, although I have not spent much time looking into it.
For now I'll use another integrator but I'd be happy to help to solve this issue in the future - if my hypothesis is correct.

Best
Malte

Issues of writes to unknown address in odeint

I am trying to use boost library's odeint with std::vector or std::array. I have no issues with compilation in either case, however, I have runtime issues regarding write to unknown address in the case of std::array. Given below is a MWE.

I have defined two classes: ODE1, and ODE2 with std::vector, and std::array respectively in the header file ode.hpp.

#ifndef ODE_HPP_                                                                   
#define ODE_HPP_                                                                   
                                                                        
// c++ standard library headers                                                    
#include <array>                                                                   
#include <iostream>                                                                
#include <vector>                                                                  
                                                                                     
namespace ode {                                                                    
constexpr size_t n_sections = 1800;                                                
constexpr size_t n_nodes_per_section = 6;                                          
constexpr size_t n_variables = n_sections * n_nodes_per_section;                   
}  // namespace ode                                                                
                                                                                     
class ODE1 {                                                                       
 public:                                                                           
  ODE1(const double o) : V_(ode::n_variables, o) {}  // explicit constructor       
                                                                                  
  std::vector<double>& GetVoltage() { return V_; }                              
                                                                                  
  void operator()(const std::vector<double>& Y, std::vector<double>& dYdt,      
                  const double /* time */) {                                       
    using namespace ode;                                                           
    for (size_t n = 0; n < n_sections; n++) {                                      
      dYdt[n * n_nodes_per_section + 0] = Y[n * n_nodes_per_section + 0];          
      dYdt[n * n_nodes_per_section + 1] = Y[n * n_nodes_per_section + 1];          
      dYdt[n * n_nodes_per_section + 2] = Y[n * n_nodes_per_section + 2];
      dYdt[n * n_nodes_per_section + 3] = Y[n * n_nodes_per_section + 3];          
      dYdt[n * n_nodes_per_section + 4] = Y[n * n_nodes_per_section + 4];          
      dYdt[n * n_nodes_per_section + 5] = Y[n * n_nodes_per_section + 5];          
    }                                                                              
  }                                                                                
                                                                                     
 private:                                                                          
  std::vector<double> V_{};                                                        
};                                                                                 

class ODE2 {                                                                       
 using Array = std::array<double, ode::n_variables>;                              
                                                                                    
 public:                                                                           
  // constructors                                                                  
  ODE2(const double o) { V_.fill(o); }  // explicit constructor                                                       

  Array& GetVoltage() { return V_; }                                            
                                                                                  
  void operator()(const Array& Y, Array& dYdt, const double /* time */) {       
    using namespace ode;                                                        
    for (size_t n = 0; n < n_sections; n++) {                                   
      dYdt[n * n_nodes_per_section + 0] = Y[n * n_nodes_per_section + 0];       
      dYdt[n * n_nodes_per_section + 1] = Y[n * n_nodes_per_section + 1];       
      dYdt[n * n_nodes_per_section + 2] = Y[n * n_nodes_per_section + 2];       
      dYdt[n * n_nodes_per_section + 3] = Y[n * n_nodes_per_section + 3];       
      dYdt[n * n_nodes_per_section + 4] = Y[n * n_nodes_per_section + 4];       
      dYdt[n * n_nodes_per_section + 5] = Y[n * n_nodes_per_section + 5];       
    }                                                                           
  }                                                                             
                                                                                  
 private:                                                                       
  Array V_{};                                                                   
};                                                                              
#endif  //   ODE_HPP_

And the main function is defined in the file `ode.cpp'

// related header                                                                  
#include "ode.hpp"                                                                 

// other library headers                                                           
#include <boost/numeric/odeint.hpp>                                                
                                                                                                                                                                          
int main() {                                                                       
  constexpr double abs_err = 1.0e-10;                                              
  constexpr double rel_err = 1.0e-8;                                               
  constexpr double sim_duration = 20.0;              // ms                         
  constexpr double t_start = 0;                      // ms                         
  constexpr double t_stop = t_start + sim_duration;  // ms                         
  constexpr double dt = 5e-3;                        // ms                         
  {                                                                                
    std::cout << "ODE 1: ";                                                        
    using namespace boost::numeric::odeint;                                        
    using error_stepper_type = runge_kutta_dopri5<std::vector<double>>;            
                                                                                     
    ODE1 o(-60.0);                                                                 
    integrate_adaptive(make_controlled<error_stepper_type>(abs_err, rel_err), o,
                       o.GetVoltage(), t_start, t_stop, dt);                       
  }                                                                                
  {                                                                                
    std::cout << "ODE 1: ";                                                        
    using namespace boost::numeric::odeint;                                        
    using error_stepper_type =                                                     
          runge_kutta_dopri5<std::array<double, ode::n_variables>>;                  
                                                                                     
    ODE2 o(-60.0);                                                                 
    integrate_adaptive(make_controlled<error_stepper_type>(abs_err, rel_err), o,
                       o.GetVoltage(), t_start, t_stop, dt);                       
  }                                                                                
  return 0;                                                                        
}

I am compiling using the clang compiler as follows
c++ -O0 -g -fsanitize=address -fno-omit-frame-pointer -Wall -Wextra -Wpedantic -std=c++17 ode.cpp -o ode

Compilation output (warnings)

In file included from ode.cpp:5:
In file included from /usr/local/include/boost/numeric/odeint.hpp:79:
In file included from /usr/local/include/boost/numeric/odeint/stepper/generation.hpp:32:
/usr/local/include/boost/numeric/odeint/stepper/generation/generation_controlled_adams_bashforth_moulton.hpp:43:99: warning: unused parameter 'stepper' [-Wunused-parameter]
    controller_type operator()( value_type abs_error , value_type rel_error , const stepper_type &stepper )
                                                                                                  ^
/usr/local/include/boost/numeric/odeint/stepper/generation/generation_controlled_adams_bashforth_moulton.hpp:49:71: warning: unused parameter 'stepper' [-Wunused-parameter]
                                time_type max_dt, const stepper_type &stepper )
                                                                      ^
2 warnings generated.

As I mentioned above, I do not run into issues with compilation. When I run the executable, I have the following issues

AddressSanitizer:DEADLYSIGNAL
=================================================================
==29478==ERROR: AddressSanitizer: stack-overflow on address 0x7ffee1db40e0 (pc 0x00010d66140c bp 0x7ffee1dde910 sp 0x7ffee1da5d80 T0)
    #0 0x10d66140c in void boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>::do_step_impl<ODE2, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul> const&, std::__1::array<double, 7200ul> const&, double, std::__1::array<double, 7200ul>&, std::__1::array<double, 7200ul>&, double, std::__1::array<double, 7200ul>&) runge_kutta_dopri5.hpp:165
    #1 0x10d660d9c in void boost::numeric::odeint::explicit_error_stepper_fsal_base<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, (unsigned short)5, (unsigned short)5, (unsigned short)4, std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>::do_step<ODE2, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul> const&, std::__1::array<double, 7200ul> const&, double, std::__1::array<double, 7200ul>&, std::__1::array<double, 7200ul>&, double, std::__1::array<double, 7200ul>&) explicit_error_stepper_fsal_base.hpp:272
    #2 0x10d660451 in boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step<ODE2, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul> const&, std::__1::array<double, 7200ul> const&, double&, std::__1::array<double, 7200ul>&, std::__1::array<double, 7200ul>&, double&) controlled_runge_kutta.hpp:768
    #3 0x10d65ea18 in boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step<ODE2, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul>&, std::__1::array<double, 7200ul>&, double&, double&) controlled_runge_kutta.hpp:718
    #4 0x10d65d9b1 in boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step_v1<ODE2, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul>&, double&, double&) controlled_runge_kutta.hpp:902
    #5 0x10d65d112 in boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step<ODE2, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul>&, double&, double&) controlled_runge_kutta.hpp:619
    #6 0x10d65cc3e in unsigned long boost::numeric::odeint::detail::integrate_adaptive<boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::null_observer>(boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>&, double&, double, double&, boost::numeric::odeint::null_observer, boost::numeric::odeint::controlled_stepper_tag) integrate_adaptive.hpp:103
    #7 0x10d65c401 in unsigned long boost::numeric::odeint::integrate_adaptive<boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::null_observer>(boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>&, double, double, double, boost::numeric::odeint::null_observer) integrate_adaptive.hpp:42
    #8 0x10d63acde in unsigned long boost::numeric::odeint::integrate_adaptive<boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>, double>(boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>, boost::numeric::odeint::default_error_checker<double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<double, double>, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>, ODE2, std::__1::array<double, 7200ul>&, double, double, double) integrate_adaptive.hpp:83
    #9 0x10d639fac in main ode.cpp:30
    #10 0x7fff203d0f3c in start+0x0 (libdyld.dylib:x86_64+0x15f3c)

SUMMARY: AddressSanitizer: stack-overflow runge_kutta_dopri5.hpp:165 in void boost::numeric::odeint::runge_kutta_dopri5<std::__1::array<double, 7200ul>, double, std::__1::array<double, 7200ul>, double, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations, boost::numeric::odeint::initially_resizer>::do_step_impl<ODE2, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul>, std::__1::array<double, 7200ul> >(ODE2, std::__1::array<double, 7200ul> const&, std::__1::array<double, 7200ul> const&, double, std::__1::array<double, 7200ul>&, std::__1::array<double, 7200ul>&, double, std::__1::array<double, 7200ul>&)
==29478==ABORTING
ODE 1: ODE 1: zsh: abort      ./ode

Strange enough, I do not have issues when I change the parameter, n_sections, in the header from 1800 to 1000. I could not figure out the issue, is this something to with the boost library containers? Any help would be greatly appreciated.

Boost\libs\numeric\odeint failed to build due to error C2955 and error C3203 on MSVC

Issue Description:
Boost\libs\numeric\odeint failed to build due to error C2955 and error C3203 on MSVC. Could you please take a look? The master branch commit we used is a1552ab.

Build Steps:

git clone -c core.autocrlf=true --recursive https://github.com/boostorg/boost.git boost
open a VS 2019 x86 command prompt and browse to boost
.\bootstrap
.\b2 headers variant=release --build-dir=..\out\x86rel address-model=32
.\b2 variant=release --build-dir=..\out\x86rel address-model=32
.\b2 -j16 variant=release --build-dir=..\out\x86rel libs\numeric\odeint\test

Error Info:
runge_kutta_error_concepts.cpp
.\boost/bind.hpp(41): note: The practice of declaring the Bind placeholders (1, 2, ...) in the global namespace is deprecated. Please use <boost/bind/bind.hpp> + using namespace boost::placeholders, or define BOOST_BIND_GLOBAL_PLACEHOLDERS to retain the current behavior.
.\boost/numeric/odeint/algebra/detail/extract_value_type.hpp(47): error C2955: 'boost::type': use of class template requires template argument list
.\boost/type.hpp(14): note: see declaration of 'boost::type'
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(150): note: see reference to class template instantiation 'boost::numeric::odeint::detail::extract_value_type<State,void>' being compiled
with
[
State=boost::multiprecision::number<boost::multiprecision::backends::cpp_dec_float<50,int32_t,void>,boost::multiprecision::et_on>
]
.\boost/mpl/aux/has_type.hpp(20): note: see reference to class template instantiation 'error_stepper_methods' being compiled
with
[
U1=boost::multiprecision::number<boost::multiprecision::backends::cpp_dec_float<50,int32_t,void>,boost::multiprecision::et_on>
]
.\boost/mpl/aux/has_type.hpp(20): note: while compiling class template member function 'boost::mpl::aux::yes_tag boost::mpl::aux::has_type<error_stepper_methods,boost::mpl::bool_>::gcc_3_2_wknd::test(volatile const boost::mpl::aux::type_wrapper *,boost::mpl::aux::type_wrapper<U::type> *)'
with
[
U1=boost::multiprecision::number<boost::multiprecision::backends::cpp_dec_float<50,int32_t,void>,boost::multiprecision::et_on>
]
.\boost/mpl/aux_/preprocessed/plain/quote.hpp(37): note: see reference to class template instantiation 'boost::mpl::aux::has_type<error_stepper_methods,boost::mpl::bool_>' being compiled
with
[
U1=boost::multiprecision::number<boost::multiprecision::backends::cpp_dec_float<50,int32_t,void>,boost::multiprecision::et_on>
]
.\boost/mpl/aux_/preprocessed/plain/apply_wrap.hpp(39): note: see reference to class template instantiation 'boost::mpl::quote1<F,Tag>::apply' being compiled
with
[
F=error_stepper_methods,
Tag=boost::mpl::void_,
T1=boost::multiprecision::number<boost::multiprecision::backends::cpp_dec_float<50,int32_t,void>,boost::multiprecision::et_on>
]
...........
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(171): note: see reference to class template instantiation 'boost::mpl::copy<container_types,boost::mpl::inserter<boost::mpl::vector0boost::mpl::na,boost::mpl::insert_range<boost::mpl::_1,boost::mpl::endboost::mpl::_1,error_stepper_methodsboost::mpl::_2>>>' being compiled
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(151): error C3203: 'type': unspecialized class template can't be used as a template argument for template parameter 'Value', expected a real type
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(152): error C3203: 'type': unspecialized class template can't be used as a template argument for template parameter 'Value', expected a real type
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(153): error C3203: 'type': unspecialized class template can't be used as a template argument for template parameter 'Value', expected a real type
libs\numeric\odeint\test\runge_kutta_error_concepts.cpp(154): error C3203: 'type': unspecialized class template can't be used as a template argument for template parameter 'Value', expected a real type

Detailed log file:
test.log.7.log

Event-based integration

Hi!

Does odeint support event-based integration?
Sometimes it is useful to stop not on specific finish time, but discrete event, when stop condition of the integrator is not defined by the time t but by a target condition on state y (say y[0]=1.0 for example).
Some open-source libraries provide this functionality, like hipparchus(Java): https://www.hipparchus.org/hipparchus-ode/#Discrete_Events_Handling

value_type cannot have an explicit contructor

There are many places in odeint similar to the following (taken from runge_kutta_dopri5.hpp):

        stepper_base_type::m_algebra.for_each3( m_x_tmp.m_v , in , dxdt_in ,
                typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , dt*b21 ) );

Unfortunately, if value_type's constructor is marked explicit, this will fail. I ran into this trying to use a ceres a jet type for the scalars here for the purpose of automatic differentiation.

It seems like it would be more robust to use something like

const value_type one{1.0};

in these scenarios instead of the raw 1.0. And if I make that change, the ceres jet type works correctly.

I would be willing to work on a PR if this is a change the maintainers would accept.

Visual Studio warnings (max_step_checker.hpp / sprintf)

  C:\devel\boost_1_64_0\boost/numeric/odeint/integrate/max_step_checker.hpp(72): warning C4996: 'sprintf': This function or variable may be unsafe. Consider using sprintf_s instead. To disable deprecation, use _CRT_SECURE_NO_WARNINGS. See online help for details. [C:\devel\vinecopulib\build\vinecopulib.vcxproj]
  C:\devel\boost_1_64_0\boost/numeric/odeint/integrate/max_step_checker.hpp(104): warning C4996: 'sprintf': This function or variable may be unsafe. Consider using sprintf_s instead. To disable deprecation, use _CRT_SECURE_NO_WARNINGS. See online help for details. [C:\devel\vinecopulib\build\vinecopulib.vcxproj]

HTH,
Sylwester

Uninitialized member variable boost::numeric::odeint::bulrisch_stoer<>::m_dt_last

It seems that m_dt_last is used w/o initialization in the first call to try_check:

// boots/numeric/odeint/stepper/bulrisch_stoer.hpp
206        if( dt != m_dt_last )    // <-- m_dt_last indeterminate?
207        {
208          reset(); // step size changed from outside -> reset
209        }

Having ode conditionally depend on Boost::MPI (or allow standalone for cmake)

odeint's cmake currently requires Boost::mpi (see here) and I think a few other boost libraries that much be prebuilt before hand. I think mpi can be kind of a big requirement for a lot of users. It's possible for odeint to be included as a header only library if it is in standalone mode.

Could you add a cmake flag to specify whether the other boost libraries should be included as dependencies? Something like

if (not BOOST_NUMERIC_ODEINT_STANDALONE)
  target_link_libraries(boost_numeric_odeint ...)
else()
  target_compile_definitions(boost_numeric_odeint BOOST_NUMERIC_ODEINT_STANDALONE)
endif()

Version 2 of controlled_runge_kutta::try_step for FSAL steppers does not compile

Hi!
I want to use controlled_runge_kutta with runge_kutta_dopri5 and get the output of try_step regardless of its success, but version 2 of try_step does not compile. The error seems to be a missing argument here, since in the FSAL version there is no try_step method that takes six arguments.

I think the intention was to call version 4 here, so a possible fix would be to use m_dxdt also as dxdt_out, i.e. change the line (controlled_runge_kutta.hpp:685) to

return try_step( system , in , m_dxdt.m_v , t , out , m_dxdt.m_v , dt );

This does compile, and it produces correct results for single steps. However, I did not check if overwriting m_dxdt.m_v has any unwanted side effects for repeated stepping.

Retrieve last-used `dt` value from controlled stepper?

I was wondering if it is possible to retrieve the last-used dt value when using integrate_adaptive with a controlled stepper?
Use case: Calling integrate_adaptive multiple times in succession, i.e. continuing from the time the previous integration stopped at. I want to avoid starting with a dt that is too small or too large for the successive call to integrate_adaptive.

What I mean with "last-used dt" is the last adaptively chosen dt, i.e. ignoring the dt that was chosen to end up at t1; see boost docs:

If stepper is a Controlled Stepper then dt is the initial step size. The actual step size is changed according to error control of the stepper. For the last step, the step size will be reduced to ensure we end exactly at t1. If provided, the observer is called after each time step (and before the first step at t0).

I guess it would be possible to do this via the observer? However, I feel like it would be valuable if it was possible to retrieve it some other way... Is this information stored anywhere or is it overwritten when choosing the dt value to end up at t1?
Digging into the code I only found ode_iterator_base::m_dt, of which I'm not at all sure that it is the relevant value.

Any guidance on this issue would be greatly appreciated. Thanks :)

Long paths result in nonportable tar files

This isn't really a bug, but the CRAN guidelines for R packages do not like paths with long file names since technically they are not supported by the tar format, which has a 100 byte limit. Most (or all?) modern tar implementations don't have a problem with this, so the complaint is annoying, but I was wondering if there is a sensible way to shorten some of the paths in Boost.

The two that cause a problem for me in particular are these:
"boost/numeric/odeint/stepper/generation/generation_controlled_adams_bashforth_moulton.hpp"
"boost/numeric/odeint/stepper/generation/generation_runge_kutta_cash_karp54_classic.hpp"

The "generation" portion is repetitive, so perhaps that could be removed from the file names?

Sorry for the irritating complaint.

EDIT: For clarification, the paths in boost are not longer than 100 bytes, but since they're part of an R package, the full path name is longer.

Inclusion of odeint headers leads to deprecation warnings

When including #include <boost/numeric/odeint.hpp>, I get the following compilation warnings with g++ 10.2/mingw-w64 on Windows:

In file included from C:\Libraries\boost_1_74_0/boost/config/header_deprecated.hpp:18, from C:\Libraries\boost_1_74_0/boost/range/result_iterator.hpp:18, from C:\Libraries\boost_1_74_0/boost/range/metafunctions.hpp:20, from C:\Libraries\boost_1_74_0/boost/range.hpp:19, from C:\Libraries\boost_1_74_0/boost/numeric/odeint/util/resize.hpp:22, from C:\Libraries\boost_1_74_0/boost/numeric/odeint/util/state_wrapper.hpp:26, from C:\Libraries\boost_1_74_0/boost/numeric/odeint/util/ublas_wrapper.hpp:33, from C:\Libraries\boost_1_74_0/boost/numeric/odeint.hpp:25, from ..MyOwnFile.cpp C:\Libraries\boost_1_74_0/boost/range/result_iterator.hpp:20:1: note: '#pragma message: This header is deprecated. Use <boost/range/iterator.hpp> instead.' 20 | BOOST_HEADER_DEPRECATED("<boost/range/iterator.hpp>") | ^~~~~~~~~~~~~~~~~~~~~~~ C:\Libraries\boost_1_74_0/boost/range/const_reverse_iterator.hpp:20:1: note: '#pragma message: This header is deprecated. Use <boost/range/reverse_iterator.hpp> instead.' 20 | BOOST_HEADER_DEPRECATED("<boost/range/reverse_iterator.hpp>") | ^~~~~~~~~~~~~~~~~~~~~~~ C:\Libraries\boost_1_74_0/boost/range/reverse_result_iterator.hpp:20:1: note: '#pragma message: This header is deprecated. Use <boost/range/reverse_iterator.hpp> instead.' 20 | BOOST_HEADER_DEPRECATED("<boost/range/reverse_iterator.hpp>") | ^~~~~~~~~~~~~~~~~~~~~~~

no matching function for call to 'adjust_size_by_resizeability' in boost 1.85

After an update to boost 1.85 we started getting errors about the adjust_size_by_resizeability function. This seems to be an internal boost problem, unrelated to how we actually use boost. We observe this bug both with clang 18.1.4 and gcc 13.2.1. If there's some additional info I could provide to help with this issue, please let me know.

In file included from /usr/include/boost/numeric/odeint/stepper/generation.hpp:28:
In file included from /usr/include/boost/numeric/odeint/stepper/generation/generation_runge_kutta_cash_karp54.hpp:22:
In file included from /usr/include/boost/numeric/odeint/stepper/runge_kutta_cash_karp54.hpp:24:
/usr/include/boost/numeric/odeint/stepper/explicit_error_generic_rk.hpp:156:24: error: no matching function for call to 'adjust_size_by_resizeability'
  156 |             resized |= adjust_size_by_resizeability( m_F[i] , x , typename is_resizeable<deriv_type>::type() );
      |                        ^~~~~~~~~~~~~~~~~~~~~~~~~~~~
/usr/include/boost/numeric/odeint/util/resizer.hpp:31:6: note: candidate function template not viable: no known conversion from 'typename is_resizeable<deriv_type>::type' (aka 'integral_constant<bool, true>') to 'boost::true_type' (aka 'integral_constant<bool, true>') for 3rd argument
   31 | bool adjust_size_by_resizeability( ResizeWrappedState &x , const State &y , boost::true_type )
      |      ^                                                                      ~~~~~~~~~~~~~~~~
/usr/include/boost/numeric/odeint/util/resizer.hpp:43:6: note: candidate function template not viable: no known conversion from 'typename is_resizeable<deriv_type>::type' (aka 'integral_constant<bool, true>') to 'boost::false_type' (aka 'integral_constant<bool, false>') for 3rd argument
   43 | bool adjust_size_by_resizeability( ResizeWrappedState & /* x */ , const State & /* y */ , boost::false_type )
      |      ^                                                                                    ~~~~~~~~~~~~~~~~~

Types error: invalid operands to binary expression ('double' and 'state_type')

I defined new state type as

struct state_type { double x; };

and vector space operators

state_type operator * (const state_type &a, double k) { return {a.x * k}; }
state_type operator * (double k, const state_type &a) { return {a.x * k}; }
state_type operator + (const state_type &a, const state_type &b) { return {a.x + b.x}; }
state_type operator - (const state_type &a, const state_type &b) { return {a.x - b.x}; }
state_type operator / (const state_type &a, const state_type &b) { return {a.x / b.x}; }
state_type abs(state_type const& v) { return {fabs(v.x)}; }

The main function is

int main()
{
    state_type x0 {1.};
    double t0 = 0, t1 = 1, dt = 1e-3;
    using Stepper = runge_kutta_dopri5<state_type, double, state_type, double, vector_space_algebra>;
    auto f = [](state_type const& x, state_type& dx, double t) {
        dx = x;
    };
    integrate_adaptive(make_controlled<Stepper>(1e-10, 1e-10), f, x0, t0, t1, dt);
}

During compilation I am getting the error

/usr/include/boost/numeric/odeint/algebra/default_operations.hpp:443:76: error: invalid operands to binary expression ('double' and 'state_type')
[build] set_unit_value( t3 , abs( get_unit_value( t3 ) ) / ( m_eps_abs + m_eps_rel * ( m_a_x * abs( get_unit_value( t1 ) ) + m_a_dxdt * abs( get_unit_value( t2 ) ) ) ) );

If I add the requested operator

state_type operator + (double k, const state_type &a) { return {a.x + k}; }

the compilation works well. It seems that the formula

set_unit_value( t3 , abs( get_unit_value( t3 ) ) / ( m_eps_abs + m_eps_rel * ( m_a_x * abs( get_unit_value( t1 ) ) + m_a_dxdt * abs( get_unit_value( t2 ) ) ) ) );

is wrong, or type of m_eps_abs is wrong. In a vector space there must not be operator + between elements of state_type and and coefficients. Moreover the operator / between state_type elements is questionable.

Is this still maintained?

There has not been any activity or response to issues or PRs for 4 years. Is this project orphaned?

Change link to Cash-Karp Solvers

This is minimal size request. In

odeint/include/boost/numeric/odeint/stepper/runge_kutta_cash_karp54.hpp

Line 203 in 2bbc186

    
            * <a href="http://en.wikipedia.org/wiki/Cash%E2%80%93Karp_methods">en.wikipedia.org/wiki/Cash-Karp_methods</a>.

the link to wikipedia should be corrected. It is: https://en.wikipedia.org/wiki/Cash%E2%80%93Karp_method

Make OpenMP range algebra buildable with MSVC

Current implementation of openmp_range_algebra uses unsigned type (size_t) as index in for loops that are parallelized by OpenMP. This is not supported by MSVC. Additionally, as MSVC is still stuck with OpenMP 2.0 support, the OpenMP reduction syntax used by norm_inf also does not compile. The following patch fixes the issue; at least of MSVC 14.
0001-Make-OpenMP-range-algebra-buildable-with-MSVC.txt

make_controlled factory bug (setting m_max_dt to eps_rel and meps_rel and m_eps_abs to eps_abs)

Bug submitted here: headmyshoulder/odeint-v2#228

bulirsch_stoer_dense_out doesn't work with multiprecision types

I want to use bulirsch_stoer with multiprecision types. However doing the "obvious" thing with the sample code elliptic_functions.cpp results in a compile-time error.

I'm using:

Linux Fedora 39
g++ (GCC) 13.2.1 20240316 (Red Hat 13.2.1-7)
boost boost-1.81.0-8.fc39.x86_64

It looks like algebra.norm_inf needs to be extended to work with multiprecision types.

Here is the test code elliptic_function.cpp

/*
 * Adapted from sample code elliptic_functions.cpp in repo
 * [email protected]:headmyshoulder/odeint-v2.git
 *
 * Compile with one of
 *   g++ -o ell -DPRECISION=2 elliptic_functions.cpp // OK
 *   g++ -o ell -DPRECISION=3 elliptic_functions.cpp // OK
 *   g++ -o ell -DPRECISION=4 elliptic_functions.cpp // compile error
 *   g++ -o ell -DPRECISION=5 elliptic_functions.cpp // compile error
 */

#include <iostream>
#include <cmath>
#include <array>

#include <boost/numeric/odeint/config.hpp>

#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/stepper/bulirsch_stoer.hpp>
#include <boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp>

#if PRECISION == 3
// This is OK
typedef long double value_type;
#elif PRECISION == 4
// This gives compile failure
#include <boost/multiprecision/float128.hpp>
typedef boost::multiprecision::float128 value_type;
#elif PRECISION == 5
// This gives compile failure
#include <boost/multiprecision/cpp_dec_float.hpp>
typedef boost::multiprecision::cpp_dec_float_50 value_type;
#else
// This is OK
typedef double value_type;
#endif

typedef std::array< value_type , 3 > state_type;

using namespace std;
using namespace boost::numeric::odeint;

/*
 * x1' = x2*x3
 * x2' = -x1*x3
 * x3' = -m*x1*x2
 */

void rhs( const state_type &x, state_type &dxdt, const value_type /*t*/ ) {
  static const value_type m = 51/value_type(100);

  dxdt[0] = x[1]*x[2];
  dxdt[1] = -x[0]*x[2];
  dxdt[2] = -m*x[0]*x[1];
}

int main() {
  bulirsch_stoer_dense_out< state_type, value_type >
    stepper( value_type(1E-9), value_type(1E-9), value_type(1), value_type(0) );

  state_type x1 = { value_type(0), value_type(1), value_type(1) };

  value_type t = 0.0;
  value_type dt = 0.01;

  integrate_adaptive( stepper, rhs, x1, t, value_type(100), dt );
  cout << x1[0] << " " << x1[1] << " " << x1[2] << endl;
}

Compiling this with, e.g.,

g++ -c -DPRECISION=4 elliptic_functions.cpp

gives

In file included from /usr/include/boost/numeric/odeint/algebra/norm_result_type.hpp:20,
                 from /usr/include/boost/numeric/odeint/algebra/range_algebra.hpp:28,
                 from /usr/include/boost/numeric/odeint/stepper/euler.hpp:25,
                 from /usr/include/boost/numeric/odeint.hpp:27,
                 from elliptic_functions.cpp:18:
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp: In instantiation of ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’:
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:47:73:   required from ‘struct boost::numeric::odeint::detail::extract_value_type<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, void>’
/usr/include/boost/numeric/odeint/algebra/norm_result_type.hpp:28:60:   required from ‘struct boost::numeric::odeint::norm_result_type<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, void>’
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64:   required by substitution of ‘template<template<class, long unsigned int <anonymous> > class Array, class T, long unsigned int dim> static typename boost::numeric::odeint::norm_result_type<Array<T, dim> >::type boost::numeric::odeint::array_algebra::norm_inf(const Array<T, dim>&) [with Array = std::array; T = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; long unsigned int dim = 3]’
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:89:32:   required from ‘boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::error(algebra_type&, const State&, const Deriv&, Err&, Time) const [with State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Err = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; value_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; algebra_type = boost::numeric::odeint::array_algebra]’
/usr/include/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp:199:63:   required from ‘boost::numeric::odeint::controlled_step_result boost::numeric::odeint::bulirsch_stoer_dense_out<State, Value, Deriv, Time, Algebra, Operations, Resizer>::try_step(System, const StateIn&, const DerivIn&, time_type&, StateOut&, DerivOut&, time_type&) [with System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); StateIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; DerivIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; StateOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; DerivOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer; time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
/usr/include/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp:339:27:   required from ‘std::pair<_Size, _Size> boost::numeric::odeint::bulirsch_stoer_dense_out<State, Value, Deriv, Time, Algebra, Operations, Resizer>::do_step(System) [with System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer]’
/usr/include/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp:140:23:   required from ‘size_t boost::numeric::odeint::detail::integrate_adaptive(Stepper, System, State&, Time, Time, Time, Observer, boost::numeric::odeint::dense_output_stepper_tag) [with Stepper = boost::numeric::odeint::bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Observer = boost::numeric::odeint::null_observer; size_t = long unsigned int]’
/usr/include/boost/numeric/odeint/integrate/integrate_adaptive.hpp:42:38:   required from ‘size_t boost::numeric::odeint::integrate_adaptive(Stepper, System, State&, Time, Time, Time, Observer) [with Stepper = bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Observer = null_observer; size_t = long unsigned int]’
/usr/include/boost/numeric/odeint/integrate/integrate_adaptive.hpp:83:30:   required from ‘size_t boost::numeric::odeint::integrate_adaptive(Stepper, System, State&, Time, Time, Time) [with Stepper = bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; size_t = long unsigned int]’
elliptic_functions.cpp:66:21:   required from here
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:47:73: error: invalid use of incomplete type ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’
   47 |     typedef typename extract_value_type< typename S::value_type >::type type;
      |                                                                         ^~~~
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:44:8: note: definition of ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’ is not complete until the closing brace
   44 | struct extract_value_type< S , typename boost::enable_if< has_value_type<S> >::type >
      |        ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
In file included from /usr/include/boost/numeric/odeint.hpp:35:
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp: In instantiation of ‘boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::error(algebra_type&, const State&, const Deriv&, Err&, Time) const [with State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Err = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; value_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; algebra_type = boost::numeric::odeint::array_algebra]’:
/usr/include/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp:199:63:   required from ‘boost::numeric::odeint::controlled_step_result boost::numeric::odeint::bulirsch_stoer_dense_out<State, Value, Deriv, Time, Algebra, Operations, Resizer>::try_step(System, const StateIn&, const DerivIn&, time_type&, StateOut&, DerivOut&, time_type&) [with System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); StateIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; DerivIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; StateOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; DerivOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer; time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
/usr/include/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp:339:27:   required from ‘std::pair<_Size, _Size> boost::numeric::odeint::bulirsch_stoer_dense_out<State, Value, Deriv, Time, Algebra, Operations, Resizer>::do_step(System) [with System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; Resizer = boost::numeric::odeint::initially_resizer]’
/usr/include/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp:140:23:   required from ‘size_t boost::numeric::odeint::detail::integrate_adaptive(Stepper, System, State&, Time, Time, Time, Observer, boost::numeric::odeint::dense_output_stepper_tag) [with Stepper = boost::numeric::odeint::bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Observer = boost::numeric::odeint::null_observer; size_t = long unsigned int]’
/usr/include/boost/numeric/odeint/integrate/integrate_adaptive.hpp:42:38:   required from ‘size_t boost::numeric::odeint::integrate_adaptive(Stepper, System, State&, Time, Time, Time, Observer) [with Stepper = bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Observer = null_observer; size_t = long unsigned int]’
/usr/include/boost/numeric/odeint/integrate/integrate_adaptive.hpp:83:30:   required from ‘size_t boost::numeric::odeint::integrate_adaptive(Stepper, System, State&, Time, Time, Time) [with Stepper = bulirsch_stoer_dense_out<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; System = void (*)(const std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>); State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; size_t = long unsigned int]’
elliptic_functions.cpp:66:21:   required from here
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:89:32: error: no matching function for call to ‘boost::numeric::odeint::array_algebra::norm_inf(std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 3>&)’
   89 |         return algebra.norm_inf( x_err );
      |                ~~~~~~~~~~~~~~~~^~~~~~~~~
In file included from /usr/include/boost/numeric/odeint/algebra/algebra_dispatcher.hpp:29,
                 from /usr/include/boost/numeric/odeint/stepper/euler.hpp:27:
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64: note: candidate: ‘template<template<class, long unsigned int <anonymous> > class Array, class T, long unsigned int dim> static typename boost::numeric::odeint::norm_result_type<Array<T, dim> >::type boost::numeric::odeint::array_algebra::norm_inf(const Array<T, dim>&)’
  276 |     static typename norm_result_type< Array< T , dim > >::type norm_inf( const Array< T , dim > &s )
      |                                                                ^~~~~~~~
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64: note:   substitution of deduced template arguments resulted in errors seen above

Compilation issue with custom state type with clang only

In the code below, I am working with a custom state-type (a std::pair in this case). The code compiles with g++, icc, msvc. With clang however, compilation fails. In order to get it to build successfully, I apparently need to define operator+() and operator*() for the state type before including the odeint header. Link to godbolt: https://godbolt.org/z/6eavzdxbf

#include <cstdlib>
#include <utility>

#include <boost/numeric/odeint.hpp> // if included here, compilations fails with clang++

using state_type = std::pair<double, double>;
using deriv_type = state_type;
using time_type = double;

state_type operator+(state_type a, state_type const& b);
state_type operator*(state_type a, time_type b);
state_type operator*(time_type b, state_type a);

//~ #include <boost/numeric/odeint.hpp> // if included here, no error!

state_type operator+(state_type a, state_type const& b) {
	a.first += b.first;
	a.second += b.second;
	return a;
}

state_type operator*(state_type a, time_type b) {
	a.first *= b;
	a.second *= b;
	return a;
}

state_type operator*(time_type b, state_type a) {
	return a * b;
}

void harmonic_oscillator(const state_type& x, deriv_type& dxdt, time_type) {
	static auto constexpr gam = 0.1;

	dxdt.first = x.second;
	dxdt.second = -x.first - gam*x.second;
}

namespace bno = boost::numeric::odeint;
using stepper_type = bno::runge_kutta_dopri5<state_type, double, deriv_type, time_type, bno::vector_space_algebra>;

int main() {
	state_type x0{1.0, 0.0};
	
	// just to see if the operators work...
	x0 * 1.;
	2. * x0;
	x0 + x0;
	
	stepper_type stepper;
	time_type const dt = 0.005;
	for (time_type t = 0; t < 10; t += dt) {
		stepper.do_step(harmonic_oscillator, x0, t, dt);
		std::cout << t << '\t' << x0.first << '\t' << x0.second << '\n';
	}
	
	return EXIT_SUCCESS;
}

Extra ; warnings

At least these files produce extra ; warnings with GCC 13 and -Wpedantic in boost 1.84:
boost/numeric/odeint/stepper/detail/adaptive_adams_bashforth_moulton.hpp
boost/numeric/odeint/stepper/detail/adaptive_adams_coefficients.hpp
boost/numeric/odeint/stepper/detail/pid_step_adjuster.hpp
boost/numeric/odeint/stepper/detail/controlled_adams_bashforth_moulton.hpp

Example:

boost/numeric/odeint/stepper/adaptive_adams_bashforth_moulton.hpp:99:6: varoitus: ylimääräinen ”;” [-Wpedantic]
   99 |     };
      |      ^
      |      -

infinite loop in integrate_adaptive/odeint

I am running into an infinite loops applying integrate_adaptive in the dense output version (odeint/integrate/detail/integrate_adaptive.hpp l. 122).

That's an easy test case replicating the issue:

auto f = []( const double& /*y*/, double& dy_dx, double /*x*/ ) {
        dy_dx = haas::math::numeric::constant<double>::two_pi() / -35.0;
    };
    typedef boost::numeric::odeint::runge_kutta_dopri5<state_type> error_stepper_type;
    auto stepper = boost::numeric::odeint::make_dense_output<error_stepper_type>( 1.0e-6, 1.0e-5);
    // Initial conditions
    state_type state1 = 0.0;
    // Solve ODE
    boost::numeric::odeint::integrate_adaptive( stepper, f, state1, 0.0, -23.8, -23.8 / 100 );

The problem seems to rise from a mismatch between the new initialization of stepper st ( l. 144) and the while loop condition (l. 135).

I would suggest to change l. 135-137 from

        while( less_eq_with_sign( static_cast<Time>(st.current_time() + st.current_time_step()) ,
               end_time ,
               st.current_time_step() ) )

into:

                    while( less_eq_with_sign(
                        st.current_time_step(), static_cast<Time>( end_time - st.current_time() ),
                        st.current_time_step() ) )

Comments very welcome!

Update Implicit-System documentation

Hi there,

looking at the code, I think, that this line in the doc is not up-to-date.

odeint/doc/concepts/implicit_system.qbk

Line 40 in 2bbc186

    
           [[Calculate ['A := df/dx (x,t)]] [`sys.second( x , jacobi , t )`] [`void`] [Calculates the Jacobian of /f/ at /x/,/t/, the result is stored into `jacobi`] ]

I compared it to

odeint/include/boost/numeric/odeint/stepper/rosenbrock4.hpp

Line 187 in 2bbc186

jacobi_func( x , m_jac.m_v , t , m_dfdt.m_v );

If this is the case, I can update it and create a pull request.

Error using odeint with float128

There's a compile time error using odeint with the float128 introduced
somewhere between boost 1.66 and boost 1.69 (or possible between g++
8.3.1 and g++ 9.2.1).

The following copy compiles and runs successfully (integrating a
simple harmonic oscillator) with PRECISION set to 1, 2, 3, and 4 (for
float, double, long double, and float128) on a Fedora Linux 29 system
with

boost 1.66
g++ 8.3.1

However after upgrading to Fedora Linux 30, i.e.,

boost 1.69
g++ 9.2.1

the float128 precision version fails to compile. Here is the test code:

// compile with
//   g++ -O3            -DPRECISION=1 -o bug1 bug.cpp
//   g++ -O3            -DPRECISION=2 -o bug2 bug.cpp
//   g++ -O3            -DPRECISION=3 -o bug3 bug.cpp
//   g++ -O3 -lquadmath -DPRECISION=4 -o bug4 bug.cpp

#include <iostream>
#include <iomanip>
#include <vector>
#include <array>
#include <boost/multiprecision/float128.hpp>
#include <boost/numeric/odeint.hpp>

#if !defined(PRECISION)
#define PRECISION 4
#endif

namespace ode = boost::numeric::odeint;
#if PRECISION == 1
typedef float real;
#elif PRECISION == 2
typedef double real;
#elif PRECISION == 3
typedef long double real;
#elif PRECISION == 4
typedef boost::multiprecision::float128 real;
#else
typedef double real
#endif

typedef std::array<real, 2> point;

class Force {
public:
  Force() {}
  void operator() (const point& p, point &d, const real /*t*/) const {
d[0] = p[1];
d[1] = -p[0];
  }
};

class ParticleFollow {
public:
  static void follow(Force& sys, const point& p0, real t0, real ds,
		 real tmax, std::vector<point>& points) {
ode::result_of::make_dense_output
  < ode::runge_kutta_dopri5< point, real > >::type integrator =
  ode::make_dense_output(real(1.0e-8), real(0.01*1.0e-8),
			 ode::runge_kutta_dopri5< point, real >() );
integrator.initialize(p0, t0, real(1e-2));
integrator.do_step(sys);
int n = 1, i = 1;
point out;
while (true) {
  real tout = i * ds;
  while (integrator.current_time() < tout) {
    integrator.do_step(sys);
    ++n;
  }
  integrator.calc_state(tout, out);
  points.push_back(out);
  ++i;
  if (tout > tmax) break;
}
  }
};

int main() {
  std::vector<point> result;
  point p0;
  p0[0] = real(1); p0[1] = real(0);
  Force sys;
  ParticleFollow::follow(sys, p0, 0, 0.1, 6.29, result);
  std::cout << std::fixed << std::setprecision(4);
  for (size_t i = 0; i < result.size(); ++i)
std::cout << result[i][0] << " " << result[i][1] << "\n";
}

and here is the output from the compiler when PRECISION=4:

In file included from /usr/include/boost/numeric/odeint/algebra/norm_result_type.hpp:20,
	     from /usr/include/boost/numeric/odeint/algebra/range_algebra.hpp:28,
	     from /usr/include/boost/numeric/odeint/stepper/euler.hpp:25,
	     from /usr/include/boost/numeric/odeint.hpp:27,
	     from bug.cpp:12:
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp: In instantiation of ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’:
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:47:73:   required from ‘struct boost::numeric::odeint::detail::extract_value_type<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, void>’
/usr/include/boost/numeric/odeint/algebra/norm_result_type.hpp:28:60:   required from ‘struct boost::numeric::odeint::norm_result_type<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, void>’
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64:   required by substitution of ‘template<template<class, long unsigned int <anonymous> > class Array, class T, long unsigned int dim> static typename boost::numeric::odeint::norm_result_type<Array<T, dim> >::type boost::numeric::odeint::array_algebra::norm_inf(const Array<T, dim>&) [with Array = std::array; T = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; long unsigned int dim = 2]’
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:89:40:   required from ‘boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::error(boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::algebra_type&, const State&, const Deriv&, Err&, Time) const [with State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Err = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::algebra_type = boost::numeric::odeint::array_algebra]’
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:768:20:   required from ‘boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step(System, const StateIn&, const DerivIn&, boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type&, StateOut&, DerivOut&, boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type&) [with System = Force; StateIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; DerivIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; StateOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; DerivOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; ErrorStepper = boost::numeric::odeint::runge_kutta_dopri5<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; ErrorChecker = boost::numeric::odeint::default_error_checker<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>; StepAdjuster = boost::numeric::odeint::default_step_adjuster<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; Resizer = boost::numeric::odeint::initially_resizer; boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
/usr/include/boost/numeric/odeint/stepper/dense_output_runge_kutta.hpp:338:17:   required from ‘std::pair<typename Stepper::stepper_type::time_type, typename Stepper::stepper_type::time_type> boost::numeric::odeint::dense_output_runge_kutta<Stepper, boost::numeric::odeint::explicit_controlled_stepper_fsal_tag>::do_step(System) [with System = Force; Stepper = boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >, boost::numeric::odeint::default_error_checker<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>; typename Stepper::stepper_type::time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
bug.cpp:51:27:   required from here
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:47:73: error: invalid use of incomplete type ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’
   47 |     typedef typename extract_value_type< typename S::value_type >::type type;
  |                                                                         ^~~~
/usr/include/boost/numeric/odeint/algebra/detail/extract_value_type.hpp:44:8: note: definition of ‘struct boost::numeric::odeint::detail::extract_value_type<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, void>’ is not complete until the closing brace
   44 | struct extract_value_type< S , typename boost::enable_if< has_value_type<S> >::type >
  |        ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
In file included from /usr/include/boost/numeric/odeint.hpp:35,
	     from bug.cpp:12:
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp: In instantiation of ‘boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::error(boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::algebra_type&, const State&, const Deriv&, Err&, Time) const [with State = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Deriv = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Err = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; Time = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Value = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; Algebra = boost::numeric::odeint::array_algebra; Operations = boost::numeric::odeint::default_operations; boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::value_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>; boost::numeric::odeint::default_error_checker<Value, Algebra, Operations>::algebra_type = boost::numeric::odeint::array_algebra]’:
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:768:20:   required from ‘boost::numeric::odeint::controlled_step_result boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::try_step(System, const StateIn&, const DerivIn&, boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type&, StateOut&, DerivOut&, boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type&) [with System = Force; StateIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; DerivIn = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; StateOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; DerivOut = std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>; ErrorStepper = boost::numeric::odeint::runge_kutta_dopri5<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; ErrorChecker = boost::numeric::odeint::default_error_checker<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>; StepAdjuster = boost::numeric::odeint::default_step_adjuster<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >; Resizer = boost::numeric::odeint::initially_resizer; boost::numeric::odeint::controlled_runge_kutta<ErrorStepper, ErrorChecker, StepAdjuster, Resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>::time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
/usr/include/boost/numeric/odeint/stepper/dense_output_runge_kutta.hpp:338:17:   required from ‘std::pair<typename Stepper::stepper_type::time_type, typename Stepper::stepper_type::time_type> boost::numeric::odeint::dense_output_runge_kutta<Stepper, boost::numeric::odeint::explicit_controlled_stepper_fsal_tag>::do_step(System) [with System = Force; Stepper = boost::numeric::odeint::controlled_runge_kutta<boost::numeric::odeint::runge_kutta_dopri5<std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >, boost::numeric::odeint::default_error_checker<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::numeric::odeint::array_algebra, boost::numeric::odeint::default_operations>, boost::numeric::odeint::default_step_adjuster<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off> >, boost::numeric::odeint::initially_resizer, boost::numeric::odeint::explicit_error_stepper_fsal_tag>; typename Stepper::stepper_type::time_type = boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>]’
bug.cpp:51:27:   required from here
/usr/include/boost/numeric/odeint/stepper/controlled_runge_kutta.hpp:89:40: error: no matching function for call to ‘boost::numeric::odeint::array_algebra::norm_inf(std::array<boost::multiprecision::number<boost::multiprecision::backends::float128_backend, boost::multiprecision::et_off>, 2>&)’
   89 |         return algebra.norm_inf( x_err );
  |                                        ^
In file included from /usr/include/boost/numeric/odeint/algebra/algebra_dispatcher.hpp:29,
	     from /usr/include/boost/numeric/odeint/stepper/euler.hpp:27,
	     from /usr/include/boost/numeric/odeint.hpp:27,
	     from bug.cpp:12:
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64: note: candidate: ‘template<template<class, long unsigned int <anonymous> > class Array, class T, long unsigned int dim> static typename boost::numeric::odeint::norm_result_type<Array<T, dim> >::type boost::numeric::odeint::array_algebra::norm_inf(const Array<T, dim>&)’
  276 |     static typename norm_result_type< Array< T , dim > >::type norm_inf( const Array< T , dim > &s )
  |                                                                ^~~~~~~~
/usr/include/boost/numeric/odeint/algebra/array_algebra.hpp:276:64: note:   substitution of deduced template arguments resulted in errors seen above

controlled_runge_kutta complex State error with Eigen.

I get an error in Line 89:
return algebra.norm_inf( x_err ); cant convert std::complex to double
i attempted to work around this by inserted a abs()
return abs(algebra.norm_inf( x_err ));
is this correct?

My state type is Eigen::VectorXcd
my Value type is double.
I'm using
#include <boost/numeric/odeint/external/eigen/eigen.hpp>

Thank you in advance and best Regards!

set up an appveyor build

Incompatibility with Thrust v1.10.0 (CUDA 11.2)

Using Thrust v1.10.0 (included in CUDA v11.2), compiling the thrust example (test_external/thrust) fails with

../../include/boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp(97): error: class template "boost::numeric::odeint::algebra_dispatcher" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp(96): error: class template "boost::numeric::odeint::operations_dispatcher" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(176): error: class template "boost::numeric::odeint::is_resizeable" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(177): error: class template "boost::numeric::odeint::resize_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(177): error: class template "boost::numeric::odeint::resize_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(178): error: class template "boost::numeric::odeint::same_size_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(178): error: class template "boost::numeric::odeint::same_size_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(179): error: class template "boost::numeric::odeint::copy_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(179): error: class template "boost::numeric::odeint::copy_impl" has already been defined

../../include/boost/numeric/odeint/external/thrust/thrust_resize.hpp(179): error: class template "boost::numeric::odeint::copy_impl" has already been defined

The problem arises, because in Thrust v1.10.0 the types thrust::cpp::vector, thrust::omp::vector, thrust::tbb::vector and thrust::cuda::vector are all aliases of the same type thrust::detail::vector_base.

Add Multi-GPU Functionality

Kindly add the functionality to do the integration using multi-GPU platforms (on a single compute node) and/or distributed across multiple nodes of a server.

Kindly also add some relevant examples.

Thank you

boostorg / odeint Goto Github PK

odeint's Introduction

odeint's People

Contributors

Stargazers

Watchers

Forkers

odeint's Issues

Recommend Projects

Recommend Topics

Recommend Org