conducting_tube.py

#!/usr/bin/env python3

# simulate the field inside a hollow conductor with N number of
# equally space thin wire at a distance R from the center.

from PIL import Image
from numpy import arange
import math

#-----------------------------------------------------------------------------    
class Vector(object):
    def __init__(self, x, y):
        self.x = x
        self.y = y
        
    def add(self, other): return Vector(self.x + other.x, self.y + other.y)
    def negate(self): return Vector(-self.x, -self.y)
    def scalar_mul(self, scalar): return Vector(scalar * self.x, scalar * self.y)
    def subtract(self, other): return self.add(other.negate())
    def mag(self): return (self.x**2 + self.y**2)**.5
    def to_unit(self): return self.scalar_mul(1/self.mag())
    
    def dist(self, other):
        dx = self.x - other.x
        dy = self.y - other.y
        return (dx*dx + dy*dy)**.5
    
    def add_update(self, other):
        # in place
        self.x += other.x
        self.y += other.y

    def cross_hatz(self):
        a1, a2, a3 = self.x, self.y, 0
        b1, b2, b3 = 0, 0, 1

        x = a2*b3 - a3*b2
        y = a3*b1 - a1*b3
        return Vector(x, y)
        
    def __repr__(self):
        return "<%f, %f>" % (self.x, self.y)

#-----------------------------------------------------------------------------    
class Field(object):    
    def __init__(self, size):
        self.center = Vector(size/2, size/2)
        self.point_vecs = {} # Map: Point -> B Vector
        
        # init vecs.
        for x in range(size):
            for y in range(size):
                # initialize these points with B vecs of <0, 0>
                self.point_vecs[Vector(x, y)] = Vector(0, 0)
                
    def query_wire(self, wire):
        print("processing wire @: ", wire.pos)
        fudge_factor = 1e12
        mu0 = 4 * math.pi * 1e-7 * fudge_factor

        # superposition!
        # calculate the associated B vector and add it to the existing
        # vector
        
        # for each field point:
        for fp in self.point_vecs:
            # find the distance between fp and the wire.       
            d = fp.dist(wire.pos)
            if d < 5: continue # avoid division by zero. and large
                                # magnitudes

            # get the B vector
            magB = (mu0 * wire.current) / (2 * math.pi * d)            
            r = fp.subtract(wire.pos)
            vecB = r.cross_hatz().to_unit().scalar_mul(magB)

            self.point_vecs[fp].add_update(vecB)
            
    def setup_color_map(self):
        # find the min,max magnitudes and determine a suitable
        # gradient for the field magnitude.
        self.min_mag = min(x.mag() for x in self.point_vecs.values())
        self.max_mag = max(x.mag() for x in self.point_vecs.values())

    def get_color(self, point):        
        x = self.point_vecs[point].mag()
        a = self.min_mag
        b = self.max_mag        
        d = 200
        color = int((x/b)*d)        
        if color %2 == 0: return (color+45, color+45, color+45)
        return (color+55, color+55, color+55)
        
    def paint(self, img):
        self.setup_color_map()
        for fp in self.point_vecs:
            #print(fp, self.get_color(fp))
            color = self.get_color(fp)
            img.putpixel((fp.x, fp.y), color)


#-----------------------------------------------------------------------------    
class Wire(object):
    def __init__(self, pos, current):
        self.pos = pos
        self.current = current

        
#-----------------------------------------------------------------------------    
def build_wires(num_wires, radius, center):
    wires = []
    # divide 2*pi radians into num_wires pieces.
    # start at theta = 0
    
    dtheta = (2 * math.pi) / num_wires
    for angle in arange(0, 2 * math.pi, dtheta):
        # these are pixel locations in screen coordinates.
        # +y points down.
        x = center.x + radius * math.cos(angle)
        y = center.y + radius * math.sin(angle)
        pos = Vector(x, y)
        current = 100.234 # Amps? 
        wire = Wire(pos, current)
        wires.append(wire)        
    return wires

def make_image(num_wires):
    # top left is <0, 0> using screen coordinates.    
    R = 100
    size = 2 * (100 + R)
    fld = Field(size)    
    center = Vector(size/2, size/2)    
    wires = build_wires(num_wires, R, center)

    for w in wires: fld.query_wire(w)
    
    img = Image.new("RGB", (size, size))
    fld.paint(img)
    img.save("hollow-tube-field-" + str(num_wires) + ".jpg", "JPEG")

def main():
    num_wires = 5
    make_image(num_wires)
    # for n in range(1, 50):
    #     make_image(n)
    
if __name__ == "__main__":
    main()