# Feil i Euler og Euler midtpunkt som funksjon av dt for dot x = x cos t

import numpy as np
import matplotlib.pyplot as plt


x0 = 1
tMax = 90


dtl = 10.**(-np.arange(0,5))
e = np.zeros(len(dtl))
em = np.zeros(len(dtl))

k = 0
for dt in dtl:
    n = int(tMax/dt)
    t = np.zeros(n)
    dx = np.zeros(n)
    x = np.zeros(n)
    xm = np.zeros(n)
    x[0] = x0
    xm[0] = x0
    dth = dt/2
    for i in range(n-1):
        t[i+1] = (i+1)*dt
        dx[i] = np.cos(t[i])*x[i]
        x[i+1] = x[i] + dx[i]*dt
        xh = xm[i] + np.cos(t[i])*xm[i]*dth
        xm[i+1] = xm[i] +  np.cos(t[i]+dth)*xh*dt
    e[k] = x[-1]-x0*np.exp(np.sin(t[-1]))
    em[k] = xm[-1]-x0*np.exp(np.sin(t[-1]))
    print(dt,k,e[k],em[k])
    k = k+1

plt.loglog(dtl,abs(e),'o', label='Eulers metode')
plt.loglog(dtl,abs(em),'o', label='Eulers midtpunktmetode')
plt.loglog(dtl,dtl*50, label='Stigningstall = 1')
plt.loglog(dtl,dtl**2*0.05, label='Stigningstall = 2')
plt.legend(loc='upper left')
plt.xlabel(r'$\Delta t$')
plt.ylabel(r'$F(\Delta t)$')
plt.savefig('eulerMidtpunktVariabelDt.pdf')
plt.show()