Converting to C++ HELP

[marcr]

I've been working on an assignment for a few days now and I'm growing more and more frustrated. The program is not finding the right answer and I can't figure out why. Here is my code in C++:

#include <iostream>
#include <cmath>
#include <iostream>
#include <fstream>
using namespace std;

int size();

void initial( double& tstop, double& t, double& dt , double y0[]);

void math( double& t, double rhs[], double y0[]) ;


int main()
{
  std :: ofstream myfile;
  myfile.open ("ass101.dat");

  double G,Msun,x,vx,vy;
  double xn,yn,vxn,vyn;
  double t, tstop, dt, thalf;
  double *y, *y1 , *rhs, *y0 ;
  double *k1, *k2, *k3, *k4 ;
  int neqns;
  neqns = size();
  y=new double[neqns];
  rhs=new double[neqns];
  y1=new double[neqns];
  k1=new double[neqns];
  k2=new double[neqns];
  k3=new double[neqns];
  k4=new double[neqns];

  initial(tstop, t, dt , y);

  while(t<tstop)
    {
      math(t, rhs, y ) ;
      for(int i=0; i < neqns; i++) k1[i] = dt * rhs[i];
      math(t, rhs, k1 ) ;
      for(int i=0; i < neqns; i++) k2[i] = (.5 * k1[i]);
      math(t, rhs, k2 ) ;
      for(int i=0; i < neqns; i++) k3[i] = (.5 * k2[i]);
      math(t, rhs, k3 ) ;
      for(int i=0; i < neqns; i++) k4[i] = ( k3[i] );

      for(int i=0; i < neqns; i++) y[i] = y[i] + (k1[i] + 2.0*k2[i] + 2.0*k3[i] + k4[i] ) / 6.0;

      myfile << y[0] << " " << y[1] << " " << t << "\n";

      t=t+dt;

    }
  cout << "\n" << " X = " << y[0] << "   " << "  Y = " << y[1] << "\n" "\n";
  delete [] y; delete [] rhs; delete [] y1;
  return(0);
}




//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

int size()
{

  return(4);
    }

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

void initial( double& tstop, double& t, double& dt , double y0[])
{
  y0[0] = 2.0;      //  x
  y0[1] = 0.0;      //  y
  y0[2] = -2.5;     // vx
  y0[3] = 3.25;     // vy
  t = 0.0;
  dt = .001/365.0;
  tstop = 2.286637;

    
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

void math(double& t, double rhs[], double y0[]) 
{
const double   Msun = 1;
const double   G = 39.47;

 double x, y, vx, vy;

  x = y0[0];
  y = y0[1];
  vx = y0[2];
  vy = y0[3];

  rhs[0] = vx;
  rhs[1] = vy;
  rhs[2] = - G*Msun*x/ pow(( pow(x,2) + pow(y,2) ),1.5);
  rhs[3] = - G*Msun*y/ pow(( pow(x,2) + pow(y,2) ),1.5);

There is another bit of code that the prof. put up online except it is in Fortran. I'm having trouble converting from Fortran to C++.

P.S. I'm new to this cite and I hope I didn't intrude or do anything wrong
Thank you for your time, I really do appreciate it.

-marc


FORTRAN CODE:

program ode_rk4 
implicit none integer :: neqns                        ! Number of equations 

integer :: i                                          ! Implicit DO-loop index 
real(kind=kind(1.0d0)), allocatable :: y(:), rhs(:)   ! y, & r.h.s 
real(kind=kind(1.0d0)), allocatable :: k1(:), k2(:), k3(:), k4(:) 
real(kind=kind(1.0d0)) :: t, tstop, dt                ! Time, Stopping time, & timestep 

call prob_size(neqns)                                 ! Get # of equations 

allocate(y(neqns))                                    ! Allocate y & k arrays to correct size 
allocate(k1(neqns), k2(neqns), k3(neqns),k4(neqns),rhs(neqns)) 

call prob_init_conds(t,tstop,dt,y)                    ! Get initial conditions 

open(unit=17,file='results.dat')                      ! Open a file for output 

do while(t < tstop)                                   ! Integrate forward in time 

 call prob_rhs(t,y,rhs)                            ! Evaluate right-hand-side of ODEs for k1
 
 k1 = dt*rhs 
 call prob_rhs(t+0.5d0*dt,y+0.5d0*k1,rhs)          ! Evaluate right-hand-side of ODEs for k2 
 k2 = dt*rhs 
 call prob_rhs(t+0.5d0*dt,y+0.5d0*k2,rhs)          ! Evaluate right-hand-side of ODEs for k3 
 k3 = dt*rhs 
 call prob_rhs(t+dt,y+k3,rhs)                      ! Evaluate right-hand-side of ODEs for k4 
 k4 = dt*rhs 
 
 y = y+(k1+2.0d0*k2+2.0d0*k3+k4)/6.0d0              ! Apply rk4 method 

 t = t+dt                                           ! Increment time 

 write(17,'(t2,20(1x,es12.5))') t,(y(i),i=1,neqns)  ! Write out results 

enddo 

close(unit=17) 
stop 
end program ode_rk4 

[admkrk]

Could you be a bit more specific as to what the right answer might be? What is it supposed to be and what is the program supposed to do? What you posted doesn't even compile for me, but that is just the missing "}" at the end and I don't have "ass101.dat" even if it did. Code from another language that doesn't explain its purpose doesn't help much.

If you could narrow it down to where you think the problem is, someone might be better able to help.