diff --git a/README.md b/README.md
deleted file mode 100644
index f1ec59983..000000000
--- a/README.md
+++ /dev/null
@@ -1,3 +0,0 @@
-# Базовый курс (моделирование на python)
-Репозиторий для прохождения базового курса по языку программирования python,
-в рамках общего цикла лекций "Моделирование на python".
diff --git a/extratasks gifs/extratask1_a.gif b/extratasks gifs/extratask1_a.gif
new file mode 100644
index 000000000..08e51f77e
Binary files /dev/null and b/extratasks gifs/extratask1_a.gif differ
diff --git a/extratasks gifs/extratask1_b.gif b/extratasks gifs/extratask1_b.gif
new file mode 100644
index 000000000..9ddf65334
Binary files /dev/null and b/extratasks gifs/extratask1_b.gif differ
diff --git a/extratasks gifs/extratask1_c.gif b/extratasks gifs/extratask1_c.gif
new file mode 100644
index 000000000..3ab1bcac7
Binary files /dev/null and b/extratasks gifs/extratask1_c.gif differ
diff --git a/extratasks gifs/extratask1_d.gif b/extratasks gifs/extratask1_d.gif
new file mode 100644
index 000000000..9f4736255
Binary files /dev/null and b/extratasks gifs/extratask1_d.gif differ
diff --git a/extratasks/lab_extratask1_a.py b/extratasks/lab_extratask1_a.py
new file mode 100644
index 000000000..fba2af0aa
--- /dev/null
+++ b/extratasks/lab_extratask1_a.py
@@ -0,0 +1,221 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 50
+t = np.linspace(10**(-7), 1.5 * 10**(-7), 50*frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3,
+     x4, v_x4, y4, v_y4) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    
+    dxdt4 = v_x4
+    dv_xdt4 = (
+      	    - G * m1 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+    dydt4 = v_y4
+    dv_ydt4 = (
+      	    - G * m1 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3,
+            dxdt4, dv_xdt4, dydt4, dv_ydt4)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+q = 1.1 *10**(-13)
+m = 1.1 * 10**(-27)
+a = 0.008
+v = 4 * 10**7
+alpha = np.pi/4
+
+m1 = m
+q1 = -q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = - q
+
+m4 = m
+q4 = q
+
+x10 = a
+v_x10 = - v * np.cos(alpha)
+y10 = a
+v_y10 = v * np.sin(alpha)
+
+x20 = a
+v_x20 = v * np.cos(alpha)
+y20 = - a
+v_y20 = v * np.sin(alpha)
+
+x30 = - a
+v_x30 = v * np.cos(alpha)
+y30 = - a
+v_y30 = - v * np.cos(alpha)
+
+x40 = - a
+v_x40 = - v * np.cos(alpha)
+y40 = a
+v_y40 = - v * np.cos(alpha)
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30,
+      x40, v_x40, y40, v_y40)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(4):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(4):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=70)
+
+plt.axis('equal')
+edge = 5 * a
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('extratask1_a.gif')
\ No newline at end of file
diff --git a/extratasks/lab_extratask1_b.py b/extratasks/lab_extratask1_b.py
new file mode 100644
index 000000000..9c3988832
--- /dev/null
+++ b/extratasks/lab_extratask1_b.py
@@ -0,0 +1,221 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 50
+t = np.linspace(10**(-7), 1.5 * 10**(-7), 50*frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3,
+     x4, v_x4, y4, v_y4) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    
+    dxdt4 = v_x4
+    dv_xdt4 = (
+      	    - G * m1 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+    dydt4 = v_y4
+    dv_ydt4 = (
+      	    - G * m1 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3,
+            dxdt4, dv_xdt4, dydt4, dv_ydt4)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+q = 1.1 *10**(-13)
+m = 1.1 * 10**(-27)
+a = 0.008
+v = 3 * 10**7
+alpha = np.pi/4
+
+m1 = m
+q1 = -q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = - q
+
+m4 = m
+q4 = q
+
+x10 = a
+v_x10 = - v * np.cos(alpha)
+y10 = a
+v_y10 = v * np.sin(alpha)
+
+x20 = a
+v_x20 = v * np.cos(alpha)
+y20 = - a
+v_y20 = v * np.sin(alpha)
+
+x30 = - a
+v_x30 = v * np.cos(alpha)
+y30 = - a
+v_y30 = - v * np.cos(alpha)
+
+x40 = - a
+v_x40 = - v * np.cos(alpha)
+y40 = a
+v_y40 = - v * np.cos(alpha)
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30,
+      x40, v_x40, y40, v_y40)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(4):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(4):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=70)
+
+plt.axis('equal')
+edge = 5 * a
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('extratask1_b.gif')
\ No newline at end of file
diff --git a/extratasks/lab_extratask1_c.py b/extratasks/lab_extratask1_c.py
new file mode 100644
index 000000000..5455238e6
--- /dev/null
+++ b/extratasks/lab_extratask1_c.py
@@ -0,0 +1,221 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 50
+t = np.linspace(10**(-7), 1.5 * 10**(-7), 50*frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3,
+     x4, v_x4, y4, v_y4) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    
+    dxdt4 = v_x4
+    dv_xdt4 = (
+      	    - G * m1 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+    dydt4 = v_y4
+    dv_ydt4 = (
+      	    - G * m1 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3,
+            dxdt4, dv_xdt4, dydt4, dv_ydt4)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+q = 1.1 *10**(-13)
+m = 1.1 * 10**(-27)
+a = 0.008
+v = 2 * 10**7
+alpha = np.pi/4
+
+m1 = m
+q1 = -q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = - q
+
+m4 = m
+q4 = q
+
+x10 = a
+v_x10 = - v * np.cos(alpha)
+y10 = a
+v_y10 = v * np.sin(alpha)
+
+x20 = a
+v_x20 = v * np.cos(alpha)
+y20 = - a
+v_y20 = v * np.sin(alpha)
+
+x30 = - a
+v_x30 = v * np.cos(alpha)
+y30 = - a
+v_y30 = - v * np.cos(alpha)
+
+x40 = - a
+v_x40 = - v * np.cos(alpha)
+y40 = a
+v_y40 = - v * np.cos(alpha)
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30,
+      x40, v_x40, y40, v_y40)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(4):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(4):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=70)
+
+plt.axis('equal')
+edge = 5 * a
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('extratask1_c.gif')
\ No newline at end of file
diff --git a/extratasks/lab_extratask1_d.py b/extratasks/lab_extratask1_d.py
new file mode 100644
index 000000000..1f94891f1
--- /dev/null
+++ b/extratasks/lab_extratask1_d.py
@@ -0,0 +1,221 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 50
+t = np.linspace(10**(-7), 1.5 * 10**(-7), 50*frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3,
+     x4, v_x4, y4, v_y4) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (x1 - x4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            - G * m4 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q4 / m1 * (y1 - y4) 
+               / ((x1 - x4)**2 + (y1 - y4)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (x2 - x4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            - G * m4 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q4 / m2 * (y2 - y4) 
+               / ((x2 - x4)**2 + (y2 - y4)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (x3 - x4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            - G * m4 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q4 / m3 * (y3 - y4) 
+               / ((x3 - x4)**2 + (y3 - y4)**2)**1.5
+              )
+    
+    dxdt4 = v_x4
+    dv_xdt4 = (
+      	    - G * m1 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (x4 - x1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (x4 - x2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (x4 - x3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+    dydt4 = v_y4
+    dv_ydt4 = (
+      	    - G * m1 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            - G * m2 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            - G * m3 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+            + k * q4 * q1 / m4 * (y4 - y1) 
+               / ((x4 - x1)**2 + (y4 - y1)**2)**1.5
+            + k * q4 * q2 / m4 * (y4 - y2) 
+               / ((x4 - x2)**2 + (y4 - y2)**2)**1.5
+            + k * q4 * q3 / m4 * (y4 - y3) 
+               / ((x4 - x3)**2 + (y4 - y3)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3,
+            dxdt4, dv_xdt4, dydt4, dv_ydt4)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+q = 1.1 *10**(-13)
+m = 1.1 * 10**(-27)
+a = 0.008
+v = 1 * 10**7
+alpha = np.pi/4
+
+m1 = m
+q1 = -q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = - q
+
+m4 = m
+q4 = q
+
+x10 = a
+v_x10 = - v * np.cos(alpha)
+y10 = a
+v_y10 = v * np.sin(alpha)
+
+x20 = a
+v_x20 = v * np.cos(alpha)
+y20 = - a
+v_y20 = v * np.sin(alpha)
+
+x30 = - a
+v_x30 = v * np.cos(alpha)
+y30 = - a
+v_y30 = - v * np.cos(alpha)
+
+x40 = - a
+v_x40 = - v * np.cos(alpha)
+y40 = a
+v_y40 = - v * np.cos(alpha)
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30,
+      x40, v_x40, y40, v_y40)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(4):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(4):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=70)
+
+plt.axis('equal')
+edge = 5 * a
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('extratask1_d.gif')
\ No newline at end of file
diff --git a/lec1.py b/lec1.py
new file mode 100644
index 000000000..d3c6244c4
--- /dev/null
+++ b/lec1.py
@@ -0,0 +1,153 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 365
+seconds_in_year = 365 * 24 * 60 * 60
+seconds_in_day = 24 * 60 * 60
+years = 0.5
+t = np.linspace(0, years*seconds_in_year, frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+x10 = 149 * 10**9
+v_x10 = 0
+y10 = 0
+v_y10 = 30000
+
+x20 = - 149 * 10**9
+v_x20 = 1
+y20 = 0
+v_y20 = - 30000
+
+x30 = 0
+v_x30 = 15000
+y30 = 149 * 10**9
+v_y30 = 0
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30)
+
+m1 = 1.1 * 10**(30)
+q1 = - 1.1 * 10**(20)
+
+m2 = 2.1 * 10**(30)
+q2 = 2.1 * 10**(20)
+
+m3 = 3.6 * 10**(30)
+q3 = - 3.1 * 10**(20)
+
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(3):
+    balls.append(plt.plot([], [], 'o', color='r'))
+    balls_lines.append(plt.plot([], [], '-', color='r'))
+
+def animate(i):
+    for j in range(3):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=30)
+
+plt.axis('equal')
+edge = 2 * x10
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('lec_13_N_body.gif')
\ No newline at end of file
diff --git a/lec_13_N_body.gif b/lec_13_N_body.gif
new file mode 100644
index 000000000..8f4913b38
Binary files /dev/null and b/lec_13_N_body.gif differ
diff --git a/tasks gifs/task1.gif b/tasks gifs/task1.gif
new file mode 100644
index 000000000..971c4bf81
Binary files /dev/null and b/tasks gifs/task1.gif differ
diff --git a/tasks gifs/task2.gif b/tasks gifs/task2.gif
new file mode 100644
index 000000000..e1e0d43ad
Binary files /dev/null and b/tasks gifs/task2.gif differ
diff --git a/tasks gifs/task3.gif b/tasks gifs/task3.gif
new file mode 100644
index 000000000..97d7417e8
Binary files /dev/null and b/tasks gifs/task3.gif differ
diff --git a/tasks/lab_task1.py b/tasks/lab_task1.py
new file mode 100644
index 000000000..359c79051
--- /dev/null
+++ b/tasks/lab_task1.py
@@ -0,0 +1,125 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 600
+seconds_in_year = 365 * 24 * 60 * 60
+years = 20
+t = np.linspace(0, years*seconds_in_year, frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+x10 = -7.45 * 149 * 10**9
+v_x10 = 0
+y10 = 0
+v_y10 = -8435
+
+x20 = 5.185 * 149 * 10**9
+v_x20 = 0
+y20 = 0
+v_y20 = 24556
+
+x30 = 4.515 * 149 * 10**9
+v_x30 = 0
+y30 = 0
+v_y30 = -21893
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30)
+
+G = 6.67 * 10**(-11)
+M_son = 1.989 * 10**30
+
+m1 = 1.06 * M_son
+
+m2 = 0.96 * M_son
+
+m3 = 0.67 * M_son
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(3):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(3):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=30)
+
+plt.axis('equal')
+edge = 1.8 * abs(x10)
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('task1.gif')
\ No newline at end of file
diff --git a/tasks/lab_task2.py b/tasks/lab_task2.py
new file mode 100644
index 000000000..c05a574d8
--- /dev/null
+++ b/tasks/lab_task2.py
@@ -0,0 +1,132 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 365
+seconds_in_year = 365 * 24 * 60 * 60
+seconds_in_day = 24 * 60 * 60
+years = 0.5
+t = np.linspace(0, years*seconds_in_year, frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+k = 8.98755 * 10**9
+q = 2 * 10**5
+m = 1
+a = 149 * 10**9
+vx = 10**4
+
+m1 = m
+q1 = q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = - q
+
+x10 = -a/2
+v_x10 = vx
+y10 = - np.sin(np.pi/6) * a / 3
+v_y10 = 0
+
+x20 = a/2
+v_x20 = -vx
+y20 = - np.sin(np.pi/6) * a / 3
+v_y20 = 0
+
+x30 = 0
+v_x30 = 0
+y30 = np.sin(np.pi/6) * a * 2 / 3
+v_y30 = 0
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(3):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(3):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=30)
+
+plt.axis('equal')
+edge = 2 * abs(x10)
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('task2.gif')
\ No newline at end of file
diff --git a/tasks/lab_task3.py b/tasks/lab_task3.py
new file mode 100644
index 000000000..2140de17a
--- /dev/null
+++ b/tasks/lab_task3.py
@@ -0,0 +1,156 @@
+import numpy as np
+from scipy.integrate import odeint
+import matplotlib.pyplot as plt
+from matplotlib.animation import FuncAnimation
+
+# Определяем переменную величину
+frames = 365
+seconds_in_year = 365 * 24 * 60 * 60
+seconds_in_day = 24 * 60 * 60
+years = 0.5
+t = np.linspace(0, years*seconds_in_year, frames)
+
+# Определяем функцию для системы диф. уравнений
+def move_func(s, t):
+    (x1, v_x1, y1, v_y1,
+     x2, v_x2, y2, v_y2,
+     x3, v_x3, y3, v_y3) = s
+
+    # Динамика первого тела под влиянием второго и третьего
+    dxdt1 = v_x1
+    dv_xdt1 = (
+      	    - G * m2 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m3 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q2 / m1 * (x1 - x2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (x1 - x3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+              )
+    dydt1 = v_y1
+    dv_ydt1 = (
+      	    - G * m2 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            - G * m2 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+            + k * q1 * q2 / m1 * (y1 - y2) 
+               / ((x1 - x2)**2 + (y1 - y2)**2)**1.5
+            + k * q1 * q3 / m1 * (y1 - y3) 
+               / ((x1 - x3)**2 + (y1 - y3)**2)**1.5
+    	      )
+
+    # Динамика второго тела под влиянием первого и третьего
+    dxdt2 = v_x2
+    dv_xdt2 = (
+      	    - G * m1 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q1 / m2 * (x2 - x1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (x2 - x3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+    	      )
+    dydt2 = v_y2
+    dv_ydt2 = (
+      	    - G * m1 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            - G * m3 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+            + k * q2 * q1 / m2 * (y2 - y1) 
+               / ((x2 - x1)**2 + (y2 - y1)**2)**1.5
+            + k * q2 * q3 / m2 * (y2 - y3) 
+               / ((x2 - x3)**2 + (y2 - y3)**2)**1.5
+              )
+
+    # Динамика третьего тела под влиянием второго и первого
+    dxdt3 = v_x3
+    dv_xdt3 = (
+      	    - G * m1 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q1 / m3 * (x3 - x1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (x3 - x2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+    dydt3 = v_y3
+    dv_ydt3 = (
+      	    - G * m1 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            - G * m2 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+            + k * q3 * q1 / m3 * (y3 - y1) 
+               / ((x3 - x1)**2 + (y3 - y1)**2)**1.5
+            + k * q3 * q2 / m3 * (y3 - y2) 
+               / ((x3 - x2)**2 + (y3 - y2)**2)**1.5
+              )
+
+    return (dxdt1, dv_xdt1, dydt1, dv_ydt1,
+            dxdt2, dv_xdt2, dydt2, dv_ydt2,
+            dxdt3, dv_xdt3, dydt3, dv_ydt3)
+
+# Определяем начальные значения и параметры, 
+# входящие в систему диф. уравнений
+G = 6.67 * 10**(-11)
+k = 8.98755 * 10**9
+q = 10**(-20)
+m = 1
+a = 100
+
+m1 = ((5*k*(q**2))/(G*m) - m)/4
+q1 = q
+
+m2 = m
+q2 = q
+
+m3 = m
+q3 = q
+
+x10 = 0
+v_x10 = 0
+y10 = 0
+v_y10 = 0
+
+x20 = - a
+v_x20 = 0
+y20 = 0
+v_y20 = 0
+
+x30 = a
+v_x30 = 0
+y30 = 0
+v_y30 = 0
+
+s0 = (x10, v_x10, y10, v_y10,
+      x20, v_x20, y20, v_y20,
+      x30, v_x30, y30, v_y30)
+
+# Решаем систему диф. уравнений
+sol = odeint(move_func, s0, t)
+
+# Строим решение в виде графика и анимируем
+fig, ax = plt.subplots()
+
+balls = []
+balls_lines = []
+
+for i in range(3):
+    balls.append(plt.plot([], [], 'o', color=f'C{i}'))
+    balls_lines.append(plt.plot([], [], '-', color=f'C{i}'))
+
+def animate(i):
+    for j in range(3):
+        balls[j][0].set_data([sol[i, 4*j]], [sol[i, 4*j+2]])
+        balls_lines[j][0].set_data([sol[:i, 4*j]], [sol[:i, 4*j+2]])
+
+ani = FuncAnimation(fig, animate, frames=frames, interval=30)
+
+plt.axis('equal')
+edge = 2 * a
+ax.set_xlim(-edge, edge)
+ax.set_ylim(-edge, edge)
+
+ani.save('task3.gif')
\ No newline at end of file