Завдання 1

import numpy as np
import pandas as pd
import matplotlib.gridspec as gridspec
import matplotlib.pyplot as plt

def dx_dt(x, y):
    return k0 - k1 * x * y

def dy_dt(x, y):
    return k1 * x * y - k2 * y


def runge_kute_system(x, y, h, a, b):
    X = [x]
    Y = [y]
    X_Y = []

    # print(f'h = {h}')
    t = a
    # print(f't = {t}  |  x = {x}  |  y = {y}')
    while round(t, 3) <= b:
        t = t + h

        k1 = h * dx_dt(x, y)
        k2 = h * dx_dt(x + h / 2, y + k1 / 2)
        k3 = h * dx_dt(x + h / 2, y + k2 / 2)
        k4 = h * dx_dt(x + h, y + k3)
        dx = (k1 + 2 * k2 + 2 * k3 + k4) / 6
        x += dx
        X.append(x)

        r1 = h * dy_dt(x, y)
        r2 = h * dy_dt(x + h / 2, y + k1 / 2)
        r3 = h * dy_dt(x + h / 2, y + k2 / 2)
        r4 = h * dy_dt(x + h, y + k3)
        dy = (r1 + 2 * r2 + 2 * r3 + r4) / 6
        y += dy
        Y.append(y)

        # print(f't = {t:.3f}  |  x = {x:.3f}  |  y = {y:.3f}')

    X_Y.append(X)
    X_Y.append(Y)

    # print('\n------------------------------------------------------------\n')
    return X_Y


def grafik(k=3):
    plt.title(f"k0={k0}; k1={k1}; k2={k2}")
    plt.xlabel("x(t)")
    plt.ylabel("y(t)")
    plt.grid()
    lines = []
    labels = []
    for i in dani:
        x_y = runge_kute_system(i[0], i[1], i[2], i[3], i[4])
        line, = plt.plot(x_y[0], x_y[1])
        lines.append(line)
        labels.append(f'x_0={i[0]}; y_0={i[1]}; h={i[2]}')
    plt.legend(lines, labels, loc='center left', bbox_to_anchor=(1, 0.5))
    plt.show()

    if k == 3:
        plt.title(f"x(t) (k0={k0}; k1={k1}; k2={k2})")
        plt.xlabel("t")
        plt.ylabel("x")
        plt.grid()
        lines = []
        labels = []
        for i in dani:
            x_y = runge_kute_system(i[0], i[1], i[2], i[3], i[4])
            plt.plot(np.arange(a, b+2*h, h), x_y[0], label=f'x(t); x_0={i[0]}; y_0={i[1]}; h={i[2]}')
        plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))
        plt.show()

        plt.title(f"y(t) (k0={k0}; k1={k1}; k2={k2})")
        plt.xlabel("t")
        plt.ylabel("y")
        plt.grid()
        lines = []
        labels = []
        for i in dani:
            x_y = runge_kute_system(i[0], i[1], i[2], i[3], i[4])
            plt.plot(np.arange(a, b + 2 * h, h), x_y[1], label=f'y(t); x_0={i[0]}; y_0={i[1]}; h={i[2]}')
        plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))
        plt.show()



def info():

    A = np.array([[-(k1 * k0) / k2, -k2],
                  [(k1 * k0) / k2, 0]])
    B = np.array([k2/k1, k0/k2])
    stat_point = np.array([k2/k1, k0/k2])
    print(f"Стаціонарна точка: ({stat_point[0]}, {stat_point[1]})")
    eigenvalues, eigenvectors = np.linalg.eig(A)
    print("lambda: ", eigenvalues)

    if eigenvalues[0] != eigenvalues[1] and np.iscomplexobj(eigenvalues[0]) == False:

        if (eigenvalues[0] < 0) and (eigenvalues[1] < 0):
            print(f'Харкатер точки спокою: стійкий вузол ; Стійкість точки спокою: асимптотично стійка.')

        if (eigenvalues[0] > 0) and (eigenvalues[1] > 0):
            print(f'Харкатер точки спокою: нестійкий вузол ; Стійкість точки спокою: нестійка.')

        if (eigenvalues[0] > 0 and eigenvalues[1] < 0) or (eigenvalues[0] < 0 and eigenvalues[1] > 0):
            print(f'Харкатер точки спокою: сідло ; Стійкість точки спокою: нестійка.')

    elif np.iscomplexobj(eigenvalues[0]) == True:

        if (eigenvalues[0] < 0) and (eigenvalues[1] < 0):
            print(f'Харкатер точки спокою: стійкий фокус ; Стійкість точки спокою: асимптотично стійка.')

        if (eigenvalues[0] > 0) and (eigenvalues[1] > 0):
            print(f'Харкатер точки спокою: нестійкий фокус ; Стійкість точки спокою: нестійка.')

        if (eigenvalues[0] > 0 and eigenvalues[1] < 0) or (eigenvalues[0] < 0 and eigenvalues[1] > 0):
            print(f'Харкатер точки спокою: центр ; Стійкість точки спокою: стійка.')

    else:
        if (eigenvalues[0] < 0) and (eigenvalues[1] < 0):
            print(f'Харкатер точки спокою: стійкий вузол ; Стійкість точки спокою: асимптотично стійка.')

        else:
            print(f'Харкатер точки спокою: нестійкий вузол ; Стійкість точки спокою: нестійка.')

    if (4 * k2 ** 2) > (k1 * k0):
        print('4*k2^2 > k1*k0, підкореневий вираз від’ємний, а особова точка – фокус.')
    else:
        print('4*k2^2 < k1*k0, особова точка – вузол.')
    return stat_point

def final_grafik():
    global k0
    global k1
    global k2
    k0 = 2
    k1 = 10
    k2 = 2
    fig = plt.figure(figsize=(12, 6))
    gs = gridspec.GridSpec(2, 2, width_ratios=[3, 1], height_ratios=[2, 1])

    ax_main = plt.subplot(gs[:, 0])
    ax_main.set_title(f"Дiаграма параметрiв для двох областей")
    ax_main.set_xlabel("k_2")
    ax_main.set_ylabel("k_1 * k_0")
    ax_main.grid()
    ax_main.plot(k_2, 4 * k_2 ** 2, label='k_0 * k_1 = 4k^2')
    ax_main.legend()

    ax_small1 = plt.subplot(gs[0, 1])
    ax_small1.set_title("k_0 * k_1 > 4*k^2")
    ax_small1.set_xlabel("x(t)")
    ax_small1.set_ylabel("y(t)")
    for i in dani:
        x_y = runge_kute_system(i[0], i[1], i[2], i[3], i[4])
        ax_small1.plot(x_y[0], x_y[1])
    ax_small1.grid()


    k0 = 2
    k1 = 10
    k2 = 10
    ax_small2 = plt.subplot(gs[1, 1])
    ax_small2.set_title("k_0 * k_1 < 4*k^2")
    ax_small2.set_xlabel("x(t)")
    ax_small2.set_ylabel("y(t)")
    ax_small2.grid()
    for i in dani:
        x_y = runge_kute_system(i[0], i[1], i[2], i[3], i[4])
        ax_small2.plot(x_y[0], x_y[1])

    plt.tight_layout()
    plt.show()


k0 = 2
k1 = 10
k2 = 2

h = 0.001
a = 0
b = 5

point = info()
dani = [[0.01, 2.0, h, a, b],
        [0.01, 0.01, h, a, b],
        [1.2, 1.0, h, a, b],
        [0.01, 0.5, h, a, b],
        [0.01, 1.0, h, a, b],
        [1.2, 0.6, h, a, b],
        [0.01, 1.5, h, a, b],
        [0.01, 0.2, h, a, b],
        [0.0, 0.0, h, a, b]]


grafik()

k0 = 2
k1 = 10
k2 = 10

point = info()
grafik()


k_2 = np.linspace(0, 10, 20)
k_1 = np.linspace(0, 10, 20)

final_grafik()