forked from lgcrego/Dynemol
-
Notifications
You must be signed in to change notification settings - Fork 0
/
ga_driver.f
199 lines (147 loc) · 5.71 KB
/
ga_driver.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
module GA_driver_m
use type_m
use constants_m
use MPI_definitions_m , only : slave
use parameters_m , only : spectrum , DP_Moment , GaussianCube , Alpha_Tensor , OPT_parms
use Solvated_m , only : DeAllocate_TDOS , DeAllocate_PDOS , DeAllocate_SPEC
use GA_m , only : Genetic_Algorithm , Dump_OPT_parameters
use GA_QCModel_m , only : GA_eigen , Mulliken , GA_DP_Analysis , AlphaPolar
use cost_EH , only : evaluate_cost , REF_DP , REF_Alpha
use Semi_Empirical_Parms , only : EH_atom
use Multipole_Routines_m , only : Util_multipoles
use Structure_Builder , only : Generate_Structure , Extended_Cell , Unit_Cell , Basis_Builder , ExCell_basis
use Oscillator_m , only : Optical_Transitions
use Data_Output , only : Dump_stuff
use Psi_squared_cube_format , only : Gaussian_Cube_Format
use DOS_m
public :: GA_driver
private
contains
!
!
!
!====================
subroutine GA_driver
!====================
implicit none
! local variables ...
integer :: i , nr , N_of_residues , info = 0
integer , allocatable :: MOnum(:)
real*8 :: first_cost , DP(3) , Alpha_ii(3)
logical :: DIPOLE_
type(R_eigen) :: UNI
type(f_grid) :: TDOS , SPEC
type(f_grid) , allocatable :: PDOS(:)
type(STO_basis) , allocatable :: OPT_basis(:)
! preprocessing stuff ...................................
DIPOLE_ = ( spectrum .OR. DP_Moment )
CALL DeAllocate_TDOS( TDOS , flag="alloc" )
CALL DeAllocate_PDOS( PDOS , flag="alloc" )
CALL DeAllocate_SPEC( SPEC , flag="alloc" )
N_of_residues = size( Unit_Cell%list_of_residues )
! reading command line arguments for plotting MO cube files ...
CALL Read_Command_Lines_Arguments( MOnum )
!.........................................................
! setting up the system ...
CALL Generate_Structure(1)
CALL Basis_Builder( Extended_Cell, ExCell_basis )
! setting up constraints ...
CALL GA_eigen( Extended_Cell, ExCell_basis, UNI )
! calculates the cost on input, for future comparison ...
first_cost = evaluate_cost(Extended_Cell, UNI, ExCell_basis)
If( DIPOLE_ ) CALL Util_multipoles
! Optimization of Huckel parameters ...
CALL Genetic_Algorithm( ExCell_basis, OPT_basis )
! calculations with new parameters ...
CALL GA_eigen( Extended_Cell, OPT_basis, UNI , info )
CALL Total_DOS( UNI%erg, TDOS )
do nr = 1 , N_of_residues
CALL Partial_DOS( Extended_Cell , UNI , PDOS , nr )
end do
If( DIPOLE_ ) CALL GA_DP_Analysis( Extended_Cell, OPT_basis, UNI%L, UNI%R, DP )
If( Alpha_Tensor .AND. DP_Moment ) CALL AlphaPolar( Extended_Cell, OPT_basis , Alpha_ii )
If( spectrum ) CALL Optical_Transitions( Extended_Cell, OPT_basis, UNI , SPEC )
!----------------------------------------------
! printing zone ...
! compare costs to evalualte otimization ...
Print*, " "
Print 210 , evaluate_cost( Extended_Cell, UNI, OPT_basis, ShowCost=.true. ) , first_cost
!Print 154, DP, sqrt( dot_product(DP,DP) )
!Print 189 , Alpha_ii , sum( Alpha_ii ) / three
Print*, " "
Print*, "dE1 = ",UNI%erg(50) - UNI%erg(49) , 6.4937d0
Print*, "dE2 = ",UNI%erg(51) - UNI%erg(49) , 7.9044d0
Print*, "dE3 = ",UNI%erg(50) - UNI%erg(48) , 8.3523d0
Print*, "dE4 = ",UNI%erg(49) - UNI%erg(48) , 1.8585d0
Print*, "dE5 = ",UNI%erg(51) - UNI%erg(50) , 1.4106d0
Print*, "dE5 = ",UNI%erg(48) - UNI%erg(47) , 0.0552d0
Print*, "dE7 = ",UNI%erg(47) - UNI%erg(46) , 0.1978d0
Print*, "dE8 = ",UNI%erg(52) - UNI%erg(51) , 0.1260d0
Print*, "dE9 = ",UNI%erg(53) - UNI%erg(51) , 0.1540d0
Print*, "dE10 = ",UNI%erg(48) - UNI%erg(46) , 0.2531d0
CALL Dump_OPT_parameters( OPT_basis , output='STDOUT' )
! Population analysis ...
If( GaussianCube ) then
do i = 1 , size(MOnum)
CALL Gaussian_Cube_Format( UNI%L(MOnum(i),:) , UNI%R(:,MOnum(i)) , MOnum(i) , 0.d0 )
end do
Print 220 , MOnum(:)
end if
If( OPT_parms ) then
Print 445
Print 45 , ( EH_atom(i)% EHSymbol , EH_atom(i)% residue , i = 1,size(EH_atom) )
else
Print*, ">> OPT_parms were not used <<"
end if
!----------------------------------------------
CALL Dump_stuff( TDOS , PDOS , SPEC )
CALL DeAllocate_TDOS( TDOS , flag="dealloc" )
CALL DeAllocate_PDOS( PDOS , flag="dealloc" )
CALL DeAllocate_SPEC( SPEC , flag="dealloc" )
include 'formats.h'
end subroutine GA_driver
!
!
!
!===============================================
subroutine Read_Command_Lines_Arguments( MOnum )
!===============================================
implicit none
integer , allocatable , intent(out) :: MOnum(:)
! local variables ...
integer :: i , MO_total , MO_first , MO_last
character(6) :: MOstr
MO_total = COMMAND_ARGUMENT_COUNT()
IF( MO_total == 3) then
call get_command_argument(2,MOstr)
select case (MOstr)
case( ":" ) ! <== orbitals within a range ...
CALL GET_COMMAND_ARGUMENT(1, MOstr)
read( MOstr,*) MO_first
CALL GET_COMMAND_ARGUMENT(3, MOstr)
read( MOstr,*) MO_last
MO_total = MO_last - MO_first + 1
allocate( MOnum(MO_total) )
do i = 1 , MO_total
MOnum(i) = MO_first + (i-1)
end do
case default ! <== list of 3 orbitals ...
allocate( MOnum(MO_total) )
do i = 1 , MO_total
CALL GET_COMMAND_ARGUMENT(i, MOstr)
read( MOstr,*) MOnum(i)
end do
end select
ELSE
! arbitrary (/= 3) list of orbitals ...
allocate( MOnum(MO_total) )
do i = 1 , MO_total
CALL GET_COMMAND_ARGUMENT(i, MOstr)
read( MOstr,*) MOnum(i)
end do
end IF
end subroutine Read_Command_Lines_Arguments
!
!
!
end module GA_driver_m