Sì scusate, posto tutto il programma e stavolta sotto forma di codice. Ho corretto anche quelle cosette che mi hai suggerito.

codice:
#include <iostream>

#include <cmath>

#include <fstream>

#include <iomanip>

#include <mpi.h>


#define G 6.67428e-11

#define M 1.9891e30

#define MIN 0.

#define MAX 1000000000.




void orb_to_car(double a, double e, double AM, double n, double v[],int i);

void ele_pia(int i);

void runge5 (double k[], double h, double m[], int N, int mype);

double f(double g[], int i, double m[], int N);




int x,y,vx,vy;



int main(int argc,char** argv)

{

using namespace std;

int nprocs, mype;

MPI_Status status;

MPI_Init(&argc,&argv);

MPI_Comm_size(MPI_COMM_WORLD, &nprocs);

MPI_Comm_rank(MPI_COMM_WORLD, &mype);



int N;    


if (mype == 0)

  {

  cout << "Inserisci il numero di pianeti: ";

  cin >> N;

  MPI_Bcast(&N,1,MPI_DOUBLE,0,MPI_COMM_WORLD);

  }



int i,n=4*N;


double h=10000.;

double *v=new double[n];

double *m=new double[N];

double a,e,AM,nn;




if (mype == 0)

{

for (i =0; i<N ; i++)

  {

  cout << "Inserisci il semiasse maggiore a del pianeta " << i+1 << ": ";

  cin >> a;

  cout << "Inserisici l'eccentricità e del pianeta " << i+1 << ": ";

  cin >> e;

  cout << "Inserisci l'anomalia media M del pianeta " << i+1 << ": ";

  cin >> AM;

  cout << "Inserisci la massa del pianeta " << i+1 << ": ";

  cin >> m[i];

  nn=sqrt(G*(M+m[i])/(pow(a,3)));

  orb_to_car(a,e,AM,nn,v,i);

  ele_pia(i);

  cout << "x: " << v[x] << "y: " << v[y] << "vx: " << v[vx] << "vy: " << v[vy] << endl;

  }

}

MPI_Barrier(MPI_COMM_WORLD);

MPI_Bcast(v,n,MPI_DOUBLE,0,MPI_COMM_WORLD);

MPI_Bcast(m,N,MPI_DOUBLE,0,MPI_COMM_WORLD);



ofstream out;


if (mype == 0) out.open("dati.dat");



for (double t=MIN; t<=MAX; t+=h)
             //<---------------------------------------------------------------------
  {

  runge5(v,h,m,N,mype);

  MPI_Barrier(MPI_COMM_WORLD);

  if (mype==0)

    {

    for (i=0; i<N;i++)

      {

      ele_pia(i);

      out << v[x] << " " << v[y] << " ";

      }

    out << endl;

    }

  }

if (mype == 0) out.close();


MPI_Finalize();

return 0;

}




void ele_pia(int i)

{

x=4*i;

y=x+1;

vx=y+1;

vy=vx+1;

}




void orb_to_car(double a, double e, double AM, double n, double v[],int i)

{

double U,Up,f,r;

ele_pia(i);

U=f=r=0.;

Up=1.;

for (int l=0; fabs(U-Up)>0.000001; l++)

  {

  Up=U;

  U=Up-(Up-(e*sin(Up*(M_PI/180)))-AM)/(1-e*cos(Up*(M_PI/180)));

  }

r=a*(1-e*cos(U*(M_PI/180)));

f=2*atan(sqrt((1+e)/(1-e))*tan(U*(M_PI/180)*0.5))*(180/M_PI);

v[x]=r*cos(f*(M_PI/180));

v[y]=r*sin(f*(M_PI/180));

v[vx]=-(n*a*sin(U*(M_PI/180)))/(1-e*cos(U*(M_PI/180)));

v[vy]=(n*a*cos(U*(M_PI/180))*sqrt(1-e*e))/(1-e*cos(U*(M_PI/180)));

}




void runge5 (double k[], double h, double m[], int N, int mype)

{

double *F1=new double[4*N];

double *F2=new double[4*N];

double *F3=new double[4*N];

double *F4=new double[4*N];

double *F5=new double[4*N];

double *F6=new double[4*N];

double *t2=new double[4*N];

double *t3=new double[4*N];

double *t4=new double[4*N];

double *t5=new double[4*N];

double *t6=new double[4*N];


int i,q;

for (int j=0; j<4; j++)

  {

  for(int s=0;s< N/4.;s++)

    {

	q=4*j+16*s;

	if((mype==j) && (q<4*N))

	  {

      for (i=q; i<q+4; i++) t2[i]= k[i]+h*(2./9.)*(F1[i]=f(k,i,m,N));

      for (i=q; i<q+4; i++) t3[i]= k[i]+h*((1./12.)*F1[i]+(1./4.)*(F2[i]=f(t2,i,m,N)));

      for (i=q; i<q+4; i++) t4[i]= k[i]+h*((69./128.)*F1[i]-(243./128.)*F2[i]+(135./64.)*(F3[i]=f(t3,i,m,N)));

      for (i=q; i<q+4; i++) t5[i]= k[i]+h*(-(17./12.)*F1[i]+(27./4.)*F2[i]-(27./5.)*F3[i]+(16./15.)*(F4[i]=f(t4,i,m,N)));

      for (i=q; i<q+4; i++) t6[i]= k[i]+h*((65./432.)*F1[i]-(5./16.)*F2[i]+(13./16.)*F3[i]+(4./27.)*F4[i]+(5./144.)*(F5[i]=f(t5,i,m,N)));

      for (i=q; i<q+4; i++) F6[i]= f(t6,i,m,N);

    

      for (i=q; i<q+4; i++) k[i]+=h*((47./450.)*F1[i]+(12./25.)*F3[i] +(32./225.)*F4[i]+(1./30.)*F5[i]+(6./25.)*F6[i]);

	  MPI_Bcast(&k[q],4,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD);

      }

    }

  }

}




double f(double g[], int i,double m[], int N)

{

int p=i/4;

ele_pia(p);

if (i%4==0) return (g[vx]);

if (i%4==1) return (g[vy]);

double sum=0.;

if (i%4==2)

  {

  for (int j=0; j<4*N; j+=4)

    {

    if (j/4!=p)

    sum+=G*m[j/4]*(((g[j]-g[x])/(pow(pow(g[j]-g[x],2)+pow(g[j+1]-g[y],2),1.5)))+(g[j]/(pow(g[j]*g[j]+g[j+1]*g[j+1],1.5))));

    }

  return ((-G*(M+m[p])*g[x]/(pow(g[x]*g[x]+g[y]*g[y],1.5)))+sum);

  }

if (i%4==3)

    {

  for (int j=1; j<4*N; j+=4)

    {

    if (j/4!=p)

    sum+=G*m[j/4]*(((g[j]-g[y])/(pow(pow(g[j-1]-g[x],2)+pow(g[j]-g[y],2),1.5)))+(g[j]/(pow(g[j-1]*g[j-1]+g[j]*g[j],1.5))));

    }

  return ((-G*(M+m[p])*g[y]/(pow(g[x]*g[x]+g[y]*g[y],1.5)))+sum);

  }

else return 0;

}
Ho evidenziato il ciclo for in questione. Il programma così scritto può essere compilato senza problemi, anche se forse non darà in output ciò che effettivamente è stato progettato di dare. Tuttavia è sufficiente per poter mostrare il mio problema. Provate a inserire il codice
codice:
if (mype==0) cout << "-------->" << v[0];
fuori e dentro il ciclo for e dovrebbe darvi il primo un numero, il secondo "nan".