以下の様なコードですと、p2 が毎回変更されてしまいますので、p2 は最後に生成した ImplicitSeries インスタンスだけになってしまいます。

    p2 = plot_implicit(...)
    p2.extend(p2)

ですので、以下の様に書き変えてみると良いかと思います。

    for i in range(len(points)):
      G = (x-points[i][0])**2 + (y-points[i][1])**2 - weights[i]**2
      p = plot_implicit(G, (x,0,100), (y,0,100), show=False, line_color='r')
      if i == 0:
        p2 = p
      else:
        p2.extend(p)

**追記**

> pointsの座標点も同時に描画したいのですが、ご存知でしょうか。

> sympy だけで処理しようとすると、座標点を線で結んでしまいますね。。。matplotlib で sympy の式をプロットする様にすれば上手く行くかもしれませんが、少し調べてみます。

sympy の `plot_implicit()` の内部処理を眺めると、方程式から決定される矩形領域を塗り潰しています(fill)。その矩形領域の 4 隅の座標は `ImplicitSeries.get_points()` メソッドで取得できるので、取り出した座標データを matplotlib で処理できる様な形式に再構成します。

    import numpy as np
    import matplotlib.pyplot as plt
    from sympy import *
    
    points = np.array([[10,20], [20, 75], [45, 15], [70, 53], [90, 30]])
    weights = np.array([5, 5.1, 5.2, 4.9, 4.8])
    rangex, rangey = (0, 100), (0, 100)
    
    fig = plt.figure()
    ax = fig.add_subplot(1, 1, 1)
    ax.set_xlim(*rangex)
    ax.set_ylim(*rangey)
    ax.set_aspect('equal')

    x, y = symbols("x y")
    for i in range(len(points)):
      px, py = points[i]
      ax.plot(px, py, marker='.', markersize=3, color='blue')
      plt.annotate('(%d, %d)' % (px, py), 
         xy=(px, py), xytext=(-15, 15),
         textcoords='offset points', ha='right', va='bottom',
         bbox=dict(boxstyle='round,pad=0.2', fc='yellow', alpha=0.5),
         arrowprops=dict(arrowstyle = '->', connectionstyle='arc3,rad=0'))
      G = (x-px)**2 + (y-py)**2 - weights[i]**2
      p = plot_implicit(G, (x,)+rangex, (y,)+rangey, show=False)[0].get_points()[0]
      px, py = [], []
      for (ix, iy) in p:
        px.extend([ix.start, ix.start, ix.end, ix.end, None])
        py.extend([iy.start, iy.end, iy.end, iy.start, None])
      ax.fill(px, py, facecolor='red')
    
    plt.show()

[![sympy の実行例][1]][1]


  [1]: https://i.sstatic.net/hotzP.png