Cant exit while loop on Simpson's Rule

Question

I am trying to calculate an integral using Simpson's Rule formula.The catch is that the value of the integral is the one that satisfies the following condition:You find the Ih and Ih/2.If the absolute of (Ih-Ih/2)<error the loop is complete.Otherwise you repeat the process with half the h,which means it calculates the absolute of (Ih/2-Ih/4) and so on and so on.

while True:
    ###Ih part
    h = (b - a) / N
    y1 = np.linspace(a, b, N)
    Ez11 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y1 - fa))) / (
                (aa ** 2 - 2 * ra * R * np.cos(y1 - fa)) ** (3 / 2))
    I11 = (h/3) * (Ez11[0] + 2*sum(Ez11[:N-2:2]) \
            + 4*sum(Ez11[1:N-1:2]) + Ez11[N-1])
    #####Ih/2 part
    h = (b-a)/(2*N)
    y2 = np.linspace(a, b, 2*N)
    Ez22 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y2 - fa))) / (
            (aa ** 2 - 2 * ra * R * np.cos(y2 - fa)) ** (3 / 2))
    print(Ez22)
    I22 = (h/ 3) * (Ez22[0] + 2 * sum(Ez22[:N - 2:2]) \
                     + 4 * sum(Ez22[1:N - 1:2]) + Ez22[N - 1])
    # error condition I1=Ih I2=Ih/2
    if np.abs(I11 - I22) < error:
        break
    else:
        N = 2*N # h/2
        print(np.abs(I11 - I22))

As far as I can tell,my approach should be correct.However the loop goes on and on,never to stop.
My code is as follows:

import numpy as np
from scipy.integrate import simps
import scipy.integrate as integrate
import scipy.special as special


# variables
a = 0
b = np.pi * 2
N = 100
ra = 0.1  # ρα
R = 0.05
fa = 35 * (np.pi / 180)  # φα
za = 0.4
Q = 10 ** (-6)
k = 9 * 10 ** 9
aa = np.sqrt(ra ** 2 + R ** 2 + za ** 2)
error = 0.55 * 10 ** (-8)
h=(b-a)/N
I1 = np.nan
I11 = np.nan

#Simpsons section

############   Ez


#automated Simpson

while True:
    ###Ih part
    y1 = np.linspace(a, b, N)
    Ez1 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y1 - fa))) / (
                (aa ** 2 - 2 * ra * R * np.cos(y1 - fa)) ** (3 / 2))
    print(len(Ez1))
    I1 = simps(Ez1, y1)
    #####Ih/2 part
    y2 = np.linspace(a, b, 2*N)
    Ez2 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y2 - fa))) / (
            (aa ** 2 - 2 * ra * R * np.cos(y2 - fa)) ** (3 / 2))
    I2 = simps(Ez2, y2)
    # error condition I1=Ih I2=Ih/2
    if np.abs(I1 - I2) < error:
        break
    else:
        N *= 2  # h/2

#custom-made Simpson

N = 100


while True:
    ###Ih part
    h = (b - a) / N
    y1 = np.linspace(a, b, N)
    Ez11 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y1 - fa))) / (
                (aa ** 2 - 2 * ra * R * np.cos(y1 - fa)) ** (3 / 2))
    I11 = (h/3) * (Ez11[0] + 2*sum(Ez11[:N-2:2]) \
            + 4*sum(Ez11[1:N-1:2]) + Ez11[N-1])
    #####Ih/2 part
    h = (b-a)/(2*N)
    y2 = np.linspace(a, b, 2*N)
    Ez22 = (np.sqrt(ra ** 2 + R ** 2 - 2 * ra * R * np.cos(y2 - fa))) / (
            (aa ** 2 - 2 * ra * R * np.cos(y2 - fa)) ** (3 / 2))
    print(Ez22)
    I22 = (h/ 3) * (Ez22[0] + 2 * sum(Ez22[:N - 2:2]) \
                     + 4 * sum(Ez22[1:N - 1:2]) + Ez22[N - 1])
    # error condition I1=Ih I2=Ih/2
    if np.abs(I11 - I22) < error:
        break
    else:
        N = 2*N # h/2
        print(np.abs(I11 - I22))

print(I1)
print(I11)

Simpson's Rule is as follows:

After a while it's stuck in this situation

The 5.23 part is the absolute diff of those 2 which shouldnt be that high.