Octopus
xc.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, 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 xc_oct_m
22 use comm_oct_m
23 use debug_oct_m
28 use global_oct_m
29 use grid_oct_m
30 use io_oct_m
32 use iso_c_binding
33 use, intrinsic :: iso_fortran_env
38 use math_oct_m
39 use mesh_oct_m
42 use mpi_oct_m
44 use parser_oct_m
47 use space_oct_m
51 use xc_cam_oct_m
52 use xc_f03_lib_m
53#ifdef HAVE_LIBXC_FUNCS
54 use xc_f03_funcs_m
55#endif
56 use xc_fbe_oct_m
59
60 implicit none
61
62 private
63 public :: &
64 xc_t, &
65 xc_init, &
66 xc_end, &
68 xc_get_vxc, &
70 xc_get_fxc, &
71 xc_get_kxc, &
79
80 ! A Structure that contains the quantities needed to compute the functionals
82 real(real64), pointer :: rho(:,:) ! A pointer to the full density
83
84 real(real64), allocatable :: dens(:,:) ! Density (in the local frame of the magnetization)
85 real(real64), allocatable :: gdens(:,:,:) ! Gradient of the density
86 real(real64), allocatable :: ldens(:,:) ! Laplacian of the density
87 real(real64), allocatable :: tau(:,:) ! Kinetic energy density
89
90
91 type xc_t
92 private
93 integer, public :: family
94 integer, public :: flags
95 integer, public :: kernel_family
96 type(xc_functional_t), public :: functional(2,2)
98
99 type(xc_functional_t), public :: kernel(2,2)
100 real(real64), public :: kernel_lrc_alpha
101 real(real64), public :: kernel_proca_a_zero
102 real(real64), public :: kernel_proca_a_one
103 type(xc_cam_t), public :: cam
104
105 logical :: use_gi_ked
106 integer :: xc_density_correction
107 logical :: xcd_optimize_cutoff
108 real(real64) :: xcd_ncutoff
109 logical :: xcd_minimum
110 logical :: xcd_normalize
111 logical :: parallel
112
113 type(internal_quantities_t) :: quantities
114
115 contains
116 procedure :: compute_exchange => xc_compute_exchange
117
118 end type xc_t
119
121 real(real64), public, parameter :: xc_tiny = 1.0e-12_real64
122
123 integer, parameter :: &
124 LR_NONE = 0, &
125 lr_x = 1
126
127contains
128
129 ! ---------------------------------------------------------
130 subroutine xc_write_info(xcs, iunit, namespace)
131 type(xc_t), intent(in) :: xcs
132 integer, optional, intent(in) :: iunit
133 type(namespace_t), optional, intent(in) :: namespace
134
135 integer :: ifunc
136
137 push_sub(xc_write_info)
138
139 write(message(1), '(a)') "Exchange-correlation:"
140 call messages_info(1, iunit=iunit, namespace=namespace)
141
142 do ifunc = func_x, func_c
143 call xc_functional_write_info(xcs%functional(ifunc, 1), iunit, namespace)
144 end do
145
146 if (abs(xcs%cam%alpha + xcs%cam%beta) > m_epsilon) then
147 write(message(1), '(1x)')
148 write(message(2), '(a,f8.5)') "Exact exchange mixing = ", xcs%cam%alpha
149 write(message(3), '(a,f8.5)') "Exact exchange for short-range beta = ", xcs%cam%beta
150 write(message(4), '(a,f8.5)') "Exact exchange range-separate omega = ", xcs%cam%omega
151 call messages_info(4, iunit=iunit, namespace=namespace)
152 end if
153
154
155 pop_sub(xc_write_info)
156 end subroutine xc_write_info
157
158
159 ! ---------------------------------------------------------
160 subroutine xc_init(xcs, namespace, ndim, periodic_dim, nel, x_id, c_id, xk_id, ck_id, hartree_fock, ispin)
161 type(xc_t), intent(out) :: xcs
162 type(namespace_t), intent(in) :: namespace
163 integer, intent(in) :: ndim
164 integer, intent(in) :: periodic_dim
165 real(real64), intent(in) :: nel
166 integer, intent(in) :: x_id
167 integer, intent(in) :: c_id
168 integer, intent(in) :: xk_id
169 integer, intent(in) :: ck_id
170 logical, intent(in) :: hartree_fock
171 integer, intent(in) :: ispin
172
173 integer :: isp, xc_major, xc_minor, xc_micro
174 logical :: ll
175 type(block_t) :: blk
176 type(xc_cam_t) :: cam_ext
178 push_sub(xc_init)
180 call xc_f03_version(xc_major, xc_minor, xc_micro)
182 xcs%family = 0
183 xcs%flags = 0
184 xcs%kernel_family = 0
185
186 call parse()
187
188 !we also need XC functionals that do not depend on the current
189 !get both spin-polarized and unpolarized
191 ! TODO: check that nel should not be spin polarized here
192 do isp = 1, 2
193
194 call xc_functional_init(xcs%functional(func_x, isp), namespace, x_id, ndim, nel, isp)
195 call xc_functional_init(xcs%functional(func_c, isp), namespace, c_id, ndim, nel, isp)
197 call xc_functional_init(xcs%kernel(func_x, isp), namespace, xk_id, ndim, nel, isp)
198 call xc_functional_init(xcs%kernel(func_c, isp), namespace, ck_id, ndim, nel, isp)
199
200 end do
202 xcs%family = ior(xcs%family, xcs%functional(func_x,1)%family)
203 xcs%family = ior(xcs%family, xcs%functional(func_c,1)%family)
205 xcs%flags = ior(xcs%flags, xcs%functional(func_x,1)%flags)
206 xcs%flags = ior(xcs%flags, xcs%functional(func_c,1)%flags)
207
208 xcs%kernel_family = ior(xcs%kernel_family, xcs%kernel(func_x,1)%family)
209 xcs%kernel_family = ior(xcs%kernel_family, xcs%kernel(func_c,1)%family)
210
212 if (xc_is_not_size_consistent(xcs, namespace) .and. periodic_dim > 0) then
213 message(1) = "Cannot perform a periodic calculation with a functional"
214 message(2) = "that depends on the number of electrons."
215 call messages_fatal(2, namespace=namespace)
216 end if
217
218 ! Take care of hybrid functionals (they appear in the correlation functional)
219 xcs%cam = cam_null
220
221 ll = (hartree_fock) &
222 .or.(xcs%functional(func_x,1)%id == xc_oep_x) &
223 .or. family_is_hybrid(xcs)
224 if (ll) then
225 if ((xcs%functional(func_x,1)%id /= 0).and.(xcs%functional(func_x,1)%id /= xc_oep_x)) then
226 message(1) = "You cannot use an exchange functional when performing"
227 message(2) = "a Hartree-Fock calculation or using a hybrid functional."
228 call messages_fatal(2, namespace=namespace)
229 end if
230
231 if (periodic_dim == ndim) then
232 call messages_experimental("Fock operator (Hartree-Fock, OEP, hybrids) in fully periodic systems", namespace=namespace)
233 end if
234
235 ! get the mixing coefficient for hybrids
236 if (family_is_hybrid(xcs)) then
237 if( .not. cam_ext%is_null() ) call set_hybrid_params(xcs, namespace, cam_ext)
238 call xc_f03_hyb_cam_coef(xcs%functional(func_c,1)%conf, xcs%cam%omega, &
239 xcs%cam%alpha, xcs%cam%beta)
240 call xc_f03_hyb_cam_coef(xcs%functional(func_c,2)%conf, xcs%cam%omega, &
241 xcs%cam%alpha, xcs%cam%beta)
242 else
243 ! we are doing Hartree-Fock plus possibly a correlation functional
244 xcs%cam = cam_exact_exchange
245 end if
246
247 ! reset certain variables
248 xcs%functional(func_x,1)%family = xc_family_oep
249 xcs%functional(func_x,1)%id = xc_oep_x
250 xcs%functional(func_x,2)%family = xc_family_oep
251 xcs%functional(func_x,2)%id = xc_oep_x
252 if (.not. hartree_fock) then
253 xcs%family = ior(xcs%family, xc_family_oep)
254 end if
255 end if
256
257 if (in_family(xcs%family, [xc_family_lca])) then
258 call messages_not_implemented("LCA current functionals", namespace) ! not even in libxc!
259 end if
260
261 call messages_obsolete_variable(namespace, 'MGGAimplementation')
262 call messages_obsolete_variable(namespace, 'CurrentInTau', 'XCUseGaugeIndependentKED')
263
264 if (xcs%functional(func_x, 1)%id == xc_mgga_x_tb09 .and. periodic_dim /= 3) then
265 message(1) = "mgga_x_tb09 functional can only be used for 3D periodic systems"
266 call messages_fatal(1, namespace=namespace)
267 end if
268
269 if (family_is_mgga(xcs%family) .or. family_is_nc_mgga(xcs%family)) then
270 !%Variable XCUseGaugeIndependentKED
271 !%Type logical
272 !%Default yes
273 !%Section Hamiltonian::XC
274 !%Description
275 !% If true, when evaluating the XC functional, a term including the (paramagnetic or total) current
276 !% is added to the kinetic-energy density such as to make it gauge-independent.
277 !% Applies only to meta-GGA (and hybrid meta-GGA) functionals.
278 !%End
279 call parse_variable(namespace, 'XCUseGaugeIndependentKED', .true., xcs%use_gi_ked)
280 end if
281
282 pop_sub(xc_init)
283
284 contains
285
286 subroutine parse()
287
288 push_sub(xc_init.parse)
289
290 ! the values of x_id, c_id, xk_id, and c_id are read outside the routine
291
292 !%Variable XCKernelLRCAlpha
293 !%Type float
294 !%Default 0.0
295 !%Section Hamiltonian::XC
296 !%Description
297 !% Set to a non-zero value to add a long-range correction for solids to the kernel.
298 !% This is the <math>\alpha</math> parameter defined in S. Botti <i>et al.</i>, <i>Phys. Rev. B</i>
299 !% 69, 155112 (2004). The <math>\Gamma = \Gamma` = 0</math> term <math>-\alpha/q^2</math> is taken
300 !% into account by introducing an additional pole to the polarizability (see R. Stubner
301 !% <i>et al.</i>, <i>Phys. Rev. B</i> 70, 245119 (2004)). The rest of the terms are included by
302 !% multiplying the Hartree term by <math>1 - \alpha / 4 \pi</math>. The use of non-zero
303 !% <math>\alpha</math> in combination with <tt>HamiltonianVariation</tt> = <tt>V_ext_only</tt>
304 !% corresponds to account of only the <math>\Gamma = \Gamma` = 0</math> term.
305 !% Applicable only to isotropic systems. (Experimental)
306 !%End
307
308 call parse_variable(namespace, 'XCKernelLRCAlpha', m_zero, xcs%kernel_lrc_alpha)
309 if (abs(xcs%kernel_lrc_alpha) > m_epsilon) then
310 call messages_experimental("Long-range correction to kernel", namespace=namespace)
311 end if
312
313 !%Variable XCProcaResonantTerm
314 !%Type float
315 !%Default 0.0
316 !%Section Hamiltonian :: XC
317 !%Description
318 !% This variable is the a_0 parameter in the Proca equation
319 !% Set to a nonzero value to give to the Maxwell equation a linear term correction as suggested in
320 !% J. K. Dewhurst <i>et al.</i>, <i>Phys. Rev. B</i> 111, L060302 (2025).
321 !%End
322 call parse_variable(namespace, 'XCProcaResonantTerm', m_zero, xcs%kernel_proca_a_zero)
323 if (abs(xcs%kernel_proca_a_zero) > m_epsilon) then
324 call messages_experimental("Proca equation a0", namespace=namespace)
325 end if
326
327 !%Variable XCProcaDamping
328 !%Type float
329 !%Default 0.0
330 !%Section Hamiltonian :: XC
331 !%Description
332 !% This variable is the a_1 parameter in the Proca equation
333 !% Set to a nonzero value to give to the Maxwell equation a first-order derivative term correction as suggested in
334 !% J. K. Dewhurst <i>et al.</i>, <i>Phys. Rev. B</i> 111, L060302 (2025).
335 !%End
336 call parse_variable(namespace, 'XCProcaDamping', m_zero, xcs%kernel_proca_a_one)
337 if (abs(xcs%kernel_proca_a_one) > m_epsilon) then
338 call messages_experimental("Proca equation a1", namespace=namespace)
339 end if
340
341 !%Variable XCDensityCorrection
342 !%Type integer
343 !%Default none
344 !%Section Hamiltonian::XC::DensityCorrection
345 !%Description
346 !% This variable controls the long-range correction of the XC
347 !% potential using the <a href=http://arxiv.org/abs/1107.4339>XC density representation</a>.
348 !%Option none 0
349 !% No correction is applied.
350 !%Option long_range_x 1
351 !% The correction is applied to the exchange potential.
352 !%End
353 call parse_variable(namespace, 'XCDensityCorrection', lr_none, xcs%xc_density_correction)
354
355 if (xcs%xc_density_correction /= lr_none) then
356 call messages_experimental('XC density correction', namespace=namespace)
357
358 if(ispin /= unpolarized) then
359 call messages_not_implemented('XCDensityCorrection with SpinComponents /= unpolarized', namespace)
360 end if
361
362 !%Variable XCDensityCorrectionOptimize
363 !%Type logical
364 !%Default true
365 !%Section Hamiltonian::XC::DensityCorrection
366 !%Description
367 !% When enabled, the density cutoff will be
368 !% optimized to replicate the boundary conditions of the exact
369 !% XC potential. If the variable is set to no, the value of
370 !% the cutoff must be given by the <tt>XCDensityCorrectionCutoff</tt>
371 !% variable.
372 !%End
373 call parse_variable(namespace, 'XCDensityCorrectionOptimize', .true., xcs%xcd_optimize_cutoff)
374
375 !%Variable XCDensityCorrectionCutoff
376 !%Type float
377 !%Default 0.0
378 !%Section Hamiltonian::XC::DensityCorrection
379 !%Description
380 !% The value of the cutoff applied to the XC density.
381 !%End
382 call parse_variable(namespace, 'XCDensityCorrectionCutoff', m_zero, xcs%xcd_ncutoff)
383
384 !%Variable XCDensityCorrectionMinimum
385 !%Type logical
386 !%Default true
387 !%Section Hamiltonian::XC::DensityCorrection
388 !%Description
389 !% When enabled, the cutoff optimization will
390 !% return the first minimum of the <math>q_{xc}</math> function if it does
391 !% not find a value of -1 (<a href=http://arxiv.org/abs/1107.4339>details</a>).
392 !% This is required for atoms or small
393 !% molecules, but may cause numerical problems.
394 !%End
395 call parse_variable(namespace, 'XCDensityCorrectionMinimum', .true., xcs%xcd_minimum)
396
397 !%Variable XCDensityCorrectionNormalize
398 !%Type logical
399 !%Default true
400 !%Section Hamiltonian::XC::DensityCorrection
401 !%Description
402 !% When enabled, the correction will be
403 !% normalized to reproduce the exact boundary conditions of
404 !% the XC potential.
405 !%End
406 call parse_variable(namespace, 'XCDensityCorrectionNormalize', .true., xcs%xcd_normalize)
407
408 end if
409
410 !%Variable ParallelXC
411 !%Type logical
412 !%Default true
413 !%Section Execution::Parallelization
414 !%Description
415 !% When enabled, additional parallelization
416 !% will be used for the calculation of the XC functional.
417 !%End
418 call messages_obsolete_variable(namespace, 'XCParallel', 'ParallelXC')
419 call parse_variable(namespace, 'ParallelXC', .true., xcs%parallel)
420
421
422 !%Variable HybridCAMParameters
423 !%Type block
424 !%Section Hamiltonian::XC
425 !%Description
426 !% This variable specifies the <math>\alpha, \beta, \omega<math> for CAM-type
427 !% hybrid functionals. Defaults are zero.
428 !%End
429 cam_ext = cam_null
430
431 if(parse_block(namespace, 'HybridCamParameters', blk) == 0) then
432 call parse_block_float(blk, 0, 0, cam_ext%alpha)
433 call parse_block_float(blk, 0, 1, cam_ext%beta)
434 call parse_block_float(blk, 0, 2, cam_ext%omega)
435 call parse_block_end(blk)
436 end if
437
438 if(.not. cam_ext%is_null()) then
439 call cam_ext%print(namespace, msg="Info: Setting external CAM parameters")
440 endif
441
442 pop_sub(xc_init.parse)
443 end subroutine parse
444
445 end subroutine xc_init
446
447
448 ! ---------------------------------------------------------
449 subroutine xc_end(xcs)
450 type(xc_t), intent(inout) :: xcs
451
452 integer :: isp
453
454 push_sub(xc_end)
455
456 do isp = 1, 2
457 call xc_functional_end(xcs%functional(func_x, isp))
458 call xc_functional_end(xcs%functional(func_c, isp))
459 call xc_functional_end(xcs%kernel(func_x, isp))
460 call xc_functional_end(xcs%kernel(func_c, isp))
461 end do
462 xcs%family = 0
463 xcs%flags = 0
464
465 pop_sub(xc_end)
466 end subroutine xc_end
467
468 ! ---------------------------------------------------------
473 logical pure function xc_is_orbital_dependent(xcs)
474 type(xc_t), intent(in) :: xcs
475
476 xc_is_orbital_dependent = family_is_hybrid(xcs) .or. &
477 in_family(xcs%functional(func_x,1)%family, [xc_family_oep]) .or. &
478 in_family(xcs%family, [xc_family_mgga, xc_family_nc_mgga])
479
480 end function xc_is_orbital_dependent
481
482 ! ---------------------------------------------------------
484 pure logical function family_is_gga(family, only_collinear)
485 integer, intent(in) :: family
486 logical, optional, intent(in) :: only_collinear
487
488 if(optional_default(only_collinear, .false.)) then
489 family_is_gga = in_family(family, [xc_family_gga, xc_family_hyb_gga, &
490 xc_family_mgga, xc_family_hyb_mgga, xc_family_libvdwxc])
491 else
492 family_is_gga = in_family(family, [xc_family_gga, xc_family_hyb_gga, &
493 xc_family_mgga, xc_family_hyb_mgga, xc_family_libvdwxc, xc_family_nc_mgga])
494 end if
495 end function family_is_gga
496
497 !----------------------------------------------------------------------
502 pure logical function family_is_supported(family)
503 integer, intent(in) :: family
504
505 family_is_supported = in_family(family, [xc_family_lda, xc_family_hyb_lda, xc_family_gga, xc_family_hyb_gga, &
506 xc_family_mgga, xc_family_hyb_mgga, xc_family_libvdwxc])
507 end function family_is_supported
508
509 ! ---------------------------------------------------------
511 pure logical function family_is_mgga(family, only_collinear)
512 integer, optional, intent(in) :: family
513 logical, optional, intent(in) :: only_collinear
514
515 if(optional_default(only_collinear, .false.)) then
516 family_is_mgga = in_family(family, [xc_family_mgga, xc_family_hyb_mgga])
517 else
518 family_is_mgga = in_family(family, [xc_family_mgga, xc_family_hyb_mgga, xc_family_nc_mgga])
519 end if
520 end function family_is_mgga
521
522 ! ---------------------------------------------------------
524 logical pure function family_is_mgga_with_exc(xcs)
525 type(xc_t), intent(in) :: xcs
526
527 integer :: ixc
528
530 do ixc = 1, 2
531 if (in_family(xcs%functional(ixc, 1)%family, [xc_family_mgga, xc_family_hyb_mgga, xc_family_nc_mgga]) &
532 .and. xc_functional_is_energy_functional(xcs%functional(ixc, 1))) then
534 end if
535 end do
536 end function family_is_mgga_with_exc
537
539 logical pure function family_is_hybrid(xcs)
540 type(xc_t), intent(in) :: xcs
541
542 integer :: ixc
543
545 do ixc = 1, 2
546 if (in_family(xcs%functional(ixc, 1)%family, [xc_family_hyb_lda, xc_family_hyb_gga, xc_family_hyb_mgga])) then
548 end if
549 end do
550 end function family_is_hybrid
551
552 pure logical function in_family(family, xc_families)
553 integer, intent(in) :: family
554 integer, intent(in) :: xc_families(:)
555
556 in_family = bitand(family, sum(xc_families)) /= 0
557 end function in_family
558
560 logical function xc_compute_exchange(xc, theory_level)
561 class(xc_t), intent(in) :: xc
562 integer, intent(in) :: theory_level
563
564 integer, parameter :: exchange_theory_level(3) = [hartree, hartree_fock, rdmft]
565
566 push_sub(xc_compute_exchange)
567
568 xc_compute_exchange = any(exchange_theory_level == theory_level) .or. &
569 (theory_level == generalized_kohn_sham_dft .and. family_is_hybrid(xc)) .or. &
570 (xc%functional(func_x, 1)%id == xc_oep_x_slater) .or. &
571 (bitand(xc%family, xc_family_oep) /= 0)
572
573 pop_sub(xc_compute_exchange)
574
575 end function xc_compute_exchange
576
577 ! ---------------------------------------------------------
579 subroutine copy_global_to_local(global, local, n_block, nspin, ip)
580 real(real64), intent(in) :: global(:,:)
581 real(real64), intent(out) :: local(:,:)
582 integer, intent(in) :: n_block
583 integer, intent(in) :: nspin
584 integer, intent(in) :: ip
585
586 integer :: ib, is
587
588 push_sub(copy_global_to_local)
589
590 do is = 1, nspin
591 !$omp parallel do
592 do ib = 1, n_block
593 local(is, ib) = global(ib + ip - 1, is)
594 end do
595 end do
596
598 end subroutine copy_global_to_local
599
600 ! ---------------------------------------------------------
601 subroutine copy_local_to_global(local, global, n_block, spin_channels, ip)
602 real(real64), intent(in) :: local(:,:)
603 real(real64), intent(inout) :: global(:,:)
604 integer, intent(in) :: n_block
605 integer, intent(in) :: spin_channels
606 integer, intent(in) :: ip
607
608 integer :: ib, is
609
611
612 do is = 1, spin_channels
613 !$omp parallel do
614 do ib = 1, n_block
615 global(ib + ip - 1, is) = global(ib + ip - 1, is) + local(is, ib)
616 end do
617 end do
618
620 end subroutine copy_local_to_global
621
622 ! ---------------------------------------------------------
624 subroutine set_hybrid_params(xcs, namespace, cam_ext)
625 type(xc_t), intent(inout) :: xcs
626 type(namespace_t), intent(in) :: namespace
627 type(xc_cam_t), intent(in) :: cam_ext
628
629 real(real64), parameter :: default_alpha_pbe0 = 0.25_real64
630 real(real64) :: parameters(3)
631
632 push_sub(set_hybrid_params)
633
634 ! LibXC expects an array of reals, ordered [alpha, beta, omega]
635 parameters = cam_ext%as_array()
636 xcs%cam%alpha = parameters(1)
637
638 select case(xcs%functional(func_c, 1)%id)
639
640 case(xc_hyb_gga_xc_pbeh, xc_hyb_lda_xc_lda0) ! original PBE0/LDA0 in libxc
641 if(parameters(1) < m_zero) parameters(1) = default_alpha_pbe0
642 call xc_f03_func_set_ext_params(xcs%functional(func_c, 1)%conf, parameters)
643 call xc_f03_func_set_ext_params(xcs%functional(func_c, 2)%conf, parameters)
644 write(message(1), '(a,f6.3,a)') 'Info: Setting mixing parameter (' , parameters(1) ,').'
645 call messages_info(1)
646
647 case(xc_hyb_gga_xc_cam_pbeh, xc_hyb_lda_xc_cam_lda0)
648 xcs%cam%beta = parameters(2)
649 xcs%cam%omega = parameters(3)
650 ! check parameters
651 call xc_f03_hyb_cam_coef(xcs%functional(func_c,1)%conf, xcs%cam%omega, &
652 xcs%cam%alpha, xcs%cam%beta)
653 call xc_f03_hyb_cam_coef(xcs%functional(func_c,2)%conf, xcs%cam%omega, &
654 xcs%cam%alpha, xcs%cam%beta)
655 call xcs%cam%print(namespace, msg="Setting CAM parameters:")
656
657 case default
658 assert(.false.)
659
660 end select
661
662 pop_sub(set_hybrid_params)
663 end subroutine set_hybrid_params
664
665 ! ---------------------------------------------------------
667 logical function xc_is_not_size_consistent(xcs, namespace)
668 type(xc_t), intent(in) :: xcs
669 type(namespace_t), intent(in) :: namespace
670
671 xc_is_not_size_consistent = xc_functional_is_not_size_consistent(xcs%functional(func_x,1), namespace) &
672 .or. xc_functional_is_not_size_consistent(xcs%functional(func_c,1), namespace)
673 end function xc_is_not_size_consistent
675 ! ---------------------------------------------------------
677 logical pure function xc_is_energy_functional(xcs)
678 type(xc_t), intent(in) :: xcs
679 xc_is_energy_functional = xc_functional_is_energy_functional(xcs%functional(func_x,1)) &
680 .or. xc_functional_is_energy_functional(xcs%functional(func_c,1))
681 end function xc_is_energy_functional
682
683#include "xc_vxc_inc.F90"
684#include "xc_fxc_inc.F90"
685#include "xc_kxc_inc.F90"
686
687#include "xc_vxc_nc_inc.F90"
688
689end module xc_oct_m
690
691
692!! Local Variables:
693!! mode: f90
694!! coding: utf-8
695!! End:
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
integer, parameter, public unpolarized
Parameters...
real(real64), parameter, public m_zero
Definition: global.F90:190
real(real64), parameter, public m_epsilon
Definition: global.F90:206
This module implements the underlying real-space grid.
Definition: grid.F90:119
Definition: io.F90:116
This module is intended to contain "only mathematical" functions and procedures.
Definition: math.F90:117
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:120
subroutine, public messages_not_implemented(feature, namespace)
Definition: messages.F90:1097
character(len=512), private msg
Definition: messages.F90:167
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1029
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:416
subroutine, public messages_experimental(name, namespace)
Definition: messages.F90:1069
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:600
integer function, public parse_block(namespace, name, blk, check_varinfo_)
Definition: parser.F90:621
This module handles spin dimensions of the states and the k-point distribution.
This module defines the unit system, used for input and output.
type(xc_cam_t), parameter, public cam_null
All CAM parameters set to zero.
Definition: xc_cam.F90:152
type(xc_cam_t), parameter, public cam_exact_exchange
Use only Hartree Fock exact exchange.
Definition: xc_cam.F90:155
subroutine, public xc_functional_write_info(functl, iunit, namespace)
Write functional information.
subroutine, public xc_functional_init(functl, namespace, id, ndim, nel, spin_channels)
integer, parameter, public xc_family_nc_mgga
integer, parameter, public xc_oep_x
Exact exchange.
subroutine, public xc_functional_end(functl)
integer, parameter, public func_c
integer, parameter, public func_x
Definition: xc.F90:116
subroutine, public xc_write_info(xcs, iunit, namespace)
Definition: xc.F90:226
subroutine xc_compute_vxc(der, xcs, st, psolver, namespace, space, quantities, ispin, vxc, ex, ec, deltaxc, vtau, ex_density, ec_density, stress_xc)
Definition: xc.F90:1024
subroutine, public xc_init(xcs, namespace, ndim, periodic_dim, nel, x_id, c_id, xk_id, ck_id, hartree_fock, ispin)
Definition: xc.F90:256
subroutine, public xc_get_kxc(xcs, mesh, namespace, rho, ispin, kxc)
Definition: xc.F90:2371
pure logical function family_is_supported(family)
Is the xc family internally supported by Octopus.
Definition: xc.F90:598
pure logical function, public family_is_mgga(family, only_collinear)
Is the xc function part of the mGGA family.
Definition: xc.F90:607
subroutine set_hybrid_params(xcs, namespace, cam_ext)
Sets external parameters for some hybrid functionals.
Definition: xc.F90:720
logical pure function, public family_is_mgga_with_exc(xcs)
Is the xc function part of the mGGA family with an energy functional.
Definition: xc.F90:620
integer, parameter lr_x
Definition: xc.F90:218
subroutine, public xc_end(xcs)
Definition: xc.F90:545
logical pure function, public xc_is_energy_functional(xcs)
Is one of the x or c functional is not an energy functional.
Definition: xc.F90:773
subroutine, public xc_get_vxc(gr, xcs, st, kpoints, psolver, namespace, space, rho, ispin, rcell_volume, vxc, ex, ec, deltaxc, vtau, ex_density, ec_density, stress_xc, force_orbitalfree)
Definition: xc.F90:799
logical pure function, public family_is_hybrid(xcs)
Returns true if the functional is an hybrid functional.
Definition: xc.F90:635
logical function, public xc_is_not_size_consistent(xcs, namespace)
Is one of the x or c functional is not size consistent.
Definition: xc.F90:763
subroutine copy_local_to_global(local, global, n_block, spin_channels, ip)
Definition: xc.F90:697
subroutine, public xc_get_nc_vxc(gr, xcs, st, kpoints, space, namespace, rho, vxc, ex, ec, vtau, ex_density, ec_density)
This routines is similar to xc_get_vxc but for noncollinear functionals, which are not implemented in...
Definition: xc.F90:2557
pure logical function family_is_nc_mgga(family)
Returns true is the functional is a noncollinear functional.
Definition: xc.F90:3130
subroutine, public xc_get_fxc(xcs, gr, namespace, rho, ispin, fxc, fxc_grad, fxc_grad_spin)
Returns the exchange-correlation kernel.
Definition: xc.F90:2029
pure logical function family_is_gga(family, only_collinear)
Is the xc function part of the GGA family.
Definition: xc.F90:580
logical function xc_compute_exchange(xc, theory_level)
Theory levels and functionals for which exacy exchange is required.
Definition: xc.F90:656
pure logical function, public in_family(family, xc_families)
Definition: xc.F90:648
subroutine copy_global_to_local(global, local, n_block, nspin, ip)
make a local copy with the correct memory order for libxc
Definition: xc.F90:675
logical pure function, public xc_is_orbital_dependent(xcs)
Is the xc family orbital dependent.
Definition: xc.F90:569
int true(void)
subroutine parse()
Definition: xc.F90:382