[table of contents] [master index] [comments] [modules] [programs] [variables] [types] [procedures]
MODULE
basis_mod - create symmetry-adapted basis functions
PURPOSE
Define data structures and procedures involving reciprocal
lattice vectors, stars, and symmetry-adapted basis functions
SOURCE
module basis_mod
use const_mod ! Defines constant long, and variable dim=1,2, or 3
use io_mod, only : output
use string_mod, only : int_string
use version_mod, only : version_type, output_version
use unit_cell_mod, only : G_basis, make_G_basis, output_unit_cell
use group_mod ! Defines types and operations for crystal groups
use space_groups_mod, only : space_groups
use grid_mod ! Utilities for fft grids
implicit none
private
! Public procedures
public:: make_basis, make_dGsq ! needed by scf
public:: make_waves, make_stars, release_basis
public:: f_basis
public:: reorder_basis
public:: scattering_intensity
public:: output_waves
! Public variables
public:: group
public:: N_star
public:: N_wave, wave, star_of_wave, coeff
public:: G_max, which_wave
public:: star_begin, star_end, star_count
public:: star_cancel, star_invert, basis_sign
public:: wave_of_star
public:: sort
public:: valid_wave
! Declaration of public variables
type(group_type) :: group ! space group
integer :: N_wave ! # of reciprocal lattice G vectors
integer :: N_star ! # of stars
integer :: G_max(3) ! maximum of values G(1), G(2), G(3)
integer :: wave(:,:) ! (dim,N_wave) list of G vectors
real(long) :: Gsq(:) ! (N_wave) list of values of |G|^{2}
complex(long):: coeff(:) ! (N_wave) coefficient of wave in star
integer :: star_of_wave(:) ! (N_wave) index of star containing wave
integer :: which_wave(:,:,:)! (-G_max(1):G_max(1), ..., ... )
! indices of waves listed by components
integer :: star_begin(:) ! (N_star) index of first wave in star
integer :: star_end(:) ! (N_star) index of last wave in star
integer :: star_count(:) ! (N_star) # of G vectors in star
integer :: wave_of_star(:,:)! (dim,N_star) components of one wave from
! the star, serves as a label
logical :: star_cancel(:) ! (N_star) see wave_format
integer :: star_invert(:) ! (N_star) see wave_format
integer :: basis_sign(:) ! (N_star) see wave_format
VARIABLE
integer N_wave = # of G vectors (reciprocal lattice vectors)
VARIABLE
integer N_star = # of stars
VARIABLE
integer G_max(3) = Max values of components of integer waves
VARIABLE
integer wave(:,:) - allocatable, dimension(dim,N_wave)
wave(1:dim,i) = integer wavevectors # i, Bravais basis
VARIABLE
real(long) Gsq(:) - allocatable, dimension(N_wave)
Gsq(i) = value of |G|^{2} for wavevector i
VARIABLE
complex(long) coeff(:) - allocatable, dimension(N_wave)
coeff(i) = complex coefficients of wavevector i
in symmetry adapted function for star
VARIABLE
integer star_of_wave(:) - allocatable, dimension(N_wave)
star_of_wave(i) = index of star containing wave i
VARIABLE
integer which_wave(G1,G2,G3) = index of wavevector G=(G1,G2,G3) If wave(i) = (G1,G2,G3), then which_wave(G1,G2,G3) = i allocatable, dimension( -G_max(1):G_max(1) , ... , ... )
VARIABLE
integer star_begin(:) - allocatable, dimension(N_star)
star_begin(i) = index of first G in star i
VARIABLE
integer star_end(:) - allocatable, dimension(N_star)
star_end(i) = index of last G in star i
VARIABLE
integer star_count(:) - allocatable, dimension(N_star)
star_count(i) = # of G vectors in star i
VARIABLE
integer star_invert(:) - allocatable, dimension(N_star)
star_invert(i) = invert flag for star i, values (0,1,-1)
See discussion of wave_format
VARIABLE
integer basis_sign(:) - allocatable, dimension(N_star)
basis_sign(i) = sign of basis function i under inversion
+1 -> even under inversion
-1 -> odd under inversion
See discussion of wave_format
VARIABLE
integer wave_of_star(:,:) - dimension(dim, N_star)
wave_of_star(:,i) = G vector from star i, serves as a label
See discussion of wave_format
VARIABLE
logical star_cancel(:) - allocatable, dimension(N_star)
star_cancel(i) = true if star i is cancelled
COMMENT
i) G vectors are listed in the array wave in non-decreasing order
of eigen (for some G_basis passed to make_waves), and (after
output from make_stars) are grouped by stars, each of which
forms contiguous block of vectors.
ii) Cancelled waves and stars are included in the lists of waves
and stars only if logical variable keep_cancel is passed to
make_waves with a value .true. If keep_cancel = true,
cancelled waves are asigned values coeff()=0.0 and cancelled
starts are assigned star_cancel = true. If keep_cancel = false,
then star_cancel = false for all stars.
iii) Within each star, vectors are listed in "descending" order,
if the vectors in wave are read as dim=2, or 3 digit numbers,
with more significant digits on the left. For example the
[111] star of a cubic structure is listed in the order:
1 1 1 (first)
1 1 -1
1 -1 1
1 -1 -1
-1 1 1
-1 1 -1
-1 -1 1
-1 -1 -1 (last)
iv) If star number i is closed under inversion, so that the
negative -G of each wave G in the star is also an element
of that star, then we assign star_invert(i) = 0. For
space groups that contain inversion symmetry (i.e., for
centro-symmetric groups), all stars must be closed under
inversion.
v) For non-centrosymmetric space groups, pairs of stars may be
related by inversion, so that the -G of each vector in star
number i is found not in star i but in a second star. Such
pairs of stars that are related by inversion are always listed
consecutively. The first star in such a pair is assigned
star_invert(i) = 1, and the second has star_invert(i+1) = -1.
vi) The combination of conventions (iii) and (v) implies that the
inverse of the first vector in a star with star_invert = 1 is
the last vector in the next star, which is the partner star
with star_invert = -1. More generally, the inverse of the ith
member of the first star in such a pair is the ith to last
vector of the second star in the pair.
vi) wave_of_star is first wave in star for star_invert=0 or 1, and
last wave in star for star_invert=-1
vii) The absolute magnitudes of the coefficients in array coeff()
for (un-cancelled) stars are chosen so that each stars yields
an orthonormal but generally complex function, i.e., values
of coeff() have absolute magnitude 1.0/sqrt(N_in_star).
viii) The phases of the coefficients in the array coeff() are chosen
such that the coefficient of wave_of_star is positive and real.
This is sufficient to guarantee that a vector -G has value of
coeff that is the complex conjugate of the coeff for G, whether
G and -G are in the same star star (for star_invert = 0) or in
two consecutive stars (for star_invert = +-1). It also
guarantees that the star basis function generated by the first
star in an inversion-related pair is the complex conjugate of
that generated by the second star in the pair. As a result, if
the first star in a pair with star_invert = +-1 is multiplied
by a coefficient and the second by the complex conjugate of
this coefficient, the resulting sum of plane waves will yield
a real basis function.
A star with star_invert = 0 will yield a normalized real basis
function, while two linearly independent real basis functions
may be constructed out of a consecutive pair of stars with
star_invert=1 and -1, by multiplying both by 1/sqrt(2) to
obtain a cosine-like function, and by multiplying the first
star by i/sqrt(2) and the second by -i/sqrt(2), to obtain a
sine-like real function. The number of real basis functions
will thus be equal to the number of stars N_star.
ix) The stars are used to construct N_star real basis functions
f_1(r), ... , f_{N_star}(r), which are defined and indexed
as follows:
For stars that are closed under inversion, with star_invert=0,
and a real star basis function phi_j(r), we create a single
real basis function that is is equal to the normalized star
basis function
f_j(r) = phi_j(r) ,
and that is assigned a basis function index equal to star
index of the star.
For pairs of stars that are related by inversion, with star
indices j and j+1, with associated star basis functions
phi_{j}(r) and phi_{j+1}(r) we construct one cosine-like
basis function, with basis i, given by the sum
f_j(r) = [ phi_{j}(r) + phi_{j+1}(r) ]/sqrt(2) ,
which is even under inversion, and a second sine-like basis,
with a basis function index j+1, given by
f_{j+1}(r) = [ phi_{j}(r) - phi_{j+1}(r) ]/sqrt(-2)
which is odd under inversion. Note that sqrt(-2) is a pure
imaginary number.
xi) For non-centrosymmetric groups, basis functions may be
either either even or odd under inversion. These two cases
are distinguished by the value of variable basis_sign:
basis_sign(i) = 1 implies f_{i}(r) = + f_{i}(-r)
basis_sign(i) = -1 implies f_{i}(r) = - f_{i}(-r)
In a centro-symmetric group basis_sign(i) = 1 for all
basis functions.
In a non-centrosymmetric group, a basis function that is
closed under inversion, so that star_invert(i) = 0, may
be either even or odd or under inversion, and may thus
have either value of basis_sign(i). This is why we need
to distinguish basis_sign from star_invert.
In a non-centrosymmetric group, pairs of sine-like and
cosine-like basis functions that are made constructed
from a pairs of inversion-related stars are indexed such
that:
If star_invert(i) = 1, then basis_sign(i) = 1
If star_invert(i) = -1, then basis_sign(i) = 1
SUBROUTINE
release_basis()
PURPOSE
release the memory allocated by constructing basis information
SOURCE
subroutine release_basis
SUBROUTINE
make_basis(R_basis,G_basis,group_name,Gabs_max)
PURPOSE
Generates information about basis functions needed by scf code:
i) Generates reciprocal lattice basis vectors G_basis, using
input R_basis (Bravais basis vectors)
ii) Reads group generators from file named by the character
string group_name, and generates remainder of the group
iii) Calls make_waves and make_stars to generate all waves
and stars with |G| < Gabs_max, excluding cancelled stars
SOURCE
subroutine make_basis( &
R_basis, &! real lattice basis vectors
G_basis, &! reciprocal lattice basis vectors
group_name, &! file for generators of space group
N_grids, &! maximum value of |G|
grid_flag &! true if PS method used
)
logical, intent(IN) :: grid_flag
real(long) , intent(IN) :: R_basis(:,:) ! real lattice basis
real(long) , intent(OUT) :: G_basis(:,:) ! Reciprocal basis
character(*), intent(INOUT) :: group_name ! file for group
integer , intent(IN) :: N_grids(3) ! # of grid points
FUNCTION
real(long) f_basis(i_basis,G_basis,R)
RETURN
Value of basis function i_basis at Cartesian position R(3)
a crystal with reciprocal basis vector array G_basis(3,3)
SOURCE
real(long) function f_basis( &
i_basis, &! index for basis function
G_basis, &! reciprocal basis vectors
R &! Cartesian position
)
! Dummy arguments
integer, intent(IN) :: i_basis
real(long), intent(IN) :: G_basis(3,3), R(3)
SUBROUTINE
output_waves(wave_port)
ARGUMENTS
wave_port = unit number for output file
PURPOSE
Output list of waves to file unit io_port
Format: i, wave, star_of_wave, |coeff|, coeff/|coeff|, |G|
SOURCE
subroutine output_waves(wave_port, group_name) integer :: wave_port character(*) :: group_name
SUBROUTINE
output_stars(o_port,keep_cancel)
ARGUMENTS
o_port = io-port of the wave file
PURPOSE
Output list of stars to supplied file o_port
Format: i, wave_of_star, star_count, star_invert
If keep_cancel = .true., then also print cancel
SOURCE
subroutine output_stars(o_unit,keep_cancel) use io_mod, only : output use string_mod, only : int_string integer, optional :: o_unit logical, optional :: keep_cancel
SUBROUTINE
reorder_basis(new_basis,old_wave,old_basis,swap,new_stars)
PURPOSE
Reorders basis functions in new_wave to match ordering in
old_wave.
COMMENT
Assumes each wave in new_wave and old_wave corresponds to
a single basis function, and that there are no repeats in
old_wave or new_wave
SOURCE
subroutine reorder_basis(new_basis,old_wave,old_basis,swap,new_stars) integer, intent(IN) :: old_wave(:,:) ! (dim,N_old) real(long), intent(IN) :: old_basis(:,:) ! (dim,N_new) real(long), intent(OUT) :: new_basis(:,:) ! (dim,N_max) integer, intent(OUT), optional :: swap(:) ! (N_max) logical, intent(OUT), optional :: new_stars(:) ! (N_max) ! N_max = max(N_old,N_new)
SUBROUTINE
make_dGsq(dGsq,dGG)
PURPOSE
Calculate dGsq given dGG_basis
SOURCE
subroutine make_dGsq( dGsq, dGG) real(long), intent(OUT) :: dGsq(:) real(long), intent(IN) :: dGG(:,:)
SUBROUTINE
scattering_intensity( G_basis, field, contrast, q, scat )
PURPOSE
Calculate the scattering intensities give a set of fields
(e.g. density fields of SCFT solution ) and some scattering
densities (e.g. electron densities)
SOURCE
subroutine scattering_intensity( G_basis, field, contrast, q, scat ) implicit none real(long), intent(IN), dimension(:,:) :: G_basis real(long), intent(IN), dimension(:,:) :: field real(long), intent(IN), dimension(:) :: contrast real(long), intent(OUT), dimension(:) :: q real(long), intent(OUT), dimension(:) :: scat