forked from lgcrego/Dynemol
-
Notifications
You must be signed in to change notification settings - Fork 0
/
diagnostic.f
183 lines (137 loc) · 5.34 KB
/
diagnostic.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
module diagnostic_m
use type_m
use omp_lib
use constants_m
use DOS_m
use parameters_m , only : spectrum , DP_Moment , &
survival , EnvField_ , &
Alpha_Tensor , &
GaussianCube , &
HFP_Forces
use Solvated_M , only : DeAllocate_TDOS , &
DeAllocate_PDOS , &
DeAllocate_SPEC
use Structure_Builder , only : Unit_Cell , &
Extended_Cell , &
Generate_Structure , &
Basis_Builder , &
ExCell_basis
use GA_QCModel_m , only : Mulliken , GA_Eigen
use DP_main_m , only : Dipole_Matrix
use Dielectric_Potential , only : Environment_SetUp
use Oscillator_m , only : Optical_Transitions
use Psi_squared_cube_format , only : Gaussian_Cube_Format
use Data_Output , only : Dump_stuff
use Embedded_FF_Alpha , only : AlphaPolar
use HuckelForces_m , only : HuckelForces
public :: diagnostic
private
contains
!
!
!
!====================
subroutine diagnostic
!====================
implicit none
! local variables ...
integer :: i , nr , N_of_residues
integer , allocatable :: MOnum(:)
real*8 :: DP(3)
type(R_eigen) :: UNI
type(f_grid) :: TDOS , SPEC
type(f_grid) , allocatable :: PDOS(:)
! preprocessing stuff ...................................
IF ( survival ) pause " >>> quit: diagnostic driver does not carry q_dynamics calculations <<< "
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 )
!.........................................................
CALL Generate_Structure(1)
CALL Basis_Builder( Extended_Cell, ExCell_basis )
If( any([DP_Moment,Spectrum,EnvField_]) ) CALL Dipole_Matrix( Extended_Cell, ExCell_basis, UNI%L, UNI%R , DP )
If( EnvField_ ) CALL Environment_SetUp( Extended_Cell )
If( Alpha_Tensor .AND. DP_Moment ) CALL AlphaPolar( Extended_Cell, ExCell_basis )
! this Eigen is MPI free ...
CALL GA_Eigen( Extended_Cell, ExCell_basis, UNI )
! save energies of the TOTAL system
OPEN(unit=9,file='system-ergs.dat',status='unknown')
do i = 1 , size(ExCell_basis)
write(9,*) i , UNI%erg(i)
end do
CLOSE(9)
CALL Total_DOS( UNI%erg , TDOS )
do nr = 1 , N_of_residues
CALL Partial_DOS( Extended_Cell , UNI , PDOS , nr )
end do
If( Spectrum ) CALL Optical_Transitions( Extended_Cell, ExCell_basis, UNI , SPEC )
If( HFP_Forces ) CALL HuckelForces( Extended_Cell, ExCell_basis, UNI )
Print*, " "
Print*, "dE1 = ", UNI%erg(115) - UNI%erg(114) , 2.6470
Print*, "dE2 = ", UNI%erg(114) - UNI%erg(113) , 0.3040
Print*, "dE3 = ", UNI%erg(115) - UNI%erg(113) , 2.9510
Print*, "dE4 = ", UNI%erg(113) - UNI%erg(112) , 0.8950
Print*, "dE5 = ", UNI%erg(112) - UNI%erg(111) , 0.4360
Print*, "dE6 = ", UNI%erg(117) - UNI%erg(116) , 1.6000
If( GaussianCube .AND. (size(MOnum) > 0) ) then
! use this to check the orbitals separately ...
! Extended_Cell% coord = Extended_Cell% coord * TWO
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
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 diagnostic
!
!
!
!===============================================
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 diagnostic_m