Eu tenho um sistema de 4 equações acopladas para resolver e um parâmetro Gamma [i] para iterar. Desde que eu sou novo em C ++, meu código é muito rudimentar. Se parece sofisticado e elegante em certas partes, é só porque adaptei o código do autor do odeint. :)

Esta questão está relacionada com (http://stackoverflow.com/questions/12060111/using-odeint-function-definition/12066958#comment16253600_12066958) mas não exatamente o mesmo. Por favor, não apague isso. :(

Perguntas foram inseridas entre as linhas de código.

#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;

}
};

Questão 1:

No meu programa original, eu tinha algo parecido para canalizar a saída para um arquivo 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;
}
   }

Como eu adapto o streaming_observer para fazer o que meu save () faz? Basicamente, eu quero gerar arquivos .csv para cada iteração i. Neste ponto, estou fazendo do jeito feio, ou seja, compilando tudo, abrindo um prompt de comando do Windows e, em seguida, canalizando a saída do exe para um arquivo de texto. Isso gera um arquivo grande com todas as iterações jogadas lá.

Isso se torna muito doloroso para analisar um grande número de iterações.

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");

Pergunta # 2: Os resultados que eu obtenho aqui não concordam com o que eu recebo com o que eu recebo usando um método simples explícito de Euler. Por isso, gostaria de ver os coeficientes RK4 (de preferência despejá-los em um arquivo) ou os passos intermediários. Como posso obter essa informação?

//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 ) );
}
    }

Obrigado por todos aqueles que ajudaram!

Edit: Existe alguma maneira de obter uma estimativa de erro em execução? Note que u [0] * u [0] + u [1] * u [1] = 1 em todos os momentos.