-
Notifications
You must be signed in to change notification settings - Fork 6
/
Copy pathga_driver.f
206 lines (153 loc) · 5.7 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
200
201
202
203
204
205
206
module GA_driver_m
use type_m
use constants_m
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 , parse_EH_cost_function
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 , Cube_Coef , Cube_Zeta
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
integer , allocatable :: MOnum(:)
real*8 :: first_cost , 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(:)
real*8 :: DP(3) = d_zero
! 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 )
If( DIPOLE_ ) &
then
CALL Util_multipoles
CALL GA_DP_Analysis( Extended_Cell, ExCell_basis, UNI%L, UNI%R, DP )
end if
! parse directives in cost_tuning.inpt ...
call parse_EH_cost_function
! calculates the cost on input, for future comparison ...
first_cost = evaluate_cost(Extended_Cell, UNI, ExCell_basis, DP )
! Optimization of Huckel parameters ...
CALL Genetic_Algorithm( ExCell_basis, OPT_basis )
! calculations with new parameters ...
CALL GA_eigen( Extended_Cell, OPT_basis, UNI , flag=2 )
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*, " "
If( DIPOLE_ ) &
then
Print 210 , evaluate_cost( Extended_Cell, UNI, OPT_basis, DP ,ShowCost=.true. ) , first_cost
Print 154, DP, sqrt( dot_product(DP,DP) )
!Print 189 , Alpha_ii , sum( Alpha_ii ) / three
else
Print 210 , evaluate_cost( Extended_Cell, UNI, OPT_basis, ShowCost=.true. ) , first_cost
end if
CALL Dump_OPT_parameters( OPT_basis , output='STDOUT' )
! Population analysis ...
If( GaussianCube ) &
then
! do it again because the STO have been optimized
do i = 1 , size(OPT_basis)
Cube_Coef(i,:) = OPT_basis(i) % coef(:)
Cube_Zeta(i,:) = OPT_basis(i) % zeta(:) * a_Bohr
end do
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 == 4) then
call get_command_argument(3,MOstr)
select case (MOstr)
case( "-" ) ! <== orbitals within a range ...
CALL GET_COMMAND_ARGUMENT(2, MOstr)
read( MOstr,*) MO_first
CALL GET_COMMAND_ARGUMENT(4, 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