Octopus
pseudo.F90
Go to the documentation of this file.
1!! Copyright (C) 2018 X. Andrade
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the Lesser GNU General Public License as published by
5!! the Free Software Foundation; either version 3, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19#include "global.h"
20
21module pseudo_oct_m
22 use, intrinsic :: iso_fortran_env
23
24 implicit none
25
26 private
27
28 public :: &
29 pseudo_t, &
32 pseudo_end, &
53 pseudo_dij, &
63
64 !the following sets of values have to match with those on base.hpp
65 integer, parameter, public :: &
66 PSEUDO_TYPE_ULTRASOFT = 30, &
70
71 integer, parameter, public :: &
72 PSEUDO_STATUS_SUCCESS = 0, &
79
80 integer, parameter, public :: &
81 PSEUDO_FORMAT_FILE_NOT_FOUND = 773, &
83 pseudo_format_upf1 = 775, &
84 pseudo_format_upf2 = 776, &
85 pseudo_format_qso = 777, &
86 pseudo_format_psml = 778, &
87 pseudo_format_psf = 779, &
88 pseudo_format_cpi = 780, &
89 pseudo_format_fhi = 781, &
90 pseudo_format_hgh = 782, &
92
93 ! we only define these values here, the specific functionals are
94 ! obtained from libxc
95 integer, parameter, public :: &
96 PSEUDO_EXCHANGE_UNKNOWN = -2, &
98
99 integer, parameter, public :: &
100 PSEUDO_CORRELATION_UNKNOWN = -2, &
102
103 type pseudo_t
104 private
105 integer(int64) :: dummy
106 end type pseudo_t
107
108 interface
109
110 ! -------------------------------------------------
111
112 integer function pseudo_detect_format(filename) bind(c)
113 use iso_c_binding
114 import :: pseudo_t
115 implicit none
117 character(kind=c_char), intent(in) :: filename(*)
118 end function pseudo_detect_format
119
120 ! -------------------------------------------------
121
122 subroutine pseudo_init(pseudo, filename, fmt, ierr)
123 use iso_c_binding
124 import :: pseudo_t
125 implicit none
126
127 type(pseudo_t), intent(out) :: pseudo
128 character(kind=c_char), intent(in) :: filename(*)
129 integer, intent(in) :: fmt
130 integer, intent(out) :: ierr
131 end subroutine pseudo_init
132
133 ! -------------------------------------------------
134
135 subroutine pseudo_end(pseudo)
136 import :: pseudo_t
137 implicit none
138
139 type(pseudo_t), intent(inout) :: pseudo
140 end subroutine pseudo_end
141
142 ! -------------------------------------------------
143
144 integer function pseudo_type(pseudo)
145 import :: pseudo_t
146 implicit none
147
148 type(pseudo_t), intent(in) :: pseudo
149 end function pseudo_type
150
151 ! -------------------------------------------------
152
153 integer function pseudo_format(pseudo)
154 import :: pseudo_t
155 implicit none
156
157 type(pseudo_t), intent(in) :: pseudo
158 end function pseudo_format
159
160 ! -------------------------------------------------
161
162 real(real64) function pseudo_valence_charge(pseudo)
163 import :: pseudo_t
164 import :: real64
165 implicit none
167 type(pseudo_t), intent(in) :: pseudo
168 end function pseudo_valence_charge
169
170 ! -------------------------------------------------
171
172 real(real64) function pseudo_mesh_spacing(pseudo)
173 import :: pseudo_t
174 import :: real64
175 implicit none
176
177 type(pseudo_t), intent(in) :: pseudo
178 end function pseudo_mesh_spacing
179
180 ! -------------------------------------------------
181
182 integer function pseudo_mesh_size(pseudo)
183 import :: pseudo_t
184 implicit none
185
186 type(pseudo_t), intent(in) :: pseudo
187 end function pseudo_mesh_size
188
189 ! -------------------------------------------------
191 real(real64) function pseudo_mass(pseudo)
192 import :: pseudo_t
193 import :: real64
194 implicit none
195
196 type(pseudo_t), intent(in) :: pseudo
197 end function pseudo_mass
199 ! -------------------------------------------------
201 integer function pseudo_lmax(pseudo)
202 import :: pseudo_t
203 implicit none
204
205 type(pseudo_t), intent(in) :: pseudo
206 end function pseudo_lmax
208 ! -------------------------------------------------
209
210 integer function pseudo_llocal(pseudo)
211 import :: pseudo_t
212 implicit none
213
214 type(pseudo_t), intent(in) :: pseudo
215 end function pseudo_llocal
216
217 ! -------------------------------------------------
218
219 integer function pseudo_nchannels(pseudo)
220 import :: pseudo_t
221 implicit none
222
223 type(pseudo_t), intent(in) :: pseudo
224 end function pseudo_nchannels
225
226 ! -------------------------------------------------
227
228 integer function pseudo_nprojectors(pseudo)
229 import :: pseudo_t
230 implicit none
231
232 type(pseudo_t), intent(in) :: pseudo
233 end function pseudo_nprojectors
234
235 ! -------------------------------------------------
236
237 integer function pseudo_nprojectors_per_l(pseudo, l)
238 import :: pseudo_t
239 implicit none
240
241 type(pseudo_t), intent(in) :: pseudo
242 integer, intent(in) :: l
243 end function pseudo_nprojectors_per_l
244
245 ! -------------------------------------------------
246
247 subroutine pseudo_grid(pseudo, grid)
248 import :: pseudo_t
249 import :: real64
250 implicit none
251
252 type(pseudo_t), intent(in) :: pseudo
253 real(real64), intent(out) :: grid
254 end subroutine pseudo_grid
255
256 ! -------------------------------------------------
258 subroutine pseudo_grid_weights(pseudo, weight)
259 import :: pseudo_t
260 import :: real64
261 implicit none
262
263 type(pseudo_t), intent(in) :: pseudo
264 real(real64), intent(out) :: weight
265 end subroutine pseudo_grid_weights
266
267 ! -------------------------------------------------
268
269 subroutine pseudo_local_potential(pseudo, local_potential)
270 import :: pseudo_t
271 import :: real64
272 implicit none
273
274 type(pseudo_t), intent(in) :: pseudo
275 real(real64), intent(out) :: local_potential
276 end subroutine pseudo_local_potential
278 ! -------------------------------------------------
279
280 subroutine pseudo_projector(pseudo, l, ic, projector)
281 import :: pseudo_t
282 import :: real64
283 implicit none
284
285 type(pseudo_t), intent(in) :: pseudo
286 integer, intent(in) :: l
287 integer, intent(in) :: ic
288 real(real64), intent(out) :: projector
289 end subroutine pseudo_projector
290
291 ! -------------------------------------------------
292
293 real(real64) function pseudo_dij(pseudo, l, ic, jc)
294 import :: pseudo_t
295 import :: real64
296 implicit none
297
298 type(pseudo_t), intent(in) :: pseudo
299 integer, intent(in) :: l
300 integer, intent(in) :: ic
301 integer, intent(in) :: jc
302 end function pseudo_dij
303
304 ! -------------------------------------------------
306 subroutine pseudo_radial_potential(pseudo, l, radial_potential)
307 import :: pseudo_t
308 import :: real64
309 implicit none
310
311 type(pseudo_t), intent(in) :: pseudo
312 integer, intent(in) :: l
313 real(real64), intent(inout) :: radial_potential
315
316 ! -------------------------------------------------
317
318 subroutine pseudo_radial_function(pseudo, l, radial_function)
319 import :: pseudo_t
320 import :: real64
321 implicit none
322
323 type(pseudo_t), intent(in) :: pseudo
324 integer, intent(in) :: l
325 real(real64), intent(out) :: radial_function
326 end subroutine pseudo_radial_function
327
328 ! -------------------------------------------------
329
330 subroutine pseudo_nlcc_density(pseudo, nlcc_density)
331 import :: pseudo_t
332 import :: real64
333 implicit none
334
335 type(pseudo_t), intent(in) :: pseudo
336 real(real64), intent(out) :: nlcc_density
337 end subroutine pseudo_nlcc_density
338
339 ! -------------------------------------------------
340
341 subroutine pseudo_density(pseudo, density)
342 import :: pseudo_t
343 import :: real64
344 implicit none
345
346 type(pseudo_t), intent(in) :: pseudo
347 real(real64), intent(out) :: density
348 end subroutine pseudo_density
349
350 ! -------------------------------------------------
351
352 integer function pseudo_nwavefunctions(pseudo)
353 import :: pseudo_t
354 implicit none
355
356 type(pseudo_t), intent(in) :: pseudo
357 end function pseudo_nwavefunctions
358
359 ! -------------------------------------------------
360
361 subroutine pseudo_wavefunction(pseudo, index, n, l, occ, wf)
362 import :: pseudo_t
363 import :: real64
364 implicit none
365
366 type(pseudo_t), intent(in) :: pseudo
367 integer, intent(in) :: index
368 integer, intent(out) :: n
369 integer, intent(out) :: l
370 real(real64), intent(out) :: occ
371 real(real64), intent(out) :: wf
372 end subroutine pseudo_wavefunction
373
374 ! -------------------------------------------------
376 integer function pseudo_exchange(pseudo)
377 import :: pseudo_t
378 implicit none
379
380 type(pseudo_t), intent(in) :: pseudo
381 end function pseudo_exchange
382
383 ! -------------------------------------------------
384
385 integer function pseudo_correlation(pseudo)
386 import :: pseudo_t
387 implicit none
389 type(pseudo_t), intent(in) :: pseudo
390 end function pseudo_correlation
391
392 ! -------------------------------------------------
393
394 integer function pseudo_projector_2j(pseudo, l, ic)
395 import :: pseudo_t
396 implicit none
397
398 type(pseudo_t), intent(in) :: pseudo
399 integer, intent(in) :: l
400 integer, intent(in) :: ic
402
403 ! -------------------------------------------------
404
405 integer function pseudo_wavefunction_2j(pseudo, ii)
406 import :: pseudo_t
407 implicit none
408
409 type(pseudo_t), intent(in) :: pseudo
410 integer, intent(in) :: ii
411 end function pseudo_wavefunction_2j
412
413 end interface
414
415contains
416
417 ! -------------------------------------------------
418
419 logical function pseudo_has_nlcc(pseudo)
420 type(pseudo_t), intent(in) :: pseudo
421
422 interface
423 integer function pseudo_has_nlcc_low(pseudo)
424 import :: pseudo_t
425 implicit none
426
427 type(pseudo_t), intent(in) :: pseudo
428 end function pseudo_has_nlcc_low
429 end interface
430
431 pseudo_has_nlcc = (pseudo_has_nlcc_low(pseudo) /= 0)
432
433 end function pseudo_has_nlcc
434
435 ! -------------------------------------------------
437 logical function pseudo_has_density(pseudo)
438 type(pseudo_t), intent(in) :: pseudo
439
440 interface
441 integer function pseudo_has_density_low(pseudo)
442 import :: pseudo_t
443 implicit none
444
445 type(pseudo_t), intent(in) :: pseudo
446 end function pseudo_has_density_low
447 end interface
448
449 pseudo_has_density = (pseudo_has_density_low(pseudo) /= 0)
450
451 end function pseudo_has_density
452
453
454 ! -------------------------------------------------
455
456 logical pure function pseudo_has_total_angular_momentum(pseudo)
457 type(pseudo_t), intent(in) :: pseudo
458
459 interface
460 integer pure function pseudo_has_total_angular_momentum_low(pseudo)
461 import :: pseudo_t
462 implicit none
463
464 type(pseudo_t), intent(in) :: pseudo
465 end function pseudo_has_total_angular_momentum_low
466 end interface
467
468 pseudo_has_total_angular_momentum = (pseudo_has_total_angular_momentum_low(pseudo) /= 0)
469
472 ! -------------------------------------------------
473
474 logical function pseudo_has_radial_function(pseudo, l)
475 type(pseudo_t), intent(in) :: pseudo
476 integer, intent(in) :: l
477
478 interface
479 integer function pseudo_has_radial_function_low(pseudo, l)
480 import :: pseudo_t
481 implicit none
482
483 type(pseudo_t), intent(in) :: pseudo
484 integer, intent(in) :: l
485 end function pseudo_has_radial_function_low
486 end interface
487
488 pseudo_has_radial_function = (pseudo_has_radial_function_low(pseudo, l) /= 0)
490 end function pseudo_has_radial_function
491
492end module pseudo_oct_m
493
494!! Local Variables:
495!! mode: f90
496!! coding: utf-8
497!! End:
integer, parameter, public pseudo_type_kleinman_bylander
Definition: pseudo.F90:160
logical function, public pseudo_has_radial_function(pseudo, l)
Definition: pseudo.F90:570
integer, parameter, public pseudo_format_qso
Definition: pseudo.F90:175
logical function, public pseudo_has_density(pseudo)
Definition: pseudo.F90:533
integer, parameter, public pseudo_format_hgh
Definition: pseudo.F90:175
integer, parameter, public pseudo_format_upf2
Definition: pseudo.F90:175
integer, parameter, public pseudo_type_paw
Definition: pseudo.F90:160
logical function, public pseudo_has_nlcc(pseudo)
Definition: pseudo.F90:515
integer, parameter, public pseudo_format_cpi
Definition: pseudo.F90:175
integer, parameter, public pseudo_status_unsupported_type
Definition: pseudo.F90:166
integer, parameter, public pseudo_status_file_not_found
Definition: pseudo.F90:166
integer, parameter, public pseudo_type_semilocal
Definition: pseudo.F90:160
integer, parameter, public pseudo_format_psp8
Definition: pseudo.F90:175
integer, parameter, public pseudo_status_unsupported_type_ultrasoft
Definition: pseudo.F90:166
logical pure function, public pseudo_has_total_angular_momentum(pseudo)
Definition: pseudo.F90:552
integer, parameter, public pseudo_format_psf
Definition: pseudo.F90:175
integer, parameter, public pseudo_status_unknown_format
Definition: pseudo.F90:166
integer, parameter, public pseudo_status_unsupported_type_paw
Definition: pseudo.F90:166
integer, parameter, public pseudo_format_fhi
Definition: pseudo.F90:175
integer, parameter, public pseudo_correlation_any
Definition: pseudo.F90:194
integer, parameter, public pseudo_format_upf1
Definition: pseudo.F90:175
integer, parameter, public pseudo_format_unknown
Definition: pseudo.F90:175
integer, parameter, public pseudo_status_format_not_supported
Definition: pseudo.F90:166
integer, parameter, public pseudo_exchange_any
Definition: pseudo.F90:190
integer, parameter, public pseudo_format_psml
Definition: pseudo.F90:175