3

我试图实现一个简单的 ode 来测试 boost.odeint 是否支持 boost.units 的使用。但是，我的示例在编译时失败。是我的代码，还是 boost.odeint 不支持 boost.Units 进行维度分析？

#include<boost/units/systems/si.hpp>
#include<boost/numeric/odeint.hpp>

typedef boost::units::quantity<boost::units::si::length,double> state_type;
typedef boost::units::quantity<velocity,double> dev_type;
typedef boost::units::quantity<boost::units::si::time,double> time_type;


using namespace boost::units::si;

void exponential_decay(const state_type &x, dev_type &dxdt, time_type t){
    dxdt = -0.0001*x/second;
}

int main(){
    state_type startValue = 10*meter;

    using namespace boost::numeric::odeint;
    auto steps = integrate(exponential_decay, startValue, 0.0*second,
            10.0*second, 0.1*second);

    return 0;
}
4

2 回答 2

2

因此，在最初未能找到 state_type、代数和操作的工作集¹之后，@Arash 提供了缺失的链接。

但是，为了完整起见，您无需携带重型机械 ( fusion_algebra.hpp)。这些适用于更一般的情况，例如状态类型不是单个值的情况。

在这个简单的例子中，真正需要的只是指定代数，而不是使用默认值。这涉及声明 a stepper_type，例如：

using stepper_type = runge_kutta4<length_type, double, velocity_type,
                            time_type, vector_space_algebra>;

接下来，我们选择一个允许我们提供该步进器的迭代算法重载：

Live On Coliru

#include <boost/numeric/odeint.hpp>
#include <boost/units/systems/si.hpp>

typedef boost::units::quantity<boost::units::si::length, double> length_type;
typedef boost::units::quantity<boost::units::si::velocity, double> velocity_type;
typedef boost::units::quantity<boost::units::si::time, double> time_type;

namespace si = boost::units::si;

void exponential_decay(const length_type &x, velocity_type &dxdt, time_type /*t*/) { dxdt = -0.0001 * x / si::second; }

int main() {
    using stepper_type = boost::numeric::odeint::runge_kutta4<
        length_type, double, velocity_type, time_type, boost::numeric::odeint::vector_space_algebra>;

    length_type startValue = 10 * si::meter;
    auto steps = integrate_const(
            stepper_type{}, exponential_decay,
            startValue, 0.0 * si::second, 10.0 * si::second,
            0.1 * si::second);

    std::cout << "Steps: " << steps << "\n";
}

印刷

Steps: 100

¹ 仅查看http://www.boost.org/doc/libs/1_66_0/libs/numeric/odeint/doc/html/boost_numeric_odeint/odeint_in_detail/state_types__algebras_and_operations.html和http://www.boost.org/doc /libs/1_66_0/doc/html/boost_units/Quantities.html

于 2018-02-28T00:32:41.793 回答
1

用于boost::fusion解决问题。

#include <iostream>
#include <vector>

#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/algebra/fusion_algebra.hpp>
#include <boost/numeric/odeint/algebra/fusion_algebra_dispatcher.hpp>

#include <boost/units/systems/si/length.hpp>
#include <boost/units/systems/si/time.hpp>
#include <boost/units/systems/si/velocity.hpp>
#include <boost/units/systems/si/acceleration.hpp>
#include <boost/units/systems/si/io.hpp>

#include <boost/fusion/container.hpp>

using namespace std;
using namespace boost::numeric::odeint;
namespace fusion = boost::fusion;
namespace units = boost::units;
namespace si = boost::units::si;

typedef units::quantity< si::time , double > time_type;
typedef units::quantity< si::length , double > length_type;
typedef units::quantity< si::velocity , double > velocity_type;
typedef units::quantity< si::acceleration , double > acceleration_type;
typedef units::quantity< si::frequency , double > frequency_type;

typedef fusion::vector< length_type  > state_type;
typedef fusion::vector< velocity_type > deriv_type;

void exponential_decay( const state_type &x , deriv_type &dxdt , time_type t )
{
    fusion::at_c< 0 >( dxdt ) = (-0.0001/si::second)* fusion::at_c< 0 >( x );
}

void observer(const state_type &x,const time_type &t )
{
}

int main( int argc , char**argv )
{
    typedef runge_kutta_dopri5< state_type , double , deriv_type , time_type > stepper_type;

    state_type startValue( 1.0 * si::meter );

    auto steps=integrate_adaptive(
        make_dense_output( 1.0e-6 , 1.0e-6 , stepper_type()),
        exponential_decay,
        startValue , 0.0 * si::second , 100.0 * si::second , 0.1 * si::second , observer);

    std::cout<<"steps: "<<steps<<std::endl;

    return 0;
}
于 2018-02-28T00:53:58.750 回答