2012-08-27 72 views
2

我有一個系統的4個耦合方程來解決和參數Gamma [i]來迭代。由於我對C++相當陌生,因此我的代碼是非常基本的。如果它在某些部分看起來複雜和優雅,那只是因爲我修改了odeint作者的代碼。 :)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做我救什麼()不?基本上,我要生成的.csv每個迭代i文件。在這一點上,我做它的醜陋的方式,即編譯所有內容,打開一個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: 結果我到這裏不就是我與我開始使用一個簡單明確的歐拉法得到同意。因此,我希望看到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。

回答

2

問題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。

+0

它應該是count_中的sprintf(fn,「%s_%04lu.csv」,count);如果它是count_,則會出現訪問衝突錯誤。如果沒有,count似乎沒有被定義。對不起,我是一個總noob。 –

+0

你是否正確初始化狀態類型?它是一個向量,需要手動定義大小(在構造函數中或通過push_back()的resize())。 – headmyshoulder

+0

是的,我做到了。編譯器抱怨「無法確定函數模板std :: count的目標是哪個實例」。 –