Python, численные методы, дифференциальные уравнения

Помогите, пожалуйста, понять в чем ошибка.

Необходимо реализовать метод Коуэла. Проблема в сравнении узлов на 2-ой итерации с 1-ой. Т.е.: на первой итерации сравнивались узлы полученные методом Эйлера с методом Коуэла (узлы сравнивались не с 1-го, т.к сам метод подразумевает уточнение узлов с 3-го до n-1), после уточнения узлов снова идем на метод Коуэла, который теперь строится на узлах предыдущей итерации.

Проблема в сравнении 1-ой итерации со второй и так далее, пока разность между ними не станет меньше eps. Код:

import numpy as np

x_k=2
eps=0.001
n_1 = 10
y_0=1
x_0=1
def du (x_0,y_0):
        return  1/(x_0+1)-5/(y_0**2+1)
h = (x_k-x_0) /n_1
#Метод Эйлера
f_arr = [0]*n_1
f_arr[0]=y_0
for i in range (0, n_1, 1):
       f_arr[i]=y_0+h*du (x_0,y_0)
       x1=x_0+h
       x_0=x1
       y_0=f_arr[i]
       f_arr[i] = round(f_arr[i], 6)
       print('el=',f_arr[i])
#метод Коуэла 
x=[0]*n_1
f=[0]*n_1
for i in range (0, n_1, 1):
        yarray =[0]*n_1
        x[i]=x_0+h*i
        f[i]=1/(x1+1)-5/(f_arr[i]**2+1)
yarray =[0]*n_1
for i in range (3, n_1-1):
        yarray[1]=f_arr[1]
        yarray[0]=f_arr[0]
        yarray[i]=yarray[i-1]+h*(-f[i+1]+13*f[i]+13*f[i-1]-f[i-2])/24
        ghgh=0
        while abs(yarray[i]-f_arr[i]) >= eps:
                f[i]=1/(x[i]+1)-5/(yarray[i]**2+1)
                ghgh+=1  
                yarray[i]=yarray[i-1]+h*(-f[i+1]+13*f[i]+13*f[i-1]-f[i-2])/24
                if ghgh>100: break
        yarray[i] = round(yarray[i], 6)
        print('kol=', yarray[i])