Задачу решает код с циклом:
x_old=ones((m,2),dtype='d')
x_new=ndarray((m,2),dtype='d')
...
for i in range(m):
if (x_new[i][0]-xc)**2/a/a+(x_new[i][1]-yc)**2/b/b>(1.+eps):
x_old[i]=x_new[i]
x_old=ones((m,2),dtype='d')
x_new=ndarray((m,2),dtype='d')
...
for i in range(m):
if (x_new[i][0]-xc)**2/a/a+(x_new[i][1]-yc)**2/b/b>(1.+eps):
x_old[i]=x_new[i]
n_a = 1/(a**2)
n_b = 1/(b**2)
n_c = 1.+eps
func = lambda n_new, n_old: (n_new[0]-xc)**2*n_a + (n_new[1]-yc)**2*n_b > n_c ? n_new : n_old
func = lambda n_new, n_old: n_new if (n_new[0]-xc)**2*n_a + (n_new[1]-yc)**2*n_b > n_c else n_old