Octopus
ks_potential.F90
Go to the documentation of this file.
1!! Copyright (C) 2024 N. Tancogne-Dejean
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
27 use accel_oct_m
28 use atom_oct_m
29 use debug_oct_m
31 use global_oct_m
34 use, intrinsic :: iso_fortran_env
36 use math_oct_m
37 use mesh_oct_m
40 use mpi_oct_m
41 use mix_oct_m
47 use space_oct_m
49 use types_oct_m
50 use unit_oct_m
53
54 implicit none
55
56 private
57 public :: &
64
66 private
67 integer :: np, nspin
68 integer :: theory_level = independent_particles
69 logical :: needs_vtau = .false.
70
71 real(real64), public, allocatable :: vhartree(:)
72 real(real64), public, allocatable :: vxc(:,:)
73 real(real64), public, allocatable :: vhxc(:,:)
74 ! !! other contributions: photon xc potential, constrained DFT, ...
75 ! !! Check v_ks_calc for more details
76 real(real64), allocatable :: vtau(:,:)
77 real(real64), allocatable :: lapl_vtau(:,:)
78
79 type(accel_mem_t) :: vtau_accel
80 type(accel_mem_t) :: lapl_vtau_accel
81
82 type(derivatives_t), pointer :: der => null()
83
84 contains
85 procedure :: init => ks_potential_init
86 procedure :: end => ks_potential_end
87 procedure :: set_vtau => ks_potential_set_vtau
88 procedure :: add_vhxc => ks_potential_add_vhxc
89 procedure :: dmult_vhxc => dks_potential_mult_vhxc
90 procedure :: zmult_vhxc => zks_potential_mult_vhxc
91 procedure :: load => ks_potential_load_vhxc
92 procedure :: dump => ks_potential_dump_vhxc
93 procedure :: init_interpolation => ks_potential_init_interpolation
94 procedure :: run_zero_iter => ks_potential_run_zero_iter
95 procedure :: interpolation_new => ks_potential_interpolation_new
96 procedure :: get_interpolated_potentials => ks_potential_get_interpolated_potentials
97 procedure :: set_interpolated_potentials => ks_potential_set_interpolated_potentials
98 procedure :: interpolate_potentials => ks_potential_interpolate_potentials
99 procedure :: dapply_vtau_psi => dks_potential_apply_vtau_psi
100 procedure :: zapply_vtau_psi => zks_potential_apply_vtau_psi
101 procedure :: dapply_lapl_vtau_psi => dks_potential_apply_lapl_vtau_psi
102 procedure :: zapply_lapl_vtau_psi => zks_potential_apply_lapl_vtau_psi
103 procedure :: dcurrent_mass_renormalization => dks_potential_current_mass_renormalization
104 procedure :: zcurrent_mass_renormalization => zks_potential_current_mass_renormalization
105 procedure :: output_potentials => ks_potential_output_potentials
106 procedure :: store_potentials => ks_potential_store_copy
107 procedure :: restore_potentials => ks_potential_restore_copy
108 procedure :: check_convergence => ks_potential_check_convergence
109 procedure :: perform_interpolation => ks_potential_perform_interpolation
110 procedure :: mix_potentials => ks_potential_mix_potentials
111 end type ks_potential_t
112
113 ! We sometimes need to store a copy of the potentials, say for interpolation
114 ! This is used to store the corresponding data, still keeping them private
116 private
117 real(real64), public, allocatable :: vhxc(:,:)
118 real(real64), allocatable :: vtau(:,:)
119 real(real64), allocatable :: lapl_vtau(:,:)
120 contains
121 procedure :: copy_vhxc_to_buffer => xc_copied_potentials_copy_vhxc_to_buffer
124
125contains
126
128 subroutine ks_potential_init(this, der, np, np_part, nspin, theory_level, needs_vtau)
129 class(ks_potential_t), intent(inout) :: this
130 type(derivatives_t), target, intent(in) :: der
131 integer, intent(in) :: np, np_part
132 integer, intent(in) :: nspin
133 integer, intent(in) :: theory_level
134 logical, intent(in) :: needs_vtau
135
136 push_sub(ks_potential_init)
137
138 this%theory_level = theory_level
139 this%needs_vtau = needs_vtau
140 this%np = np
141 this%nspin = nspin
142 this%der => der
143
144 ! In the case of spinors, vxc_11 = hm%vxc(:, 1), vxc_22 = hm%vxc(:, 2), Re(vxc_12) = hm%vxc(:. 3);
145 ! Im(vxc_12) = hm%vxc(:, 4)
146 safe_allocate(this%vhxc(1:np, 1:nspin))
147 this%vhxc(1:np, 1:nspin) = m_zero
148
149 if (theory_level /= independent_particles) then
150
151 safe_allocate(this%vhartree(1:np_part))
152 this%vhartree=m_zero
153
154 safe_allocate(this%vxc(1:np, 1:nspin))
155 this%vxc=m_zero
156
157 if (needs_vtau) then
158 safe_allocate(this%vtau(1:np_part, 1:nspin))
159 this%vtau = m_zero
160 safe_allocate(this%lapl_vtau(1:np, 1:nspin))
161 this%lapl_vtau = m_zero
165 call accel_create_buffer(this%lapl_vtau_accel, accel_mem_read_only, type_float, &
167 end if
169 end if
170 end if
173 end subroutine ks_potential_init
176 subroutine ks_potential_end(this)
177 class(ks_potential_t), intent(inout) :: this
178
179 push_sub(ks_potential_end)
181 safe_deallocate_a(this%vhxc)
182 safe_deallocate_a(this%vhartree)
183 safe_deallocate_a(this%vxc)
184 safe_deallocate_a(this%vtau)
185 safe_deallocate_a(this%lapl_vtau)
187 if (this%needs_vtau .and. accel_is_enabled()) then
188 call accel_release_buffer(this%vtau_accel)
189 call accel_release_buffer(this%lapl_vtau_accel)
190 end if
193 end subroutine ks_potential_end
196 subroutine ks_potential_set_vtau(this, vtau)
197 class(ks_potential_t), intent(inout) :: this
198 real(real64), contiguous, intent(in) :: vtau(:,:)
202 assert(this%needs_vtau)
204 call lalg_copy(this%np, this%nspin, vtau, this%vtau)
207
208 pop_sub(ks_potential_set_vtau)
209 end subroutine ks_potential_set_vtau
213 class(ks_potential_t), intent(inout) :: this
215 integer :: offset, ispin
218
219 assert(this%needs_vtau)
220
221 do ispin = 1, this%nspin
222 call dderivatives_lapl(this%der, this%vtau(:, ispin), this%lapl_vtau(:, ispin))
223 end do
224
225 if (accel_is_enabled()) then
226 offset = 0
227 do ispin = 1, this%nspin
228 call accel_write_buffer(this%vtau_accel, this%np, this%vtau(:, ispin), offset = offset)
229 call accel_write_buffer(this%lapl_vtau_accel, this%np, this%lapl_vtau(:, ispin), offset = offset)
230 offset = offset + accel_padded_size(this%np)
231 end do
232 end if
233
236
238 subroutine ks_potential_add_vhxc(this, pot, nspin)
239 class(ks_potential_t), intent(in) :: this
240 real(real64), contiguous, intent(inout) :: pot(:,:)
241 integer, optional, intent(in) :: nspin
242
243 integer :: nspin_
244
245 push_sub(ks_potential_add_vhxc)
246
247 assert(not_in_openmp())
248 nspin_ = optional_default(nspin, this%nspin)
249 assert(size(pot, dim=1) >= this%np)
250 assert(size(pot, dim=2) >= nspin_)
251
252 !TODO: We could skip this for IPA
253
254 call lalg_axpy(this%np, nspin_, m_one, this%vhxc, pot)
255
256 pop_sub(ks_potential_add_vhxc)
257 end subroutine ks_potential_add_vhxc
258
259 ! -----------------------------------------------------------------
261 subroutine ks_potential_dump_vhxc(this, restart, mesh, ierr)
262 class(ks_potential_t), intent(in) :: this
263 type(restart_t), intent(in) :: restart
264 class(mesh_t), intent(in) :: mesh
265 integer, intent(out) :: ierr
266
267 integer :: iunit, err, err2(2), isp
268 character(len=MAX_PATH_LEN) :: filename
269 character(len=100) :: lines(2)
270
272
273 ierr = 0
274
275 if (restart%skip() .or. this%theory_level == independent_particles) then
277 return
278 end if
279
280 if (debug%info) then
281 message(1) = "Debug: Writing Vhxc restart."
282 call messages_info(1)
283 end if
284
285 !write the different components of the Hartree+XC potential
286 iunit = restart%open('vhxc')
287 lines(1) = '# #spin #nspin filename'
288 lines(2) = '%vhxc'
289 call restart%write(iunit, lines, 2, err)
290 if (err /= 0) ierr = ierr + 1
292 err2 = 0
293 do isp = 1, this%nspin
294 filename = get_filename_with_spin('vhxc', this%nspin, isp)
295 write(lines(1), '(i8,a,i8,a)') isp, ' | ', this%nspin, ' | "'//trim(adjustl(filename))//'"'
296 call restart%write(iunit, lines, 1, err)
297 if (err /= 0) err2(1) = err2(1) + 1
298
299 call restart%write_mesh_function(filename, mesh, this%vhxc(:,isp), err)
300 if (err /= 0) err2(2) = err2(2) + 1
301
302 end do
303 if (err2(1) /= 0) ierr = ierr + 2
304 if (err2(2) /= 0) ierr = ierr + 4
305
306 lines(1) = '%'
307 call restart%write(iunit, lines, 1, err)
308 if (err /= 0) ierr = ierr + 4
309
310 ! MGGAs and hybrid MGGAs have an extra term that also needs to be dumped
311 if (this%needs_vtau) then
312 lines(1) = '# #spin #nspin filename'
313 lines(2) = '%vtau'
314 call restart%write(iunit, lines, 2, err)
315 if (err /= 0) ierr = ierr + 8
316
317 err2 = 0
318 do isp = 1, this%nspin
319 filename = get_filename_with_spin('vtau', this%nspin, isp)
320 write(lines(1), '(i8,a,i8,a)') isp, ' | ', this%nspin, ' | "'//trim(adjustl(filename))//'"'
321 call restart%write(iunit, lines, 1, err)
322 if (err /= 0) err2(1) = err2(1) + 16
323
324 call restart%write_mesh_function(filename, mesh, this%vtau(:,isp), err)
325 if (err /= 0) err2(1) = err2(1) + 1
326
327 end do
328 if (err2(1) /= 0) ierr = ierr + 32
329 if (err2(2) /= 0) ierr = ierr + 64
330
331 lines(1) = '%'
332 call restart%write(iunit, lines, 1, err)
333 if (err /= 0) ierr = ierr + 128
334 end if
335
336 call restart%close(iunit)
337
338 if (debug%info) then
339 message(1) = "Debug: Writing Vhxc restart done."
340 call messages_info(1)
341 end if
342
344 end subroutine ks_potential_dump_vhxc
345
346
347 ! ---------------------------------------------------------
349 subroutine ks_potential_load_vhxc(this, restart, mesh, ierr)
350 class(ks_potential_t), intent(inout) :: this
351 type(restart_t), intent(in) :: restart
352 class(mesh_t), intent(in) :: mesh
353 integer, intent(out) :: ierr
354
355 integer :: err, err2, isp
356 character(len=MAX_PATH_LEN) :: filename
357
358 push_sub(ks_potential_load_vhxc)
359
360 ierr = 0
361
362 if (restart%skip() .or. this%theory_level == independent_particles) then
363 ierr = -1
365 return
366 end if
367
368 if (debug%info) then
369 message(1) = "Debug: Reading Vhxc restart."
370 call messages_info(1)
371 end if
372
373 err2 = 0
374 do isp = 1, this%nspin
375 filename = get_filename_with_spin('vhxc', this%nspin, isp)
376
377 call restart%read_mesh_function(filename, mesh, this%vhxc(:,isp), err)
378 if (err /= 0) err2 = err2 + 1
379
380 end do
381 if (err2 /= 0) ierr = ierr + 1
382
383 ! MGGAs and hybrid MGGAs have an extra term that also needs to be read
384 err2 = 0
385 if (this%needs_vtau) then
386 do isp = 1, this%nspin
387 filename = get_filename_with_spin('vtau', this%nspin, isp)
388
389 call restart%read_mesh_function(filename, mesh, this%vtau(:,isp), err)
390 if (err /= 0) err2 = err2 + 1
391
392 end do
393
394 if (err2 /= 0) ierr = ierr + 2
395 end if
396
397 if (debug%info) then
398 message(1) = "Debug: Reading Vhxc restart done."
399 call messages_info(1)
400 end if
401
403 end subroutine ks_potential_load_vhxc
404
406 subroutine ks_potential_init_interpolation(this, vksold, order)
407 class(ks_potential_t), intent(in) :: this
408 type(potential_interpolation_t), intent(inout) :: vksold
409 integer, optional, intent(in) :: order
410
412
413 call potential_interpolation_init(vksold, this%np, this%nspin, this%needs_vtau, order=order)
414
417
419 subroutine ks_potential_run_zero_iter(this, vksold)
420 class(ks_potential_t), intent(in) :: this
421 type(potential_interpolation_t), intent(inout) :: vksold
422
424
425 call potential_interpolation_run_zero_iter(vksold, this%np, this%nspin, this%vhxc, vtau=this%vtau)
426
428 end subroutine ks_potential_run_zero_iter
429
431 subroutine ks_potential_interpolation_new(this, vksold, current_time, dt)
432 class(ks_potential_t), intent(in) :: this
433 type(potential_interpolation_t), intent(inout) :: vksold
434 real(real64), intent(in) :: current_time
435 real(real64), intent(in) :: dt
436
438
439 call potential_interpolation_new(vksold, this%np, this%nspin, current_time, dt, this%vhxc, vtau=this%vtau)
440
442 end subroutine ks_potential_interpolation_new
443
447 subroutine ks_potential_get_interpolated_potentials(this, vksold, history, storage)
448 class(ks_potential_t), intent(inout) :: this
449 type(potential_interpolation_t), intent(inout) :: vksold
450 integer, intent(in) :: history
451 type(xc_copied_potentials_t), optional, intent(inout) :: storage
452
453 if (this%theory_level == independent_particles) return
454
456
457 if (.not. present(storage)) then
458 if (this%needs_vtau) then
459 call potential_interpolation_get(vksold, this%np, this%nspin, history, this%vhxc, vtau = this%vtau)
461 else
462 call potential_interpolation_get(vksold, this%np, this%nspin, history, this%vhxc)
463 end if
464 else
465 call ks_potential_storage_allocate(this, storage)
466 if (this%needs_vtau) then
467 call potential_interpolation_get(vksold, this%np, this%nspin, history, storage%vhxc, vtau = storage%vtau)
468 else
469 call potential_interpolation_get(vksold, this%np, this%nspin, history, storage%vhxc)
470 end if
471 end if
472
475
477 subroutine ks_potential_set_interpolated_potentials(this, vksold, history)
478 class(ks_potential_t), intent(in) :: this
479 type(potential_interpolation_t), intent(inout) :: vksold
480 integer, intent(in) :: history
481
482 if (this%theory_level == independent_particles) return
483
485
486 if (this%needs_vtau) then
487 call potential_interpolation_set(vksold, this%np, this%nspin, history, this%vhxc, vtau = this%vtau)
488 else
489 call potential_interpolation_set(vksold, this%np, this%nspin, history, this%vhxc)
490 end if
491
494
495
497 subroutine ks_potential_interpolate_potentials(this, vksold, order, current_time, dt, interpolation_time)
498 class(ks_potential_t), intent(inout) :: this
499 type(potential_interpolation_t), intent(inout) :: vksold
500 integer, intent(in) :: order
501 real(real64), intent(in) :: current_time
502 real(real64), intent(in) :: dt
503 real(real64), intent(in) :: interpolation_time
504
505 if (this%theory_level == independent_particles) return
506
508
509 if (this%needs_vtau) then
510 call potential_interpolation_interpolate(vksold, order, current_time, dt, interpolation_time, &
511 this%vhxc, vtau = this%vtau)
513 else
514 call potential_interpolation_interpolate(vksold, order, current_time, dt, interpolation_time, &
515 this%vhxc)
516 end if
517
520
521 ! ---------------------------------------------------------
522 subroutine vtau_set_vout(field, this)
523 type(mixfield_t), intent(inout) :: field
524 type(ks_potential_t), intent(in) :: this
525
527
528 call mixfield_set_vout(field, this%vtau)
529
530 pop_sub(vtau_set_vout)
531 end subroutine vtau_set_vout
532
533 ! ---------------------------------------------------------
534 subroutine vtau_set_vin(field, this)
535 type(mixfield_t), intent(inout) :: field
536 type(ks_potential_t), intent(in) :: this
537
538 push_sub(vtau_set_vin)
539
540 call mixfield_set_vin(field, this%vtau)
541
543 end subroutine vtau_set_vin
544
545 ! ---------------------------------------------------------
546 subroutine vtau_get_vnew(field, this)
547 type(mixfield_t), intent(in) :: field
548 type(ks_potential_t), intent(inout) :: this
549
550 push_sub(vtau_get_vnew)
551
552 call mixfield_get_vnew(field, this%vtau)
553
555
556 pop_sub(vtau_get_vnew)
557 end subroutine vtau_get_vnew
558
560 subroutine ks_potential_output_potentials(this, namespace, how, dir, space, mesh, pos, atoms, grp)
561 class(ks_potential_t), intent(in) :: this
562 type(namespace_t), intent(in) :: namespace
563 integer(int64), intent(in) :: how
564 character(len=*), intent(in) :: dir
565 class(space_t), intent(in) :: space
566 class(mesh_t), intent(in) :: mesh
567 real(real64), optional, intent(in) :: pos(:,:)
568 type(atom_t), optional, intent(in) :: atoms(:)
569 type(mpi_grp_t), optional, intent(in) :: grp
570
571 character(len=MAX_PATH_LEN) :: fname
572 integer :: is, err
573
575
576 call dio_function_output(how, dir, 'vh', namespace, &
577 space, mesh, this%vhartree, units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
578
579 do is = 1, this%nspin
580 fname = get_filename_with_spin('vxc', this%nspin, is)
581 call dio_function_output(how, dir, fname, namespace, space, &
582 mesh, this%vxc(:, is), units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
583
584 if (this%needs_vtau) then
585 fname = get_filename_with_spin('vtau', this%nspin, is)
586 call dio_function_output(how, dir, fname, namespace, space, &
587 mesh, this%vtau(:, is), units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
588 end if
589 end do
590
593
595 subroutine ks_potential_storage_allocate(this, copy)
596 class(ks_potential_t), intent(in) :: this
597 type(xc_copied_potentials_t), intent(inout) :: copy
598
600
601 if (.not. allocated(copy%vhxc)) then
602 safe_allocate(copy%vhxc(1:this%np, 1:this%nspin))
603 end if
604 if (this%needs_vtau) then
605 if (.not. allocated(copy%vtau)) then
606 safe_allocate(copy%vtau(1:this%np, 1:this%nspin))
607 end if
608 end if
609
611 end subroutine ks_potential_storage_allocate
612
614 subroutine ks_potential_store_copy(this, copy)
615 class(ks_potential_t), intent(in) :: this
616 type(xc_copied_potentials_t), intent(inout) :: copy
618 if (this%theory_level == independent_particles) return
619
621
622 call ks_potential_storage_allocate(this, copy)
623 call lalg_copy(this%np, this%nspin, this%vhxc, copy%vhxc)
624 if (this%needs_vtau) then
625 call lalg_copy(this%np, this%nspin, this%vtau, copy%vtau)
626 end if
627
630
632 subroutine ks_potential_restore_copy(this, copy)
633 class(ks_potential_t), intent(inout) :: this
634 type(xc_copied_potentials_t), intent(in) :: copy
635
636 if (this%theory_level == independent_particles) return
637
639
640 assert(allocated(copy%vhxc))
641 call lalg_copy(this%np, this%nspin, copy%vhxc, this%vhxc)
642 if (this%needs_vtau) then
643 assert(allocated(copy%vtau))
644 call this%set_vtau(copy%vtau)
645 end if
646
648 end subroutine ks_potential_restore_copy
649
651 subroutine xc_copied_potentials_copy_vhxc_to_buffer(this, np, nspin, pnp, buffer)
652 class(xc_copied_potentials_t), intent(in) :: this
653 integer(int64), intent(in) :: np
654 integer, intent(in) :: nspin
655 integer(int64), intent(in) :: pnp
656 type(accel_mem_t), intent(inout) :: buffer
657
658 integer :: ispin
659
661
662 do ispin = 1, nspin
663 call accel_write_buffer(buffer, np, this%vhxc(:, ispin), offset = (ispin - 1)*pnp)
664 end do
665
668
670 subroutine xc_copied_potentials_end(this)
671 type(xc_copied_potentials_t), intent(inout) :: this
672
673 safe_deallocate_a(this%vhxc)
674 safe_deallocate_a(this%vtau)
675 end subroutine xc_copied_potentials_end
676
683 real(real64) function ks_potential_check_convergence(this, copy, mesh, rho, qtot) result(diff)
684 class(ks_potential_t), intent(in) :: this
685 type(xc_copied_potentials_t), intent(in) :: copy
686 class(mesh_t), intent(in) :: mesh
687 real(real64), intent(in) :: rho(:,:)
688 real(real64), intent(in) :: qtot
689
690 real(real64), allocatable :: diff_pot(:)
691 integer :: ip, is
692
693 push_sub(ks_potential_check_convergence)
694
695 safe_allocate(diff_pot(1:this%np))
696
697 !$omp parallel
698 !$omp do
699 do ip = 1, this%np
700 diff_pot(ip) = abs(copy%vhxc(ip, 1) - this%vhxc(ip, 1))*rho(ip, 1)
701 end do
702 do is = 2, this%nspin
703 !$omp do
704 do ip = 1, this%np
705 diff_pot(ip) = diff_pot(ip) + abs(copy%vhxc(ip, is) - this%vhxc(ip, is))*rho(ip, is)
706 end do
707 end do
708 !$omp end parallel
709 diff = dmf_integrate(mesh, diff_pot) / qtot
710
711 safe_deallocate_a(diff_pot)
712
713 pop_sub(ks_potential_check_convergence)
715
717 subroutine ks_potential_perform_interpolation(this, vksold, times, current_time)
718 class(ks_potential_t), intent(inout) :: this
719 type(potential_interpolation_t), intent(inout) :: vksold
720 real(real64), intent(in) :: times(:)
721 real(real64), intent(in) :: current_time
722
724
725 assert(size(times) == 3)
726
727 call interpolate(times, vksold%v_old(:, :, 1:3), current_time, vksold%v_old(:, :, 0))
728 if (this%needs_vtau) then
729 call interpolate(times, vksold%vtau_old(:, :, 1:3), current_time, vksold%vtau_old(:, :, 0))
730 end if
731
734
736 subroutine ks_potential_mix_potentials(this, vold, dt)
737 class(ks_potential_t), intent(in) :: this
738 type(xc_copied_potentials_t), intent(inout) :: vold
739 real(real64), intent(in) :: dt
740
741 integer :: ispin, ip
742
744
745 assert(not_in_openmp())
747 !$omp parallel private(ip)
748 do ispin = 1, this%nspin
749 !$omp do simd
750 do ip = 1, this%np
751 vold%vhxc(ip, ispin) = m_half*dt*(this%vhxc(ip, ispin) - vold%vhxc(ip, ispin))
752 end do
753
754 if (this%needs_vtau) then
755 !$omp do simd
756 do ip = 1, this%np
757 vold%vtau(ip, ispin) = m_half*dt*(this%vtau(ip, ispin) - vold%vtau(ip, ispin))
758 end do
759 end if
760 end do
761 !$omp end parallel
762
764 end subroutine ks_potential_mix_potentials
766
767#include "undef.F90"
768#include "complex.F90"
769#include "ks_potential_inc.F90"
770
771#include "undef.F90"
772#include "real.F90"
773#include "ks_potential_inc.F90"
774
775end module ks_potential_oct_m
776
777!! Local Variables:
778!! mode: f90
779!! coding: utf-8
780!! End:
constant times a vector plus a vector
Definition: lalg_basic.F90:173
Copies a vector x, to a vector y.
Definition: lalg_basic.F90:188
subroutine, public accel_release_buffer(this, async)
Definition: accel.F90:919
pure logical function, public accel_is_enabled()
Definition: accel.F90:419
integer, parameter, public accel_mem_read_only
Definition: accel.F90:196
type(debug_t), save, public debug
Definition: debug.F90:158
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
subroutine, public dderivatives_lapl(der, ff, op_ff, ghost_update, set_bc, factor)
apply the Laplacian to a mesh function
real(real64), parameter, public m_zero
Definition: global.F90:191
integer, parameter, public independent_particles
Theory level.
Definition: global.F90:237
logical pure function, public not_in_openmp()
Definition: global.F90:493
real(real64), parameter, public m_one
Definition: global.F90:192
subroutine, public dio_function_output(how, dir, fname, namespace, space, mesh, ff, unit, ierr, pos, atoms, grp, root)
Top-level IO routine for functions defined on the mesh.
A module to handle KS potential, without the external potential.
subroutine ks_potential_set_interpolated_potentials(this, vksold, history)
Set the interpolated potentials to history.
subroutine zks_potential_current_mass_renormalization(this, gpsi, space_dim, ndim, ispin)
Nonlocal contribution of vtau for the current.
subroutine dks_potential_apply_lapl_vtau_psi(this, mesh, d, ispin, psib, vpsib)
Wrapper to hamiltonian_elec_base_local_sub to hide the data of lapl_vtau.
subroutine ks_potential_add_vhxc(this, pot, nspin)
Adds vHxc to the potential.
subroutine dks_potential_current_mass_renormalization(this, gpsi, space_dim, ndim, ispin)
Nonlocal contribution of vtau for the current.
subroutine, public vtau_set_vin(field, this)
subroutine ks_potential_get_interpolated_potentials(this, vksold, history, storage)
Get the interpolated potentials from history.
subroutine zks_potential_apply_vtau_psi(this, mesh, d, ispin, psib, vpsib)
Wrapper to hamiltonian_elec_base_local_sub to hide the data of vtau.
subroutine, public vtau_set_vout(field, this)
subroutine zks_potential_mult_vhxc(this, mf, ispin)
Multiply a mesh function by vHxc.
subroutine ks_potential_init(this, der, np, np_part, nspin, theory_level, needs_vtau)
Allocate the memory for the KS potentials.
subroutine ks_potential_load_vhxc(this, restart, mesh, ierr)
Loads the vhxc potentials.
subroutine ks_potential_interpolate_potentials(this, vksold, order, current_time, dt, interpolation_time)
Interpolate potentials to a new time.
subroutine ks_potential_update_vtau_buffer(this)
Update vtau GPU buffer.
subroutine ks_potential_run_zero_iter(this, vksold)
Run zero iter for the interpolation.
subroutine ks_potential_store_copy(this, copy)
Copy the potentials to a storage object.
subroutine ks_potential_mix_potentials(this, vold, dt)
Replace vold potentials by 0.5*dt(vold + vhxc)
subroutine ks_potential_perform_interpolation(this, vksold, times, current_time)
Perform a time interpolation of the potentials.
subroutine, public vtau_get_vnew(field, this)
subroutine ks_potential_init_interpolation(this, vksold, order)
Initialize the potential interpolation.
real(real64) function ks_potential_check_convergence(this, copy, mesh, rho, qtot)
Check the convergence of a vhxc for predictor-corrector.
subroutine ks_potential_restore_copy(this, copy)
Copy the potentials from a storage object.
subroutine ks_potential_storage_allocate(this, copy)
Copy the potentials to a storage object.
subroutine ks_potential_dump_vhxc(this, restart, mesh, ierr)
Dumps the vhxc potentials.
subroutine ks_potential_interpolation_new(this, vksold, current_time, dt)
New interpolation point for the interpolation.
subroutine ks_potential_end(this)
Releases the memory for the KS potentials.
subroutine dks_potential_mult_vhxc(this, mf, ispin)
Multiply a mesh function by vHxc.
subroutine, public xc_copied_potentials_end(this)
Finalizer for the copied potentials.
subroutine dks_potential_apply_vtau_psi(this, mesh, d, ispin, psib, vpsib)
Wrapper to hamiltonian_elec_base_local_sub to hide the data of vtau.
subroutine ks_potential_set_vtau(this, vtau)
Set vtau and update the corresponding GPU buffer.
subroutine xc_copied_potentials_copy_vhxc_to_buffer(this, np, nspin, pnp, buffer)
Copy the vhxc potential to a gpu buffer.
subroutine ks_potential_output_potentials(this, namespace, how, dir, space, mesh, pos, atoms, grp)
Outputs vh, vxc, and vtau potentials.
subroutine zks_potential_apply_lapl_vtau_psi(this, mesh, d, ispin, psib, vpsib)
Wrapper to hamiltonian_elec_base_local_sub to hide the data of lapl_vtau.
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
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:594
this module contains the low-level part of the output system
Definition: output_low.F90:117
character(len=max_path_len) function, public get_filename_with_spin(output, nspin, spin_index)
Returns the filame as output, or output-spX is spin polarized.
Definition: output_low.F90:234
subroutine, public potential_interpolation_interpolate(potential_interpolation, order, time, dt, t, vhxc, vtau)
subroutine, public potential_interpolation_set(potential_interpolation, np, nspin, i, vhxc, vtau)
subroutine, public potential_interpolation_new(potential_interpolation, np, nspin, time, dt, vhxc, vtau)
subroutine, public potential_interpolation_run_zero_iter(potential_interpolation, np, nspin, vhxc, vtau)
subroutine, public potential_interpolation_get(potential_interpolation, np, nspin, i, vhxc, vtau)
subroutine, public potential_interpolation_init(potential_interpolation, np, nspin, mgga_with_exc, order)
This module handles spin dimensions of the states and the k-point distribution.
type(type_t), public type_float
Definition: types.F90:135
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:134
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
Describes mesh distribution to nodes.
Definition: mesh.F90:187
Quantities used in mixing: Input, output and new potentials, and the residuals.
Definition: mix.F90:173
This is defined even when running serial.
Definition: mpi.F90:144