Calculating Air Resistance

Written in C++ this is a simple console application which calculates the trajectory of a projectile with air resistance applied. The coordinates are printed to a text file.

This was created for the first Maths & Physics Assignment of year 1 as an added extra.

#include "stdafx.h"
#include <iostream>
#include <math.h>
#include <fstream>
#include <vector>
#include <string>
using namespace std;
const double PI = 3.14159265359;

struct XYCoords
{
 double x;
 double y;
};

vector<XYCoords> plotPoints;
XYCoords temp;
double x = 0; // x co-ordinate
double y = 0; // y co-ordinate

void fillVector()
{
 temp.x = x;
 temp.y = y;
 plotPoints.push_back(temp);
}

void openToFile(string filename, const vector<XYCoords>& v)
{
 ofstream openstream{ filename };
 if (!openstream) { cout << "Error, can't open output file" + filename << endl; }

 for (int i = 0; i < v.size(); ++i)
 {
  openstream << "x " << v[i].x << " y " << v[i].y << endl;
 }

}

int main()
{
 //constants of grenade
 double m = 0.39; //mass kg
 double radius = 0.03175; // radius of grenade in meters
 double ar = PI*(radius*radius); //area of grenade
 double cD = 0.5; // drag co-efficient
 double p = 1.225; //air density
 double d = (p*cD*ar) / 2; //drag
 double g = 9.81; //gravity

 //initial conditions
 double deltaT = 0.1; //time interval
 double deltaT2 = deltaT*deltaT; //delta time squared
 double v = 35; // velocity
 double theta = 60; // angle at launch
 double vx = v*cos(theta*(PI / 180)); // initial x velocity
 double vy = v*sin(theta*(PI / 180)); // initial y velocity
 double t = 0; // time
 double ax = -(d / m)*v*vx; // initial acceleration on x
 double ay = -g - (d / m)*v*vy; // initial acceleration on y

 //start loop
 int loop = 1; // counts the number of plot points
 while (y > -0.001) // start a loop and continue looping until the y co-ordinate goes below -0.001
 {
  //Step 5 - calculate acceleration components
  ax = -(d / m)*v*vx;
  ay = -g - (d / m)*v*vy;

  //Step 6 - Push Co-ordinates back into Vector
  fillVector();
  ++loop;

  //Step 7 - calculate new velocity
  vx = vx + ax*deltaT;
  vy = vy + ay*deltaT;

  //Step Eight - calculate new co-ordinates
  x = x + (vx*deltaT) + (0.5*ax*deltaT2);
  y = y + (vy*deltaT) + (0.5*ay*deltaT2);

  //Step 9 - increment time by delta time
  t = t + deltaT;

  //Step 9 - recalculate velocity
  v = sqrt((vx*vx) + (vy*vy));

 }

 //write to a text file
 openToFile("airResistanceCoords.txt", plotPoints);

    return 0;
}