[table of contents] [master index] [comments] [modules] [programs] [variables] [types] [procedures]
MODULE
response_pd_mod
PURPOSE
Calculate response of periodic structure to periodic perturbation, for use in construction of approximate Jacobian in iterate_mod. Calculate approximate functional derivatives of monomer concentration field rho with respect to changes in the periodic potential field omega, i.e., derivatives of the coefficients in an expansion of rho with respect to coefficients in the corresponding expansion of omega. Also calculate derivatives of rho with respect to changes in unit cell dimensions. These derivatives can be to construct an approximate Jacobian for use in quasi-Newton iteration schemes.
SOURCE
module response_pd_mod use const_mod use io_mod implicit none PRIVATE PUBLIC :: init_response_pd ! allocate memory to array corrlt PUBLIC :: response_pd_omega ! calculate d_rho/d_omega PUBLIC :: response_pd_cell ! calculate d_rho/d_cell_param PUBLIC :: make_correlation ! block-block correlation of ideal chains PUBLIC :: corrlt ! used by spinodal_mod PUBLIC :: response_dstress_dcell ! calculate dstress/dcell only real(long), allocatable :: corrlt(:,:,:) ! correlation function
SUBROUTINE
subroutine init_response_pd(N_monomer,N)
PURPOSE
Allocate memory for module array corrlt. This routine must be called before make_correlation is called
ARGUMENTS
N_monomer = number monomers N = number of basis functions
SOURCE
subroutine init_response_pd(N_monomer,N) integer :: N_monomer, N
SUBROUTINE
response_pd_omega( N, cut, omega, drho_domega )
PURPOSE
Calculates matrix representation of response function d rho/d omega
ARGUMENTS
N - total number of symmetry-adapted basis functions cut - number of basis functions treated numerically omega - cofficients of monomer chemical potential field drho_domega - response matrix drho_domega(alpha,i,beta,j) = d rho(alpha,i) / d omega(beta,j)
COMMENT
Long wavelength entries, with i,j in [1, cut], are calculated numerically, by numerical differentiation. Short wavelength entries, i,j in (cut, N], are approximated by the response of a homogeneous ideal gas, which is diagonal in k-space.
SOURCE
subroutine response_pd_omega( N, cut, omega, drho_domega ) use chemistry_mod, only : ensemble, N_monomer use unit_cell_mod, only : N_cell_param use scf_mod, only : density integer, intent(IN) :: N integer, intent(IN) :: cut real(long), intent(IN) :: omega(:,:) !(N_monomer, N) real(long), intent(OUT) :: drho_domega(:,:,:,:) !(N_monomer,N,N_monomer,N)
SUBROUTINE
response_pd_cell(N,omega,drho_dcell,dstress_dcell)
PURPOSE
drho_dcell(alpha,i,j) = d rho(alpha,i) / d cell_param(j) dstress_dcell(i,j) = d stress(i) / d cell_param(j)
SOURCE
subroutine response_pd_cell(N,omega,drho_dcell,dstress_dcell) use chemistry_mod, only : ensemble, N_monomer use scf_mod, only : density, scf_stress use grid_mod, only : make_ksq use basis_mod, only : make_dGsq use unit_cell_mod, only : N_cell_param, cell_param, &! G_basis, dGG_basis, make_unit_cell integer, intent(IN) :: N real(long), intent(IN) :: omega(:,:) real(long), intent(OUT) :: drho_dcell(:,:,:) real(long), intent(OUT) :: dstress_dcell(:,:)
SUBROUTINE
response_dstress_dcell(N,omega,drho_dcell,dstress_dcell)
PURPOSE
dstress_dcell(i,j) = d stress(i) / d cell_param(j)
SOURCE
subroutine response_dstress_dcell(N,omega,dstress_dcell) use chemistry_mod, only : ensemble, N_monomer use scf_mod, only : density, scf_stress use grid_mod, only : make_ksq use basis_mod, only : make_dGsq use unit_cell_mod, only : N_cell_param, cell_param, &! G_basis, dGG_basis, make_unit_cell integer, intent(IN) :: N real(long), intent(IN) :: omega(:,:) real(long), intent(OUT) :: dstress_dcell(:,:)
SUBROUTINE
make_correlation(N,[lambda])
PURPOSE
calculate correlation function of ideal chains melt
SOURCE
subroutine make_correlation(N, lambda) use chemistry_mod use group_mod use basis_mod, only : wave_of_star use scf_mod, only : density use unit_cell_mod, only : G_basis integer :: N ! # of stars real(long), optional :: lambda(:) ! $k^2$