antenna_field.py

import numpy as np
 
mu_0 = 4 * np.pi * 10e-7
 
def field_contribution(x, y, pos, I, r):
    v = np.asarray([x,y])
    v = v - pos
    if np.abs(np.sum(v*v)) < r*r:
        return np.asarray([0,0])
    factor = mu_0 / (2 * np.pi * np.dot(v,v))
    return np.array([ I * v[1] * factor, -I * v[0] * factor ])
        
def get_x(val):
    return val[0]
 
def get_y(val):
    return val[1]
 
def limit_vec(c, limit):
    abs = np.sqrt(np.sum(c*c))
    if abs > limit:
        return np.zeros(c.shape)
    else:
        return c
 
def limit_val(s, limit):
    if s > limit:
        return s-s
    else:
        return s
    
vgx = np.vectorize(get_x)
vgy = np.vectorize(get_y)
vsum = np.vectorize(np.sum)
field = np.vectorize(field_contribution, otypes=[np.ndarray], excluded=['pos','I', 'r'])
limit = np.vectorize(limit_vec, otypes=[np.ndarray], excluded=['limit'])
limits = np.vectorize(limit_val, excluded=['limit'])

def calc_field(d, I, c_w, c_h, npc, r):
    xs = np.linspace(-c_w, c_w, num=2 * c_w * npc + 1)
    ys = np.linspace(-c_h, c_h, num=2 * c_h * npc + 1)
 
    X,Y = np.meshgrid(xs, ys)
 
    c1 = field(X, Y, pos=[d, 0], I=I, r=r)
    c2 = field(X, Y, pos=[-d, 0], I=-I, r=r)
    c = c1 + c2
    
    return X,Y,c

def coilInField(c, width):
    f = np.zeros(c.shape)
    yf = vgy(c)
    w = width // 2
    for x in range(w, Y.shape[1] - w):
        for y in range(0, Y.shape[0]):
            f[y][x] = np.abs(np.sum(yf[y][x-w : x+w]))
    return f / width

X,Y,c = calc_field(0.8, 1, 2, 1, 300, 0.025)
#np.save('C:/work/antenna_field/field.npy', c)
#del c
#c = np.load('C:/work/antenna_field/field.npy', 'r')

f = coilInField(c, 5)
#np.save('C:/work/antenna_field/coil_strength.npy', f)
#del f
#f = np.load('C:/work/antenna_field/coil_strength.npy', 'r')

#lc = limit(c, 1e-5)
#lf = limit(f, 10)

import matplotlib
matplotlib.use('SVG')
import matplotlib.pyplot as plt

plt.figure(figsize=[10.4, 22.22])
plt.axis('off')
plt.contour(X, Y, f, [12e-6,4.1e-6,3.9e-6,6e-6,3e-6,1e-6])
plt.gca().set_aspect('equal', adjustable='box')

plt.savefig('C:/work/antenna_field/xxx.svg', bbox_inches='tight')#, dpi=300)

#unused plots
#plt.figure()
#QP = plt.quiver(X, Y, vgx(lc), vgy(lc))
#plt.quiver(X, Y, np.zeros(lc.shape), vgy(lc))
#plt.contour(X,Y,abs,[30e-6,20e-6,15e-6,10e-6,8e-6,6e-6,5e-6,4.25e-6,4.5e-6,4e-6,3e-6,2e-6,1e-6])   
#plt.gca().set_aspect('equal', adjustable='box')