Merged tests with tokamak branch

.github/workflows/main.yaml

@@ -0,0 +1,43 @@
+name: CI
+
+# We can specify which Github events will trigger a CI build
+on: push
+
+# now define a single job 'build' (but could define more)
+jobs:
+
+  build:
+
+    strategy:
+      matrix:
+        os: [ubuntu-latest, macos-latest, windows-latest]
+        python-version: ["3.9", "3.10"]
+
+    runs-on: ${{ matrix.os }}
+
+    # a job is a seq of steps
+    steps:
+
+    # Next we need to checkout out repository, and set up Python
+    # A 'name' is just an optional label shown in the log - helpful to clarify progress - and can be anything
+    - name: Checkout repository
+      uses: actions/checkout@v2
+
+    - name: Set up Python 3.9
+      uses: actions/setup-python@v2
+      with:
+        python-version: ${{ matrix.python-version }}
+
+    - name: Install Python dependencies
+      run: |
+        python3 -m pip install --upgrade pip
+        pip3 install -r requirements.txt
+
+
+    - name: Test with PyTest
+      run: |
+        python -m pytest --cov=fusion_toolbox tests/
+
+    - name: Check style with Pylint
+      run: |
+        python3 -m pylint --fail-under=0 --reports=y fusion_toolbox

fusion_toolbox/Bfield.py

@@ -38,11 +38,12 @@ def B(self,xyz):
         for i,pt in enumerate(self.pts[:-1]): #Skip last point because it repeats
             B_out += self.BGreen(xyz,pt,drs[i])
 
-        return B_out
+        return B_out*self.I
 
     def BGreen(self, xyz_samples, xyz_center, dl):
         """Evaluates the B field at all sample locations due to a wire segment. Uses the formula:
         dB = mu_0/4pi * I * dl x (r-r_0)/|r-r_0|^3
+        The current I 
 
         Parameters
         ----------
@@ -58,8 +59,9 @@ def BGreen(self, xyz_samples, xyz_center, dl):
         """
         _r = xyz_samples - xyz_center[None,:]
         Bvecs = np.cross(dl,_r)/np.linalg.norm(_r,axis=1)[:,None]**3
-        return 1e-7 * self.I * Bvecs #mu_0/4pi = 1e-7 H/m
-    def plot(self):
+        return 1e-7 * Bvecs #mu_0/4pi = 1e-7 H/m
+
+    def plot_contour(self):
         """
         Displays a 3D plot of the coil
         """
@@ -80,6 +82,25 @@ def plot(self):
         ax.set_box_aspect(np.ptp(limits, axis = 1))
 
         plt.show()
+    def plot_B_slice(self,axis,r,w1,w2,n1,n2):
+        ax = plt.figure().add_subplot()
+        #ax = plt.figure().add_subplot(projection='3d')
+        X,Y,Z = gen_plane_pts(axis,r,w1,w2,n1,n2)
+        if axis==0:
+            C1,C2 = Y,Z
+        elif axis==1:
+            C1,C2 = X,Z
+        elif axis==2:
+            C1,C2 = X,Y
+        Bpts = np.vstack((X.ravel(),Y.ravel(),Z.ravel())).T
+        B_samples = self.B(Bpts).T
+        B_clean = np.delete(B_samples,axis,axis=0).T
+        B_2D = np.zeros((n2,n1,2))
+        for i,B in enumerate(B_clean):
+            B_2D[(i//n1)%n2,i%n1] = B
+        ax.streamplot(C1,C2,B_2D[:,:,0],B_2D[:,:,1],density=1.5)
+        ax.set_aspect(1)
+        plt.show()
 
 class PFCoil(Coil):
     def __init__(self,R,Z0,I,npts=100):
@@ -181,3 +202,38 @@ def make_TFset(self,phi,I):
 
 #coils = CMod.make_PFset(4,1,1)
 #coils = CMod.make_PFset(4,-1,1)
+def gen_plane_pts(axis,r,w1,w2,n1,n2):
+    """
+    Generates a 2D grid of points on a plane in 3D which is normal to one of the three axes.
+
+    Parameters
+    ----------
+    axis : int
+        Which axis the plane should be normal to. Valid values are 0, 1, or 2 for x, y, or z.
+    r : np.ndarray of shape (3)
+        Displacement from origin of the center of the plane
+    w1 : float
+        Physical width of first dimension of grid
+    w2 : float
+        Physical width of second dimension of grid
+    n1 : int
+        Number of grid points along side 1 of the grid.
+    n2 : int
+        Number of grid points along side 2 of the grid.
+    """
+    if axis==0:
+        yspan = np.linspace(r[1]-w1/2,r[1]+w1/2,n1)
+        zspan = np.linspace(r[2]-w2/2,r[2]+w2/2,n2)
+        Y,Z = np.meshgrid(yspan,zspan)
+        X = np.zeros_like(Y)
+    elif axis==1:
+        xspan = np.linspace(r[1]-w1/2,r[1]+w1/2,n1)
+        zspan = np.linspace(r[2]-w2/2,r[2]+w2/2,n2)
+        X,Z = np.meshgrid(xspan,zspan)
+        Y = np.zeros_like(Z)
+    elif axis==2:
+        xspan = np.linspace(r[1]-w1/2,r[1]+w1/2,n1)
+        yspan = np.linspace(r[2]-w2/2,r[2]+w2/2,n2)
+        X,Y = np.meshgrid(xspan,yspan)
+        Z = np.zeros_like(X)
+    return X,Y,Z

requirements.txt

@@ -0,0 +1,5 @@
+matplotlib>=3.0.0
+numpy>=1.0.0
+pytest>=7.0.0
+pytest-cov>=4.0.0
+pylint>=2.0.0

tests/test_empty.py

@@ -0,0 +1,4 @@
+
+
+def test_empty():
+    pass

tests/test_magnetics.py

@@ -0,0 +1,45 @@
+import numpy as np
+import numpy.testing as npt
+from math import pi
+from fusion_toolbox.Bfield import *
+import pytest
+
+MU_0 = 4*pi*1E-7 # T * m / A
+
+def generate_pf_coil_xyz(R, N = 100, z = 0): 
+    '''
+    Get xyz coordinates of points along a pf coil 
+    R: radius [m]
+    N: number of points 
+    z: pf coil height [m]
+    Returns: xyz - np.ndarray [N,3]
+    '''
+    thetas =  np.linspace(0, 2*pi, N)
+    xyz = np.zeros((N,3))
+    xyz[:,0] = R * np.cos(thetas)
+    xyz[:,1] = R * np.sin(thetas)
+    xyz[:,2] = z
+    return xyz
+
+def analytic_B_center_of_pf_coil(I, R): 
+    return np.array([0,0,MU_0 * I / (2 * R)])
+
+@pytest.mark.parametrize(
+    "I, R",
+    [[1E3,0.3], [1E5, 0.02], [1E7, 0.004], [1E2, 1]
+    ])
+def test_B_center_of_pf_circular_coil(I, R): 
+    '''
+    I: current [A]
+    R: coil radius [m]
+    Checks the calculated field at the center of a circular coil against the analytic solution
+    '''
+
+    # Generate test coil and calculate the field at the center
+    xyz = generate_pf_coil_xyz(R, N = int(1E4))
+    test_coil = Coil(xyz, I)
+    B_center = test_coil.B([np.array([0,0,0])])
+
+    # Compare to analytic calculation
+    B_analytic = [analytic_B_center_of_pf_coil(I,R)]
+    npt.assert_almost_equal(B_center, B_analytic, decimal = 4)