import numpy
import matplotlib.pyplot as plt

a = 0
b = 1

c = -1
d = 1

def f(x):
    return (3 * x ** 3 * numpy.arcsin(x) + x**2 * numpy.sqrt(1 - x**2) + 2 * numpy.sqrt(1 - x**2) - 2) / (9 * x ** 3)

def integrate(x, n):
    sum = 0
    for t in numpy.linspace((b-a)/n, b, n-1):
        sum += t**2 * numpy.arcsin(x * t)
    return sum * (b-a)/n

def calculate_max_diff(Y1, Y2):
    return max(abs(Y1-Y2))

points_num = 20
N = [5, 10, 20, 50, 100, 200]
Y_approx = []
max_diff = []

X = numpy.linspace(c, d, points_num).reshape(-1, 1)
Y = f(X)

for i in range(len(N)):
    Y_approx.append(integrate(X, N[i]))
    max_diff.append(calculate_max_diff(Y, Y_approx[i]))

fig = plt.figure()
ax1 = fig.add_subplot(111)

ax1.scatter(X, Y, label='original')
for i in range(len(N)):
    ax1.plot(X, Y_approx[i], label=str(N[i]) + ' diff = ' + str(max_diff[i]))
plt.legend(loc='upper left')
plt.show()