eclipse.py

d = np.linspace(0, 2, num=10000)
A = 1 - 2 /np.pi * ((1 - d/2)*np.sqrt(d - d**2/4) + np.arctan((1-d/2)/(np.sqrt(d-d**2/4))))

fig, axes = plt.subplots(1, 3, figsize = (16, 3))
axes[0].plot(d/2, A)
axes[0].grid()
axes[1].plot(d/2, np.gradient(A, d/2))
axes[2].plot(d/2, np.gradient(np.gradient(A, d/2), d/2))
for ax in axes:
    ax.set_xlim(0, 1)
axes[0].set_ylim(0, 1)
axes[1].set_ylim(ymin=0)
axes[2].set_ylim(ymin=0)
plt.show()