[table of contents] [master index] [comments] [modules] [programs] [variables] [types] [procedures]
PURPOSE
Define crystal unit cell and lattice basis vectors
SOURCE
module unit_cell_mod use const_mod use io_mod implicit none private ! Public procedures public :: input_unit_cell ! read dim, crystal_system, cell_param public :: output_unit_cell ! write dim, crystal_system, cell_param public :: standard_cell_param ! return parameters a,b,c, alpha,beta,gamma public :: make_unit_cell ! create lattice basis vectors, etc. public :: define_unit_cell ! reset cell parameters public :: make_G_basis ! make G_basis from R_basis ! Public variables public :: crystal_system, N_cell_param, cell_param public :: R_basis, G_basis, dGG_basis character(30) :: crystal_system ! type of crystal cell (cubic, etc.) integer :: N_cell_param ! # of unit cell parameters real(long) :: cell_param(6) ! unit cell parameters real(long) :: R_basis(:,:) ! (dim,dim) lattice bases a_i real(long) :: G_basis(:,:) ! (dim,dim) reciprocal bases b_i real(long) :: dGG_basis(:,:,:) ! (dim,dim,6) derivatives of b_i.b_j
VARIABLE
character(30) crystal_system = string identifying crystal system Allowed values: 3D crystal systems (dim = 3): 'cubic', 'tetragonal', 'orthorhombic', 'monoclinic', 'hexagonal', 'triclinic' 2D crystal systems (dim = 2) 'square', 'rectangular', 'rhombus', 'hexagonal', 'oblique' 1D crystal system (dim = 1): 'lamellar'
VARIABLE
integer N_cell_param = # of unit cell parameters Different values needed for different crystal systems (e.g., 1 for cubic, 6 for triclinic)
VARIABLE
real(long) cell_param(6) - array of cell parameters Only elements 1:N_cell_param are used
VARIABLE
real(long) R_basis(:,:) - dimension(dim,dim) R_basis(i,:) = a_i = Bravais lattice basis vector i
VARIABLE
real(long) G_basis(:,:) - dimension(dim,dim) G_basis(i,:) = b_i = reciprocal lattice basis vector i
VARIABLE
real(long) dGG_basis(:,:,:) - dimension(dim,dim,6) dGG_basis(i,j,k) = d(b_i.dot.b_j)/d(cell_param(k)) Needed in calculation of stress by perturbation theory
SUBROUTINE
input_unit_cell(i_unit, fmt_flag)
PURPOSE
Read data needed to construct unit cell from file i_unit. Inputs dim, crystal_system, N_cell_param, and cell_param If necessary, allocates R_basis, G_basis, related arrays
ARGUMENTS
integer i_unit - unit # of input file character(1) fmt_flag - flag specifying input format
COMMENT
Allowed values of fmt_flag: F -> formatted ascii input U -> unformatted input
SOURCE
subroutine input_unit_cell(i_unit,fmt_flag) integer, intent(IN) :: i_unit character(len = 1), intent(IN) :: fmt_flag
SUBROUTINE
output_unit_cell(o_unit,fmt_flag)
PURPOSE
Write crystal_system, N_cell_param, and cell_param to file
ARGUMENTS
integer o_unit - unit # of output file character(1) fmt_flag - flag specifying output format
COMMENT
Allowed values of fmt_flag: F -> formatted output U -> unformatted output
SOURCE
subroutine output_unit_cell(o_unit,fmt_flag) integer, intent(IN) :: o_unit character(len = 1), intent(IN) :: fmt_flag
FUNCTION
standard_cell_param(cell_param)
PURPOSE
Returns array (a, b, c, alpha, beta, gamma) for 3-d systems a, b, c are lengths of the three Bravais basis vectors alpha is the angle beween b, c beta is the angle between a, c gamma is the angle between a, b
RETURN
standard_cell_param(1:3) = (a,b,c) standard_cell_param(4:6) = (alpha,beta,gamma)
SOURCE
function standard_cell_param(cell_param) real(long), dimension(6), intent(IN) :: cell_param real(long), dimension(6) :: standard_cell_param
SUBROUTINE
make_unit_cell
PURPOSE
Constructs Bravais and reciprocal lattice vectors, and related arrays, from knowledge of module input variables.
COMMENT
All inputs and outputs are module variables, rather than explicit parameters. Inputs: crystal_system, N_cell_param, and cell_param Outputs: R_basis, G_basis, dRR_basis, dGG_basis
SOURCE
subroutine make_unit_cell
SUBROUTINE
define_unit_cell( mylattice, my_N, my_param )
PURPOSE
Modify crystal system and/or unit cell parameters
ARGUMENTS
mylattice - crystal system my_N - number of cell parameters my_param - array of cell parameters
SOURCE
subroutine define_unit_cell( mylattice, my_N, my_param ) character(*), intent(IN) :: mylattice integer, intent(IN) :: my_N real(long), intent(IN) :: my_param(:)
SUBROUTINE
make_G_basis(R_basis,G_basis)
PURPOSE
Construct array G_basis of reciprocal lattice basis vectors from input array R_basis of Bravais lattice basis vectors
SOURCE
subroutine make_G_basis(R_basis,G_basis) use group_mod, only : Inverse real(long), intent(IN) :: R_basis(:,:) ! (dim,dim) Bravais real(long), intent(OUT) :: G_basis(:,:) ! (dim,dim) reciprocal