Repository files navigation

A 3D electromagnetic FDTD simulator written in Python. The FDTD simulator hasan optional PyTorch backend, enabling FDTD simulations on a GPU.

Thefdtd-library can be installed withpip:

pip install fdtd

The development version can be installed by cloning the repository

git clone http://github.com/flaport/fdtd

and linking it with pip

pip install -e fdtd

Development dependencies can be installed with

pip install -e fdtd[dev]

• python 3.6+
• numpy
• scipy
• matplotlib
• tqdm
• pytorch (optional)

All improvements or additions (for example new objects, sources or detectors) arewelcome. Please make a pull-request 😊.

read the documentation here:https://fdtd.readthedocs.org

Thefdtd library is simply imported as follows:

importfdtd

Thefdtd library allows to choose a backend. The"numpy" backend is thedefault one, but there are also several additional PyTorch backends:

• "numpy" (defaults to float64 arrays)
• "torch" (defaults to float64 tensors)
• "torch.float32"
• "torch.float64"
• "torch.cuda" (defaults to float64 tensors)
• "torch.cuda.float32"
• "torch.cuda.float64"

For example, this is how to choose the"torch" backend:

fdtd.set_backend("torch")

In general, the"numpy" backend is preferred for standard CPU calculationswith"float64" precision. In general,"float64" precision is alwayspreferred over"float32" for FDTD simulations, however,"float32" mightgive a significant performance boost.

The"cuda" backends are only available for computers with a GPU.

The FDTD grid defines the simulation region.

# signaturefdtd.Grid(shape:Tuple[Number,Number,Number],grid_spacing:float=155e-9,permittivity:float=1.0,permeability:float=1.0,courant_number:float=None,)

A grid is defined by itsshape, which is just a 3D tuple ofNumber-types(integers or floats). If the shape is given in floats, it denotes the width,height and length of the grid in meters. If the shape is given in integers, itdenotes the width, height and length of the grid in terms of thegrid_spacing. Internally, these numbers will be translated to three integers:grid.Nx,grid.Ny andgrid.Nz.

Agrid_spacing can be given. For stability reasons, it is recommended tochoose a grid spacing that is at least 10 times smaller than thesmallestwavelength in the grid. This means that for a grid containing a source withwavelength1550nm and a material with refractive index of3.1, therecommended minimumgrid_spacing turns out to be50nm

For thepermittivity andpermeability floats or arrays with the followingshapes

• (grid.Nx, grid.Ny, grid.Nz)
• or(grid.Nx, grid.Ny, grid.Nz, 1)
• or(grid.Nx, grid.Ny, grid.Nz, 3)

are expected. In the last case, the shape implies the possibility for differentpermittivity for each of the major axes (so-calleduniaxial orbiaxialmaterials). Internally, these variables will be converted (for performancereasons) to their inversesgrid.inverse_permittivity array and agrid.inverse_permeability array of shape(grid.Nx, grid.Ny, grid.Nz, 3). Itis possible to change those arrays after making the grid.

Finally, thecourant_number of the grid determines the relation between thetime_step of the simulation and thegrid_spacing of the grid. If not given,it is chosen to be the maximum number allowed by theCourant-Friedrichs-LewyCondition:1 for1D simulations,1/√2 for2D simulations and1/√3 for3Dsimulations (the dimensionality will be derived by the shape of the grid). Forstability reasons, it is recommended not to change this value.

grid=fdtd.Grid(shape= (25e-6,15e-6,1),# 25um x 15um x 1 (grid_spacing) --> 2D FDTD)print(grid)

Grid(shape=(161,97,1), grid_spacing=1.55e-07, courant_number=0.70)

An other option to locally change thepermittivity orpermeability in thegrid is to add anObject to the grid.

# signaturefdtd.Object(permittivity:Tensorlike,name:str=None)

An object defines a part of the grid with modified update equations, allowingto introduce for example absorbing materials or biaxial materials for whichmixing between the axes are present throughPockels coefficients or manymore. In this case we'll make an object with a differentpermittivity thanthe grid it is in.

Just like for the grid, theObject expects apermittivity to be a floats oran array of the following possible shapes

• (obj.Nx, obj.Ny, obj.Nz)
• or(obj.Nx, obj.Ny, obj.Nz, 1)
• or(obj.Nx, obj.Ny, obj.Nz, 3)

Note that the valuesobj.Nx,obj.Ny andobj.Nz are not given to theobject constructor. They are in stead derived from its placing in the grid:

grid[11:32,30:84,0]=fdtd.Object(permittivity=1.7**2,name="object")

Several things happen here. First of all, the object is given the space[11:32, 30:84, 0] in the grid. Because it is given this space, the object'sNx,Ny andNz are automatically set. Furthermore, by supplying a name tothe object, this name will become available in the grid:

print(grid.object)

    Object(name='object')        @ x=11:32, y=30:84, z=0:1

A second object can be added to the grid:

grid[13e-6:18e-6,5e-6:8e-6,0]=fdtd.Object(permittivity=1.5**2)

Here, a slice with floating point numbers was chosen. These floats will bereplaced by integerNx,Ny andNz during the registration of the object.Since the object did not receive a name, the object won't be available as anattribute of the grid. However, it is still available via thegrid.objectslist:

print(grid.objects)

[Object(name='object'), Object(name=None)]

This list stores all objects (i.e. of typefdtd.Object) in the order thatthey were added to the grid.

Similarly as to adding an object to the grid, anfdtd.LineSource can also beadded:

# signaturefdtd.LineSource(period:Number=15,# timesteps or secondsamplitude:float=1.0,phase_shift:float=0.0,name:str=None,)

And also just like anfdtd.Object, anfdtd.LineSource size is defined by itsplacement on the grid:

grid[7.5e-6:8.0e-6,11.8e-6:13.0e-6,0]=fdtd.LineSource(period=1550e-9/ (3e8),name="source")

However, it is important to note that in this case aLineSource is added tothe grid, i.e. the source spans the diagonal of the cube defined by the slices.Internally, these slices will be converted into lists to ensure this behavior:

print(grid.source)

    LineSource(period=14, amplitude=1.0, phase_shift=0.0, name='source')        @ x=[48, ... , 51], y=[76, ... , 83], z=[0, ... , 0]

Note that one could also have supplied lists to index the grid in the firstplace. This feature could be useful to create aLineSource of arbitraryshape.

# signaturefdtd.LineDetector(name=None)

Adding a detector to the grid works the same as adding a source

grid[12e-6, :,0]=fdtd.LineDetector(name="detector")

print(grid.detector)

    LineDetector(name='detector')        @ x=[77, ... , 77], y=[0, ... , 96], z=[0, ... , 0]

# signaturefdtd.PML(a:float=1e-8,# stability factorname:str=None)

Although, having an object, source and detector to simulate is in principleenough to perform an FDTD simulation, One also needs to define a grid boundaryto prevent the fields to be reflected. One of those boundaries that can beadded to the grid is aPerfectly MatchedLayer orPML. Theseare basically absorbing boundaries.

# x boundariesgrid[0:10, :, :]=fdtd.PML(name="pml_xlow")grid[-10:, :, :]=fdtd.PML(name="pml_xhigh")# y boundariesgrid[:,0:10, :]=fdtd.PML(name="pml_ylow")grid[:,-10:, :]=fdtd.PML(name="pml_yhigh")

A simple summary of the grid can be shown by printing out the grid:

print(grid)

Grid(shape=(161,97,1), grid_spacing=1.55e-07, courant_number=0.70)sources:    LineSource(period=14, amplitude=1.0, phase_shift=0.0, name='source')        @ x=[48, ... , 51], y=[76, ... , 83], z=[0, ... , 0]detectors:    LineDetector(name='detector')        @ x=[77, ... , 77], y=[0, ... , 96], z=[0, ... , 0]boundaries:    PML(name='pml_xlow')        @ x=0:10, y=:, z=:    PML(name='pml_xhigh')        @ x=-10:, y=:, z=:    PML(name='pml_ylow')        @ x=:, y=0:10, z=:    PML(name='pml_yhigh')        @ x=:, y=-10:, z=:objects:    Object(name='object')        @ x=11:32, y=30:84, z=0:1    Object(name=None)        @ x=84:116, y=32:52, z=0:1

Running a simulation is as simple as using thegrid.run method.

grid.run(total_time:Number,progress_bar:bool=True)

Just like for the lengths in the grid, thetotal_time of the simulationcan be specified as an integer (number oftime_steps) or as a float (inseconds).

grid.run(total_time=100)

Let's visualize the grid. This can be done with thegrid.visualize method:

# signaturegrid.visualize(grid,x=None,y=None,z=None,cmap="Blues",pbcolor="C3",pmlcolor=(0,0,0,0.1),objcolor=(1,0,0,0.1),srccolor="C0",detcolor="C2",show=True,)

This method will by default visualize all objects in the grid, as well as thefield intensity at the currenttime_step at a certainx,yORz-plane. Bysettingshow=False, one can disable the immediate visualization of thematplotlib image.

grid.visualize(z=0)

An as quick as possible explanation of the FDTD discretization of the Maxwellequations.

An electromagnetic FDTD solver solves the time-dependent Maxwell Equations

curl(H)=ε*ε0*dE/dtcurl(E)=-µ*µ0*dH/dt

These two equations are calledAmpere's Law andFaraday's Law respectively.

In these equations, ε and µ are the relative permittivity and permeabilitytensors respectively. ε0 and µ0 are the vacuum permittivity and permeabilityand their square root can be absorbed into E and H respectively, such thatE := √ε0*E andH := õ0*H.

Doing this, the Maxwell equations can be written as update equations:

E+=c*dt*inv(ε)*curl(H)H-=c*dt*inv(µ)*curl(E)

The electric and magnetic field can then be discretized on a grid withinterlaced Yee-coordinates, which in 3D looks like this:

According to the Yee discretization algorithm, there are inherently two typesof fields on the grid:E-type fields on integer grid locations andH-typefields on half-integer grid locations.

The beauty of these interlaced coordinates is that they enable a very naturalway of writing the curl of the electric and magnetic fields: the curl of anH-type field will be an E-type field and vice versa.

This way, the curl of E can be written as

curl(E)[m,n,p]= (dEz/dy-dEy/dz,dEx/dz-dEz/dx,dEy/dx-dEx/dy)[m,n,p]=( ((Ez[m,n+1,p]-Ez[m,n,p])/dy- (Ey[m,n,p+1]-Ey[m,n,p])/dz),                      ((Ex[m,n,p+1]-Ex[m,n,p])/dz- (Ez[m+1,n,p]-Ez[m,n,p])/dx),                      ((Ey[m+1,n,p]-Ey[m,n,p])/dx- (Ex[m,n+1,p]-Ex[m,n,p])/dy) )=(1/du)*( ((Ez[m,n+1,p]-Ez[m,n,p])- (Ey[m,n,p+1]-Ey[m,n,p])), [assumedx=dy=dz=du]                             ((Ex[m,n,p+1]-Ex[m,n,p])- (Ez[m+1,n,p]-Ez[m,n,p])),                             ((Ey[m+1,n,p]-Ey[m,n,p])- (Ex[m,n+1,p]-Ex[m,n,p])) )

this can be written efficiently with array slices (note that the factor(1/du) was left out):

defcurl_E(E):curl_E=np.zeros(E.shape)curl_E[:,:-1,:,0]+=E[:,1:,:,2]-E[:,:-1,:,2]curl_E[:,:,:-1,0]-=E[:,:,1:,1]-E[:,:,:-1,1]curl_E[:,:,:-1,1]+=E[:,:,1:,0]-E[:,:,:-1,0]curl_E[:-1,:,:,1]-=E[1:,:,:,2]-E[:-1,:,:,2]curl_E[:-1,:,:,2]+=E[1:,:,:,1]-E[:-1,:,:,1]curl_E[:,:-1,:,2]-=E[:,1:,:,0]-E[:,:-1,:,0]returncurl_E

The curl for H can be obtained in a similar way (note again that the factor(1/du) was left out):

defcurl_H(H):curl_H=np.zeros(H.shape)curl_H[:,1:,:,0]+=H[:,1:,:,2]-H[:,:-1,:,2]curl_H[:,:,1:,0]-=H[:,:,1:,1]-H[:,:,:-1,1]curl_H[:,:,1:,1]+=H[:,:,1:,0]-H[:,:,:-1,0]curl_H[1:,:,:,1]-=H[1:,:,:,2]-H[:-1,:,:,2]curl_H[1:,:,:,2]+=H[1:,:,:,1]-H[:-1,:,:,1]curl_H[:,1:,:,2]-=H[:,1:,:,0]-H[:,:-1,:,0]returncurl_H

The update equations can now be rewritten as

E+= (c*dt/du)*inv(ε)*curl_HH-= (c*dt/du)*inv(µ)*curl_E

The number(c*dt/du) is a dimensionless parameter called theCourant numbersc. For stability reasons, the Courant number should always be smaller than1/√D, withD the dimension of the simulation. This can be intuitively beunderstood as the condition that information should always travel slower thanthe speed of light through the grid. In the FDTD method described here,information can only travel to the neighboring grid cells (through applicationof the curl). It would therefore takeD time steps to travel over thediagonal of aD-dimensional cube (square in2D, cube in3D), the Courantcondition follows then automatically from the fact that the length of thisdiagonal is1/√D.

This yields the final update equations for the FDTD algorithm:

E+=sc*inv(ε)*curl_HH-=sc*inv(µ)*curl_E

This is also how it is implemented:

classGrid:# ... [initialization]defstep(self):self.update_E()self.update_H()defupdate_E(self):self.E+=self.courant_number*self.inverse_permittivity*curl_H(self.H)defupdate_H(self):self.H-=self.courant_number*self.inverse_permeability*curl_E(self.E)

Ampere's Law can be updated to incorporate a current density:

curl(H)=J+ε*ε0*dE/dt

Making again the usual substitutionssc := c*dt/du,E := √ε0*E andH := õ0*H, the update equations can be modified to include the current density:

E+=sc*inv(ε)*curl_H-dt*inv(ε)*J/√ε0

Making one final substitutionEs := -dt*inv(ε)*J/√ε0 allows us to write thisin a very clean way:

E+=sc*inv(ε)*curl_H+Es

Where we defined Es as theelectric field source term.

It is often useful to also define amagnetic field source termHs, which would bederived from themagnetic current density if it were to exist. In the same way,Faraday's update equation can be rewritten as

H-=sc*inv(µ)*curl_E+Hs

classSource:# ... [initialization]defupdate_E(self):# electric source function heredefupdate_H(self):# magnetic source function hereclassGrid:# ... [initialization]defupdate_E(self):# ... [electric field update equation]forsourceinself.sources:source.update_E()defupdate_H(self):# ... [magnetic field update equation]forsourceinself.sources:source.update_H()

When a material has aelectric conductivity σ, a conduction-current willensure that the medium is lossy. Ampere's law with a conduction current becomes

curl(H)=σ*E+ε*ε0*dE/dt

Making the usual substitutions, this becomes:

E(t+dt)-E(t)=sc*inv(ε)*curl_H(t+dt/2)-dt*inv(ε)*σ*E(t+dt/2)/ε0

This update equation depends on the electric field on a half-integer time step (amagnetic field time step). We need to substituteE(t+dt/2)=(E(t)+E(t+dt))/2 tointerpolate the electric field to the correct time step.

    (1+0.5*dt*inv(ε)*σ/√ε0)*E(t+dt)=sc*inv(ε)*curl_H(t+dt/2)+ (1-0.5*dt*inv(ε)*σ/ε0)*E(t)

Which, yield the new update equations:

f=0.5*inv(ε)*σ*sc*du/(ε0*c)E*=inv(1+f)* (1-f)E+=inv(1+f)*sc*inv(ε)*curl_H

Note that the more complicated the permittivity tensor ε is, the more timeconsuming this algorithm will be. It is therefore sometimes a nice hack totransfer the absorption to the magnetic domain by introducing a(nonphysical) magnetic conductivity, because the permeability tensor µ isusually just equal to one:

f=0.5*inv(μ)*σm*sc*du/(μ0*c)H*=inv(1+f)* (1-f)H+=inv(1+f)*sc*inv(µ)*curl_E

The electromagnetic energy density can be given by

e= (1/2)*ε*ε0*E**2+ (1/2)*µ*µ0*H**2

making the above substitutions, this becomes in simulation units:

e= (1/2)*ε*E**2+ (1/2)*µ*H**2

The Poynting vector is given by

P=E×H

Which in simulation units becomes

P=c*E×H

The energy introduced by a sourceEs can be derived from tracking the changein energy density

de=ε*Es·E+ (1/2)*ε*Es**2

This could also be derived from Poyntings energy conservation law:

de/dt=-grad(S)-J·E

where the first term just describes the redistribution of energy in a volumeand the second term describes the energy introduced by a current density.

Note: although it is unphysical, one could also have introduced a magneticsource. This source would have introduced the following energy:

de=ε*Hs·H+ (1/2)*µ*Hs**2

Since the µ-tensor is usually just equal to one, using a magnetic source termis often more efficient.

Similarly, one can also keep track of the absorbed energy due to an electricconductivity in the following way:

f=0.5*inv(ε)*σ*sc*du/(ε0*c)Enoabs=E+sc*inv(ε)*curl_HE*=inv(1+f)* (1-f)E+=inv(1+f)*sc*inv(ε)*curl_HdE=Enoabs-Ee_abs+=ε*E*dE+0.5*ε*dE**2

or if we want to keep track of the absorbed energy by magnetic a magneticconductivity:

f=0.5*inv(μ)*σm*sc*du/(μ0*c)Hnoabs=E+sc*inv(µ)*curl_EH*=inv(1+f)* (1-f)H+=inv(1+f)*sc*inv(µ)*curl_EdH=Hnoabs-He_abs+=µ*H*dH+0.5*µ*dH**2

The electric term and magnetic term in the energy density are usually of thesame size. Therefore, the same amount of energy will be absorbed by introducingamagnetic conductivity σm as by introducing aelectric conductivity σ if:

inv(µ)*σm/µ0=inv(ε)*σ/ε0

Assuming we want periodic boundary conditions along theX-direction, then wehave to make sure that the fields atXlow andXhigh are the same. This hasto be enforced after performing the update equations:

Note that the electric fieldE is dependent oncurl_H, which means that thefirst indices ofE will not be updated through the update equations. It'sthose indices that need to be set through the periodic boundary condition.Concretely:E[0] needs to be set to equalE[-1]. For the magnetic field,the inverse is true:H is dependent oncurl_E, which means that its lastindices will not be set. This has to be done by the boundary condition:H[-1]needs to be set equal toH[0]:

classPeriodicBoundaryX:# ... [initialization]defupdate_E(self):self.grid.E[0, :, :, :]=self.grid.E[-1, :, :, :]defupdate_H(self):self.grid.H[-1, :, :, :]=self.grid.H[0, :, :, :]classGrid:# ... [initialization]defupdate_E(self):# ... [electric field update equation]# ... [electric field source update equations]forboundaryinself.boundaries:boundary.update_E()defupdate_H(self):# ... [magnetic field update equation]# ... [magnetic field source update equations]forboundaryinself.boundaries:boundary.update_H()

a Perfectly Matched Layer (PML) is the state of the art forintroducing absorbing boundary conditions in an FDTD grid.A PML is an impedance-matched absorbing area in the grid. It turns out thatfor a impedance-matching condition to hold, the PML can only be absorbing ina single direction. This is what makes a PML in fact a nonphysical material.

Consider Ampere's law for theEz component, where we use the following substitutions:E := √ε0*E,H := õ0*H andσ := inv(ε)*σ/ε0 arealready introduced:

ε*dEz/dt+ε*σ*Ez=c*dHy/dx-c*dHx/dy

This becomes in the frequency domain:

iω*ε*Ez+ε*σ*Ez=c*dHy/dx-c*dHx/dy

We can split this equation in a x-propagating wave and a y-propagating wave:

iω*ε*Ezx+ε*σx*Ezx=iω*ε*(1+σx/iω)*Ezx=c*dHy/dxiω*ε*Ezy+ε*σy*Ezy=iω*ε*(1+σy/iω)*Ezy=-c*dHx/dy

We can define theS-operators as follows

Su=1+σu/iωwithuin {x,y,z}

In general, we prefer to add a stability factorau and a scaling factorku toSu:

Su=ku+σu/(iω+au)withuin {x,y,z}

Summing the two equations forEz back together after dividing by the respectiveS-operator gives

iω*ε*Ez= (c/Sx)*dHy/dx- (c/Sy)*dHx/dy

Converting this back to the time domain gives

ε*dEz/dt=c*sx[*]dHy/dx-c*sx[*]dHx/dy

wheresx denotes the inverse Fourier transform of(1/Sx) and[*] denotes a convolution.The expression forsu can be proven [after some derivation] to look as follows:

su= (1/ku)*δ(t)+Cu(t)withuin {x,y,z}

whereδ(t) denotes the Dirac delta function andC(t) an exponentiallydecaying function given by:

Cu(t)=-(σu/ku**2)*exp(-(au+σu/ku)*t)forallt>0anduin {x,y,z}

Plugging this in gives:

dEz/dt= (c/kx)*inv(ε)*dHy/dx- (c/ky)*inv(ε)*dHx/dy+c*inv(ε)*Cx[*]dHy/dx-c*inv(ε)*Cx[*]dHx/dy= (c/kx)*inv(ε)*dHy/dx- (c/ky)*inv(ε)*dHx/dy+c*inv(ε)*Фez/duwithdu=dx=dy=dz

This can be written as an update equation:

Ez+= (1/kx)*sc*inv(ε)*dHy- (1/ky)*sc*inv(ε)*dHx+sc*inv(ε)*Фez

Where we definedФeu as

Фeu=Ψeuv-Ψezwwithu,v,win {x,y,z}

andΨeuv as the convolution updating the componentEu by taking the derivative ofHw in thev direction:

Ψeuv=dv*Cv[*]dHw/dvwithu,v,win {x,y,z}

This can be rewritten [after some derivation] as an update equation in itself:

Ψeuv=bv*Ψeuv+cv*dv*(dHw/dv)=bv*Ψeuv+cv*dHwwithu,v,win {x,y,z}

Where the constantsbu andcu are derived to be:

bu=exp(-(au+σu/ku)*dt)withuin {x,y,z}cu=σu*(bu-1)/(σu*ku+au*ku**2)withuin {x,y,z}

The final PML algorithm for the electric field now becomes:

1. UpdateФe=[Фex, Фey, Фez] by using the update equation for theΨ-components.
2. Update the electric fields the normal way
3. AddФe to the electric fields.

or as python code:

classPML(Boundary):# ... [initialization]defupdate_phi_E(self):# update convolutionself.psi_Ex*=self.bEself.psi_Ey*=self.bEself.psi_Ez*=self.bEc=self.cEHx=self.grid.H[self.locx]Hy=self.grid.H[self.locy]Hz=self.grid.H[self.locz]self.psi_Ex[:,1:, :,1]+= (Hz[:,1:, :]-Hz[:, :-1, :])*c[:,1:, :,1]self.psi_Ex[:, :,1:,2]+= (Hy[:, :,1:]-Hy[:, :, :-1])*c[:, :,1:,2]self.psi_Ey[:, :,1:,2]+= (Hx[:, :,1:]-Hx[:, :, :-1])*c[:, :,1:,2]self.psi_Ey[1:, :, :,0]+= (Hz[1:, :, :]-Hz[:-1, :, :])*c[1:, :, :,0]self.psi_Ez[1:, :, :,0]+= (Hy[1:, :, :]-Hy[:-1, :, :])*c[1:, :, :,0]self.psi_Ez[:,1:, :,1]+= (Hx[:,1:, :]-Hx[:, :-1, :])*c[:,1:, :,1]self.phi_E[...,0]=self.psi_Ex[...,1]-self.psi_Ex[...,2]self.phi_E[...,1]=self.psi_Ey[...,2]-self.psi_Ey[...,0]self.phi_E[...,2]=self.psi_Ez[...,0]-self.psi_Ez[...,1]defupdate_E(self):# update PML located at self.locself.grid.E[self.loc]+= (self.grid.courant_number*self.grid.inverse_permittivity[self.loc]*self.phi_E        )classGrid:# ... [initialization]defupdate_E(self):forboundaryinself.boundaries:boundary.update_phi_E()# ... [electric field update equation]# ... [electric field source update equations]forboundaryinself.boundaries:boundary.update_E()

The same has to be applied for the magnetic field.

These update equations for the PML were based onSchneider, Chap. 11.

As a bare FDTD library, this is dimensionally agnostic for any unit system you may choose.No conversion factors are applied within the library API; this is left to the user.(The code used to calculate the Courant limit may be a sticking point depending on the time scale involved).

However, as noted above (H := õ0*H), it is generally good numerical practice to scale all values toget the maximum precision from floating-point types.

In particular, a scaling scheme detailed in"Novel architectures for brain-inspired photonic computers", Chapters 4.1.2 and 4.1.6, is highly recommended.

A set of conversion functions to and from reduced units are available for users in conversions.py.

You can run a linter in the root usingpylint fdtd.

© Floris laporte -MIT License