import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm

N = 1
a1 = 0.5 + 1/N
b1 = 1.5 + 1/(N + 1)
r1 = 1 - 1/(N + 2)
a2 = 1.5 + 1/N
b2 = 0.5 + 1/(N + 1)
r2 = 0.5 + 1/(N + 2)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x2 = np.linspace(-1, 1, 100)
y2 = np.linspace(-1, 1, 100)
x2, y2 = np.meshgrid(x2, y2)
z2 = np.zeros(x2.shape)

ax.plot_surface(x2, y2, z2, rstride=4, cstride=4, cmap=cm.ocean) #плоскость


u = np.linspace(0, 2 * np.pi, 100)
v = np.linspace(0, np.pi, 100)

r=0.5

x = r * np.outer(np.cos(u), np.sin(v))
y = r * np.outer(np.sin(u), np.sin(v))
z = 0.5 + r * np.outer(np.ones(np.size(u)), np.cos(v))
ax.plot_surface(x, y, z,  rstride=4, cstride=4, color='r',) #сфера

ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlim3d(-1,1)

plt.show()