diff --git "a/\320\2403212/shamatulskiy_409870/lab1/lab1.docx" "b/\320\2403212/shamatulskiy_409870/lab1/lab1.docx"
new file mode 100644
index 0000000..bff1c8e
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab1/lab1.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab1/lab1.py" "b/\320\2403212/shamatulskiy_409870/lab1/lab1.py"
new file mode 100644
index 0000000..ed65db4
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab1/lab1.py"
@@ -0,0 +1,69 @@
+import numpy as np
+
+def gauss_seidel(A, b, tol, max_iter=10000):
+    n = len(b)
+    x = np.zeros(n)
+    for k in range(1, max_iter+1):
+        x_old = x.copy()
+        for i in range(n):
+            s = np.dot(A[i, :i], x[:i]) + np.dot(A[i, i+1:], x_old[i+1:])
+            x[i] = (b[i] - s)/A[i, i]
+        err = np.abs(x - x_old)
+        if np.all(err < tol):
+            return x, k, err
+    return x, k, err
+
+def make_diagonally_dominant(A, b):
+    n = A.shape[0]
+    for i in range(n):
+        for j in range(i, n):
+            if abs(A[j, i]) >= np.sum(np.abs(A[j, :])) - abs(A[j, i]):
+                if j != i:
+                    A[[i, j]], b[[i, j]] = A[[j, i]].copy(), b[[j, i]].copy()
+                break
+    for i in range(n):
+        if abs(A[i, i]) < np.sum(np.abs(A[i, :])) - abs(A[i, i]):
+            return False
+    return True
+
+def read_matrix():
+    mode = input("Чтение матрицы из файла? (y/n): ")
+    if mode.lower() == "y":
+        fn = input("Имя файла: ")
+        data = np.loadtxt(fn)
+        A = data[:, :-1]
+        b = data[:, -1]
+    else:
+        n = int(input("n: "))
+        A = np.zeros((n, n))
+        b = np.zeros(n)
+        for i in range(n):
+            row = list(map(float, input(f"Строка {i+1} (n+1 чисел): ").split()))
+            A[i, :] = row[:-1]
+            b[i] = row[-1]
+    return A, b
+
+def read_tol():
+    mode = input("Чтение точности из файла? (y/n): ")
+    if mode.lower() == "y":
+        fn = input("Имя файла с точностью: ")
+        tol = float(open(fn).read().strip())
+    else:
+        tol = float(input("Точность: "))
+    return tol
+
+def main():
+    A, b = read_matrix()
+    tol = read_tol()
+    if not make_diagonally_dominant(A, b):
+        print("Невозможно добиться диагонального преобладания")
+        return
+    x, iters, err = gauss_seidel(A, b, tol)
+    norm = np.linalg.norm(A, ord=np.inf)
+    print("Норма матрицы (сумма по строкам, бесконечная):", norm)
+    print("Решение x:", x)
+    print("Число итераций:", iters)
+    print("Вектор погрешностей:", err)
+
+if __name__ == "__main__":
+    main()
diff --git "a/\320\2403212/shamatulskiy_409870/lab2/lab2.py" "b/\320\2403212/shamatulskiy_409870/lab2/lab2.py"
new file mode 100644
index 0000000..be62797
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab2/lab2.py"
@@ -0,0 +1,125 @@
+import math
+import sys
+import matplotlib.pyplot as plt
+
+def f1(x):
+    return x**3 - 2*x - 5
+
+def f2(x):
+    return math.sin(x) - 0.5*x
+
+def f3(x):
+    return math.exp(-x) - x
+
+functions = {
+    '1': ('x^3 - 2x - 5', f1),
+    '2': ('sin(x) - 0.5x', f2),
+    '3': ('e^(-x) - x', f3),
+}
+
+def phi3(x):
+    return math.exp(-x)
+
+phis = {
+    '3': phi3,
+}
+
+def bisection(f,a,b,eps):
+    fa,fb=f(a),f(b)
+    if fa*fb>0: return None,None,None
+    n=0
+    while (b-a)/2>eps:
+        c=(a+b)/2
+        if f(a)*f(c)<=0:
+            b=c
+        else:
+            a=c
+        n+=1
+    return (a+b)/2,f((a+b)/2),n
+
+def secant(f,a,b,eps):
+    fa,fb=f(a),f(b)
+    if fa==fb: return None,None,None
+    n=0
+    while abs(b-a)>eps:
+        c=b - fb*(b-a)/(fb-fa)
+        a,b,fa,fb=b,c,fb,f(c)
+        n+=1
+    return b,f(b),n
+
+def simple_iter(f,phi,a,b,eps):
+    mp = max(abs((phi(a+1e-6)-phi(a))/1e-6),abs((phi(b)-phi(b-1e-6))/1e-6))
+    if mp>=1: return 'bad',None,None
+    x0=(a+b)/2
+    n=0
+    while True:
+        x1=phi(x0)
+        n+=1
+        if abs(x1-x0)<eps:
+            return x1,f(x1),n
+        x0=x1
+
+def get_input():
+    if input('Ввод из файла? (y/n) ').lower()=='y':
+        fn=input('Имя файла: ')
+        with open(fn) as g:
+            data=g.read().split()
+        it=iter(data)
+        eq=next(it)
+        method=next(it)
+        a,b,eps=map(float,(next(it),next(it),next(it)))
+    else:
+        print('Уравнения:')
+        for k,v in functions.items():
+            print(k,v[0])
+        eq=input('Выберите номер функции: ')
+        print('Методы: 1-бисекция, 2-секущие, 3-простая итерация')
+        method=input('Выберите метод: ')
+        a=float(input('a= '))
+        b=float(input('b= '))
+        eps=float(input('ε= '))
+    return eq,method,a,b,eps
+
+def output(root,val,iter_count):
+    if input('Вывод в файл? (y/n) ').lower()=='y':
+        fn=input('Имя файла: ')
+        with open(fn,'w') as g:
+            g.write(f'Корень: {root}\n')
+            g.write(f'f(root): {val}\n')
+            g.write(f'Итераций: {iter_count}\n')
+    else:
+        print('Корень:',root)
+        print('f(root)=',val)
+        print('Итераций:',iter_count)
+
+def plot(f,a,b):
+    xs=[a+i*(b-a)/1000 for i in range(1001)]
+    ys=[f(x) for x in xs]
+    plt.plot(xs,ys)
+    plt.axhline(0, color='black')
+    plt.show()
+
+def main():
+    eq,method,a,b,eps=get_input()
+    desc,f=functions.get(eq,(None,None))
+    if f is None:
+        print('Неизвестная функция'); sys.exit()
+    fa,fb=f(a),f(b)
+    if method in ('1','2') and fa*fb>0:
+        print('На интервале нет единственного корня'); sys.exit()
+    if method=='1':
+        root,val,it=bisection(f,a,b,eps)
+    elif method=='2':
+        root,val,it=secant(f,a,b,eps)
+    else:
+        phi=phis.get(eq)
+        if not phi:
+            print('Итерационная форма не задана'); sys.exit()
+        root,val,it=simple_iter(f,phi,a,b,eps)
+        if root=='bad':
+            print('Условие сходимости не выполнено'); sys.exit()
+    output(root,val,it)
+    plot(f,a,b)
+
+if __name__=='__main__':
+    main()
diff --git "a/\320\2403212/shamatulskiy_409870/lab2/lab2_system.py" "b/\320\2403212/shamatulskiy_409870/lab2/lab2_system.py"
new file mode 100644
index 0000000..71c5067
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab2/lab2_system.py"
@@ -0,0 +1,93 @@
+import math
+import sys
+import numpy as np
+import matplotlib.pyplot as plt
+
+def sys1(x):
+    u,v = x
+    return np.array([u**2 + v**2 - 4, u*v - 1])
+
+def phi1(x):
+    u,v = x
+    return np.array([math.sqrt(max(0,4 - v**2)), 1/u])
+
+def sys2(x):
+    u,v = x
+    return np.array([math.sin(u) + v - 1, u + math.cos(v) - 1])
+
+def phi2(x):
+    u,v = x
+    return np.array([math.asin(1 - v), math.acos(1 - u)])
+
+systems = {
+    '1': (sys1, phi1),
+    '2': (sys2, phi2),
+}
+
+def jacobian(phi, x, h=1e-6):
+    n = len(x)
+    J = np.zeros((n,n))
+    for i in range(n):
+        dx = np.zeros(n); dx[i] = h
+        J[:,i] = (phi(x+dx) - phi(x-dx))/(2*h)
+    return J
+
+def simple_iter(sysf, phi, x0, eps, maxit=1000):
+    J = jacobian(phi, x0)
+    if max(abs(np.linalg.eigvals(J))) >= 1:
+        print('Условие сходимости не выполнено'); sys.exit()
+    x_prev = np.array(x0)
+    errors = []
+    for k in range(1, maxit+1):
+        x_new = phi(x_prev)
+        err = np.abs(x_new - x_prev)
+        errors.append(err)
+        if np.all(err < eps):
+            return x_new, k, errors
+        x_prev = x_new
+    print('Не сошлось за max итераций'); sys.exit()
+
+def get_input():
+    print('Системы:')
+    print('1: x^2 + y^2 = 4; x*y = 1')
+    print('2: sin(x) + y = 1; x + cos(y) = 1')
+    sysn = input('Выберите номер системы: ')
+    x0 = float(input('x0= ')), float(input('y0= '))
+    eps = float(input('ε= '))
+    return sysn, x0, eps
+
+def output(x, it, err):
+    if input('Вывод в файл? (y/n) ').lower()=='y':
+        fn = input('Имя файла: ')
+        with open(fn, 'w') as g:
+            g.write(f'x = {x}\n')
+            g.write(f'Итераций = {it}\n')
+            g.write(f'Ошибка последней итерации = {err}\n')
+    else:
+        print('x =', x)
+        print('Итераций =', it)
+        print('Ошибка последней итерации =', err)
+
+def plot(sysf, a=-2, b=2):
+    xs = np.linspace(a, b, 400)
+    ys = np.linspace(a, b, 400)
+    X, Y = np.meshgrid(xs, ys)
+    Z1 = sysf([X, Y])[0]
+    Z2 = sysf([X, Y])[1]
+    plt.contour(X, Y, Z1, levels=[0])
+    plt.contour(X, Y, Z2, levels=[0])
+    plt.axhline(0, color='black')
+    plt.axvline(0, color='black')
+    plt.show()
+
+def main():
+    sysn, x0, eps = get_input()
+    sysf, phi = systems.get(sysn, (None,None))
+    if not sysf:
+        print('Неизвестная система'); sys.exit()
+    x, it, errors = simple_iter(sysf, phi, x0, eps)
+    output(x, it, list(errors[-1]))
+    plot(sysf)
+
+if __name__=='__main__':
+    main()
diff --git "a/\320\2403212/shamatulskiy_409870/lab2/report.docx" "b/\320\2403212/shamatulskiy_409870/lab2/report.docx"
new file mode 100644
index 0000000..f6418bd
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab2/report.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab3/report.docx" "b/\320\2403212/shamatulskiy_409870/lab3/report.docx"
new file mode 100644
index 0000000..3394bd5
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab3/report.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab3/task.py" "b/\320\2403212/shamatulskiy_409870/lab3/task.py"
new file mode 100644
index 0000000..147e89e
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab3/task.py"
@@ -0,0 +1,168 @@
+import math
+import sys
+
+
+def f1(x):
+    """f(x) = x^2"""
+    return x**2
+
+def f2(x):
+    """f(x) = sin(x)"""
+    return math.sin(x)
+
+def f3(x):
+    """f(x) = e^x"""
+    return math.exp(x)
+
+def f4(x):
+    """f(x) = 1/x (при x ≠ 0)"""
+    return 1.0 / x
+
+def f5(x):
+    """f(x) = x^3 + 2x"""
+    return x**3 + 2*x
+
+functions = [f1, f2, f3, f4, f5]
+
+def rectangle_method(f, a, b, n, mode="left"):
+    h = (b - a) / n
+    total = 0.0
+
+    if mode == "left":
+        for i in range(n):
+            x = a + i * h
+            total += f(x)
+
+    elif mode == "right":
+        for i in range(1, n+1):
+            x = a + i * h
+            total += f(x)
+
+    elif mode == "middle":
+        for i in range(n):
+            x = a + (i + 0.5) * h
+            total += f(x)
+    else:
+        raise ValueError("Неизвестный режим прямоугольников.")
+
+    return total * h
+
+def trapezoid_method(f, a, b, n):
+    h = (b - a) / n
+    total = 0.5 * (f(a) + f(b))
+    for i in range(1, n):
+        x = a + i * h
+        total += f(x)
+    return total * h
+
+def simpson_method(f, a, b, n):
+    if n % 2 != 0:
+        n += 1
+
+    h = (b - a) / n
+    total = f(a) + f(b)
+    odd_sum = 0.0
+    even_sum = 0.0
+    for k in range(1, n):
+        x = a + k * h
+        if k % 2 == 0:
+            even_sum += f(x)
+        else:
+            odd_sum += f(x)
+    total += 4 * odd_sum + 2 * even_sum
+    return total * h / 3.0
+
+def integrate_with_precision(f, a, b, eps, method_name="rectangle", mode="left", n_init=4):
+    if method_name == "rectangle":
+        method_func = lambda f,a,b,n: rectangle_method(f,a,b,n,mode=mode)
+        p = 2 if mode == "middle" else 1
+        if p == 2:
+            runge_factor = 1.0/3.0
+        else:
+            runge_factor = 1.0  # Если p=1
+
+    elif method_name == "trapezoid":
+        method_func = trapezoid_method
+        runge_factor = 1.0/3.0
+
+    elif method_name == "simpson":
+        method_func = simpson_method
+        runge_factor = 1.0/15.0
+
+    else:
+        raise ValueError("Неизвестный метод: " + method_name)
+
+    n = n_init
+    I_old = method_func(f, a, b, n)
+    while True:
+        n *= 2
+        I_new = method_func(f, a, b, n)
+        # оценка погрешности по правилу Рунге:
+        error_est = abs(I_new - I_old) * runge_factor
+
+        if error_est < eps:
+            return I_new, n
+        I_old = I_new
+
+def main():
+    print("Выберите функцию для интегрирования:")
+    for idx, func in enumerate(functions, start=1):
+        print(f"{idx}) {func.__doc__}")
+
+    choice = int(input("Введите номер функции (1..5): "))
+    if choice < 1 or choice > len(functions):
+        print("Некорректный выбор функции. Завершение.")
+        return
+    f = functions[choice - 1]
+
+    a = float(input("Введите левую границу интегрирования a = "))
+    b = float(input("Введите правую границу интегрирования b = "))
+    if a == b:
+        print("Границы интегрирования совпадают, интеграл будет 0.")
+        return
+
+    eps = float(input("Введите требуемую точность (например, 1e-6): "))
+
+    n_init = int(input("Введите начальное число разбиений (n >= 2), по умолчанию 4: ") or 4)
+    if n_init < 2:
+        n_init = 4
+
+    print("\nВыберите метод интегрирования:")
+    print("1) Прямоугольники (левые)")
+    print("2) Прямоугольники (правые)")
+    print("3) Прямоугольники (средние)")
+    print("4) Трапеции")
+    print("5) Симпсон")
+
+    method_choice = int(input("Номер метода: "))
+
+    if method_choice == 1:
+        method_name = "rectangle"
+        mode = "left"
+    elif method_choice == 2:
+        method_name = "rectangle"
+        mode = "right"
+    elif method_choice == 3:
+        method_name = "rectangle"
+        mode = "middle"
+    elif method_choice == 4:
+        method_name = "trapezoid"
+        mode = ""
+    elif method_choice == 5:
+        method_name = "simpson"
+        mode = ""
+    else:
+        print("Некорректный выбор метода. Завершение.", file=sys.stderr)
+        return
+
+    I, N = integrate_with_precision(f, min(a, b), max(a, b), eps, method_name, mode, n_init)
+
+    # вывод рез-тов
+    print(f"\nРезультаты:")
+    print(f"Метод: {method_name}, вариант = {mode if mode else '—'}")
+    print(f"Число разбиений, достигнутое для заданной точности: {N}")
+    print(f"Приблизительное значение интеграла: {I:.8f}")
+    print(f"Точность (заданная): {eps}")
+
+if __name__ == "__main__":
+    main()
diff --git "a/\320\2403212/shamatulskiy_409870/lab4/report.docx" "b/\320\2403212/shamatulskiy_409870/lab4/report.docx"
new file mode 100644
index 0000000..38bc118
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab4/report.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab4/task.py" "b/\320\2403212/shamatulskiy_409870/lab4/task.py"
new file mode 100644
index 0000000..569da87
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab4/task.py"
@@ -0,0 +1,179 @@
+import sys
+import math
+from pathlib import Path
+
+import numpy as np
+import matplotlib.pyplot as plt
+
+# -----------------------------------------------------
+# helpers ------------------------------------------------
+
+def read_points_interactive():
+    print("Введите пары x y (пустая строка — конец ввода):")
+    pts = []
+    while True:
+        line = input().strip()
+        if not line:
+            break
+        parts = line.replace(",", ".").split()
+        if len(parts) != 2:
+            print("Ожидалось два числа через пробел")
+            continue
+        pts.append((float(parts[0]), float(parts[1])))
+    return pts
+
+
+def read_points_from_file() -> list:
+    fname = input("Имя файла с данными: ").strip()
+    path = Path(fname)
+    if not path.exists():
+        print("Файл не найден"); sys.exit(1)
+    pts = []
+    for line in path.read_text().splitlines():
+        line = line.strip()
+        if not line:
+            continue
+        parts = line.replace(",", ".").split()
+        if len(parts) != 2:
+            continue
+        pts.append((float(parts[0]), float(parts[1])))
+    return pts
+
+
+def request_points():
+    mode = input("Читать точки из файла? (y/n): ").strip().lower()
+    pts = read_points_from_file() if mode == "y" else read_points_interactive()
+    if not (8 <= len(pts) <= 12):
+        print("Нужно 8–12 точек"); sys.exit(1)
+    xs, ys = map(np.array, zip(*pts))
+    return xs, ys
+
+
+# -----------------------------------------------------
+# метрики ------------------------------------------------
+
+# Сигма - среднеквадратичное отклонение
+def rms(pred, real):
+    return math.sqrt(np.mean((pred - real) ** 2))
+# R^2 -
+def r2_score(pred, real):
+    ss_res = np.sum((real - pred) ** 2)
+    ss_tot = np.sum((real - real.mean()) ** 2)
+    return 1 - ss_res / ss_tot if ss_tot else 0.0
+
+
+# -----------------------------------------------------
+# модели ------------------------------------------------
+
+def model_linear(xs, ys):
+    c = np.polyfit(xs, ys, 1)
+    return np.polyval(c, xs), c
+
+def model_poly2(xs, ys):
+    c = np.polyfit(xs, ys, 2)
+    return np.polyval(c, xs), c
+
+def model_poly3(xs, ys):
+    c = np.polyfit(xs, ys, 3)
+    return np.polyval(c, xs), c
+
+def model_expo(xs, ys):
+    if np.any(ys <= 0):
+        return None, None
+    b, ln_a = np.polyfit(xs, np.log(ys), 1)
+    a = math.exp(ln_a)
+    return a * np.exp(b * xs), (a, b)
+
+def model_log(xs, ys):
+    if np.any(xs <= 0):
+        return None, None
+    b, a = np.polyfit(np.log(xs), ys, 1)
+    return a + b * np.log(xs), (a, b)
+
+def model_power(xs, ys):
+    if np.any(xs <= 0) or np.any(ys <= 0):
+        return None, None
+    b, ln_a = np.polyfit(np.log(xs), np.log(ys), 1)
+    a = math.exp(ln_a)
+    return a * xs ** b, (a, b)
+
+MODELS = {
+    "linear":   (model_linear,    "y = a₀ + a₁ x"),
+    "poly2":    (model_poly2,     "y = a₂ x² + a₁ x + a₀"),
+    "poly3":    (model_poly3,     "y = a₃ x³ + a₂ x² + a₁ x + a₀"),
+    "expo":     (model_expo,      "y = A · e^{Bx}"),
+    "log":      (model_log,       "y = A + B · ln x"),
+    "power":    (model_power,     "y = A · x^{B}"),
+}
+
+# -----------------------------------------------------
+# main --------------------------------------
+
+def main():
+    print("Будут исследоваться следующие формы зависимостей:")
+    for name, (_, formula) in MODELS.items():
+        print(f"  {name:6s}: {formula}")
+    print()
+
+    xs, ys = request_points()
+
+    results = {}
+    for name, (fn, _) in MODELS.items():
+        pred, coeffs = fn(xs, ys)
+        if pred is None:
+            continue
+        results[name] = {
+            "sigma": rms(pred, ys),
+            "coeffs": coeffs,
+            "pred": pred,
+            "r2": r2_score(pred, ys),
+        }
+
+    if not results:
+        print("Не удалось построить ни одной модели"); sys.exit(1)
+
+    best_name = min(results, key=lambda n: results[n]["sigma"])
+
+    out_mode = input("Сохранить результаты в файл? (y/n): ").strip().lower()
+    if out_mode == "y":
+        out_path = Path(input("Имя выходного файла: ").strip())
+        with out_path.open("w", encoding="utf-8") as f:
+            for name, res in results.items():
+                f.write(f"{name}: σ={res['sigma']:.3f}  R²={res['r2']:.3f}  coeffs={res['coeffs']}\n")
+            f.write(f"Лучшее приближение: {best_name} (σ={results[best_name]['sigma']:.3f})\n")
+            if "linear" in results:
+                r = np.corrcoef(xs, ys)[0, 1]
+                f.write(f"Корреляция Пирсона (linear): r={r:.3f}\n")
+        print(f"Результаты сохранены в {out_path}")
+    else:
+        for name, res in results.items():
+            print(f"{name:6s}: σ={res['sigma']:.3f}  R²={res['r2']:.3f}  coeffs={res['coeffs']}")
+        print(f"Лучшее приближение => {best_name} (σ={results[best_name]['sigma']:.3f})")
+        if "linear" in results:
+            r = np.corrcoef(xs, ys)[0, 1]
+            print(f"Корреляция Пирсона (linear): r={r:.3f}")
+
+    plt.scatter(xs, ys, label="data", zorder=3)
+    pu = np.linspace(xs.min(), xs.max(), 400)
+    for name, res in results.items():
+        if name.startswith("poly"):
+            y_plot = np.polyval(res["coeffs"], pu)
+        elif name == "linear":
+            y_plot = res["coeffs"][0] * pu + res["coeffs"][1]
+        elif name == "expo":
+            y_plot = res["coeffs"][0] * np.exp(res["coeffs"][1] * pu)
+        elif name == "log":
+            y_plot = res["coeffs"][0] + res["coeffs"][1] * np.log(pu)
+        elif name == "power":
+            y_plot = res["coeffs"][0] * pu ** res["coeffs"][1]
+        else:
+            continue
+        plt.plot(pu, y_plot, label=name)
+    plt.title("Аппроксимация точек")
+    plt.legend()
+    plt.grid(True)
+    plt.show()
+
+
+if __name__ == "__main__":
+    main()
diff --git "a/\320\2403212/shamatulskiy_409870/lab5/report.docx" "b/\320\2403212/shamatulskiy_409870/lab5/report.docx"
new file mode 100644
index 0000000..2ee7d8e
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab5/report.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab5/task.py" "b/\320\2403212/shamatulskiy_409870/lab5/task.py"
new file mode 100644
index 0000000..94fef04
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab5/task.py"
@@ -0,0 +1,129 @@
+import sys, math, numpy as np
+import matplotlib.pyplot as plt
+
+def finite_diffs(y):
+    table=[y.copy()]
+    while len(table[-1])>1:
+        table.append([table[-1][i+1]-table[-1][i] for i in range(len(table[-1])-1)])
+    return table
+
+def lagrange(xv,yv,x):
+    s=0.0
+    for i in range(len(xv)):
+        p=1.0
+        for j in range(len(xv)):
+            if i!=j:
+                p*=(x-xv[j])/(xv[i]-xv[j])
+        s+=yv[i]*p
+    return s
+
+def newton_divided(xv,yv,x):
+    div=yv.copy()
+    n=len(xv)
+    for k in range(1,n):
+        for i in range(n-1,k-1,-1):
+            div[i]=(div[i]-div[i-1])/(xv[i]-xv[i-k])
+    s=div[-1]
+    for k in range(n-2,-1,-1):
+        s=s*(x-xv[k])+div[k]
+    return s
+
+def newton_fd(xv,yv,x):
+    h=xv[1]-xv[0]
+    table=finite_diffs(yv)
+    n=len(xv)
+    if abs(x-xv[0])<abs(x-xv[-1]):
+        t=(x-xv[0])/h
+        s,fact=table[0][0],1.0
+        for k in range(1,n):
+            fact*=(t-(k-1))/k
+            s+=fact*table[k][0]
+    else:
+        t=(x-xv[-1])/h
+        s,fact=table[0][-1],1.0
+        for k in range(1,n):
+            fact*=(t+(k-1))/k
+            s+=fact*table[k][-1]
+    return s
+
+def read_points():
+    m=int(input("Количество точек: "))
+    pts=[tuple(map(float,input("x y: ").split())) for _ in range(m)]
+    pts.sort()
+    xv,yv=zip(*pts)
+    return list(xv),list(yv)
+
+def func_choice():
+    print("Доступные функции:\ny = sin(x) -- sin\ny = cos(x) -- cos\ny = exp(x) -- exp")
+    fx=input("Выберите функцию: ").strip()
+    a,b=map(float,input("Интервал a b: ").split())
+    n=int(input("Число точек: "))
+    xv=[a+i*(b-a)/(n-1) for i in range(n)]
+    if fx=="sin":
+        yv=[math.sin(v) for v in xv]
+    elif fx=="cos":
+        yv=[math.cos(v) for v in xv]
+    else:
+        yv=[math.exp(v) for v in xv]
+    return xv,yv
+
+def load_file(path):
+    with open(path,"r",encoding="utf-8") as f:
+        lines=[l.strip() for l in f if l.strip()]
+    try:
+        count=int(lines[0])
+        data=lines[1:]
+    except ValueError:
+        data=lines
+    xv,yv=[],[]
+    for l in data:
+        parts=l.split()
+        if len(parts)>=2:
+            a,b=map(float,parts[:2])
+            xv.append(a)
+            yv.append(b)
+    return xv,yv
+
+def plot_method(title,xv,yv,ys_func):
+    xs=np.linspace(min(xv),max(xv),400)
+    plt.figure()
+    plt.plot(xs,[ys_func(xx) for xx in xs])
+    plt.scatter(xv,yv,color="black",zorder=5)
+    plt.title(title)
+    plt.xlabel("x")
+    plt.ylabel("y")
+    plt.grid(True)
+
+
+def main():
+    if input("Ввод из файла? (y/n): ").lower().startswith("y"):
+        xv,yv=load_file(input("Файл: "))
+    else:
+        if input("Сгенерировать по функции? (y/n): ").lower().startswith("y"):
+            xv,yv=func_choice()
+        else:
+            xv,yv=read_points()
+    steps=[round(xv[i+1]-xv[i],10) for i in range(len(xv)-1)]
+    if len(set(steps))!=1:
+        print("Шаг по x должен быть равномерным")
+        sys.exit()
+    table=finite_diffs(yv)
+    print("Таблица конечных разностей:")
+    for i in range(len(xv)):
+        row=[f"{xv[i]:g}",f"{yv[i]:g}"]
+        for k in range(1,len(table)):
+            if i<len(table[k]):
+                row.append(f"{table[k][i]:g}")
+        print("\t".join(row))
+    xq=float(input("Введите X для интерполяции: "))
+    print("Лагранж:",lagrange(xv,yv,xq))
+    print("Ньютон разделённые:",newton_divided(xv,yv,xq))
+    print("Ньютон конечные:",newton_fd(xv,yv,xq))
+
+    plot_method("Lagrange interpolation",xv,yv,lambda x: lagrange(xv,yv,x))
+    plot_method("Newton (divided) interpolation",xv,yv,lambda x: newton_divided(xv,yv,x))
+    plot_method("Newton (finite diffs) interpolation",xv,yv,lambda x: newton_fd(xv,yv,x))
+    plt.show()
+
+if __name__=="__main__":
+    main()
\ No newline at end of file
diff --git "a/\320\2403212/shamatulskiy_409870/lab6/report.docx" "b/\320\2403212/shamatulskiy_409870/lab6/report.docx"
new file mode 100644
index 0000000..eed48e1
Binary files /dev/null and "b/\320\2403212/shamatulskiy_409870/lab6/report.docx" differ
diff --git "a/\320\2403212/shamatulskiy_409870/lab6/task.py" "b/\320\2403212/shamatulskiy_409870/lab6/task.py"
new file mode 100644
index 0000000..ed8adba
--- /dev/null
+++ "b/\320\2403212/shamatulskiy_409870/lab6/task.py"
@@ -0,0 +1,129 @@
+import math
+import matplotlib.pyplot as plt
+
+def eq1_f(x, y):
+    return x + y
+
+def eq1_exact(x, x0, y0):
+    return (y0 + 1) * math.exp(x - x0) - x - 1 + x0
+
+def eq2_f(x, y):
+    return y * x * x
+
+def eq2_exact(x, x0, y0):
+    return y0 * math.exp((x**3 - x0**3) / 3)
+
+def eq3_f(x, y):
+    return y * x
+
+def eq3_exact(x, x0, y0):
+    return y0 * math.exp((x**2 - x0**2) / 2)
+
+def euler(f, x0, y0, xn, h):
+    xs = []
+    ys = []
+    x = x0
+    y = y0
+    while x <= xn + 1e-12:
+        xs.append(x)
+        ys.append(y)
+        y = y + h * f(x, y)
+        x = x + h
+    return xs, ys
+
+def rk4(f, x0, y0, xn, h):
+    xs = []
+    ys = []
+    x = x0
+    y = y0
+    while x <= xn + 1e-12:
+        xs.append(x)
+        ys.append(y)
+        k1 = f(x, y)
+        k2 = f(x + h/2, y + h * k1 / 2)
+        k3 = f(x + h/2, y + h * k2 / 2)
+        k4 = f(x + h, y + h * k3)
+        y = y + h * (k1 + 2*k2 + 2*k3 + k4) / 6
+        x = x + h
+    return xs, ys
+
+def milne(f, x0, y0, xn, h):
+    xs_rk, ys_rk = rk4(f, x0, y0, x0 + 3*h, h)
+    xs = xs_rk[:]
+    ys = ys_rk[:]
+    fs = [f(xs[i], ys[i]) for i in range(4)]
+    i = 3
+    while xs[i] + h <= xn + 1e-12:
+        x_next = xs[i] + h
+        y_pred = ys[i-3] + (4*h/3) * (2*fs[i-2] - fs[i-1] + 2*fs[i])
+        f_pred = f(x_next, y_pred)
+        y_corr = ys[i-1] + (h/3) * (fs[i-1] + 4*fs[i] + f_pred)
+        xs.append(x_next)
+        ys.append(y_corr)
+        fs.append(f(x_next, y_corr))
+        i += 1
+    return xs, ys
+
+def runge_error(method, f, x0, y0, xn, h, p):
+    _, y_h = method(f, x0, y0, xn, h)
+    _, y_h2 = method(f, x0, y0, xn, h/2)
+    y_end_h = y_h[-1]
+    y_end_h2 = y_h2[-1]
+    return abs(y_end_h2 - y_end_h) / (2**p - 1)
+
+def main():
+    print("Выберите уравнение:")
+    print("1) y' = x + y")
+    print("2) y' = y * x^2")
+    print("3) y' = y * x")
+    choice = int(input())
+    if choice == 1:
+        f = eq1_f
+        exact = eq1_exact
+    elif choice == 2:
+        f = eq2_f
+        exact = eq2_exact
+    else:
+        f = eq3_f
+        exact = eq3_exact
+    x0 = float(input("x0 = "))
+    y0 = float(input("y0 = "))
+    xn = float(input("xn = "))
+    h = float(input("h = "))
+    eps = float(input("epsilon = "))
+    xs_e, ys_e = euler(f, x0, y0, xn, h)
+    xs_rk, ys_rk = rk4(f, x0, y0, xn, h)
+    xs_m, ys_m = milne(f, x0, y0, xn, h)
+    print("\n i\t   x\t   Euler\t\t   RK4\t\t   Milne\t\t   Exact")
+    N = len(xs_e)
+    for i in range(N):
+        x_val = xs_e[i]
+        y_ex = exact(x_val, x0, y0)
+        y_e = ys_e[i] if i < len(ys_e) else float('nan')
+        y_r = ys_rk[i] if i < len(ys_rk) else float('nan')
+        y_m = ys_m[i] if i < len(ys_m) else float('nan')
+        print(f"{i:2d}\t{x_val:8.4f}\t{y_e:10.6f}\t{y_r:10.6f}\t{y_m:10.6f}\t{y_ex:10.6f}")
+    err_euler = runge_error(euler, f, x0, y0, xn, h, 1)
+    err_rk4 = runge_error(rk4, f, x0, y0, xn, h, 4)
+    max_err_m = 0.0
+    for i, xm in enumerate(xs_m):
+        ye = exact(xm, x0, y0)
+        ym = ys_m[i]
+        max_err_m = max(max_err_m, abs(ye - ym))
+    print(f"\nОценка погрешности Эйлера (правило Рунге): {err_euler:.6e}")
+    print(f"Оценка погрешности RK4 (правило Рунге): {err_rk4:.6e}")
+    print(f"Максимальная погрешность Милна: {max_err_m:.6e}")
+    xs_exact = [x0 + i*h/10 for i in range(int((xn - x0)/(h/10)) + 1)]
+    ys_exact = [exact(x, x0, y0) for x in xs_exact]
+    plt.plot(xs_exact, ys_exact, label="Exact", color="black")
+    plt.plot(xs_e, ys_e, label="Euler", linestyle="--")
+    plt.plot(xs_rk, ys_rk, label="RK4", linestyle="-.")
+    plt.plot(xs_m, ys_m, label="Milne", linestyle=":")
+    plt.legend()
+    plt.xlabel("x")
+    plt.ylabel("y")
+    plt.title("Numerical vs Exact Solution")
+    plt.show()
+
+if __name__ == "__main__":
+    main()