ti_loh.py

import numpy as np

a = 1
b = 2
h = 0.001
x0 = 1
y0 = 2


def dif_function(x, y):
    return (6 - x ** 2 * y ** 2) / (-x ** 2)
    #return 3 * x ** 2


def first_method(x, y0):
    y = np.ones(len(x))

    for i in range(len(x)):
        if i == 0:
            y[i] = y0
        else:
            y[i] = y[i - 1] + h * dif_function(x[i - 1], y[i - 1])
    return y


def second_method(x, y0):
    y = np.ones(len(x))

    for i in range(len(x)):
        if i == 0:
            y[i] = y0
        else:
            y[i] = y[i - 1] + h * dif_function(x[i - 1] + h / 2, y[i - 1] + h / 2 * dif_function(x[i - 1], y[i - 1]))
    return y


def third_method(x, y0):
    y = np.ones(len(x))

    for i in range(len(x)):
        if i == 0:
            y[i] = y0
        else:
            y_1 = y[i - 1] + h * dif_function(x[i - 1], y[i - 1])
            y[i] = y[i - 1] + (h / 2) * (dif_function(x[i - 1], y[i - 1]) + dif_function(x[i], y_1))
    return y


def rungecut(x, y0):
    y = np.ones(len(x))

    for i in range(len(x)):
        if i == 0:
            y[i] = y0
        else:
            k1 = h * dif_function(x[i - 1], y[i - 1])
            k2 = h * dif_function(x[i - 1] + h / 2, y[i - 1] + k1 / 2)
            k3 = h * dif_function(x[i - 1] + h / 2, y[i - 1] + k2 / 2)
            k4 = h * dif_function(x[i - 1] + h, y[i - 1] + k3)

            y[i] = y[i - 1] + (1.0 / 6) * (k1 + 2 * k2 + 2 * k3 + k4)
    return y


def adams(x, y0):
    y = []
    q = []
    for i in range(5):
        if i == 0:
            y.append(y0)
            q.append(h * dif_function(x[i], y[i]))
        else:
            k1 = h * dif_function(x[i - 1], y[i - 1])
            k2 = h * dif_function(x[i - 1] + h / 2, y[i - 1] + k1 / 2)
            k3 = h * dif_function(x[i - 1] + h / 2, y[i - 1] + k2 / 2)
            k4 = h * dif_function(x[i - 1] + h, y[i - 1] + k3)

            y.append(y[i - 1] + (1.0 / 6) * (k1 + 2 * k2 + 2 * k3 + k4))
            q.append(h * dif_function(x[i], y[i]))

    for i in range(5, len(x)):
        # qi = q[i - 1]
        # qi_1 = q[i - 1] - q[i - 2]
        # qi_2 = q[i - 1] - 2 * q[i - 2] + q[i - 3]
        # qi_3 = q[i - 1] - 3 * q[i - 2] + 3 * q[i - 3] - q[i - 4]
        # qi_4 = q[i - 1] - 4 * q[i - 2] + 6 * q[i - 3] - 4 * q[i - 4] + q[i - 5]

        # y_new = y[i - 1] + qi + (1.0 / 2) * qi_1 + (5.0 / 12) * qi_2 + (3.0 / 8) * qi_3 + (251.0 / 720) * qi_4
        y_new = y[i-1]+ (1901*q[i-1] -2774*q[i-2] + 2616 * q[i-3] - 1274*q[i-4] +251*q[i-5]) / 720
        y.append(y_new)
        q.append(h * dif_function(x[i], y[i]))
    return y


def diff(x, y):
    dif = []
    for i in range(len(x)):
        dif.append(abs(2 / x[i] - y[i]))
        #dif.append(abs(x[i] ** 3 + 1 - y[i]))
    return dif


x = np.linspace(1, 2, 1.0 / h + 1)
print('First method:')
y = first_method(x, y0)
difference = diff(x, y)
print(difference[len(x) - 1])
print('Second method:')
y = second_method(x, y0)
difference = diff(x, y)
print(difference[len(x) - 1])

print('Third method:')
y = third_method(x, y0)
difference = diff(x, y)
print(difference[len(x) - 1])

print('Fourth method:')
y = rungecut(x, y0)
difference = diff(x, y)
print(difference[len(x) - 1])

print('Fifth method:')
y = adams(x, y0)
difference = diff(x, y)
# for xi in difference:
#     print(xi)
print(difference[len(x) - 1])

print('End')