python绘制空间自由电荷电场

            
              import numpy as np
from mayavi import mlab
from sympy import *
from mayavi.mlab import *

q = 1 #电荷量
e0 = 4*3.14*8.85 #4×pi×ε0 放大10的12次方

def cacualte():
    x, y, z = np.mgrid[-1:1:20j, -1:1:20j,-1:1:20j]#x y z网格化，形成填充三维坐标数组

    r=np.sqrt(x**2+y**2+z**2)#与原点距离
   
    E_scale = q/(e0*(r**2))#电场大小
    Exoy = E_scale*(np.sqrt(x**2+y**2)/r)#电场在XOY投影大小
    Ex = Exoy*(x/np.sqrt(x**2+y**2))#电场在X分量大小
    Ey = Exoy*(y/np.sqrt(x**2+y**2))#电场在Y分量大小
    Ez = E_scale*(z/r)#电场在Z分量大小
    
    u = Ex #三坐标分量赋值
    v = Ey
    w = Ez
    
    scalars = x * x * 0.5 + y * y + z * z * 2.0#切片范围
    ball = mlab.points3d(np.array(0), np.array(0.), np.array(0),color=(1,0,0),extent = [-0.01,0.01,-0.01,0.01,-0.01,0.01])#插入小球代表电荷
    obj1 = volume_slice(x,y,z,scalars, plane_orientation='x_axes',extent = [-1,1,-1,1,1,1])#插入切片以观察矢量场
    obj = quiver3d(x, y, z, u, v, w, mask_points=1)#绘制3d矢量场
    return obj
    
cacualte()