BOOST:ODEINT 迭代突然停止

BOOST:ODEINT Sudden Iteration stop

我是 C++ 世界的新手,我在使用 boost 库时遇到了一些问题。在我的问题中,我想用 5 个方程求解 ODE 系统。这不是一个棘手的问题。作为迭代方法,我同时使用了 integreate(rhs, x0, t0, tf, size_step, write_output)integreate_adaptive(stepper, sys, x0, t0, tf, size_step, write_output)。这两种方法实际上都对方程进行了积分,但给出了无意义的结果,几乎随机地将步长从 0.001 更改为 5。方程式和数据是正确的。我该怎么做才能解决这个问题?这是代码:

#include <iostream>
#include <vector>
#include <boost/numeric/odeint.hpp>
#include <fstream>
#include <boost/array.hpp>

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

//DATA
double Lin = 20000;                             // kg/h
double Gdry = 15000;                            // kg/h
double P = 760;                                 // mmHg
double TinH2O = 50;                             // °C
double ToutH2O = 25;                            // °C
double Tinair = 20;                             // °C
double Z = 0.5;                                 // relative humidity
double Cu = 0.26;                               // kcal/kg*K
double CpL = 1;                                 // kcal/kg*K
double DHev = 580;                              // kcal/kg
double hga = 4000;                              // kcal/m/h/K
double hla = 30000;                             // kcal/m/h/K
double A = -49.705;                             // Pev 1st coeff    mmHg vs °C
double B = 2.71;                                // Pev 2nd coeff    mmHg vs °C
double Usair = 0.62*(A + B*Tinair) / P;
double Uair = Z*Usair;
double Kua = hga / Cu;
double L0 = 19292;                              // kg/h


typedef vector< double > state_type;
vector <double> pack_height;
vector <double> Umidity;
vector <double> T_liquid;
vector <double> T_gas;
vector <double> Liquid_flow;
vector <double> Gas_flow;

void rhs(const state_type& x , state_type& dxdt , const double z )
{// U Tl Tg L G

double Ti = (hla*x[1] + hga*x[2] + Kua*DHev*(x[0] - 0.62*A / P)) / (hla + hga + Kua*DHev*0.62*B / P);
double Ui = 0.62*(A + B*Ti) / P;

dxdt[0] = Kua*(Ui - x[0]) / Gdry / 100;
dxdt[1] = hla*(x[1] - Ti) / x[3] / CpL / 100;
dxdt[2] = hga*(Ti - x[2]) / Gdry / Cu / 100;
dxdt[3] = Kua*(Ui - x[0]) / 100;
dxdt[4] = Kua*(Ui - x[0]) / 100;
}

void write_output(const state_type& x, const double z)
{
pack_height.push_back(z);
Umidity.push_back(x[0]);
T_liquid.push_back(x[1]);
T_gas.push_back(x[2]);
Liquid_flow.push_back(x[3]);
Gas_flow.push_back(x[4]);

cout << z << "      " << x[0] << "      " << x[1] << "      " << x[2] << "      " << x[3] << "      " << x[4] << endl;
}

int main()
{
state_type x(5);
x[0] = Uair;
x[1] = ToutH2O;
x[2] = Tinair;
x[3] = L0;
x[4] = Gdry;

double z0 = 0.0;
double zf = 5.5;
double stepsize = 0.001;
integrate( rhs , x , z0 , zf , stepsize , write_output );

return 0;
}

这是我从提示中得到的最终结果:

    0         0.00183349      25           20           19292      15000
    0.001     0.00183356      25           20           19292      15000
    0.0055    0.0018339       25.0002      20.0001      19292      15000
    0.02575   0.00183542      25.001       20.0007      19292      15000
    0.116875  0.00184228      25.0046      20.003       19292.1    15000.1
    0.526938  0.00187312      25.0206      20.0135      19292.6    15000.6
    2.37222   0.00201203      25.0928      20.0608      19294.7    15002.7
    5.5       0.00224788      25.2155      20.142       19298.2    15006.2

只有第一次迭代具有正确要求的步长..显然解决方案不正确..我该怎么办?先感谢您。 :)

如果你read the documentation,那么你会发现常数step-size例程是integrate_constintegrate_n_steps,或者可能是integrate_adaptive加上non-controlled步进器。

integrate 的简短调用使用具有自适应步长的标准 dopri5 步进器,因此改变步长并不奇怪。您可以使用步进器的密集输出在等距时间插值。