Fandom

Stanford University Wiki

FDTD Plus

229pages on
this wiki
Add New Page
Talk0 Share
Header master


This is a Wiki for a finite-difference time-domain (FDTD) simulation software package developed by Sunil Sandhu of Fan Group at Stanford University.

The purpose of this Wiki is to

  • document the software package.
  • allow group members to post feedback on the software.

Overview Edit

The name of the FDTD software package is FDTD Plus. A parallel version of FDTD Plus executable is available to Fan group members at several servers,

   fdtd_plus_mpi <directory location> <input file name> <x-proc> <y-proc> <z-proc>
  • <directory location> is the directory where the input file is located and is the location where all FDTD Plus outputs are performed.
  • <x-proc>, <y-proc> and <z-proc> are the number of processors in the 3 directions.
    • In the next release, the allocation of processors in the 3 directions will be done internally.
  • During the run, errors, warnings, status of the run etcetera are output to a file called output.dat.
  • Email centaur(at)stanford.edu if you can't access any of the above sites or for any other questions.

Convenience script (not implemented in Hera) Edit

There is a convenience script available for submitting FDTD Plus jobs into any queue in the above platforms using a minimum set of arguments. The usage of the script at the command line is

             qsub_fdtd_plus <options>
  • To obtain information on the available options of qsub_fdtd_plus type
             qsub_fdtd_plus --help

Example usageEdit

  • An example use is
             qsub_fdtd_plus -o=$SCRATCH/run01 -f=$HOME/input/input01.txt -t=00:30:00 -x=2

which will submit a FDTD Plus run with the following settings

  • user input file path $HOME/input/input01.txt (this file will be copied to the working directory indicated in the -o option),
  • working directory (i.e. where all outputs are done) at $SCRATCH/run01 (note: run01 directory is created if it does not exist),
  • for a wall time of 30min,
  • using 2 processors in the x direction (default value for -x,-y,-z directions is 1 processor).

Access to scriptEdit

For Shanhuiserve, you can use it right away at the command line. As for other platforms, the following steps have to be performed before qsub_fdtd_plus can be used at the command line.

  • For SDSC, add the following at the end of the .soft file (.soft is located in your $HOME directory),
         PATH += /users/ux455370/share

User input fileEdit

The FDTD Plus simulation is set up based on instructions written in a user input file.

Input file syntaxEdit

The input file consists of instructions enclosed in parentheses (round brackets).

            syntax of an instruction:
                          (instruction)
            error! no closing parenthesis:
                          (instruction_1)
                          (instructio_2     <-instruction_2 missing closing parenthesis
                          (instruction_3)

Certain instructions introduced in later sections have nested instructions.

           nested instruction example:
                           (instruction_1 (instruction_2) )


The main rule in the user input file is never have two parentheses together.

     this will result in an error: 
                           (instruction_1 (instruction_2))  <-closed parentheses together
            
     this fixes the above error: 
                           (instruction_1 (instruction_2) )


An instruction can be disabled by adding a slash "\" to the open parenthesis.

           disable instruction: 
                            (\instruction)
           only instruction 2 is disabled: 
                            (instruction_1 (\instruction_2) )

Having "(\" together with another parenthesis (e.g. "(\(" ) will result in an error too.

A block of instructions can be disabled by enclosing them in "(\" and ")".

          instructions 3-5 are disabled:
                            (instruction_1)   
                            (instruction_2) 
                       (\            
                            (instruction_3)   
                            (instruction_4)   
                            (instruction_5) 
                        )
                            (instruction_6)   
                            (instruction_7)

Comments can be added anywhere outside instructions. If parentheses (i.e. round brackets) need to be used for comments, the open parenthesis must have a slash attached to it, "(\".

  • Having a single "(\" without an enclosing ")" will result in an error too.
         this will result in an error:
                       This comment is fine.
                       #This is another example of a comment that is fine too
                       (\This is a comment that is fine.. I guess you get the idea!)
                       (\This comment has an error..  <-no closing parenthesis 
                                   (instruction)

Available instructionsEdit

There are five groups of instructions, namely Options, Materials, Sources, Objects, and Outputs .

  • All instructions except for Objects can be placed in any order. Source Objects and Material Objects are processed separately according to the order they are written.
    • The order of Source Objects have no relation with the order of Material Objects (see comment in Sources).


The following sections describe the different possible user input instructions available.

Example fileEdit

Click here for an example input file.

Options Edit

This user input sets up the dimensions, boundaries and cell type for the simulation.

        syntax:
                (Option OptionType <parameters> )
  • In the definition of the different options below, any option that is not marked "(optional)" is a required option.

CellParamEdit

Defines the type of cell used in the FDTD computation including the field type.

         (Option CellParam (CellType      value) 
                           (FieldType     value)
                           (PrecisionType value) 
         )

CellTypeEdit

The type of cells used to represent the materials in the computation region,
value\in\{linear\},

  • linear is a field independent material.
    • A linear cell only communicates E- and H-fields to neighboring cells.

FieldTypeEdit

The mathematical field type of the cell physical fields that are communicated to other cells,
value\in\{complex, real\}.

  • Most simulations can be performed using FieldType = real. If a material type or output type requires FieldType = complex, this will be highlighted in its description below.

PrecisionType (optional)Edit

The floating point precision of fields, coefficients etc in the cell,
value\in\{double, float\}.

  • Default value = double.
  • To include the PrecisionType = float during compilation of FDTD Plus, use the FDTD_PLUS_PRECISION_TYPE_FLOAT flag.
    • The FDTD Plus executable available in the servers above do not have the PrecisionType = float option.

DimensionParamEdit

Defines the spatial and time dimensions of the compute region.

     (Option DimensionParam (LengthScale   value)         
                            (Center        x_value y_value z_value) 
                            (Size          x_value y_value z_value) 
                            (Resolution    x_value y_value z_value)
                            (TimeParam     (Start value) (End value) ) 
                            (CourantNumber value)  
     )

LengthScale (optional)Edit

The length scale of the compute region,
value\in\mathbb{R}_{>0}.

  • Unit = meters.
  • LengthScale does not need to be defined for lossless, non-dispersive, linear materials.
    • Since the Maxwell equations are scale invariant in the case of lossless, non-dispersive, linear materials, a length scale (i.e. LengthScale) can be used to represent a unit of distance in these cases.
    • In the description of material types below, materials requiring the definition of LengthScale will be highlighted.

CenterEdit

The center coordinates of the compute region,
x\_value,y\_value,z\_value\in\mathbb{R}

  • Unit = LengthScale.

SizeEdit

The size of the compute region,
x\_value,y\_value,z\_value\in\mathbb{R}_{>0}

  • Unit = LengthScale.

ResolutionEdit

The spatial increment in the compute region,
x\_value,y\_value,z\_value\in\mathbb{R}_{>0}

  • Unit = LengthScale.

TimeParamEdit

The time parameters for the FDTD run.

Start (optional)Edit

The start time of computation
value\in\mathbb{R}

  • Unit = LengthScale / c, where c = speed of light in vacuum.
  • Default value = 0.
EndEdit

The end time of computation
value\in\mathbb{R}_{\geq{Start}}

  • Unit = LengthScale / c.

CourantNumber (optional)Edit

The courant number,
value\in\mathbb{R}_{>0}.

  • Default value = 0.5\cdot\sqrt{3}\approx{0.866}.
  • Using CourantNumber, the FDTD time increment is calculated as
                  \Delta{t}=\frac{S}{\sqrt{\frac{1}{\Delta{x}^{2}}+\frac{1}{\Delta{y}^{2}}+\frac{1}{\Delta{z}^{2}}}}
  • where,
  • In order to ensure numerical stability, CourantNumber should be slightly less than n_{min},
    • n_{min} = minimum refractive index in the compute region.

BoundaryParamEdit

Defines the boundaries of the compute region.

     (Option BoundaryParam (Type        x_value y_value z_value) 
                           (PmlParam   (Size (X value1 value2) (Y value1 value2) (Z value1 value2) )
                                       (SigmaFactor (X value1 value2) (Y value1 value2) (Z value1 value2) ) 
                                       (KappaFactor (X value1 value2) (Y value1 value2) (Z value1 value2) )    
                                       (Type value)
                                       (DecayConstant value)
                           )    
                           (BlochK     (X value) (Y value) (Z value) )
                           (Symmetry   (X value) (Y value) (Z value) )
     )

TypeEdit

The type of boundaries in the compute region,
x\_value,y\_value,z\_value\in\{pml, bloch, none\},

  • pml is a Perfectly Matched Layer (PML) boundary,
  • bloch is a Bloch boundary,
  • none results in the boundary that includes the origin of the compute region to be of type PEC (E-fields parallel to boundary set to zero), and the opposite boundary of the compute region to be of type PMC (H-fields parallel to boundary set to zero).

PmlParam (optional)Edit

The parameters for PML boundaries.

  • PmlParam is only required to be defined if any BoundaryParam::Type == pml.
SizeEdit

The size of PML boundaries.

  • Size is included in DimensionParam::Size.
  • Size is only required to be defined for directions with BoundaryParam::Type == pml.
value1Edit

The PMLboundary size for boundaries that include the origin of the compute region,
value1\in\mathbb{R}_{\geq{0}}

  • Unit = DimensionParam::LengthScale.
value2 (optional)Edit

The PMLboundary size for boundaries that do not include the origin of the compute region,
value2\in\mathbb{R}_{\geq{0}}

  • Unit = DimensionParam::LengthScale.
  • Default value = value1.
SigmaFactor (optional)Edit

This controls the imaginary part of the sparameter of UPML (see Taflove). Controls the normal incidence reflection R(0) at the PMLboundaries. Generally, R(0)=Exp\left[-\gamma\right] where \gamma=SigmaFactor\times{1.6}\times\frac{PmlParams::Size}{DimensionParams::Resolution}

  • This parameter adds flexibility to designing your PML boundaries.
  • By default SigmaFactor = 0.8. This results in \gamma=12.8 for 10 PML cells. But for 20 PML cells, \gamma=25.6 (which can cause numerical issues). Be careful or you may get undesirable reflections or numerical issues. Read Taflove for more information.
  • You should read Taflove's Computational Electrodynamics (within chapter on PML, section on Theoretical Performance of the PML) to get a better idea of this parameter.
value1Edit

The SigmaFactor for PMLboundaries that include the origin of the compute region,
value1\in\mathbb{R}_{>0}

  • Unit = Dimensionless.
  • Default value = 0.8.
value2 (optional)Edit

The SigmaFactor for PMLboundaries that do not include the origin of the compute region,
value2\in\mathbb{R}_{>0}

  • Unit = Dimensionless.
  • Default value = value1.
KappaFactor (optional)Edit

This controls the real part of the s parameter of UPML (see Taflove).

  • KappaFactor is \kappa_{max} in Taflove.
  • This parameter adds flexibility to designing your PML boundaries.
  • By default KappaFactor = 1.
  • Only use this parameter if you REALLY know what you are doing.
value1Edit

The KappaFactor for PMLboundaries that include the origin of the compute region,
value1\in\mathbb{R}_{>0}

  • Unit = Dimensionless.
  • Default value = 1.
value2 (optional)Edit

The KappaFactor for PMLboundaries that do not include the origin of the compute region,
value2\in\mathbb{R}_{>0}

  • Unit = Dimensionless.
  • Default value = value1.
TypeEdit

The type of PMLboundary,
value\in\{upml\},

DecayConstant (optional)Edit

The PMLdecay constant,
value\in\mathbb{R}_{>0}

  • Default value = e^{-16}.

BlochK (optional)Edit

The Bloch k-vector for Bloch boundaries,
value\in\mathbb{R}

  • Unit = 2\pi / DimensionParam::Size.
  • BlochK is only required to be defined for directions with BoundaryParam::Type == bloch.
  • When BoundaryParam::Type == bloch, only Bloch k-vector value = 0 allows CellParam::FieldType = real.

Symmetry (optional)Edit

For a compute region that has mirror symmetry in the x-, y- or z-direction, this option enforces the symmetry at the center plane normal to that direction (i.e. each symmetry results in a reduction by half of the corresponding DimensionParam::Size used for FDTD computation),
value\in\{pec, pmc\},

  • pec results in a Perfect Electric Conductor (PEC) center plane where the E-fields parallel to the plane are set to zero.
  • pmc results in a Perfect Magnetic Conductor (PMC) center plane where the H-fields parallel to the plane are set to zero.

OutputStructure (optional)Edit

Prints the structure within the compute region into a hdf5 file,

     (Option OutputStructure (FileName value) )
  • This option is disabled by default.

FileName (optional)Edit

The name of the output file,
value\in{s}tring.

  • Default value = "structure".
  • The suffix ".h5" is appended to FileName internally.

DefaultMaterial (optional)Edit

Sets the default material in the compute region,

     (Option DefaultMaterial value)


value\in{s}tring.

  • value must be a defined MaterialName (see Material below).
  • Default value = "vacuum". The "vacuum" material type is generated internally.

CheckInput (optional)Edit

Option to check the user input file. This option runs the simulation up to the point before cell construction.

     (Option CheckInput)
  • This option is disabled by default.
  • If OutputStructure is enabled, CheckInput will run OutputStructure.
  • Read comment in Sources for a bug in CheckInput.

Materials Edit

This user input defines a material type.

   syntax:
          (Material MaterialType (MaterialName value) 
                                 (MaterialToken value)
                                  <other parameters> 
           )

Common material parametersEdit

MaterialNameEdit

The name of the material,
value\in{s}tring.

  • Duplicate MaterialName values are not allowed.
  • "vacuum" is a reserved MaterialName.

MaterialTokenEdit

The index of the material (used by Option::OutputStructure to print material),
value\in\mathbb{I}_{>0}.

  • Duplicate MaterialToken values are not allowed.
  • 0 is a reserved MaterialToken for "vacuum".

Linear materialsEdit

These are the material types defined for CellParam::CellType = linear.

LosslessEdit

Defines a lossless MaterialType.

         (Material Lossless (MaterialName    value) 
                            (MaterialToken   value)
                            (RelPermittivity value)
                            (RelPermeability value) 
         )

RelPermittivityEdit

The relative permittivity of the material,
value\in\mathbb{R}.

RelPermeabilityEdit

The relative permeability of the material,
value\in\mathbb{R}.

LossyEdit

Defines a lossy MaterialType.

         (Material Lossy    (MaterialName    value) 
                            (MaterialToken   value)
                            (RelPermittivity value)
                            (RelPermeability value) 
                            (Conductivity    value) 
         )
  • The parameters for Lossy are the same as Lossless except for the additional conductivity parameter.
  • Option::DimensionParam::LengthScale must be defined for Lossy.

ConductivityEdit

The conductivity of the material,
value\in\mathbb{R}.

  • Unit = siemens/meter.

PRPEdit

Defines an arbitrary dispersive MaterialType based on general complex-conjugate pole-residue pairs (PRP).

         (Material Prp      (MaterialName     value) 
                            (MaterialToken    value)
                            (RelPermeability  value)
                            (EpsilonInfinity  value) 
                            (PoleResiduePair  (Pair1 residue_value pole_value )  
                                              (Pair2 residue_value pole_value ) 
                                                              |
                                                              |
                                                              |
                                              (PairN residue_value pole_value )
                             )  
         )
  • Option::DimensionParam::LengthScale must be defined for PRP.
  • PRP can only be used with Option::CellParam::FieldType = real.

EpsilonInfinityEdit

The relative permittivity at infinite frequency,
value\in\mathbb{R}.

PoleResiduePairEdit

The N pole-residue pairs,
residue\_value,\,pole\_value\, \in \,\mathbb{C}.

  • Unit = radians/second.
  • If your Xth (residue, pole) pair is (a + j b, c + j d), the format for entry into the PoleResiduePair parameter is (PairX (a, b) (c, d) ), where residue_value = (a, b) and pole_value = (c, d) in the above syntax for PRP definition.

Example PRP parametersEdit

  • Click here for gold (Au) parameters.
  • Click here for silver (Ag) parameters.

Sources Edit

This user input defines a source type.

   syntax:
          (Source SourceType (SourceName value) 
                             (Amplitude  value)
                             (Duration   value) 
                             (Center     value) 
                             (Width      value)
                             (Direction  x_value y_value z_value) 
                              <other parameters> 
          )
  • In FDTD Plus, sources and materials are populated using separate arrays during initialization.
    • Sources and materials are treated separately, up to the point of cell construction in the compute region, where a source from the source array is attached to a material in the material array.
      • The order in which source objects are input in the user input file have no relation with the order of material objects.
  • Currently, only the Material::PRP material cannot be attached to a source.
    • If this is attempted, Option::CheckInput will not catch the error. The program will abort during cell construction. This bug will be fixed in the next release.

Common source parametersEdit

SourceNameEdit

The name of the source,
value\in{s}tring.

  • Duplicate SourceName values are not allowed.

AmplitudeEdit

The amplitude of the excitation,
value\in\mathbb{R}.

  • Unit = amperes (if Option::DimensionParam::LengthScale is defined).
    • if Option::DimensionParam::LengthScale is not defined, value must be normalized to Option::DimensionParam::LengthScale by user (i.e. if source amplitude is I_{0} (unit = amperes), it should be entered as value = I_{0} / Option::DimensionParam::LengthScale ).
    • In SPP, the unit of Amplitude was not amperes. If your Amplitude value of the source in SPP was 1, enter it as 1 for FDTD Plus when Option::DimensionParam::LengthScale is not defined; if Option::DimensionParam::LengthScale is defined, enter it as the value of Option::DimensionParam::LengthScale.

DurationEdit

The duration of the excitation,
value\in\mathbb{R}_{\geq{0}}.

  • Unit = Option::DimensionParam::LengthScale / c.

CenterEdit

The center of the excitation,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale / c.

WidthEdit

The width of the excitation,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale / c.
  • A rule of thumb is to have Width \leq Center / 4. This is to allow a smooth turn on of the source.

DirectionEdit

The direction (polarization) of the source,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
  • This direction vector is normalized internally during initialization.

GaussTimeEdit

Defines a Gaussian source.

         (Source GaussTime (SourceName value) 
                           (Duration   value)
                           (Center     value)
                           (Width      value) 
                           (Frequency  value)
                           (Phase      value)
                           (Amplitude  value)
                           (Direction  x_value y_value z_value)   
        )
  • The Gaussian excitation is defined as,
             \frac{A_{0}}{\Delta_{0}\Delta_{1}}\cdot{E}xp\left[-\frac{(t-t_{0})^{2}}{2\cdot\Delta_{t}^{2}}\right]\cdot\sin\left[2\pi{f}(t-t_{0})-\phi\right]\cdot\left(u\left[t-(t_{0}-\frac{T}{2})\right]-u\left[t-(t_{0}+\frac{T}{2})\right]\right)
  • t = time (Unit = Option::DimensionParam::LengthScale / c).
  • A_{0} = Amplitude.
  • \Delta_{0},\Delta_{1} = Option::DimensionParam::Resolution.
    • \Delta_{0}\cdot\Delta_{1} is the area of the unit cell plane perpendicular to the direction of the current source.
  • \Delta_{t} = Width.
  • t_{0} = Center.
  • f = Frequency.
  • \phi = Phase.
  • T = Duration.
  • u\left[t\right] = unit step function.

FrequencyEdit

The center frequency of the excitation,
value\in\mathbb{R}.

  • Unit = 2\pi \,c / Option::DimensionParam::LengthScale.

PhaseEdit

The phase of the excitation,
value\in\mathbb{R}.

  • Unit = radians.

RectangleEdit

Defines a sinusoidal source with a Gaussian turn on/off.

         (Source Rectangle (SourceName value) 
                           (Duration   value)
                           (Center     value)
                           (Width      value) 
                           (Frequency  value)
                           (Phase      value)
                           (Amplitude  value)
                           (Direction  x_value y_value z_value)  
                           (Start        value) 
        )
  • The Rectangular excitation pulse is defined as,
             \frac{A_{0}}{\Delta_{0}\Delta_{1}}\left[\left(1-u(t-t_{start})\right)\cdot{E}xp\left[-\frac{(t-t_{start})^{2}}{2\cdot\Delta_{t}^{2}}\right]+\left(u(t-t_{start})-u(t-t_{end})\right)+u(t-t_{end})\cdot{E}xp\left[-\frac{(t-t_{end})^{2}}{2\cdot\Delta_{t}^{2}}\right]\right]\cdot\sin\left[2\pi{f}(t-t_{0})-\phi\right]\cdot\left(u\left[t-(t_{0}-\frac{T}{2})\right]-u\left[t-(t_{0}+\frac{T}{2})\right]\right)
  • t = time (Unit = Option::DimensionParam::LengthScale / c).
  • A_{0} = Amplitude.
  • \Delta_{0},\Delta_{1} = Option::DimensionParam::Resolution.
    • \Delta_{0}\cdot\Delta_{1} is the area of the unit cell plane perpendicular to the direction of the current source.
  • \Delta_{t} = Width of the Gaussian turn on/off.
  • t_{0} = Center.
  • f = Frequency.
  • \phi = Phase.
  • T = Duration.
  • u\left[t\right] = unit step function.
  • t_{start} = Start of the constant envelope or peak/end of the Gaussian turn on.
  • t_{end} = Center + Center - t_{start} = End of the constant envelope or peak/start of the Gaussian turn off.


WidthEdit

The width of the Gaussian turn on/off pulse,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale / c.
  • A rule of thumb is to have Width \leq (Start - (Center - Duration / 2)) / 4. This is to allow a smooth turn on/off of the source.

FrequencyEdit

The center frequency of the excitation,
value\in\mathbb{R}.

  • Unit = 2\pi \,c / Option::DimensionParam::LengthScale.

PhaseEdit

The phase of the excitation,
value\in\mathbb{R}.

  • Unit = radians.

StartEdit

Determines start and end of the constant envelope (or the peaks of the Gaussian turn on/off),
value\in\mathbb{R}.
The end of the constant envelope is calculated as Center + Center - Start.

  • Unit = Option::DimensionParam::LengthScale / c.
  • Default value = Center - Duration / 2 + 4 * Width.

ObjectsEdit

This user input defines the objects that can be built in the compute region.

   syntax:
          (Object ObjectType (Center x_value y_value z_value) 
                             (MaterialName/SourceName value)  
                                   <other parameters> 
                             (Array (Axis0    x_value y_value z_value)
                                    (Axis1    x_value y_value z_value)
                                    (Axis2    x_value y_value z_value)
                                    (Period   value0  value1  value2) 
                                    (NoPeriod value0  value1  value2)
                             )
           )

Common object parametersEdit

CenterEdit

The center coordinate of the object,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.

MaterialName/SourceNameEdit

The MaterialName (SourceName) of a defined Material (Source),
value\in{s}tring.

  • Defining both MaterialName and SourceName for the same object will throw an initialization error. If you want to define a source and a material at the same point, you have to use separate objects.

Array (optional)Edit

Creates an array of the object.

  • In the description of different object types, this optional parameter will be indicated as <optional array param>.
AxisN (N = 0, 1, 2)Edit

The Nth axis along which the object is cloned,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
  • Currently, if axes other than the scalar multiple of the x-, y-, z-unit axes are chosen, the result of object array population is undefined. This will be fixed in the next release.
PeriodEdit

The spatial increment of the object center along AxisN,
value0,value1,value2\in\mathbb{R}.

  • valueN is the spatial increment along AxisN.
  • Unit = Option::DimensionParam::LengthScale.
NoPeriodEdit

The number of periods along AxisN,
value0,value1,value2\in\mathbb{I}_{>0}.

  • valueN is the number of periods along AxisN.

PointEdit

Builds a point (i.e. populates one cell).

          (Object Point (Center x_value y_value z_value)
                        (MaterialName/SourceName value) 
                            <optional array param>   
          )

PlaneEdit

Builds a plane,

          (Object Plane (Center   x_value y_value z_value)
                        (MaterialName/SourceName value)
                        (Axis0    x_value y_value z_value)
                        (Axis1    x_value y_value z_value)
                        (Size     value0  value1)
                            <optional array param>
          )
AxisN (N = 0, 1)Edit

The two axes that are parallel to any two non-parallel edges of the plane,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
  • Currently, if axes other than the scalar multiple of the x-, y-, z-unit axes are chosen, the result of the plane population is undefined. This will be fixed in the next release.
SizeEdit

The size along AxisN,
value0,value1,value2\in\mathbb{R}.

  • valueN is the size along AxisN.
  • Unit = Option::DimensionParam::LengthScale.

BlockEdit

Builds a parallelepiped,

          (Object Block (Center   x_value y_value z_value)
                        (MaterialName/SourceName value)
                        (Axis0    x_value y_value z_value)
                        (Axis1    x_value y_value z_value)
                        (Axis2    x_value y_value z_value)
                        (Size     value0  value1  value2)
                            <optional array param>
          )
AxisN (N = 0, 1, 2)Edit

The three axes that are parallel to any three non-parallel edges of the parallelepiped,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
  • Currently, if axes other than the scalar multiple of the x-, y-, z-unit axes are chosen, the result of the block population is undefined. This will be fixed in the next release.
SizeEdit

The size along AxisN,
value0,value1,value2\in\mathbb{R}.

  • valueN is the size along AxisN.
  • Unit = Option::DimensionParam::LengthScale.

CylinderEdit

Builds a cylinder,

          (Object Cylinder (Center x_value y_value z_value)
                           (MaterialName/SourceName value)
                           (Axis   x_value y_value z_value)
                           (Radius value)
                           (Height value)
                            <optional array param>
          )
AxisEdit

The axis of the cylinder,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
HeightEdit

The size along Axis,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
RadiusEdit

The radius of the cylinder,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.

SphereEdit

Builds a sphere,

          (Object Sphere (Center x_value y_value z_value)
                           (MaterialName/SourceName value)
                           (Radius value)
                            <optional array param>
          )
RadiusEdit

The radius of the sphere,
value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.

OutputsEdit

This user input defines the outputs of the FDTD simulation.

   syntax:
          (Output OutputType (StartTime value) 
                             (EndTime   value)
                             (TimeStep  value)
                             (FileName  value)   
                             <other parameters> 
           )
  • The leap frog algorithmn of FDTD results in field components being stored at different time points and different spatial points.
  • Before any output computation below, each field component is spatially averaged (with the corresponding field components in its neighboring cells) to the same spatial point wihin the cell; in addition, time averaging is performed to obtain all field components at the same time point.

Common output parametersEdit

StartTimeEdit

The start time of the output,
value\in\mathbb{R}

  • Unit = Option::DimensionParam::LengthScale / c.

EndTimeEdit

The end time of the output,
value\in\mathbb{R}_{>StartTime}

  • Unit = Option::DimensionParam::LengthScale / c.

TimeStepEdit

The time step of the output,
value\in\mathbb{R}

  • Unit = Option::DimensionParam::LengthScale / c.

FileNameEdit

The name of the file,
value\in{s}tring.

FieldTimeEdit

Output selected fields at a particular point,

          (Output FieldTime (Position    x_value y_value z_value)
                            (OutputField value)  
                            (StartTime   value) 
                            (EndTime     value)
                            (TimeStep    value)
                            (FileName    value)   
          )
  • Click here for help on post-processing of output data.

PositionEdit

The position coordinate at which output is performed,
x\_value,y\_value,z\_value\in\mathbb{R}.

OutputFieldEdit

These are the fields requested. Possible values are dependent on value of CellParam::CellType.

  • See below for material dependent output definitions.

FieldSpaceEdit

Output selected fields in the whole compute region,

          (Output FieldSpace (OutputField value)  
                             (StartTime   value) 
                             (EndTime     value)
                             (TimeStep    value)
                             (FileName    value)   
          )
  • The output is performed into a hdf5 file.
    • The suffix ".h5" is appended to FileName internally.
  • If Option::CellParam::FieldType=complex, both real and imaginary parts of field are output.
  • todo (in future release):
    • To allow the output of the field in a plane instead of the whole region.

FluxTimePlaneEdit

At each requested TimeStep, output the sum of the time average Poynting vector over the plane,

          (Output FluxTimePlane (Center          x_value y_value z_value)
                                (NormalDirection value)
                                (Size            value0 value1)  
                                (StartTime       value) 
                                (EndTime         value)
                                (TimeStep        value)
                                (FileName        value)   
          )
  • Click here for help on post-processing of output data.
  • Output as Power / \left(Option::DimensionParam::LengthScale\right)^{2} where Power unit = Watts.
  • The sum of the time average Poynting vector over the plane is defined as,
          \frac{1}{2}\Delta_{1}\Delta_{2}\sum_{plane}\hat{n}\cdot{R}e[E(t,r)^{*}\times{H(t,r)}]
  • \Delta_{1},\Delta_{2} are the Option::DimensionParam::Resolution along the two perpendicular axes parallel to the plane (see Size below).
  • \hat{n} is the unit normal axis of the plane (i.e. NormalDirection).
  • FluxTimePlane requires CellParam::FieldType = complex.
    • since the calculation of the time average Poynting vector requires E(t,r) and H(t,r) to be complex fields.

CenterEdit

The center coordinate of the plane,
x\_value,y\_value,z\_value\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.

NormalDirectionEdit

The normal direction of the plane,
value\in\{\,positive_x, positive_y,positive_z,negative_x,negative_y,negative_z\,\}.

SizeEdit

The size along the two perpendicular axis parallel to the plane,
value0,value1\in\mathbb{R}.

  • Unit = Option::DimensionParam::LengthScale.
  • For NormalDirection = positive_x, negative_x, value0 corresponds to size along y-axis and value1 corresponds to size along z-axis.
  • For NormalDirection = positive_y, negative_y, value0 corresponds to size along x-axis and value1 corresponds to size along z-axis.
  • For NormalDirection = positive_z, negative_z, value0 corresponds to size along x-axis and value1 corresponds to size along y-axis.


  • To do
    • Implement time average Poynting vector through a box.

FluxSpectraPlaneEdit

Output the transmitted power spectra through a plane,

          (Output FluxSpectraPlane (Center          x_value y_value z_value)
                                   (NormalDirection value)
                                   (Size            value0 value1)  
                                   (StartTime       value) 
                                   (EndTime         value)
                                   (TimeStep        value)   
                                   (StartFrequency  value) 
                                   (EndFrequency    value)
                                   (FrequencyStep   (Value value) (Force ) )
                                   (FileName        value)   
          )
  • Click here for help on post-processing of output data.
  • Output as Power / \left(Option::DimensionParam::LengthScale\right)^{2} where Power unit = Watts.
  • The first step in calculating the transmitted power spectra involves taking the Discrete Fourier Transform (DFT) of each field parallel to the plane from StartTime to EndTime at time increment TimeStep. The DFT for field F(t,r) at spatial point r and frequency \omega is computed as
          F(\omega,r)=\sum_{t'=StartTime}^{EndTime}Re\left[F(t',r)\right]\cdot{e}^{j \omega{t'}}
  • At EndTime, the transmitted power at each frequency \omega through the plane is obtained by summing the complex Poynting vector across the plane,
          T(\omega)=\frac{1}{2}\Delta_{0}\Delta_{1}\sum_{plane}\hat{n}\cdot\left[E(\omega,r)^{*}\times{H}(\omega,r)\right]
  • \Delta_{0},\Delta_{1} are the Option::DimensionParam::Resolution along the two perpendicular axes parallel to the plane (see Size in FluxTimePlane).
  • \hat{n} is the unit normal axis of the plane (i.e. NormalDirection).
  • FluxSpectraPlane can be performed with CellParam::FieldType = real.
    • since E(t,r) and H(t,r) in the DFT summation above can be real fields.


  • The plane parameters Center, NormalDirection and Size are the same as FluxTimePlane.

EndTime, StartTimeEdit

  • Additional restrictions (besides those specified above)
    • If EndTime > Option::DimensionParam::TimeParam::End or EndTime <= StartTime, error is thrown during initialization.
    • If StartTime < Option::DimensionParam::TimeParam::Start, error is thrown during initialization.

TimeStep (optional)Edit

  • See above for definition.
  • TimeStep is optional for FluxSpectraPlane.
  • Default value = FDTD time increment.

StartFrequencyEdit

The start frequency of the flux spectra calculation,
value\in\mathbb{R}_{>0}

  • Unit = 2\pi c / Option::DimensionParam::LengthScale.

EndFrequencyEdit

The end frequency of the flux spectra calculation,
value\in\mathbb{R}_{\leq{f}_{max}}

  • f_{max}\, =\,1 / \left(2\cdot{T}imeStep\right).
  • Unit = 2\pi\,c / Option::DimensionParam::LengthScale.
  • value is automatically adjusted to f_{max} if value\, > \,f_{max} (in order to avoid aliasing).

FrequencyStepEdit

ValueEdit

The frequency step of the flux spectra calculation,
value\in\mathbb{R}_{>0}

  • Unit = 2\pi\, c / Option::DimensionParam::LengthScale.
Force (optional)Edit
  • If Force is not specified, Value is automatically adjusted to 1 / (EndTime - StartTime).
  • Force is disabled by default.
  • Currently, you need to leave a space after Force (i.e. (Force ) is fine but (Force) will throw an error). This bug will be fixed in the next release.

Material specific output definitionsEdit

Linear material outputsEdit

This section defines the material dependent parameters and computation related to outputs for CellParam::CellType = linear.

OutputFieldEdit

The fields requested,
value\in\{all, ex,ey,ez,hx,hy,hz\}.

  • OutputField can be any combination of above values. For example, (OutputField ex ey hz) will only output Ex-, Ey- and Hz-fields.
  • H-fields are output as H / \sqrt{c\cdot\epsilon_{0}}, where H unit = amperes / meter.
  • E-fields are output as E / \sqrt{c\cdot\mu_{0}}, where E unit = volts / meter.

Example user input file 1Edit

Field imageEdit

Below is an image of the Ey field generated from running the following input file.

Ey69696

User inputEdit

  # This is a 2D simulation of a waveguide in the \Gamma{-}K direction of
  # a hexagonal lattice photonic crystal consisting of air holes in Silicon.
  # The outer holes of the waveguide have a reduced radius.
  # 2D here means only one cell is used in the z-direction with Bloch kz = 0.
   (Option CellParam      (CellType linear)
                          (FieldType complex)
   )
   (Option BoundaryParam  (Type bloch pml bloch)
                          (PmlParam (Type upml) (Size (Y 1.2) )
                           // SigmaFactor here gives R[0]=Exp[-10] 
                          (SigmaFactor (Y 0.26) ) )
                          (BlochK (X 0.3660) (Z 0) )
   )
   (Option DimensionParam (Center 0 0 0)
                          (Resolution  0.05 0.05 0.125)
                          (Size 1 20.4 0.125)
                          (TimeParam (End 3500) )
   )
   // To enable CheckInput, remove "\".
   (\Option CheckInput)
   // File will be output to structure.h5 by default.
   (Option OutputStructure)
   // RelPermittivity here is from effective index of  Si.
   (Material Lossless (MaterialName Si)
                      (MaterialToken 3)
                      (RelPermittivity 8.6696892)
                      (RelPermeability 1)
   )
   (Source GaussTime (SourceName Source1)
                     (Duration 1200)
                     (Center 600)
                     (Width 125)
                     (Frequency 0.251)
                     (Phase 0)
                     (Amplitude 1)
                     (Direction 0 1 0)
   )
   (Output FieldTime (Position 0 0 0)
                     (StartTime 0)
                     (EndTime 3500)
                     (TimeStep  0.19)
                     (FileName field_time_center.dat)
                     (OutputField ey hz)
   )
   (Output FieldSpace (StartTime 1500)
                      (EndTime 1500.3)
                      (TimeStep 0.2)
                      (FileName field_space)
                      (OutputField ey)
    )
    // Build a source.
    (Object Point (SourceName Source1) (Center 0 0 0) )

    // Build  tri-lattice air-holes in Si.
    (Object  Block (MaterialName Si)
                   (Center 0 0 0)
                   (Axis0  1 0 0)
                   (Axis1  0 1 0)
                   (Axis2  0 0 1)
                   (Size 1 20.0 0.125)
    )
    //Note : vacuum is the default material in FDTD+ (it is generated internally).
    (Object Cylinder (MaterialName vacuum)
                     (Center 0 -8.66025403784 0)
                     (Axis 0 0 1)
                     (Radius 0.3)
                     (Height 0.125)
                     (Array  (NoPeriod 1 11 1)
                             (Period 1 1.73205080757 1)
                             (PeriodAxis0 1 0 0)
                             (PeriodAxis1 0 1 0)
                             (PeriodAxis2 0 0 1)
                     )
    )
    (Object Cylinder (MaterialName vacuum)
                     (Center -0.5 -9.52627944163 0)
                     (Axis 0 0 1)
                     (Radius 0.3)
                     (Height 0.125)
                     (Array (NoPeriod 2 12 1)
                            (Period 1 1.73205080757 1)
                            (PeriodAxis0 1 0 0)
                            (PeriodAxis1 0 1 0)
                            (PeriodAxis2 0 0 1)
                      )
    )
     // Create a waveguide.
    (Object Cylinder (MaterialName Si)
                     (Center 0 0 0)
                     (Axis 0 0 1)
                     (Radius 0.35)
                     (Height 0.125)
    )
    // These two instructions reduce the hole size of the outer row of holes
    // surrounding the waveguide.
    (Object Cylinder (MaterialName Si)
                     (Center -0.5 -0.866025403784 0)
                     (Axis 0 0 1)
                     (Radius 0.35)
                     (Height 0.125)
                     (Array  (NoPeriod 2 2 1)
                             (Period 1 1.73205080757 0)
                             (PeriodAxis0 1 0 0)
                             (PeriodAxis1 0 1 0)
                             (PeriodAxis2 0 0 1)
                      )
    )
    (Object Cylinder (MaterialName vacuum)
                     (Center -0.5 -0.866025403784 0)
                     (Axis 0 0 1)
                     (Radius 0.265)
                     (Height 0.125)
                     (Array  (NoPeriod 2 2 1)
                             (Period 1 1.73205080757 0)
                             (PeriodAxis0 1 0 0)
                             (PeriodAxis1 0 1 0)
                             (PeriodAxis2 0 0 1)
                      )
    )

Example Material::Prp parametersEdit

Gold (Au) Prp parametersEdit

      (Material Prp (MaterialName Au)   // MaterialName value is user defined.
                    (MaterialToken 8)   // MaterialToken value is user defined.
                    (EpsilonInfinity 1)
                    (RelPermeability 1)
                    (PoleResiduePair
                        (Pair1 (2.2217e16, -6.3742e18) (-1.4982e14, 3.0413e11) )
                        (Pair2 (-6.4292e15, 1.3091e18) (-2.8617e13, -3.4173e13) )
                        (Pair3 (2.1046e16, 1.4907e15)  (-3.2597e15, -4.3949e15) )
                        (Pair4 (1.0807e15, 5.889e14)   (-4.4078e14, -3.9961e15) )
                        (Pair5 (-7.173e15, 2.5124e18)  (-3.3045e12, -1.0357e13) )
                        (Pair6 (-5.3709e15, 2.517e15)  (-1.9179e15, 7.2146e15) )
                        (Pair7 (-7.8496e14, 3.3901e18) (2.39e13, -5.5125e12) )
                        (Pair8 (-7.1488e15, 5.1946e18) (4.5398e13, 1.2772e12) )
                    )
      )

Silver (Ag) Prp parametersEdit

      (Material Prp (MaterialName Ag)   // MaterialName value is user defined.
                    (MaterialToken 9)   // MaterialToken value is user defined.
                    (EpsilonInfinity 1)
                    (RelPermeability 1)
                    (PoleResiduePair
                          (Pair1  (2.004847623185898e16, 6.373518292877005e18)
                                  (-1.450986957494932e13, -1.471286657124635e13)
                          )
                          (Pair2 (4.419387050893835e15, -2.455715772083395e16)
                                 (-1.617050856766036e15, -7.484019360934898e14)
                          )
                          (Pair3  (-1.253762773346362e16, 1.111152157244156e18)
                                  (-1.621233619602248e16, -2.197813817662167e15)
                          )
                          (Pair4 (2.397213142810919e15, -2.295684023806528e16)
                                 (-3.358496690509831e15, -4.204829676212065e15)
                          )
                          (Pair5 (-1.287525821691552e16, -2.726588225501935e15)
                                 (-1.826237040321377e14, -1.915550886065391e14)
                          )
                          (Pair6 (3.258726396082003e15, -1.140962611954672e17)
                                 (-9.157642928770140e15, -1.161886436310815e16)
                          )
                          (Pair7 (1.485343112935528e15, -4.055619609190444e14)
                                 (-4.591136080183809e14, -6.063844504839434e15)
                          )
                          (Pair8 (3.256535035230075e13, -8.975107704148206e13)
                                 (-2.459619405093500e14, -6.564514351528399e15)
                          )
                    )
      )


Silicon (Si) Prp parametersEdit

      // Fit for wavelengths up to ~800nm.
      (Material Prp (MaterialName Si)   // MaterialName value is user defined.
                    (MaterialToken 2)   // MaterialToken value is user defined.
                    (EpsilonInfinity 1)
                    (RelPermeability 1)
                    (PoleResiduePair
                       (Pair1 (-1.2393e17, -1.5487e17)
                               (-6.3301e14, -6.4586e15)
                       )
                       (Pair2 (-3.6772e16, 7.6757e16)
                               (-7.3714e14, -6.1149e15)
                       )
                       (Pair3 (1.7097e17, 1.0258e17)
                               (-5.666e14, -6.3804e15)
                       )
                       (Pair4 (-1.2185e16, 7.215e15)
                               (-3.3952e15, -8.0514e15)
                       )
                       (Pair5 (3.2464e15, -3.7598e15)
                               (-1.952e14, 5.126e15)
                       )
                       (Pair6 (-1.0e11, 7e13)
                               (-1e15, 3e15)
                       )
                       (Pair7 (-5000000, 2.93e12)
                               (-3e14, 1.9435e15)
                       )
                       (Pair8 (-5000, -2.5e12)
                               (-1.98e14, 1.1e15)
                       )
                       (Pair9 (-5000, -4.6e12)
                               (-2.5e14, 2.5e15)
                       )
                       (Pair10 (-5000, -1e12)
                               (-2.98e14, 4.75e14)
                       )
                    )
      )

Reading output data with MathematicaEdit

   (* Set the directory location. *)
   In[5]:= SimNo = "sim68";
   In[6]:= DirName = "C:\Research\FdtdPlusSimulations\\" <> SimNo ; SetDirectory[DirName]
   Out[6]= "C:\\Research\\FdtdPlusSimulations\\sim68"
   
   (* File name.  *)
   In[7]:= FileName = "field_time.dat";
   (* Read the first few lines (optional step; helps you see the headers). *)
   In[8]:= ReadList[FileName, Word, 4, RecordLists -> True, WordSeparators -> {"\n"}]
   Out[8]= {{"time_value	hz	ey	"}, {"0	(0,0)	(0,0)	"}, {"0.189038	\
                (-0.000608598,0.000132897)	(0.000814581,-0.000388972)	"}, {"0.378076	\
                (-0.000812379,1.68512e-05)	(0.00168293,-0.000526529)	"}}


   (* Read the file. *)
   (* Data is stored in list MyData. *)
   (* Note header above (i.e. "time value hz ey") is removed. *)
   In[9]:= MyString = Import[FileName, "Text"]; MyString =
   StringReplace[MyString, {"(" -> "", ")" -> "", "," -> "\t"}]; MyStream =
   StringToStream[MyString]; MyData = {}; Read[MyStream, String]; MyData =
   ReadList[MyStream, Real, RecordLists -> True,
   RecordSeparators -> {"\n"}]; Close[MyStream]; Clear[MyString];
   (* Output a line from MyData. *)   
   In[10]:= MyData[[2]]
   Out[10]= {0.189038, -0.000608598, 0.000132897, 0.000814581, -0.000388972}

Reading output data with MatlabEdit


   [t eyr eyi hzr hzi] = textread('field_time_center.dat',
                                  '%f%*c%f%*c%f%*c%*c%f%*c%f%*c','headerlines',1);

Ad blocker interference detected!


Wikia is a free-to-use site that makes money from advertising. We have a modified experience for viewers using ad blockers

Wikia is not accessible if you’ve made further modifications. Remove the custom ad blocker rule(s) and the page will load as expected.