Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Initialize the use of the gaussians to treat the QMMM
10 : !> coupling potential
11 : !> \par History
12 : !> 6.2004 created [tlaino]
13 : !> \author Teodoro Laino
14 : ! **************************************************************************************************
15 : MODULE qmmm_gaussian_init
16 : USE ao_util, ONLY: exp_radius
17 : USE cp_log_handling, ONLY: cp_get_default_logger,&
18 : cp_logger_type
19 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
20 : cp_print_key_unit_nr
21 : USE gaussian_gridlevels, ONLY: gaussian_gridlevel,&
22 : gridlevel_info_type
23 : USE input_constants, ONLY: do_qmmm_gauss,&
24 : do_qmmm_swave
25 : USE input_section_types, ONLY: section_vals_type,&
26 : section_vals_val_get
27 : USE kinds, ONLY: dp
28 : USE memory_utilities, ONLY: reallocate
29 : USE message_passing, ONLY: mp_para_env_type
30 : USE pw_env_types, ONLY: pw_env_get,&
31 : pw_env_type
32 : USE pw_pool_types, ONLY: pw_pool_p_type
33 : USE qmmm_gaussian_data, ONLY: max_geep_lib_gauss,&
34 : min_geep_lib_gauss
35 : USE qmmm_gaussian_input, ONLY: read_mm_potential,&
36 : set_mm_potential_erf,&
37 : set_mm_potential_swave
38 : USE qmmm_gaussian_types, ONLY: qmmm_gaussian_p_type,&
39 : qmmm_gaussian_type
40 : #include "./base/base_uses.f90"
41 :
42 : IMPLICIT NONE
43 : PRIVATE
44 :
45 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
46 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qmmm_gaussian_init'
47 : PUBLIC :: qmmm_gaussian_initialize
48 :
49 : CONTAINS
50 :
51 : ! **************************************************************************************************
52 : !> \brief Initialize the Gaussian QMMM Environment
53 : !> \param qmmm_gaussian_fns ...
54 : !> \param para_env ...
55 : !> \param pw_env ...
56 : !> \param mm_el_pot_radius ...
57 : !> \param mm_el_pot_radius_corr ...
58 : !> \param qmmm_coupl_type ...
59 : !> \param eps_mm_rspace ...
60 : !> \param maxradius ...
61 : !> \param maxchrg ...
62 : !> \param compatibility ...
63 : !> \param print_section ...
64 : !> \param qmmm_section ...
65 : !> \par History
66 : !> 06.2004 created [tlaino]
67 : !> \author Teodoro Laino
68 : ! **************************************************************************************************
69 808 : SUBROUTINE qmmm_gaussian_initialize(qmmm_gaussian_fns, para_env, pw_env, &
70 : mm_el_pot_radius, mm_el_pot_radius_corr, &
71 : qmmm_coupl_type, eps_mm_rspace, maxradius, maxchrg, compatibility, &
72 : print_section, qmmm_section)
73 : TYPE(qmmm_gaussian_p_type), DIMENSION(:), POINTER :: qmmm_gaussian_fns
74 : TYPE(mp_para_env_type), POINTER :: para_env
75 : TYPE(pw_env_type), POINTER :: pw_env
76 : REAL(KIND=dp), DIMENSION(:), POINTER :: mm_el_pot_radius, mm_el_pot_radius_corr
77 : INTEGER, INTENT(IN) :: qmmm_coupl_type
78 : REAL(KIND=dp), INTENT(IN) :: eps_mm_rspace
79 : REAL(KIND=dp), DIMENSION(:), POINTER :: maxradius
80 : REAL(KIND=dp), INTENT(IN) :: maxchrg
81 : LOGICAL, INTENT(IN) :: compatibility
82 : TYPE(section_vals_type), POINTER :: print_section, qmmm_section
83 :
84 : INTEGER :: i, ilevel, j, Ndim, num_geep_gauss, &
85 : output_unit
86 : LOGICAL :: Found, use_geep_lib
87 : REAL(KIND=dp) :: alpha, mymaxradius, Prefactor
88 404 : REAL(KIND=dp), DIMENSION(:), POINTER :: c_radius, radius
89 : TYPE(cp_logger_type), POINTER :: logger
90 : TYPE(gridlevel_info_type), POINTER :: gridlevel_info
91 404 : TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pools
92 : TYPE(qmmm_gaussian_type), POINTER :: mypgf
93 :
94 : ! Statements
95 :
96 404 : NULLIFY (mypgf, gridlevel_info, radius, c_radius, logger)
97 808 : logger => cp_get_default_logger()
98 404 : CALL section_vals_val_get(qmmm_section, "USE_GEEP_LIB", i_val=num_geep_gauss)
99 404 : IF (num_geep_gauss == 0) THEN
100 : use_geep_lib = .FALSE.
101 : ELSE
102 200 : use_geep_lib = .TRUE.
103 200 : CPASSERT(num_geep_gauss >= min_geep_lib_gauss)
104 200 : CPASSERT(num_geep_gauss <= max_geep_lib_gauss)
105 : END IF
106 662 : SELECT CASE (qmmm_coupl_type)
107 : CASE (do_qmmm_gauss, do_qmmm_swave)
108 : !
109 : ! Preprocessing...
110 : !
111 258 : ALLOCATE (radius(1))
112 258 : ALLOCATE (c_radius(1))
113 258 : Ndim = SIZE(radius)
114 4188 : Loop_on_all_values: DO I = 1, SIZE(mm_el_pot_radius)
115 3930 : Found = .FALSE.
116 5232 : Loop_on_found_values: DO J = 1, SIZE(radius) - 1
117 5232 : IF (mm_el_pot_radius(i) .EQ. radius(j)) THEN
118 : Found = .TRUE.
119 : EXIT Loop_on_found_values
120 : END IF
121 : END DO Loop_on_found_values
122 4188 : IF (.NOT. Found) THEN
123 : Ndim = SIZE(radius)
124 428 : radius(Ndim) = mm_el_pot_radius(i)
125 428 : c_radius(Ndim) = mm_el_pot_radius_corr(i)
126 428 : Ndim = Ndim + 1
127 428 : CALL REALLOCATE(radius, 1, Ndim)
128 428 : CALL REALLOCATE(c_radius, 1, Ndim)
129 : END IF
130 : END DO Loop_on_all_values
131 : !
132 258 : IF (Ndim - 1 > 0) THEN
133 256 : CALL REALLOCATE(radius, 1, Ndim - 1)
134 256 : CALL REALLOCATE(c_radius, 1, Ndim - 1)
135 2 : ELSE IF (Ndim - 1 == 0) THEN
136 2 : DEALLOCATE (radius)
137 2 : DEALLOCATE (c_radius)
138 : ELSE
139 0 : CPABORT("")
140 : END IF
141 : !
142 1200 : ALLOCATE (qmmm_gaussian_fns(Ndim - 1))
143 686 : DO I = 1, Ndim - 1
144 428 : NULLIFY (qmmm_gaussian_fns(I)%pgf)
145 428 : ALLOCATE (qmmm_gaussian_fns(I)%pgf)
146 428 : NULLIFY (qmmm_gaussian_fns(I)%pgf%Ak)
147 428 : NULLIFY (qmmm_gaussian_fns(I)%pgf%Gk)
148 428 : NULLIFY (qmmm_gaussian_fns(I)%pgf%grid_level)
149 : !
150 : ! Default Values
151 : !
152 428 : qmmm_gaussian_fns(I)%pgf%Elp_Radius = radius(I)
153 686 : qmmm_gaussian_fns(I)%pgf%Elp_Radius_corr = c_radius(I)
154 : END DO
155 258 : IF (ASSOCIATED(radius)) THEN
156 256 : DEALLOCATE (radius)
157 : END IF
158 258 : IF (ASSOCIATED(c_radius)) THEN
159 256 : DEALLOCATE (c_radius)
160 : END IF
161 : !
162 258 : IF (use_geep_lib) THEN
163 164 : IF (qmmm_coupl_type == do_qmmm_gauss) THEN
164 : CALL set_mm_potential_erf(qmmm_gaussian_fns, &
165 130 : compatibility, num_geep_gauss)
166 34 : ELSEIF (qmmm_coupl_type == do_qmmm_swave) THEN
167 : CALL set_mm_potential_swave(qmmm_gaussian_fns, &
168 34 : num_geep_gauss)
169 : END IF
170 : ELSE
171 : CALL read_mm_potential(para_env, qmmm_gaussian_fns, &
172 184 : (compatibility .AND. (qmmm_coupl_type == do_qmmm_gauss)), qmmm_section)
173 : END IF
174 : !
175 258 : CALL pw_env_get(pw_env, pw_pools=pools, gridlevel_info=gridlevel_info)
176 774 : ALLOCATE (maxradius(SIZE(pools)))
177 1306 : maxradius = 0.0_dp
178 686 : DO J = 1, SIZE(qmmm_gaussian_fns)
179 428 : mypgf => qmmm_gaussian_fns(J)%pgf
180 1284 : ALLOCATE (mypgf%grid_level(SIZE(mypgf%Ak)))
181 3672 : mypgf%grid_level = 0
182 428 : mymaxradius = 0.0_dp
183 3930 : DO I = 1, mypgf%number_of_gaussians
184 3672 : IF (mypgf%Gk(I) /= 0.0_dp) THEN
185 3244 : alpha = 1.0_dp/mypgf%Gk(I)
186 3244 : alpha = alpha*alpha
187 3244 : ilevel = gaussian_gridlevel(gridlevel_info, alpha)
188 3244 : Prefactor = mypgf%Ak(I)*maxchrg
189 3244 : mymaxradius = exp_radius(0, alpha, eps_mm_rspace, Prefactor, rlow=mymaxradius)
190 3244 : maxradius(ilevel) = MAX(maxradius(ilevel), mymaxradius)
191 3244 : mypgf%grid_level(i) = ilevel
192 : END IF
193 : END DO
194 : END DO
195 : !
196 : ! End of gaussian initialization...
197 : CASE DEFAULT
198 : output_unit = cp_print_key_unit_nr(logger, print_section, "PROGRAM_RUN_INFO", &
199 146 : extension=".qmmmLog")
200 146 : IF (output_unit > 0) WRITE (output_unit, '(A)') " QMMM Gaussian Data Not Initialized!"
201 550 : CALL cp_print_key_finished_output(output_unit, logger, print_section, "PROGRAM_RUN_INFO")
202 : END SELECT
203 404 : END SUBROUTINE qmmm_gaussian_initialize
204 :
205 : END MODULE qmmm_gaussian_init
|