<jupyter_start><jupyter_text>Two competing species <NAME><jupyter_code>import numpy as np
import matplotlib.pyplot as plt
import math
plt.style.use("ggplot")
#Runge-Kutt method
def runge_Kutt(f, g, t0, x0, y0, h, b):
    t = t0
    x = x0
    y = y0
    arrayX = np.array([])
    arrayY = np.array([])
    arrayT = np.array([])
    while t<b:
        k1 = g(t, x, y)
        q1 = f(t, x, y)
        
        k2 = g(t + h/2, x + (h*q1)/2, y + (h*k1)/2)
        q2 = f(t + h/2, x + (h*q1)/2, y + (h*k1)/2)
        
        k3 = g(t + h/2, x + (h*q2)/2, y + (h*k2)/2)
        q3 = f(t + h/2, x + (h*q2)/2, y + (h*k2)/2)
        
        k4 = g(t + h, x + h*q3, y + h*k3)
        q4 = f(t + h, x + h*q3, y + h*k3)
        
        y = y + h*(k1 + 2*k2 + 2*k3 + k4)/6
        x = x + h*(q1 + 2*q2 + 2*q3 + q4)/6
        
        arrayX = np.append(arrayX, x)
        arrayT = np.append(arrayT, t)
        arrayY = np.append(arrayY, y)
        
        t += h
    return arrayT, arrayX, arrayY
# a1 = 1
# a2 = 1
# b11 = 0.2
# b12 = 0.3
# b21 = 0.205
# b22 = 0.25
a1 = 1
a2 = 0.75
b11 = 1
b12 = 0.85
b21 = 0.55
b22 = 0.55

def fx(t, x, y):
    return a1*x - b11*x**2 - b12*x*y

def fy(t, x, y):
    return a2*y - b21*x*y - b22*y**2

def fxl(t, x, y):
    return a1*math.log(x) - b11*math.log(x)**2 - b12*math.log(x)*math.log(y)

def fyl(t, x, y):
    return a2*math.log(y) - b21*math.log(x)*math.log(y) - b22*math.log(y)**2


p1, p2, p3 = runge_Kutt(fxl, fyl, 0, 4, 5, 0.01, 100)
print(f"Solution: x = {p2[-1]}, y = {p3[-1]}")<jupyter_output>Solution: x = 1.0000001991111718, y = 3.910382811931909<jupyter_text>Visualization<jupyter_code>fig, ax = plt.subplots(figsize=(11,8))
n = len(p1)
t = np.linspace(0, 100, n)
t1 = np.linspace(0, 14, 20)
plt.plot(t, p2, label="first population")
plt.plot(t, p3, label="second population")

plt.legend(prop={'size': 12})
plt.show()
fig, ax = plt.subplots(figsize=(11, 8))
for i in range(15):
    p1, p2, p3 = runge_Kutt(fx, fy, 0, math.cos(i*np.pi/60), math.sin(i*np.pi/60), 0.01, 100)
    plt.plot(p2, p3,c='purple')

#stationary points
plt.scatter(0, 0, c='black', label="x = 0, y = 0")
plt.scatter(a1/b11, 0, c='red', label="x = a1/b11, y = 0")
plt.scatter(0, a2/b22, c='blue', label="x = 0, y = a2/b22")
plt.scatter((b12*a2 - b22*a1)/(b21*b12 - b22*b11), (a1*b21 - a2*b11)/(b12*b21 - b22*b11), c='green')
ax.set_xlim(-0.05, a1/b11+0.05)
ax.set_ylim(-0.05, a2/b22+0.05)
plt.legend()    
plt.show()
fig, ax = plt.subplots(figsize=(12,8))
X, Y = np.mgrid[0:a1/b11:16j,
                0:a2/b22:16j]

U = a1*X - b11*X**2 - b12*X*Y
V = a2*Y - b21*X*Y - b22*Y**2
#vector field
ax.quiver(X, Y, U, V, color='black')

#stationary points
plt.scatter(0, 0, c='black', label="x = 0, y = 0")
plt.scatter(a1/b11, 0, c='red', label="x = a1/b11, y = 0")
plt.scatter(0, a2/b22, c='blue', label="x = 0, y = a2/b22")
plt.scatter((b12*a2 - b22*a1)/(b21*b12 - b22*b11), (a1*b21 - a2*b11)/(b12*b21 - b22*b11), c='green')
ax.set_xlim(-0.05, a1/b11+0.05)
ax.set_ylim(-0.05, a2/b22+0.05)
plt.legend()    
plt.show()<jupyter_output><empty_output><jupyter_text>x1 = 0, x2 = 0<jupyter_code>x1 = 0
x2 = 0
jacobi_matrix = np.array([[a1-2*b11*x1-b12*x2, -b12*x1], [-b21*x2, a2-2*b22*x2-b21*x1]])
jacobi_matrix
D0 = (jacobi_matrix.trace())**2-4*(jacobi_matrix[0][0]*jacobi_matrix[1][1]-(jacobi_matrix[0][1]*jacobi_matrix[1][0]))
D0<jupyter_output><empty_output><jupyter_text>неустойчивый узел  x1 = a1/b11, x2 = 0<jupyter_code>x1 = a1/b11
x2 = 0
jacobi_matrix = np.array([[a1-2*b11*x1-b12*x2, -b12*x1], [-b21*x2, a2-2*b22*x2-b21*x1]])
jacobi_matrix
D = (jacobi_matrix.trace())**2-4*(jacobi_matrix[0][0]*jacobi_matrix[1][1]-(jacobi_matrix[0][1]*jacobi_matrix[1][0]))
D<jupyter_output><empty_output><jupyter_text>устойчивый узел  x1 = 0, x2 = a2/b22<jupyter_code>x1 = 0
x2 = a2/b22
jacobi_matrix1 = np.array([[a1-2*b11*x1-b12*x2, -b12*x1], [-b21*x2, a2-2*b22*x2-b21*x1]])
jacobi_matrix1
D1 = (jacobi_matrix1.trace())**2-4*(jacobi_matrix1[0][0]*jacobi_matrix1[1][1]-(jacobi_matrix1[0][1]*jacobi_matrix1[1][0]))
D1<jupyter_output><empty_output><jupyter_text>устойчивый узел  x1 = (a2b12-a1b22)/(b12b21-b22b11), x2 = (a1b21-a2b11)/(b12b21-b22b11)<jupyter_code>x1 = (a2*b12-a1*b22)/(b12*b21-b22*b11)
x2 = (a1*b21-a2*b11)/(b12*b21-b22*b11)
jacobi_matrix2 = np.array([[a1-2*b11*x1-b12*x2, -b12*x1], [-b21*x2, a2-2*b22*x2-b21*x1]])
jacobi_matrix2
D2 = (jacobi_matrix1.trace())**2-4*(jacobi_matrix1[0][0]*jacobi_matrix1[1][1]-(jacobi_matrix1[0][1]*jacobi_matrix1[1][0]))
D2<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.22726261952412719, 'fraction_numerical': 0.13564598621303092, 'mean_word_length': 3.6228160328879753, 'pattern_counts': {'":': 0, '<': 26, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '50732', 'license': 'MIT', 'n_tokens_mistral': 2654, 'n_tokens_neox': 2389, 'n_words': 393, 'path': 'Mathematical_modeling/two_competing_species.ipynb', 'repo_name': 'kirillzx/Math-projects', 'size': '4506'}