diff --git a/boris.py b/boris.py index 650fe2e..a9b5418 100644 --- a/boris.py +++ b/boris.py @@ -4,6 +4,7 @@ import matplotlib.pyplot as plt from matplotlib import rc rc('font',**{'family':'serif','serif':['Computer Modern']}) rc('text', usetex=True) +from scipy.special import ellipe,ellipk e = 1.602176565e-19 #osnovni naboj [C] m_pr = 1.672621777e-27 #masa protona [kg] @@ -53,6 +54,32 @@ def B_zemlja(x): k = -B0*Rz**3/r**5 return k*np.array([3*x[0]*x[2],3*x[1]*x[2],2*x[2]**2-x[0]**2-x[1]**2]) +def B_loop(p, a): + ''' polje magnetne zanke polmera a v xy ravnini, normiran tok 1A + glej: Sinigoj, ELMG polje, stran 200''' + mu = 1e-7 # mu/4/pi + x,y,z = p + r = np.sqrt(x**2+y**2) + R0 = np.sqrt((r+a)**2+z**2)) + m = 4*a*r/( (r+a)**2+z**2) + K = ellipk(m) + E = ellipe(m) + + Br = -mu * z/r/R0 * (2*K - (2-m)*E/(1-m)) + Bz = mu/R0 * (2*K - (2-m*(r+a)/r)*E/(1-m)) + + return np.array([Br*x/r, Br*y/r, Bz]) + + + +def B_bottle(p, a, h, I): + ''' polje magnetne steklenice visine 2h in polmera a + ''' + B1 = B_loop(p-np.array([0.,0.,h]),a)*I + B2 = B_loop(p+np.array([0.,0.,h]),a)*I + return B1+B2 + + def Wk(m,V): return 0.5*m*np.sum(V**2,1)