Added class and methods for upcoming transmission line.

gprMax/sources.py

@@ -18,11 +18,13 @@
 
 import numpy as np
 
+from gprMax.constants import c, floattype
+from gprMax.grid import Ix, Iy, Iz
 from gprMax.utilities import rvalue
 
 
 class VoltageSource:
-    """Voltage sources."""
+    """The voltage source can be a hard source if it's resistance is zero, i.e. the time variation of the specified electric field component is prescribed. If it's resistance is non-zero it behaves as a resistive voltage source."""
     
     def __init__(self):
         self.polarisation = None
@@ -34,12 +36,11 @@ class VoltageSource:
         self.resistance = None
         self.waveformID = None
 
-    def update_fields(self, abstime, timestep, updatecoeffsE, ID, Ex, Ey, Ez, G):
+    def update_E(self, abstime, updatecoeffsE, ID, Ex, Ey, Ez, G):
         """Updates electric field values for a voltage source.
             
         Args:
             abstime (float): Absolute time.
-            timestep (int): Iteration number.
             updatecoeffsE (memory view): numpy array of electric field update coefficients.
             ID (memory view): numpy array of numeric IDs corresponding to materials in the model.
             Ex, Ey, Ez (memory view): numpy array of electric field values.
@@ -74,7 +75,7 @@ class VoltageSource:
 
 
 class HertzianDipole:
-    """Hertzian dipoles, i.e. normal additive source (current density)."""
+    """The Hertzian dipole is an additive source (electric current density)."""
     
     def __init__(self):
         self.polarisation = None
@@ -85,12 +86,11 @@ class HertzianDipole:
         self.stop = None
         self.waveformID = None
 
-    def update_fields(self, abstime, timestep, updatecoeffsE, ID, Ex, Ey, Ez, G):
+    def update_E(self, abstime, updatecoeffsE, ID, Ex, Ey, Ez, G):
         """Updates electric field values for a Hertzian dipole.
             
         Args:
             abstime (float): Absolute time.
-            timestep (int): Iteration number.
             updatecoeffsE (memory view): numpy array of electric field update coefficients.
             ID (memory view): numpy array of numeric IDs corresponding to materials in the model.
             Ex, Ey, Ez (memory view): numpy array of electric field values.
@@ -116,7 +116,7 @@ class HertzianDipole:
 
 
 class MagneticDipole:
-    """Magnetic dipoles, i.e. current on a small loop."""
+    """The magnetic dipole is an additive source (magnetic current density)."""
     
     def __init__(self):
         self.polarisation = None
@@ -127,12 +127,11 @@ class MagneticDipole:
         self.stop = None
         self.waveformID = None
 
-    def update_fields(self, abstime, timestep, updatecoeffsH, ID, Hx, Hy, Hz, G):
-        """Updates electric field values for a magnetic dipole.
+    def update_H(self, abstime, updatecoeffsH, ID, Hx, Hy, Hz, G):
+        """Updates magnetic field values for a magnetic dipole.
             
         Args:
             abstime (float): Absolute time.
-            timestep (int): Iteration number.
             updatecoeffsH (memory view): numpy array of magnetic field update coefficients.
             ID (memory view): numpy array of numeric IDs corresponding to materials in the model.
             Hx, Hy, Hz (memory view): numpy array of magnetic field values.
@@ -156,3 +155,143 @@ class MagneticDipole:
             elif self.polarisation is 'z':
                 Hz[i, j, k] -= waveform.amp  * waveform.calculate_value(time, G.dt) * (G.dt / (G.dx * G.dy * G.dz))
 
+
+class TransmissionLine:
+    """The transmission line source is a one-dimensional transmission line which is attached virtually to a grid cell."""
+    
+    def __init__(self, G):
+        """
+        Args:
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+        
+        self.polarisation = None
+        self.positionx = None
+        self.positiony = None
+        self.positionz = None
+        self.start = None
+        self.stop = None
+        self.length = None
+        self.resistance = None
+        self.waveformID = None
+        
+        # Coefficients for ABC termination of end of the transmission line
+        self.abcv0 = 0
+        self.abcv1 = 0
+
+        # Spatial step of transmission line
+        self.dl = np.sqrt(3) * c * G.dt
+        
+        # Nodal position of one-way injector excitation in the transmission line
+        self.source = 10
+        
+        # Number of nodes in the transmission line; add nodes to the length to account for position of one-way injector
+        self.nl = rvalue(self.length/self.dl) + self.source
+
+        self.voltage = np.zeros(self.nl, dtype=floattype)
+        self.current = np.zeros(self.nl, dtype=floattype)
+    
+    def update_abc(self, G):
+        """Updates absorbing boundary condition at end of the transmission line.
+            
+        Args:
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+    
+        tmp = (c * G.dt - self.dl) / (c * G.dt + self.dl)
+    
+        self.voltage[0] = (self.voltage[1] - self.abcv0) + self.abcv1
+        self.abcv0 = self.voltage[0]
+        self.abcv1 = self.voltage[1]
+
+    def update_voltage(self, time, G):
+        """Updates voltage values along the transmission line.
+            
+        Args:
+            time (float): Absolute time.
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+            
+        waveform = next(x for x in G.waveforms if x.ID == self.waveformID)
+        
+        # Update all the voltage values along the line
+        for i in range(1, self.nl):
+            self.voltage[i] -= self.resistance * (c * G.dt / self.dl) * (self.current[i] - self.current[i - 1])
+        
+        # Update the voltage at the position of the one-way injector excitation
+        self.voltage[self.source] += (c * G.dt / self.dl) * waveform.amp * waveform.calculate_value(time - 0.5 * G.dt, G.dt)
+
+        # Update ABC before updating current
+        self.update_abc(G)
+
+    def update_current(self, time, G):
+        """Updates current values along the transmission line.
+            
+        Args:
+            time (float): Absolute time.
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+            
+        waveform = next(x for x in G.waveforms if x.ID == self.waveformID)
+        
+        # Update all the current values along the line
+        for i in range(0, self.nl - 1):
+            self.current[i] -= (1 / self.resistance) * (c * G.dt / self.dl) * (self.voltage[i + 1] - self.voltage[i])
+        
+        # Update the current one node before the position of the one-way injector excitation
+        self.current[self.source - 1] += (c * G.dt / self.dl) * waveform.amp * waveform.calculate_value(time - 0.5 * G.dt, G.dt) * (1 / self.resistance)
+
+    def update_E(self, abstime, Ex, Ey, Ez, G):
+        """Updates electric field value in the main grid from voltage value in the transmission line.
+            
+        Args:
+            abstime (float): Absolute time.
+            Ex, Ey, Ez (memory view): numpy array of electric field values.
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+        
+        if abstime >= self.start and abstime <= self.stop:
+            # Set the time of the waveform evaluation to account for any delay in the start
+            time = abstime - self.start
+            i = self.positionx
+            j = self.positiony
+            k = self.positionz
+            
+            self.update_voltage(time, G)
+            
+            if self.polarisation is 'x':
+                Ex[i, j, k] = - self.voltage[self.nl - 1] / G.dx
+
+            elif self.polarisation is 'y':
+                Ey[i, j, k] = - self.voltage[self.nl - 1] / G.dy
+
+            elif self.polarisation is 'z':
+                Ez[i, j, k] = - self.voltage[self.nl - 1] / G.dz
+
+    def update_H(self, abstime, Hx, Hy, Hz, G):
+        """Updates current value in transmission line from magnetic field values in the main grid.
+            
+        Args:
+            abstime (float): Absolute time.
+            Hx, Hy, Hz (memory view): numpy array of magnetic field values.
+            G (class): Grid class instance - holds essential parameters describing the model.
+        """
+        
+        if abstime >= self.start and abstime <= self.stop:
+            # Set the time of the waveform evaluation to account for any delay in the start
+            time = abstime - self.start
+            i = self.positionx
+            j = self.positiony
+            k = self.positionz
+            
+            if self.polarisation is 'x':
+                self.current[self.nl - 1] = Ix(i, j, k, G.Hy, G.Hz, G)
+
+            elif self.polarisation is 'y':
+                self.current[self.nl - 1] = Iy(i, j, k, G.Hx, G.Hz, G)
+
+            elif self.polarisation is 'z':
+                self.current[self.nl - 1] = Iz(i, j, k, G.Hx, G.Hy, G)
+
+            self.update_current(time, G)
+