"""
Python code for task 2 in 3-4
Author: Anders Haarberg Eriksen
"""
from matplotlib import pyplot as plt
import numpy as np
 
 
def plot_solution(ax, f, x_0: float, y_0: float, T: float, h: float = 0.1) -> None:
    """
    This is implemented with explicit euler
    """
    curve = [np.array([x_0, y_0])]
    n = int(T / h)
    for i in range(n):
        x, y = curve[i]
        curve.append(curve[i] + h * f(x, y))
    x = [a[0] for a in curve]
    y = [a[1] for a in curve]
    ax.plot(x, y)
 
 
# Vector field - function
f = lambda x, y: np.array([x * (2 - y), 
                           y *(x - 1)])
 
# Plot vector field
x = np.arange(0, 6, 0.25)
y = np.arange(0, 7, 0.25)
 
X, Y = np.meshgrid(x, y)
U, V = f(X, Y)
 
fig, ax = plt.subplots()
ax.quiver(X, Y, U, V, label='_nolegend_')
 
# Plot solutions
T = 6
h = 0.0001
num = 4
 
for i in range(1, num + 1):
    plot_solution(ax, f, i, i, T, h)
 
ax.legend([f'$x_0={i}$, $y_0={i}$' for i in range(1, num + 1)], loc='lower right')
 
plt.show()

import numpy as np
from typing import Callable
 
 
def calculate_integral(f: Callable[[float], float], theta_0: float, theta_1: float, n: int = 1000) -> float:
    """
    Implemented with trapezoid method
    """
    result = 0
    theta_values = np.linspace(theta_0, theta_1, n)
 
    h = (theta_1 - theta_0) / (n - 1)
    for i, theta in enumerate(theta_values[1:]):
        result += (f(theta_values[i]) + f(theta)) / 2 * h
    return result
 
 
if __name__ == '__main__':
    # Function variables
    a = 4
    b = 3
 
    # Norm of derivative of x-vector
    f = lambda theta: np.sqrt(a ** 2 * np.sin(theta) ** 2 + b ** 2 * np.cos(theta) ** 2)
 
    # Integral variables
    theta_0 = 0
    theta_1 = 2 * np.pi
    n = 10000
 
    print(f"The circumference of an ellipse with axis {a = } and {b = } is approximately {calculate_integral(f, theta_0, theta_1, n):.6f}.")

"""
Python code for task 3 in 3-5
Author: Anders Haarberg Eriksen
"""
import numpy as np
from typing import Callable
 
 
def calculate_integral(f: Callable[[float], float], theta_0: float, theta_1: float, n: int = 1000) -> float:
    """
    Implemented with trapezoid method
    """
    result = 0
    theta_values = np.linspace(theta_0, theta_1, n)
 
    h = (theta_1 - theta_0) / (n - 1)
    for i, theta in enumerate(theta_values[1:]):
        result += (f(theta_values[i]) + f(theta)) / 2 * h
    return result
 
 
if __name__ == '__main__':
    # Integrand
    f = lambda theta: (2 + 3 * 4 / 10 * np.cos(theta) * np.sin(theta)) * np.sqrt(16 * np.sin(theta) ** 2 + 9 * np.cos(theta) ** 2)
 
    # Integral variables
    theta_0 = 0
    theta_1 = np.pi / 2
    n = 10000
 
    print(f"The area of the wall is approximately {calculate_integral(f, theta_0, theta_1, n):.6f}.")