Magnetic field of a railgun with square rails (maybe)
Aug 10, 2017 11:50:05 GMT
apophys, matterbeam, and 4 more like this
Post by omnipotentvoid on Aug 10, 2017 11:50:05 GMT
Well, I've spent an unreasonable (unexpectadly so) amount of time over the last three days trying to accurately model the acceleration forces on railgun armatures. To do so, I assumed infinitely long rails of a square crossection of side length d through which the current I flows through with a constant distribution. I then tried to use the Biot-Savart law to calculate the magnetic field. The resulting solutions produced disproportionatly strong fields with rather odd scaling for d, despite having the expected form. Being unable to find a fault in the logic behind my formulation of the eqaution, I assumed that I was simply unable to solve the integral. To avoid the integral in the Biot-Savart law, I simplified the problem by looking at the square conductor as a collection if infinitely thin wires, for which the magnetic field is well known and then integrating over all these infentesimal wires. This integral again proved difficult for me to solve (apperantly I'm terrible at math and online integral calculators were producing conflicting results), though I was able to produce a solution that at least appeared reasonable even if the field had a somewhat unexpected form.
I'll post pictures of my calculations as soon as I write them down in a form that is legible.
Here are the results (note. cs-force means cross sectional force):
And heres the python code for it, for anyone thats interested (warning: it's a mess):
As far as the implications go: If armatures would be made to fit the force distribution in the armature, internal forces in the projectile would be severely reduced (perhaps almost totaly negated), allowing for much greater acceleration and thus higher velocities with shorter barrels, thus reducing weight and cost while increasing firepower.
I'll post pictures of my calculations as soon as I write them down in a form that is legible.
Here are the results (note. cs-force means cross sectional force):
And heres the python code for it, for anyone thats interested (warning: it's a mess):
import math
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import numpy as np
%matplotlib notebook
fig1 = plt.figure()
ax1 = fig1.add_subplot(111, projection='3d')
plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0))
plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
plt.ticklabel_format(style='sci', axis='z', scilimits=(0,0))
fig2 = plt.figure()
ax2 = fig2.add_subplot(111, projection='3d')
plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0))
plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
plt.ticklabel_format(style='sci', axis='z', scilimits=(0,0))
n = 100
d=.002
I=1700
ForceTotal = 0
BTotal = 0
Mass = .003
#def BfieldFunction_WRONG1 (xs,ys,d):
# y = ys
# x=xs
#
# A1=((y-.5*d)/(x+1.5*d))+math.sqrt((((y-.5*d)/(x+1.5*d))**2)+1)
# A2=((y+.5*d)/(x+1.5*d))+math.sqrt((((y+.5*d)/(x+1.5*d))**2)+1)
# A3=((y-.5*d)/(x+.5*d))+math.sqrt((((y-.5*d)/(x+.5*d))**2)+1)
# A4=((y+.5*d)/(x+.5*d))+math.sqrt((((y+.5*d)/(x+.5*d))**2)+1)
# B1= ((((x+1.5*d)**2)*(.5*np.log(A1)+(((A1**4)-1)/(8*(A1**2))))-(.5*np.log(A2)+(((A2**4)-1)/(8*(A2**2)))))-(((x+.5*d)**2)*(.5*np.log(A3)+(((A3**4)-1)/(8*(A3**2))))-(.5*np.log(A4)+(((A4**4)-1)/(8*(A4**2))))))
#
# C1=((y-.5*d)/(x-1.5*d))+math.sqrt((((y-.5*d)/(x-1.5*d))**2)+1)
# C2=((y+.5*d)/(x-1.5*d))+math.sqrt((((y+.5*d)/(x-1.5*d))**2)+1)
# C3=((y-.5*d)/(x-.5*d))+math.sqrt((((y-.5*d)/(x-.5*d))**2)+1)
# C4=((y+.5*d)/(x-.5*d))+math.sqrt((((y+.5*d)/(x-.5*d))**2)+1)
# B2= ((((x-1.5*d)**2)*(.5*np.log(C1)+(((C1**4)-1)/(8*(C1**2))))-(.5*np.log(C2)+(((C2**4)-1)/(8*(C2**2)))))-(((x-.5*d)**2)*(.5*np.log(C3)+(((C3**4)-1)/(8*(C3**2))))-(.5*np.log(C4)+(((C4**4)-1)/(8*(C4**2))))))
#
# return -(B1+B2)
#
#def BfieldFunction_WRONG2 (xs,ys,ds):
# a = xs
# b = ys
# d= ds
#
# B1 = (((3*d+2*a)**2*np.log(abs(math.sqrt(4*(0.5*d-b)**2+(3*d+2*a)**2)+2*0.5*d-2*b))-(d+2*a)**2*np.log(abs(math.sqrt(4*(0.5*d-b)**2+(d+2*a)**2)+2*0.5*d-2*b))+2*(0.5*d-b)*(math.sqrt(4*(0.5*d-b)**2+(3*d+2*a)**2)-math.sqrt(4*(0.5*d-b)**2+(d+2*a)**2)))/8)-(((3*d+2*a)**2*np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(3*d+2*a)**2)+2*-0.5*d-2*b))-(d+2*a)**2*np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(d+2*a)**2)+2*-0.5*d-2*b))+2*(-0.5*d-b)*(math.sqrt(4*(-0.5*d-b)**2+(3*d+2*a)**2)-math.sqrt(4*(-0.5*d-b)**2+(d+2*a)**2)))/8)
# B2 = (((3*d-2*a)**2*np.log(abs(math.sqrt(4*(0.5*d-b)**2+(3*d-2*a)**2)+2*0.5*d-2*b))-(d-2*a)**2*np.log(abs(math.sqrt(4*(0.5*d-b)**2+(d-2*a)**2)+2*0.5*d-2*b))+2*(0.5*d-b)*(math.sqrt(4*(0.5*d-b)**2+(3*d-2*a)**2)-math.sqrt(4*(0.5*d-b)**2+(d-2*a)**2)))/8)-(((3*d-2*a)**2*np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(3*d-2*a)**2)+2*-0.5*d-2*b))-(d-2*a)**2*np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(d-2*a)**2)+2*-0.5*d-2*b))+2*(-0.5*d-b)*(math.sqrt(4*(-0.5*d-b)**2+(3*d-2*a)**2)-math.sqrt(4*(-0.5*d-b)**2+(d-2*a)**2)))/8)
# return B1+B2
#candidate 2
def BfieldFunction_V1 (xs,ys,ds,I):
a = xs
b = ys
d= ds
B1 = (np.log(abs(math.sqrt(4*(0.5*d-b)**2+(d+2*a)**2)+2*0.5*d-2*b))-np.log(abs(math.sqrt(4*(0.5*d-b)**2+(3*d+2*a)**2)+2*0.5*d-2*b)))-(np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(d+2*a)**2)+2*-0.5*d-2*b))-np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(3*d+2*a)**2)+2*-0.5*d-2*b)))
B2 = (np.log(abs(math.sqrt(4*(0.5*d-b)**2+(d-2*a)**2)+2*0.5*d-2*b))-np.log(abs(math.sqrt(4*(0.5*d-b)**2+(3*d-2*a)**2)+2*0.5*d-2*b)))-(np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(d-2*a)**2)+2*-0.5*d-2*b))-np.log(abs(math.sqrt(4*(-0.5*d-b)**2+(3*d-2*a)**2)+2*-0.5*d-2*b)))
return ((I*1.2566*(10**(-6)))/(4*math.pi*(d**2)))*(B1+B2)
def BfieldFunction_V2 (xs,ys,ds,I):
a = xs
b = ys
d = ds
B1 =((a*math.sqrt(((4*(0.5*d-b)**2)/((d+2*a)**2))+1)-a*math.sqrt(((4*(0.5*d-b)**2)/((3*d+2*a)**2))+1))/(0.5*d-b))-((a*math.sqrt(((4*(-0.5*d-b)**2)/((d+2*a)**2))+1)-a*math.sqrt(((4*(-0.5*d-b)**2)/((3*d+2*a)**2))+1))/(-0.5*d-b))
B2 =(-(a*math.sqrt(((4*(0.5*d-b)**2)/((d-2*a)**2))+1)-a*math.sqrt(((4*(0.5*d-b)**2)/((3*d-2*a)**2))+1))/(0.5*d-b))-(-(a*math.sqrt(((4*(-0.5*d-b)**2)/((d-2*a)**2))+1)-a*math.sqrt(((4*(-0.5*d-b)**2)/((3*d-2*a)**2))+1))/(-0.5*d-b))
return (-(I*1.2566*(10**(-6)))/(4*math.pi*(d**2)))*(B1+B2)
def BfieldFunction_V3 (xs,ys,ds,I):
a = xs
b = ys
d = ds
B1 =(((d+2*b)*np.log(4*(-1.5*d-a)**2+(d+2*b)**2)+(d-2*b)*np.log(4*(-1.5*d-a)**2+(d-2*b)**2)+4*(-1.5*d-a)*(math.atan((d+2*b)/(2*(-1.5*d-a)))+math.atan((d-2*b)/(2*(-1.5*d-a)))))/4)-(((d+2*b)*np.log(4*(-.5*d-a)**2+(d+2*b)**2)+(d-2*b)*np.log(4*(-.5*d-a)**2+(d-2*b)**2)+4*(-.5*d-a)*(math.atan((d+2*b)/(2*(-.5*d-a)))+math.atan((d-2*b)/(2*(-.5*d-a)))))/4)
B2 =(((d+2*b)*np.log(4*(.5*d -a)**2+(d+2*b)**2)+(d-2*b)*np.log(4*( .5*d-a)**2+(d-2*b)**2)+4*( .5*d-a)*(math.atan((d+2*b)/(2*( .5*d-a)))+math.atan((d-2*b)/(2*( .5*d-a)))))/4)-(((d+2*b)*np.log(4*(1.5*d-a)**2+(d+2*b)**2)+(d-2*b)*np.log(4*(1.5*d-a)**2+(d-2*b)**2)+4*(1.5*d-a)*(math.atan((d+2*b)/(2*(1.5*d-a)))+math.atan((d-2*b)/(2*(1.5*d-a)))))/4)
return ((I*1.2566*(10**(-6)))/(4*math.pi*(d**2)))*(B1-B2)
#candidate 1
def BfieldFunction_V4 (xs,ys,ds,I):
a = xs
b = ys
d = ds
B1 =(((.5*d-b)*(np.log(4*(.5*d-b)**2+(3*d+2*a)**2)-np.log(4*(.5*d-b)**2+(d+2*a)**2))+abs(3*d+2*a)*math.atan((2*(.5*d-b))/abs(3*d+2*a))-abs(d+2*a)*math.atan((2*(.5*d-b))/abs(d+2*a)))/2)-(((-.5*d-b)*(np.log(4*(-.5*d-b)**2+(3*d+2*a)**2)-np.log(4*(-.5*d-b)**2+(d+2*a)**2))+abs(3*d+2*a)*math.atan((2*(-.5*d-b))/abs(3*d+2*a))-abs(d+2*a)*math.atan((2*(-.5*d-b))/abs(d+2*a)))/2)
B2 =(-((.5*d-b)*(np.log(4*(.5*d-b)**2+(3*d-2*a)**2)-np.log(4*(.5*d-b)**2+(d-2*a)**2))+abs(3*d-2*a)*math.atan((2*(.5*d-b))/abs(3*d-2*a))-abs(d-2*a)*math.atan((2*(.5*d-b))/abs(d-2*a)))/2)+(((-.5*d-b)*(np.log(4*(-.5*d-b)**2+(3*d-2*a)**2)-np.log(4*(-.5*d-b)**2+(d-2*a)**2))+abs(3*d-2*a)*math.atan((2*(-.5*d-b))/abs(3*d-2*a))-abs(d-2*a)*math.atan((2*(-.5*d-b))/abs(d-2*a)))/2)
return ((I*1.2566*(10**(-6)))/(4*math.pi*(d**2)))*(B1-B2)
def ForceFunction_V1 (x,y,d,Is):
I = Is
return (I/n)*BfieldFunction_V1(x,y,d,I)
def ForceFunction_V2 (x,y,d,Is):
I = Is
return (I/n)*BfieldFunction_V2(x,y,d,I)
def ForceFunction_V3 (x,y,d,Is):
I = Is
return (I/n)*BfieldFunction_V3(x,y,d,I)
def ForceFunction_V4 (x,y,d,Is):
I = Is
return (I/n)*BfieldFunction_V4(x,y,d,I)
X = np.arange(-.49*d,.5*d,(1/n)*d)
Y = np.arange(-.49*d,.5*d,(1/n)*d)
for y in Y:
for x in X:
ForceTotal = ForceTotal + ForceFunction_V4(x,y,d,I)
ForceTotal = ForceTotal/((n-2)**2)
for y in Y:
for x in X:
BTotal = BTotal + BfieldFunction_V4(x,y,d,I)
BTotal = BTotal/((n-2)**2)
X, Y = np.meshgrid(X,Y)
z = np.array([ForceFunction_V4(x,y,d,I) for x,y in zip(np.ravel(X),np.ravel(Y))])
Z = z.reshape(X.shape)
zB = np.array([BfieldFunction_V4(x,y,d,I) for x,y in zip(np.ravel(X),np.ravel(Y))])
ZB = zB.reshape(X.shape)
ax1.plot_surface(X,Y,Z,cmap=cm.coolwarm)
ax2.plot_surface(X,Y,ZB,cmap=cm.coolwarm)
ax1.set_xlabel('x in m')
ax1.set_ylabel('y in m')
ax1.set_zlabel('CS-Force in N')
ax2.set_xlabel('x in m')
ax2.set_ylabel('y in m')
ax2.set_zlabel('B-field in T')
plt.show()
print('Side length of {}m, Current of {}A, Projectile Mass of {}kg, Center Field of {}T, average B-field of {}T, Total force of {}N and Acceleration of {}m/s^2'.format(d, I, Mass, BfieldFunction_V4(0,0,d,I), BTotal, ForceTotal, ForceTotal/Mass))
As far as the implications go: If armatures would be made to fit the force distribution in the armature, internal forces in the projectile would be severely reduced (perhaps almost totaly negated), allowing for much greater acceleration and thus higher velocities with shorter barrels, thus reducing weight and cost while increasing firepower.
