I wrote a program in c++ that should solve differential equations. The problem is, it seems like it does not work well with ROOT. It compilates fine, but when I execute, this is what I get:
*** Break *** segmentation violation
===========================================================
There was a crash.
This is the entire stack trace of all threads:
===========================================================
#0 0x00007fc28193984a in __GI___waitpid (pid=7730, stat_loc=stat_loc
entry=0x7ffffe4ae000, options=options
entry=0) at ../sysdeps/unix/sysv/linux/waitpid.c:31
#1 0x00007fc2818b2ffb in do_system (line=<optimized out>) at ../sysdeps/posix/system.c:148
#2 0x00007fc2831d0954 in TUnixSystem::StackTrace() () from /usr/lib/root/libCore.so
#3 0x00007fc2831d29ec in TUnixSystem::DispatchSignals(ESignals) () from /usr/lib/root/libCore.so
#4 <signal handler called>
#5 0x0000000000405a8a in Runge_Kutta::Passo(double, VettoreLineare&, double) ()
#6 0x0000000000403b8a in main ()
===========================================================
The lines below might hint at the cause of the crash.
If they do not help you then please submit a bug report at
http://root.cern.ch/bugs. Please post the ENTIRE stack trace
from above as an attachment in addition to anything else
that might help us fixing this issue.
===========================================================
#5 0x0000000000405a8a in Runge_Kutta::Passo(double, VettoreLineare&, double) ()
#6 0x0000000000403b8a in main ()
===========================================================
This is my program
equazione_differenziale.c
#include "equazione_differenziale.h"
EqDifferenzialeBase :: EqDifferenzialeBase (FunzioneBase* f) {
_f=f;
};
Eulero :: Eulero (FunzioneBase*f) : EqDifferenzialeBase(f) {};
Runge_Kutta :: Runge_Kutta (FunzioneBase* f) : EqDifferenzialeBase(f) {};
VettoreLineare Eulero :: Passo (double t, VettoreLineare& x, double h) {
VettoreLineare vec(x.GetN());
vec=x+_f->Eval(t,x)*h;
return vec;
};
Protone :: Protone (double m, double q, double E0, double f, double lambda){
_m=m;
_q=q;
_E0=E0;
_f=f;
_lambda=lambda;
};
VettoreLineare Protone::Eval(double t, const VettoreLineare& v) const{
VettoreLineare y(v.GetN());
for(int i=0; i<v.GetN()/2; i++){
y.SetComponent(i, v.GetComponent(v.GetN()/2+i));
y.SetComponent(i+v.GetN()/2, (-1.)*(_q/_m)*_E0*sin(2*M_PI*(v.GetComponent(i)/_lambda)-2*M_PI*_f*t));
};
return y;
};
VettoreLineare Runge_Kutta::Passo(double t, VettoreLineare& v, double h){
VettoreLineare k1=_f->Eval(t,v);
VettoreLineare k2=_f->Eval(t+h/2.,v+k1*(h/2.));
VettoreLineare k3=_f->Eval(t+h/2.,v+k2*(h/2.));
VettoreLineare k4=_f->Eval(t+h,v+k3*h);
VettoreLineare y=v+(k1+k2*2.+k3*2.+k4)*(h/6.);
return y;
};
equazione_differenziale.h
#ifndef equazione_differenziale_h_
#define equazione_differenziale_h_
#include "Vettore.h"
#include <iostream>
#include <cmath>
class FunzioneBase {
public:
virtual VettoreLineare Eval(double t, const VettoreLineare& v) const=0;
};
class Protone: public FunzioneBase {
private:
double _m,_q,_E0,_f,_lambda;
public:
Protone(double m, double q, double E0, double f, double lambda);
virtual VettoreLineare Eval(double t, const VettoreLineare& v) const;
};
class EqDifferenzialeBase {
protected:
FunzioneBase* _f;
public:
EqDifferenzialeBase (FunzioneBase* f);
virtual VettoreLineare Passo (double t, VettoreLineare& x, double h)=0;
};
class Eulero : public EqDifferenzialeBase {
public:
Eulero (FunzioneBase* f);
virtual VettoreLineare Passo (double t, VettoreLineare& x, double h);
};
class Runge_Kutta: public EqDifferenzialeBase {
protected:
FunzioneBase* _f;
public:
Runge_Kutta (FunzioneBase* f);
virtual VettoreLineare Passo(double t, VettoreLineare& v, double h);
};
#endif
Vettore.h
#ifndef vettore_h_
#define vettore_h_
#include <iostream>
#include <cmath>
#include <fstream>
class Vettore {
protected:
unsigned int _N;
double * _v;
void Quicksort(unsigned int primo, unsigned int ultimo);
void Scambia (int a, int b);
public:
Vettore ();
Vettore (int N);
Vettore (int N, char* nomefile);
Vettore (const Vettore& v);
virtual void SetComponent (int i, double x);
void AddComponent (double x);
double GetComponent (int i) const;
void Print () const;
void Print (char* nomefile) const;
void Sort();
virtual int GetN() const;
Vettore& operator=(const Vettore & vetty);
~Vettore();
};
class VettoreLineare : public Vettore {
protected:
public:
VettoreLineare () : Vettore() {};
VettoreLineare (int N) : Vettore(N) {};
VettoreLineare (int N, char* nomefile) : Vettore(N, nomefile) {};
VettoreLineare (const Vettore& v) : Vettore(v) {};
VettoreLineare operator+(const VettoreLineare& v);
VettoreLineare operator*(double lambda);
VettoreLineare& operator=(const VettoreLineare& v);
virtual int GetN() const;
virtual void SetComponent(int i, double x);
};
Vettore.c
#include "Vettore.h"
//Default Constructor
Vettore :: Vettore () {
_N=0;
_v=NULL;
};
//N Constructor
Vettore :: Vettore (int N) {
_N=N;
_v=new double [_N];
for (int i=0; i<_N; i++)
_v[i]=0;
};
//N file-taken constructor
Vettore :: Vettore (int N, char* nomefile) {
_N=N;
_v=new double [_N];
std::ifstream input;
input.open(nomefile);
double dato;
input>>dato;
for(int i=0; i<N; i++){
_v[i]=dato;
input>>dato;
};
input.close();
};
//Copyconstructor
Vettore :: Vettore (const Vettore& V) {
_N=V.GetN();
_v=new double [_N];
for(int i=0; i<_N; i++)
_v[i]=V.GetComponent(i);
};
//Destructor
Vettore::~Vettore(){
delete[] _v;
};
//Set Component
void Vettore :: SetComponent (int i, double x) {
if (i>_N) {
std::cout<<"errore!"<<std::endl;
return ;
};
_v[i]=x;
};
//Get Component
double Vettore :: GetComponent (int i) const {
if (i>_N){
std::cout<<"errore!"<<std::endl;
return 0;
};
return _v[i];
};
//Add Component (aggiunge il valore desiderato nella coda del vettore)
void Vettore :: AddComponent (double x) {
double* a=new double [_N+1];
for(int i=0; i<_N; i++)
a[i]=_v[i];
a[_N]=x;
delete [] _v;
_v=a;
_N=_N+1;
};
//Print
void Vettore :: Print () const {
std::cout<<"Il vettore ha: "<<_N<<" componenti."<<std::endl;
for(int i=0; i<_N; i++)
std::cout<<_v[i]<<std::endl;
};
//Stampa su file
void Vettore :: Print (char* nomefile) const {
std::ofstream output;
output.open(nomefile);
output<<_N;
for(int i=0; i<_N; i++)
output<<_v[i]<<std::endl;
output.close();
};
//Get _N
int Vettore :: GetN () const {
return _N;
};
//Operatore di Assegnazione
Vettore & Vettore::operator =(const Vettore& vetty){
if (_v) delete [] _v;
_N=vetty.GetN();
_v=new double [_N];
for(int n; n<_N; n++)
_v[n]=vetty._v[n];
return *this;
};
//Algoritmo Quicksort
void Vettore :: Sort (){
Quicksort(0,GetN()-1);
};
void Vettore :: Quicksort (unsigned int primo, unsigned int ultimo) {
if(ultimo-primo<=1){
if (GetComponent(primo)>GetComponent(ultimo)) Scambia(primo, ultimo);
return;
}
double pivot= GetComponent(int((primo+ultimo)/2));
unsigned int basso= primo, alto=ultimo;
while(basso < alto) {
while (GetComponent(basso)<pivot) basso++;
while (GetComponent(alto)>pivot) alto--;
if(basso<alto) { Scambia(basso,alto); basso++;};
};
Quicksort(primo, basso-1);
Quicksort(basso, ultimo);
};
void Vettore :: Scambia(int a, int b){
double k;
k=_v[a];
_v[a]=_v[b];
_v[b]=k;
};
//Operatore somma fra vettori
VettoreLineare VettoreLineare::operator+ (const VettoreLineare& v){
VettoreLineare sum(v.GetN());
for(int i=0; i<_N; i++)
sum.SetComponent(i, _v[i]+v.GetComponent(i));
return sum;
};
//Operatore Moltiplicazione scalare
VettoreLineare VettoreLineare::operator* (double lambda){
for(int i=0; i<_N; i++)
_v[i]=_v[i]*lambda;
return *this;
};
//Operatore Assegnazione
VettoreLineare& VettoreLineare::operator= (const VettoreLineare& k){
if (_v) delete [] _v;
_N=k.GetN();
_v=new double [k.GetN()];
for (int i=0; i<_N; i++)
_v[i]=k.GetComponent(i);
return *this;
};
int VettoreLineare :: GetN() const {
return _N;
};
void VettoreLineare :: SetComponent(int i, double x) {
_v[i]=x;
};
main.c
#include "equazione_differenziale.h"
#include "Vettore.h"
#include "iostream"
#include "TGraph.h"
#include "TApplication.h"
#include "TCanvas.h"
#include "TAxis.h"
using namespace std;
int main () {
//PRIMO PUNTO
//Dichiarazione equazione
Protone* myProt=new Protone (1.67E-27, 1.60E-19, 1.E7, 0.2, 5.E8);
Runge_Kutta myKutta(myProt);
//Dichiarazione DatiIniziali
VettoreLineare DatiIniziali (2);
DatiIniziali.SetComponent(0, 0);
DatiIniziali.SetComponent(1,1E8);
//dichiarazione tempo
double t_min=0;
double h=1E-12;
//definizione variabili root
TApplication myApp("myApp",0,0);
TGraph* g = new TGraph();
//ciclo
for(int i=0;i<(1E-7-t_min)/h;i++){
double x,y;
x=t_min+i*h;
DatiIniziali= myKutta.Passo(x, DatiIniziali ,h);
y=DatiIniziali.GetComponent(0);
g->SetPoint(i,x,y);
};
//Run grafico
TCanvas *c=new TCanvas("C1","C1",1);
c->cd();
g->GetXaxis()->SetTitle("t[s]");
g->GetYaxis()->SetTitle("x[m]");
g->Draw("AL");
myApp.Run();
return 0;
};
The strange thing is that this program works on university computers, but it doesn't on neither of my two computers. I'm thinking this means I installed ROOT badly on both computers, but I sincerely wouldn't know how to proove it.
