4. Calculating and plotting band diagrams
2. Writing a control file
The best way to learn to do this is to examine a file, with a line by line explanation. The following control file specifies a lattice, geometric objects to occupy the lattice points, and a set of k-points for which optical modes will be calculated.
2. Writing a control file
The best way to learn to do this is to examine a file, with a line by line explanation. The following control file specifies a lattice, geometric objects to occupy the lattice points, and a set of k-points for which optical modes will be calculated.
Let’s take this line by line.
;ctl for the calculation of Inverse Opals
Everything after a semicolon is ignored by the program so that user notes can be made within a file.
(define-param r (/ 1 (* 2 (sqrt 2)))) ; radius of the spheres
Here, the function is “define-param” which has the format (define-param parameter_name parameter_value). The parameter name is something the user chooses, I chose r for the radius. The parameter value is (/ 1 (* 2 (sqrt 2))), which is how 1/(2√2) is written using Schemes syntax for mathematical expressions. The functions “/” (divide) and “*” (multiply) are called the same way as before, in parenthesis with the function name followed by arguments, for example (/ numerator denominator). Functions are evaluated starting from the innermost parenthesis.
(define-param n 2.2); refractive index
(define-param eps (* n n)) ; the dielectric constant of structure
It is smart to define the parameters that are used, and then use their names later in the control file. By doing this if something is changed, like the dielectric constant of the structure, it only needs to be changed in one place. Also, it is clearer later in the file if named parameters are used rather than numbers.
(set! geometry-lattice (make lattice
(basis-size (sqrt 0.5) (sqrt 0.5) (sqrt 0.5))
(basis1 0 1 1)
(basis2 1 0 1)
(basis3 1 1 0)))
The ctl file is written in Scheme syntax. A function or command is called in parenthesis, followed by what is passed on to that function. Here, the function set! is called and passed the parameter to be set, geometry-lattice, followed by the definition of the lattice (make lattice …).
(set! geometry-lattice (make lattice
(basis-size (sqrt 0.5) (sqrt 0.5) (sqrt 0.5))
(basis1 0 1 1)
(basis2 1 0 1)
(basis3 1 1 0)))
These few lines of code define a FCC lattice, irrespective of what occupies the lattice points. The function set! is part of MPB, and is used whenever a major program setting is defined. geometry-lattice specifies the basis vectors and lattice size of the computational cell (centered at the origin of the Cartesian coordinate system in real space). The units are in a, the lattice constant, but remember that when describing vectors using a basis that those vectors are in units of the basis which is in turn in units of a. basis-size defines the length of the lattice vectors. The closest neighbors to a given lattice point in an FCC lattice is sqrt(0.5) or (1/√2) away in all three directions (this is the distance between a corner and the center of a face of a cube). The three basis vectors (primitive lattice vectors), basis1, basis2, and basis3 point to the three nearest neighbors (in this case the center of the three faces intersecting on a corner of a cube). The three lattice vectors in this case point in each direction basis1, basis2, and basis3 with a length of (sqrt 0.5). In other words, taken from the first part of the output of this example:
init-params: initializing eigensolver data
Computing 16 bands with 1.000000e-07 tolerance.
Working in 3 dimensions.
Grid size is 64 x 64 x 64.
Solving for 8 bands at a time.
Creating Maxwell data...
Allocating fields...
Mesh size is 8.
Lattice vectors:
(0, 0.5, 0.5)
(0.5, 0, 0.5)
(0.5, 0.5, 0)
Cell volume = 0.25
Checking that everything is right, the length of the 1st lattice vector is √((1/2)^2+(1/2)^2 )=√(1/2) which is what we specified in the control file. Moving onto the next line of input:
(define diel (make dielectric (epsilon eps)))
This line “makes” or defines a dielectric material with epsilon “eps” and calls that material “diel.” In other words, the program is being told “Get ready! We are going to use a material with a certain dielectric soon, we call this material ‘diel’ and it has an epsilon of ‘eps.’”
(set! geometry (list (make block (center 0 0 0) (size 2 2 2) (material diel))
(make sphere (center 0 0 0) (radius r)(material air))))
set! geometry is used to define or set the geometry of what occupies each individual lattice point. The geometry starts blank, then each entry in the list overwrites what was there before. Anything that extends outside of the primitive unit cell volume is ignored. For example, we start with a block of material with a particular dielectric constant centered at the origin with size 2×2×2 (in units of a) which is WAY too big to fit in the primitive cell, but we want to be sure to start with something uniform that fills up the whole space. Then, we overwrite the dielectric with air spheres. Complex geometries can be made this way.
(set-param! resolution 64) ; use a 64x64x64 grid
(set-param! resolution 64) is how the three dimensional grid describing the real-space geometry is broken up within the primitive unit cell. Basically, the calculations work by specifying a matrix of matrixes. Each point in the primitive unit cell is specified by a matrix, the values of that matrix come from the geometry specified earlier and the mesh size. The resolution is the number of points per dimension that go into this matrix. Depending on the geometry, it may be fine to go with a resolution of 32 or 16. As the resolution goes up, the calculation time rises exponentially.
(set-param! mesh-size 8) sets the dimensions of the sub-matrixes that make up each point. It is how the average dielectric of a point is computed. The higher the mesh is, the more sensitive the calculation is to fine changes in the dielectric structure.
(set-param! num-bands 16) is the number of bands each k-point will be computed for. Each higher order band has a k-point that is an integral multiple of the first band.
(set-param! eigensolver-block-size 4) is the number of bands which the program solves at a time. A small number requires less memory, but requires more iteration. A large number may require too much memory, making it slower in spite of the fact that less iteration is needed. If calculation time is a critical issue (like when density of optical states is calculated) then it is wise to run several test calculations with different eigensolver bloack sizes to see which is fastest.
(set! k-points (interpolate 16 (list
(vector3 0.00000 0.50000 0.00000); L
(vector3 0.00000 0.00000 0.00000); GAMMA
(vector3 0.00000 0.50000 0.50000); X
(vector3 0.00000 0.62500 0.37500); U
(vector3 0.25000 0.75000 0.50000); W
(vector3 0.37500 0.75000 0.37500)))); K
set! k-points is used to specify the wavelength and direction of light within the lattice. For technical reasons it is more convenient to calculate optical modes when both the lattice and optical modes are expressed in what is called reciprocal space. It is difficult to visualize what a lattice looks like in reciprocal space, so the geometry is defined in real space in the control file and converted to reciprocal space later. The direction and wavelength of light is expressed in reciprocal space by convention in the literature, and so is used here. What follows set! k-points is a list of k-points that span the most important directions and wavelengths in the FCC lattice. These directions are from, in reciprocal space, the center of the first irreducible brillouin zone (a representation of the lattice in reciprocal space) called Γ, to each of the high symmetry points of the brillouin zone and along the edges.
![]() |
The first Brillouin zone of the FCC lattice, with high symmetry points labeled. |
(set! k-points (interpolate 16 (list
(vector3 0.00000 0.50000 0.00000); L
(vector3 0.00000 0.00000 0.00000); GAMMA
(vector3 0.00000 0.50000 0.50000); X
(vector3 0.00000 0.62500 0.37500); U
(vector3 0.25000 0.75000 0.50000); W
(vector3 0.37500 0.75000 0.37500); K
)))
(interpolate 16 (list... creates a list consisting of 8 points interpolated between each of the items of the original list, roughly multiplying the total number by 16. The values of the high symmetry k-points were found in the literature and in MPB’s online documentation. The order of the k-points in the list was chosen so that points would be interpolated along an edge of the Brillouin zone. We start at L, then travel along the line LΓ to Γ and then along ΓX to X and across to U. It would be a mistake to go from X directly to L because the line made by XL does not lay on the surface of the Brillouin zone.
(run)
This final command starts the calculation now that all the settings and definitions have been taken in by the program. There are different versions of the run function that do different things, see the MPB documentation for examples. Before the optical modes are examined, it is wise to take a look at the dielectric function defined in the control file to be sure it is what was intended.
No comments:
Post a Comment