import csv
import numpy as np
from math import exp, sqrt, sin, cos, pi

"""
  Solve simulataneous second order diffrential equations by Euler method
"""

# d2x/dt2 = force(t, x, y), d2y/dt2 = force(t, x, y)
# define function to be integrated
def force(t, x, y):
    return -x, -y

# solution: x(0.0) = 0.0, y(0.0) = 2.0, vx(0.0) = 1.0, vy(0.0) = 0.0
def xysolution(t):
    return sin(t), 2.0 * cos(t)

#===================
# parameters
#===================
outfile = 'diffeq2nd_2D_Euler.csv'
x0, vx0 = 0.0, 1.0
y0, vy0 = 2.0, 0.0
dt = 0.01
nt = 501
iprint_interval = 20

#===================
# main routine
#===================
def main(x0, y0, vx0, vy0, dt, nt):
    print("Solve simulataneous second order diffrential equations by Euler method")

    print("Write to [{}]".format(outfile))

# open outfile to write a csv file
    f = open(outfile, 'w')
    fout = csv.writer(f, lineterminator='\n')
    fout.writerow([
        't', 'x(cal)', 'x(exact)', 'y(cal)', 'y(exact)'
        ])

    print("{:^5}  {:^12}  {:^12}  {:^12}  {:^12}".format('t', 'x(cal)', 'x(exact)', 'y(cal)', 'y(exact)'))
    for i in range(0, nt):
        t0 = i * dt
        fx0, fy0 = force(t0, x0, y0)
        vx1 = vx0 + dt * fx0
        vy1 = vy0 + dt * fy0
        x1  = x0  + dt * vx0
        y1  = y0  + dt * vy0
        xexact, yexact = xysolution(t0)
        if i == 1 or i % iprint_interval == 0:
            print("t={:5.2f}  {:12.6f}  {:12.6f}  {:12.6f}  {:12.6f}".format(t0, x0, xexact, y0, yexact))
        fout.writerow([t0, x0, xexact, y0, yexact])
        x0  = x1
        vx0 = vx1
        y0  = y1
        vy0 = vy1

    f.close()

if __name__ == '__main__':
    main(x0, y0, vx0, vy0, dt, nt)



