Magnetic Field and Vector Potential from Solenoid
2025-4-28 (R.2)
Kiyoshi Yoshida
kiyoshi.yoshida@eagle.ocn.ne.jp
Radius a
(m)
Radial Width Wr
(m)
Axial Length Wz
(m)
Current
(MA)
Radial Position
(m)
Axial Position
(m)
Calculate!
from pyscript import document, display import js from numpy import abs, pi, sqrt, log, arctan, sin, cos from scipy import integrate small_angl = 1.00e-6 small_dist = 1.00e-6 small_value = 1.00e-12 def sgn(x): if x>=0.0: return 1.0 else: return -1.0 def ba0zs(r1, r2, lg, rp, zp): r0 = 0.5*(r1 + r2) s= 0.0 for i in range (4): if i==0: r = r1; z = - 0.5*lg elif i==1: r = r1; z = 0.5*lg elif i==2: r = r2; z = 0.5*lg elif i==3: r = r2; z = - 0.5*lg zj = zp - z rj = sqrt(rp**2 + r**2 - 2.0*rp*r + zj**2) sg = sgn(z*(r0 - r)) s += sg*zj*log(rj + r) b= 2.0*pi*s return b def bar(t,r1,r2,lg,xp,zp): r0 = 0.5*(r1 + r2) co = cos(t) s = 0.0 for i in range(4): if i==0: r = r1; z = -0.5*lg elif i==1: r = r1; z = 0.5*lg elif i==2: r = r2; z = 0.5*lg elif i==3: r = r2; z = -0.5*lg if abs(z) < small_dist: z = small_dist # to avoid divergence zj = zp - z rj = sqrt(xp**2 + r**2 - 2.0*xp*r*co + zj**2) sg = sgn(z*(r0 - r) ) s += sg*(rj + xp*co*log(rj + r - xp*co))*co s = -s return s def baz(t, r1, r2, lg, xp, zp): r0 = 0.50*(r1 + r2) co = cos(t) si = sin(t) s = 0.0 for i in range(4): if i==0: r = r1; z = -0.5*lg elif i==1: r = r1; z = 0.5*lg elif i==2: r = r2; z = 0.5*lg elif i==3: r = r2; z = -0.5*lg if abs(z) < small_dist: z = small_dist # to avoid divergence zj = zp - z rj = sqrt( xp**2 + r**2 - 2.0*xp*r*co + zj**2) a = zj*(r - xp*co)/(rj*xp*si) sg = sgn(z*(r0 - r) ) s += sg*(xp*si*arctan(a) - zj*log(rj + r - xp*co) + xp*co*log(rj + zj)) s = -s return s def aat(t, r1, r2, lg, xp, zp): r0 = 0.50*(r1 + r2) co = cos(t) si = sin(t) s = 0.00 for i in range(4): if i==0: r = r1; z = -0.5*lg elif i==1: r = r1; z = 0.5*lg elif i==2: r = r2; z = 0.5*lg elif i==3: r = r2; z = -0.5*lg zj = zp - z rj = sqrt(xp**2 + r**2 - 2.0*xp*r*co + zj**2) b = ((xp*co - r)*rj - xp**2 - r**2 + 2.0*xp*r*co)/(zj*xp*si) sg = sgn(z*(r0 - r)) s += sg*(zj*rj*0.50 + zj*xp*co*log(rj + r - xp*co) \ + (r**2*0.50 + xp**2*0.50 - xp**2*co**2) * log(rj + zj) + \ + xp**2*si*co*arctan(b))*co return s def ba_sole(rp, zp, ra, wr, lz, cd): r1 = ra - wr*0.5 r2 = ra + wr*0.5 cc = 1.0e-7*cd t1 = small_angl t2 = 2.0*pi - small_angl e = 0.0 if r1 > r2: print('**ERROR at basole(), r1 < r2: ', r1, r2) br = bz = 0.0 return br, bz if(abs(rp) < small_dist) : at = 0.0 br = 0.0 bz = cc * ba0zs(r1, r2, lz, rp, zp) elif abs(zp) < small_dist: a, e = integrate.quad(aat, t1, t2, args=(r1, r2, lz, rp, zp)) at = cc * a br = 0.0 b, e = integrate.quad(baz, t1, t2, args=(r1, r2, lz, rp, zp)) bz = cc * b else: # at any place a, e = integrate.quad(aat, t1, t2, args=(r1, r2, lz, rp, zp)) at = cc * a b, e = integrate.quad(bar, t1, t2, args=(r1, r2, lz, rp, zp)) br = cc * b b, e = integrate.quad(baz, t1, t2, args=(r1, r2, lz, rp, zp)) bz = cc * b return br, bz, at def bacalc(event): ra = float(js.document.getElementById("radius").value) wr = float(js.document.getElementById("wr").value) wz = float(js.document.getElementById("wz").value) ia = float(js.document.getElementById("current").value)*1.0e6 rp = float(js.document.getElementById("radial").value) zp = float(js.document.getElementById("axial").value) ja = ia/(wr*wz) br, bz, at = ba_sole(rp, zp, ra, wr, wz, ja) js.document.getElementById("msg1").innerHTML = "Br (T): "+ str(br) js.document.getElementById("msg2").innerHTML = "Bz (T): "+ str(bz) js.document.getElementById("msg3").innerHTML = "At (Tm): "+ str(at)
