[table of contents] [master index] [comments] [programs] [procedures] [variables] [types] [modules]
PURPOSE
Calculate monomer concentrations, free energies, and stresses.
Solve the modified diffusion equation for chains by a pseudo-spectral
algorithm. Use a simple Boltzmann weight on a grid for solvent.
SOURCE
module scf_mod
use const_mod
use chemistry_mod
use fft_mod
use grid_mod
use grid_basis_mod
use chain_mod
use step_mod
implicit none
private
! public procedures
public:: density_startup ! allocates arrays needed by density
public:: density ! scf calculation of rho & q
public:: scf_stress ! calculates d(free energy)/d(cell_param)
public:: mu_phi ! calculates mu from phi (canonical)
! or phi from mu (grand canonical)
public:: free_energy ! calculates helmholtz free energy
! (optionally calculates pressure)
public:: free_energy_FH ! Flory-Huggins helmholtz free energy
public:: set_omega_uniform ! sets k=0 component of omega (canonical)
public:: divide_energy ! calculates different components of free energy
! public module variable
public:: plan ! module variable, used in iterate_mod
VARIABLE
type(fft_plan) plan - Plan of grid sizes etc. used for FFTs
(Public because its used in iterate_mod)
SUBROUTINE
subroutine density_startup(N_grids,extr_order,chain_step,update_chain)
PURPOSE
Initialize FFT_plan, grid_data.
Allocate or update/re-allocate memory for chains
ARGUMENTS
N_grids = grid dimensions
extr_order = Richardson extrapolation order
chain_step = the discretized chain segment length
update_chain = true if simply the chain memory need to be re-allocated
SOURCE
subroutine density_startup(N_grids, extr_order, chain_step, update_chain) implicit none integer, intent(IN) :: N_grids(3) ! # of grid points in each direction integer, intent(IN) :: extr_order real(long), intent(IN) :: chain_step logical, intent(IN) :: update_chain
SUBROUTINE
density(N,omega,rho,qout)
PURPOSE
Main SCFT calculation. Solve the modified diffusion equation
for all polymer species, and calculate monomer density field
for all monomer types.
ARGUMENTS
N = # of basis functions
omega(N_monomer,N) = chemical potential
rho(N_monomer,N) = monomer density fields
qout(N_chain) = partition functions
q_solvent(N_solvent) = partition functions of solvent molecules
COMMENT
density_startup should be called prior to density to
allocate arrays used by density and scf_stress.
SOURCE
subroutine density( &
N, & ! # of basis functions
omega, & ! (N_monomer,N)chemical potential field
rho, & ! (N_monomer,N) monomer density field
qout, & ! (N_chain) 1-chain partition functions
q_solvent & ! (N_solvent) solvent partition functions
)
implicit none
integer, intent(IN) :: N
real(long), intent(IN) :: omega(:,:)
real(long), intent(OUT) :: rho(:,:)
real(long), intent(OUT), optional :: qout(N_chain)
real(long), intent(OUT), optional :: q_solvent(N_solvent)
SUBROUTINE
solvent_density(monomer,s_size,omega,rho_grid,bigQ_solvent)
PURPOSE
to calculate the density profile of a solvent specie
ARGUMENTS
monomer - monomer type of the solvent species
s_size - volume occupied by solvent molecule / reference volume
(volume in units where reference volume = 1)
omega - omega fields on grid, per reference volume
rho_grid - density fields on grid
bigQ_solvent - spatial average of Boltzmann factor exp(-s_size*omega)
SOURCE
subroutine solvent_density(monomer,s_size,omega,rho_grid,bigQ_solvent) implicit none real(long),intent(IN) :: s_size real(long),intent(IN) :: omega(0:,0:,0:,:) integer,intent(IN) :: monomer real(long),intent(INOUT) :: rho_grid(0:,0:,0:,:) real(long),intent(OUT) :: bigQ_solvent
SUBROUTINE
chain_density(i_chain, chain, ksq, omega)
PURPOSE
solve the PDE for a single chain
evaluate the density for each block
ARGUMENTS
i_chain - index to the chain
chain - chain_grid_type, see chain_mod
ksq - k^2 on grid, initialized in grid_mod
omega - omega fields on grid
SOURCE
subroutine chain_density(i_chain, chain, ksq, omega) implicit none integer,intent(IN) :: i_chain type(chain_grid_type),intent(INOUT) :: chain real(long),intent(IN) :: ksq(0:,0:,0:) real(long),intent(IN) :: omega(0:,0:,0:,:)
FUNCTION
scf_stress(N, size_dGsq, dGsq )
RETURN
real(long) array of dimension(size_dGsq) containing
derivatives of free energy with respect to size_dGsq
cell parameters or deformations
ARGUMENTS
N = number of basis functions
size_dGsq = number of cell parameters or deformations
dGsq = derivatives of |G|^2 w.r.t. cell parameters
dGsq(i,j) = d |G(i)|**2 / d cell_param(j)
COMMENT
Requires previous call to density, because scf_stress
uses module variables computed in density.
SOURCE
function scf_stress(N, size_dGsq, dGsq ) implicit none integer, intent(IN) :: N integer, intent(IN) :: size_dGsq real(long), intent(IN) :: dGsq(:,:)
SUBROUTINE
divide_energy(rho, omega, phi_chain, phi_solvent, Q, f_comp, ovlap)
PURPOSE
Divide free energy into components arising from binary
interaction free energy and from chain entropy
ARGUMENTS
rho = density fields
omega = potential fields
phi_chain = volume fraction of species (chain)
phi_solvent = volume fraction of species (solvent)
Q = partion function of species (chain)
f_tot = total free energy
f_comp = components of free energy (see below)
ovlap = overlap integrals
COMMENT
a) Components of f_comp array:
f_comp(1) = overall interaction energy
f_comp(2) = conformational energy of first block
f_comp(3) = conformational energy of last block
f_comp(4) = junction translational energy (diblock)
b) Calculation of junction translational entropy is correct
only for diblocks, for which there is only one junction
c) Components of overlap integral array ovlap can be used
to divide interaction energy into components arising from
interactions between specific pairs of monomer types.
SOURCE
subroutine divide_energy(rho, omega, phi_chain, phi_solvent, Q, f_tot, f_comp, ovlap) implicit none real(long), intent(IN) :: rho(:,:) ! monomer vol. frac fields real(long), intent(IN) :: omega(:,:) ! chemical potential field real(long), intent(IN) :: phi_chain(:) ! molecule vol. frac of chain mol real(long), intent(IN) :: phi_solvent(:) ! molecule vol. frac of solvent mol real(long), intent(IN) :: Q(:) ! chain partition functions real(long), intent(IN) :: f_tot ! components of free energy real(long), intent(OUT) :: f_comp(:) ! components of free energy real(long), intent(OUT) :: ovlap(:,:) ! overlap integrals,N_monomer**2
SUBROUTINE
set_omega_uniform(omega)
PURPOSE
Sets uniform (k=0) component of field omega to convention
omega(:,1) = chi(:,:) .dot. phi_mon(:)
corresponding to vanishing Lagrange multiplier field
SOURCE
subroutine set_omega_uniform(omega) real(long), intent(INOUT) :: omega(:,:)
SUBROUTINE
mu_phi(mu_chain,phi_chain,q,mu_solvent,phi_solvent,q_solvent)
PURPOSE
If ensemble = 0 (canonical), calculate mu from phi
If ensemble = 1 (grand can), calculate phi from mu
ARGUMENTS
mu_chain(N_chain) = chemical potentials of chain(units kT=1)
phi_chain(N_chain) = molecular volume fractions of chain
q(N_chain) = single chain partition functions
mu_solvent(N_solvent) = chemical potentials of solvent species
phi_solvent(N_solvent) = molecular volume fractions of solvent species
q_solvent(N_solvent) = partition functions of solvent species
SOURCE
subroutine mu_phi(mu_chain,phi_chain,q,mu_solvent,phi_solvent,q_solvent) real(long), intent(INOUT) :: mu_chain(N_chain) real(long), intent(INOUT) :: phi_chain(N_chain) real(long), intent(IN) :: q(N_chain) real(long), intent(INOUT) :: mu_solvent(N_solvent) real(long), intent(INOUT) :: phi_solvent(N_solvent) real(long), intent(IN) :: q_solvent(N_solvent)
SUBROUTINE
free_energy( N, rho, omega, phi_chain, mu_chain, phi_solvent,
mu_solvent, f_Helmholtz, [pressure] )
PURPOSE
Calculates Helmholtz free energy / monomer and (optionally)
the pressure, given phi, mu, and omega and rho fields
SOURCE
subroutine free_energy(N, rho, omega, phi_chain, mu_chain, &
phi_solvent, mu_solvent, f_Helmholtz, pressure )
integer, intent(IN) :: N ! # of basis functions
real(long), intent(IN) :: rho(:,:) ! monomer vol. frac fields
real(long), intent(IN) :: omega(:,:) ! chemical potential field
real(long), intent(IN) :: phi_chain(:) ! molecule vol. frac of chain species
real(long), intent(IN) :: mu_chain(:) ! chemical potential of chain species
real(long), intent(IN) :: phi_solvent(:) ! molecule vol. fraction of solvent species
real(long), intent(IN) :: mu_solvent(:) ! chemical potential of solvent species
real(long), intent(OUT):: f_Helmholtz ! free energy/monomer
real(long), intent(OUT), optional :: pressure
FUNCTION
real(long) function free_energy_FH(phi_chain,phi_solvent)
RETURN
Flory-Huggins Helmholtz free energy per monomer, in units
such that kT =1, for a homogeneous mixture of the specified
composition.
ARGUMENTS
phi_chain(N_chain) = molecular volume fractions of chains
phi_solvent(N_solvent) = molecular volume fractions of solvents
SOURCE
real(long) function free_energy_FH(phi_chain,phi_solvent) real(long), intent(IN) :: phi_chain(N_chain) real(long), intent(IN), optional :: phi_solvent(N_solvent) real(long) :: rho(N_monomer)