! ***************************************************************************
! *   Copyright (C) 2007 - 2015 by J.A.Purton                               *
! *   john.purton[@]stfc.ac.uk                                              *
! *                                                                         *
! *   Contributors:                                                         *
! *   -------------                                                         *
! *   A.V.Brukhno - Free Energy Diff, planar pore & algorithm optimizations *
! *   andrey.brukhno[@]stfc.ac.uk abrukhno[@]gmail.com                      *
! *                                                                         *
! *   T.L.Underwood - Lattice/Phase-Switch MC method & optimizations        *
! *   t.l.Underwood[@]bath.ac.uk                                            *
! *                                                                         *
! *   J.Grant - Complexi constant for Ewald refactor and general use        *
! *   j.grant[@]bath.ac.uk                                                  *
! ***************************************************************************

module constants_module

    use kinds_f90

    implicit none

    real(kind = wp), parameter :: BOLTZMAN = 0.831451115_wp ! boltzman constant
    real(kind = wp), parameter ::  PRSFACT = 0.163882576_wp ! conversion factor from katm to internal units

    integer, parameter :: MAXMESH = 10004 ! number of points along mesh
    integer, parameter :: MAXINTH = huge(MAXMESH)
    integer, parameter :: MAXINT1 = huge(MAXMESH)-1
    integer, parameter :: MAXINT2 = huge(MAXMESH)-2
    integer, parameter :: MAXINT3 = huge(MAXMESH)-3
    integer, parameter :: MAXINT4 = huge(MAXMESH)-4
    integer, parameter :: MAXINTL = 2147483647
    integer, parameter :: MAXINTS = 32767
    
    !TU: Starting point for range of integers signifying a 1-4 dihedral pair
    !TU: in exclusion lists
    integer, parameter :: EXCDIH = huge(MAXMESH)/2

    real(kind = wp), parameter ::  PI = 3.14159265358979323846264338328_wp
    real(kind = wp), parameter ::  TWOPI  = 2.0_wp * PI
    real(kind = wp), parameter ::  ROOTPI = 1.77245385090551602729816748334_wp
    real(kind = wp), parameter ::  RTWOPI = 1.0_wp / TWOPI
    real(kind = wp), parameter ::  TORADIANS = PI / 180.0_wp

    real(kind = wp), parameter ::  THIRD = 1.0_wp / 3.0_wp
    real(kind = wp), parameter ::  HALF  = 0.5_wp
    real(kind = wp), parameter ::  DHALF = HALF-1.e-12_wp
    real(kind = wp), parameter ::  DZERO = 1.e-12_wp
    real(kind = wp), parameter ::  FZERO = 1.e-6_wp

!the bits for fh correction - the exponentials are fiddled later
    real(kind = wp), parameter ::  PLANCK = 6.62606597_wp
    real(kind = wp), parameter ::  HBAR   = PLANCK / TWOPI
    real(kind = wp), parameter ::  AVAG   = 6.02214_wp         
    real(kind = wp), parameter ::  BOLT   = 1.3806488_wp

! complexi for use in Ewald refactor and generally
    complex(kind = wp), parameter :: COMPLEXI = (0.0_wp,1.0_wp)    
    complex(kind = wp), parameter :: COMPLEXTWOPI = TWOPI * (0.0_wp,1.0_wp)    

!STRIP_SIZE is the size of the buffer in strip mining
    integer, parameter :: STRIP_SIZE = 32
    
!internal conversion units
    real(kind = wp), parameter ::  CTOINTERNAL = 138935.4835_wp
    real(kind = wp), parameter ::  CTOEV = 14.3997584_wp
    real(kind = wp), parameter ::  EVTOINTERNAL = 9648.530821_wp
    real(kind = wp), parameter ::  KCALTOINTERNAL = 418.4_wp
    real(kind = wp), parameter ::  KJTOINTERNAL = 100.00_wp
    real(kind = wp), parameter ::  KTOINTERNAL = 0.831451115_wp
    real(kind = wp), parameter ::  INTERNALTOEV = 1.0364305e-4_wp

    real(kind = wp), parameter ::  RCORE  = 0.001_wp
    real(kind = wp), parameter ::  RCORE2 = 0.000001_wp

    integer, parameter :: uctrl = 55  ! control file
    integer, save :: uout  = 66 ! output file
    integer, parameter :: ucfg  = 77  ! config file
    integer, parameter :: uarch = 8  ! archive file
    integer, parameter :: urevc = 9  ! revcon file
    integer, parameter :: ustat = 10 ! statistics file
    integer, parameter :: urest = 11 ! restart file
    integer, parameter :: ufld  = 12 ! field file
    integer, parameter :: uplt  = 13 ! data file (PTFILE)
    integer, parameter :: utbl  = 14 ! potential table files
    integer, parameter :: urevi = 15 ! revive file
    integer, parameter :: uzdn  = 16 ! z-density
    integer, parameter :: urdf  = 17 ! rdfs
    integer, parameter :: ueng  = 18 ! energies of atoms
    integer, parameter :: udsp  = 19 ! displacements
    integer, parameter :: ufed  = 20 !AB: FED data file
    integer, parameter :: upsmc = 22 !TU: PSMC data (NB: unit 21 is taken by 'unpt', for 'VOLDAT' - see 'vol_moves.f90')
    integer, parameter :: utm   = 23 !TU: FED TM I/O
    integer, parameter :: uzch  = 26 !AB: charge z-densities
    integer, parameter :: utraj = 30 !AB: trajectory file (for conversion)
    integer, parameter :: uorient = 31 !TU: Periodic dump of molecular orientations 
    integer, parameter :: uyaml  = 32 !TU: YAML data file
    integer, parameter :: uyaml2 = 33 !TU: 2nd YAML data file - required for the Gibbs ensemble 
    integer, parameter :: uarch2 = 34 !TU: 2nd archive file - required for the Gibbs ensemble
    integer, parameter :: udcd   = 35 !TU: DCD trajectory file
    integer, parameter :: udcd2  = 36 !TU: 2nd DCD trajectory file - required for the Gibbs ensemble
    integer, parameter :: uscan = 40 !scp: output for scan

! AB: NPT move types

    integer, parameter :: NPT_OFF  = 0 ! no NPT move
    integer, parameter :: NPT_VEC  =-3 ! NPT for cell matrix - sampling cell vectors, not allowed in slit
    integer, parameter :: NPT_LOG  =-2 ! NPT for orthorhombic/cubic cell - sampling ln(V), same in slit
    integer, parameter :: NPT_INV2 =-4 ! NPT for orthorhombic       cell - sampling inverse linear dimension(s), same in slit
    integer, parameter :: NPT_INV  =-1 ! NPT for cubic              cell - sampling inverse linear dimension(s)
    integer, parameter :: NPT_LIN  = 1 ! NPT for cubic              cell - sampling linear dimension(s)
    integer, parameter :: NPT_LIN2 = 2 ! NPT for orthorhombic       cell - sampling linear dimension(s), specific in slit
    integer, parameter :: NPT_VOL  = 3 ! NPT for orthorhombic/cubic cell - sampling the volume uniformly, same in slit
    integer, parameter :: NPT_GIBBS = 10 ! Volume moves for orthorhombic/cubic cells in Gibbs ensemble (coupled boxes) - sampling volume uniformly
    integer, parameter :: NPT_ORTHOANI = 11 ! NPT for orthorhombic with anisotropic moves - sampling ln(V)

! AB: Trajectory file types

    integer, parameter :: VTK      = -1 ! Legacy VTK (for visualising atom orientations)
    integer, parameter :: DLM      = 0 ! DL_MONTE  (legacy; native)
    integer, parameter :: DLP      = 1 ! DL_POLY   (legacy; hybrid)
    integer, parameter :: DLP2     = 2 ! DL_POLY-2 (Classic)
    integer, parameter :: DLP2DCD  = 3 ! DL_POLY-2 (Classic) + DCD (Charmm; binary)
    integer, parameter :: DLP4     = 4 ! DL_POLY-4 (Latest)
    integer, parameter :: DLP4DCD  = 5 ! DL_POLY-4 (Latest)  + DCD (Charmm; binary)
    integer, parameter :: DCD      = 6 ! DCD  (only Charmm; binary)
    
! AB: below are constants added for 2D-slit extension (see slit_module.f90)

    integer, parameter :: SLIT_EXT0 = 0 ! slit with no wall having short-range interactions
    integer, parameter :: SLIT_EXT1 = 1 ! slit with one wall having short-range interactions
    integer, parameter :: SLIT_EXT2 = 2 ! slit with two walls having short-range interactions

! AB: slit cell types (must be < 0):
    integer, parameter :: SLIT_SOFT1 = -1 ! slit with one soft wall  (no Coulomb)
    integer, parameter :: SLIT_SOFT2 = -2 ! slit with two soft walls (no Coulomb)
    integer, parameter :: SLIT_HARD2 = -3 ! slit with two hard walls (no Coulomb)

! AB: slit MFA types (coultype = -2 for all):
    integer, parameter :: MFA_OFF  = 0 ! no MFA correction (NO/Explicit Coulomb)
    integer, parameter :: MFA_RECT = 1 ! MFA outside rectangular in-box cut-off (-1 from scratch; 1 from input)
    integer, parameter :: MFA_CYLD = 2 ! MFA outside cylindrical in-box cut-off (-2 from scratch; 2 from input)
    !integer, parameter :: MFA_SPHR = 3 ! MFA outside spherical   in-box cut-off (-3 from scratch; 3 from input)

! AB: collective cell types (for quicker checks):
    integer, parameter :: SLIT_EXT1_MFA0 = -10 ! slit with one soft wall, direct in-box Coulomb
    integer, parameter :: SLIT_EXT1_MFA1 = -11 ! slit with one soft wall, MFA outside rectangular in-box cut-off
    integer, parameter :: SLIT_EXT1_MFA2 = -12 ! slit with one soft wall, MFA outside cylindrical in-box cut-off
    integer, parameter :: SLIT_EXT1_MFA3 = -13 ! slit with one soft wall, MFA outside spherical   in-box cut-off

    integer, parameter :: SLIT_EXT2_MFA0 = -20 ! slit with two soft walls, direct in-box Coulomb
    integer, parameter :: SLIT_EXT2_MFA1 = -21 ! slit with two soft walls, MFA outside rectangular in-box cut-off
    integer, parameter :: SLIT_EXT2_MFA2 = -22 ! slit with two soft walls, MFA outside cylindrical in-box cut-off
    integer, parameter :: SLIT_EXT2_MFA3 = -23 ! slit with two soft walls, MFA outside spherical   in-box cut-off

    integer, parameter :: SLIT_HARD_MFA0 = -30 ! slit with two hard walls, direct in-box Coulomb
    integer, parameter :: SLIT_HARD_MFA1 = -31 ! slit with two hard walls, MFA outside rectangular in-box cut-off
    integer, parameter :: SLIT_HARD_MFA2 = -32 ! slit with two hard walls, MFA outside cylindrical in-box cut-off
    integer, parameter :: SLIT_HARD_MFA3 = -33 ! slit with two hard walls, MFA outside spherical   in-box cut-off

! AB: below are constants added for FED calculation (fed_interface_module.f90)

    ! FED FLAVOR IDENTIFIERS
    integer, parameter :: FED_NO = 0 ! NO FED 
    integer, parameter :: FED_GN = 1 ! generic
    integer, parameter :: FED_PS = 3 ! Phase Switch   / PSMC
    integer, parameter :: FED_HW = 4 ! Histogram Reweighting / HRMC
    !integer, parameter :: FED_?? = ??

    ! FED METHODOLY IDENTIFIERS
    integer, parameter :: FED_US = 1 ! Umbrella Sampling
    integer, parameter :: FED_EE = 2 ! Expanded/Extended Ensemble (iterative)
    integer, parameter :: FED_WL = 3 ! Wang-Landau on-the-fly feedback (iterative)
    integer, parameter :: FED_TM = 4 ! Transition Matrix 
    !integer, parameter :: FED_?? = ??

    ! IDENTIFIERS for NATURALLY EVOLVING order parameter (reaction coordinates)
    ! - depend on structure & need re-calculation every time!
    integer, parameter :: FED_PAR_PSMC    = -6 ! PSMC order parameter - difference in energy between passive and active boxes
    integer, parameter :: FED_PAR_PSMC_HS = -5 ! PSMC order parameter - difference in number of hard-sphere overlaps
    integer, parameter :: FED_PAR_ORDR2   = -4 ! Structure/order 2
    integer, parameter :: FED_PAR_ORDR1   = -3 ! Structure/order 1
    integer, parameter :: FED_PAR_DIST2   = -2 ! COM-distance 2
    integer, parameter :: FED_PAR_DIST1   = -1 ! COM-distance 1

    ! IDENTIFIERS for MANUALLY CONTROLLED order parameters (reaction coordinates)
    ! - a priory known/imposed within a range, or well-defined based on system macro-state
    integer, parameter :: FED_PAR_NONE  = 0 ! NO FED - NO PARAMETER
    integer, parameter :: FED_PAR_TEMP  = 1 ! Temperature
    integer, parameter :: FED_PAR_BETA  = 2 ! Beta = 1/(k_B*T)
    integer, parameter :: FED_PAR_LAMB  = 3 ! Kirkwood Lambda - Hamiltonian/coupling factor/scaling
    integer, parameter :: FED_PAR_VOLM  = 4 ! Volume
    integer, parameter :: FED_PAR_DENS  = 5 ! Density
    integer, parameter :: FED_PAR_NMOLS = 6 ! Total number of molecules in the system
    integer, parameter :: FED_PAR_LCS2  = 7 ! Liquid crystal scalar order parameter - 2nd order Legendre polynomial
    integer, parameter :: FED_PAR_LCS4  = 8 ! Liquid crystal scalar order parameter - 4th order Legendre polynomial
    integer, parameter :: FED_PAR_EXTRA = 9 ! Some other macro-state parameter - ???


!TU: IDENTIFIERS for MOLECULE ORIENTATION class

    integer, parameter :: ORIENTATION_NONE = 0 ! No orientation, e.g. for single-atom molecules
    integer, parameter :: ORIENTATION_LIN  = 1 ! Linear molecule: orientation completely determined by only one vector
    integer, parameter :: ORIENTATION_GEN  = 2 ! General molecule 

    
!TU: Atom type identifiers
        
    integer, parameter :: ATM_TYPE_CORE  = 1   ! "core"  / "c"
    integer, parameter :: ATM_TYPE_SEMI  = 2   ! "semi"  / "s"
    integer, parameter :: ATM_TYPE_METAL = 3   ! "metal" / "m"
    integer, parameter :: ATM_TYPE_SPIN  = 4   ! "spin"  / "p"

    
end module
