numerical-methods - odeint streaming observer 及相关问题

Question

我有一个由 4 个耦合方程组成的系统要求解，还有一个要迭代的参数 Gamma[i]。由于我对 C++ 很陌生，所以我的代码非常简陋。如果它在某些地方看起来精致优雅，那只是因为我改编了odeint作者的代码。 :)

这个问题与 (http://stackoverflow.com/questions/12060111/using-odeint-function-definition/12066958#comment16253600_12066958) 有关，但不完全相同。请不要删除它。 :(

代码行之间插入了问题。

#include <iostream>
#include <iterator>
#include <algorithm>
#include <boost/numeric/odeint.hpp>
#include <cmath>
#include <vector>
#include <fstream>
#include <iomanip>

using namespace std;
using namespace boost::numeric::odeint;
class NLI_class {
private:
    double gamma;
public:
 NLI_class (double r) : gamma(r) {} 

 void operator()( vector<double> &u , vector<double> &du , double z ) {
            du[0] = u[0]*u[1]*cos(u[3]); //u1
            du[1] = -u[0]*u[0]*cos(u[3]); //u2
            du[2] = gamma * (2/(u[0]*u[0]) - 1/(u[1]*u[1])); //theta
            du[3] = gamma * (1.0/(u[0]*u[0])); //phi1
            du[4] = gamma * (1.0/(u[1]*u[1])); //phi2;

}
};

问题 #1:

在我原来的程序中，我有类似这样的东西将输出通过管道传输到 csv 文件:

 inline void save(vector<double>& v, string filename)
  {
ofstream output(filename);
for(int i=0;i<v.size();++i){ 
    output << setprecision(64) << v[i] << endl;
}
   }

我如何调整 streaming_observer 来完成我的 save() 的功能？基本上，我想为每次迭代 i 生成 .csv 文件。在这一点上，我正在以丑陋的方式进行操作，即编译所有内容，打开 Windows 命令提示符，然后将 exe 输出通过管道传输到文本文件。这会生成一个大文件，其中包含所有迭代。

分析大量迭代变得非常痛苦。

struct streaming_observer {

 std::ostream &m_out;
 streaming_observer( std::ostream &out ) : m_out( out ) {}

 void operator()( const vector<double> &x , double t ) const
 {
      m_out << t;
      for( size_t i=0 ; i < x.size() ; ++i )
          m_out << "\t" << x[i];
      m_out << "\n";
 }
};




    int main(){

vector<double> x( 5 );
vector<double> Gamma;
vector<double>delta;
const double  pi=acos(-1.0); 
short delta_n=5;
const double delta_step=(2*pi)/delta_n;
const double dz = 0.01;
const double  zeta = 3.0;
const double  theta_initial=0.0;
const double  u20=tanh(zeta); 
const double  u10=sqrt(1.0-(u20*u20)); 

double d=0.0;
double G=0.0;

for(int i=0;i<=delta_n;i++){
    //When i=0, the d=0.0 and G=0.0 are pushed into the vector.
    delta.push_back(d);  
    Gamma.push_back(G);
    // Compute delta and Gamma
    d=d+delta_step; 
    G=-u10*u10*u20*sin(theta_initial+d);
}

save(delta,"delta.csv");
save(Gamma,"Gamma.csv");

问题#2:我在这里得到的结果与我使用简单显式 Euler 方法得到的结果不一致。因此，我希望看到 RK4 系数(最好将它们转储到文件中)或中间步骤。我怎样才能得到这些信息？

//Numeric Integration
    for (unsigned i = 0; i < Gamma.size(); ++i) {
        x[0] = u10;
        x[1] = u20;
        x[2] = 0.0;
        x[3] = 0.0;
        x[4] = 0.0;

        NLI_class nli_obj(Gamma[i]);
        integrate_const( runge_kutta4< vector<double > >(), nli_obj, x , 0.0 , 3.0 , dz,streaming_observer( std::cout ) );
}
    }

感谢所有帮助过的人!

编辑:有什么方法可以获得运行错误估计吗？注意 u[0]*u[0]+u[1]*u[1]=1 在任何时候。

Answer 1

问题#1:

我不明白你需要什么样的输出。但是如果你想在每次迭代后写出结果，你可以像这样实现一个输出观察器:

struct output_observer
{
    string filename_;
    size_t count_;
    output_observer( const string &filename ) : filename_( filename ) , count_( 0 ) { }
    void operator()( const state_type &x , time_type dt )
    {
        char fn[512] = "";
        sprintf( fn , "%s_%04lu.csv" , filename_.c_str() , count_ );
        ofstream fout( fn );
        for( size_t i=0 ; i<x.size() ; ++i ) fout << x[i] << "\n";
        ++count_;
    }
};

您可以简单地应用这个观察者

integrate_const( runge_kutta4< vector<double > >() , nli_obj , x ,
    0.0 , 3.0 , dz , output_observer( "filename" ) );

这是所需的功能吗？

问题 #2:

看不到runge_kutta4的中间e步。系数是经典 Runge-Kutta 方法的标准系数:http://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods

问题#3:

odeint 有几个错误步进器，它们估计一步中产生的错误。您可以使用例如 Runge_Kutta Cash Karp 算法；

runge_kutta_cash_karp54< state_type > rk;
state_type xerr;
rk.do_step( nli_obj , x , t , xerr );

一步一步估计错误并将错误结果写入xerr。