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, space, mesh, ierr)
262 class(ks_potential_t), intent(in) :: this
263 type(restart_t), intent(in) :: restart
264 class(space_t), intent(in) :: space
265 class(mesh_t), intent(in) :: mesh
266 integer, intent(out) :: ierr
267
268 integer :: iunit, err, err2(2), isp
269 character(len=MAX_PATH_LEN) :: filename
270 character(len=100) :: lines(2)
272 push_sub(ks_potential_dump_vhxc)
273
274 ierr = 0
275
276 if (restart%skip() .or. this%theory_level == independent_particles) then
278 return
279 end if
280
281 if (debug%info) then
282 message(1) = "Debug: Writing Vhxc restart."
283 call messages_info(1)
284 end if
285
286 !write the different components of the Hartree+XC potential
287 iunit = restart%open('vhxc')
288 lines(1) = '# #spin #nspin filename'
289 lines(2) = '%vhxc'
290 call restart%write(iunit, lines, 2, err)
291 if (err /= 0) ierr = ierr + 1
292
293 err2 = 0
294 do isp = 1, this%nspin
295 filename = get_filename_with_spin('vhxc', this%nspin, isp)
296 write(lines(1), '(i8,a,i8,a)') isp, ' | ', this%nspin, ' | "'//trim(adjustl(filename))//'"'
297 call restart%write(iunit, lines, 1, err)
298 if (err /= 0) err2(1) = err2(1) + 1
299
300 call restart%write_mesh_function(space, filename, mesh, this%vhxc(:,isp), err)
301 if (err /= 0) err2(2) = err2(2) + 1
302
303 end do
304 if (err2(1) /= 0) ierr = ierr + 2
305 if (err2(2) /= 0) ierr = ierr + 4
306
307 lines(1) = '%'
308 call restart%write(iunit, lines, 1, err)
309 if (err /= 0) ierr = ierr + 4
310
311 ! MGGAs and hybrid MGGAs have an extra term that also needs to be dumped
312 if (this%needs_vtau) then
313 lines(1) = '# #spin #nspin filename'
314 lines(2) = '%vtau'
315 call restart%write(iunit, lines, 2, err)
316 if (err /= 0) ierr = ierr + 8
317
318 err2 = 0
319 do isp = 1, this%nspin
320 filename = get_filename_with_spin('vtau', this%nspin, isp)
321 write(lines(1), '(i8,a,i8,a)') isp, ' | ', this%nspin, ' | "'//trim(adjustl(filename))//'"'
322 call restart%write(iunit, lines, 1, err)
323 if (err /= 0) err2(1) = err2(1) + 16
324
325 call restart%write_mesh_function(space, filename, mesh, this%vtau(:,isp), err)
326 if (err /= 0) err2(1) = err2(1) + 1
327
328 end do
329 if (err2(1) /= 0) ierr = ierr + 32
330 if (err2(2) /= 0) ierr = ierr + 64
331
332 lines(1) = '%'
333 call restart%write(iunit, lines, 1, err)
334 if (err /= 0) ierr = ierr + 128
335 end if
336
337 call restart%close(iunit)
338
339 if (debug%info) then
340 message(1) = "Debug: Writing Vhxc restart done."
341 call messages_info(1)
342 end if
343
345 end subroutine ks_potential_dump_vhxc
346
347
348 ! ---------------------------------------------------------
350 subroutine ks_potential_load_vhxc(this, restart, space, mesh, ierr)
351 class(ks_potential_t), intent(inout) :: this
352 type(restart_t), intent(in) :: restart
353 class(space_t), intent(in) :: space
354 class(mesh_t), intent(in) :: mesh
355 integer, intent(out) :: ierr
357 integer :: err, err2, isp
358 character(len=MAX_PATH_LEN) :: filename
359
360 push_sub(ks_potential_load_vhxc)
361
362 ierr = 0
363
364 if (restart%skip() .or. this%theory_level == independent_particles) then
365 ierr = -1
367 return
368 end if
369
370 if (debug%info) then
371 message(1) = "Debug: Reading Vhxc restart."
372 call messages_info(1)
373 end if
374
375 err2 = 0
376 do isp = 1, this%nspin
377 filename = get_filename_with_spin('vhxc', this%nspin, isp)
378
379 call restart%read_mesh_function(space, filename, mesh, this%vhxc(:,isp), err)
380 if (err /= 0) err2 = err2 + 1
381
382 end do
383 if (err2 /= 0) ierr = ierr + 1
384
385 ! MGGAs and hybrid MGGAs have an extra term that also needs to be read
386 err2 = 0
387 if (this%needs_vtau) then
388 do isp = 1, this%nspin
389 filename = get_filename_with_spin('vtau', this%nspin, isp)
390
391 call restart%read_mesh_function(space, filename, mesh, this%vtau(:,isp), err)
392 if (err /= 0) err2 = err2 + 1
393
394 end do
395
396 if (err2 /= 0) ierr = ierr + 2
397 end if
398
399 if (debug%info) then
400 message(1) = "Debug: Reading Vhxc restart done."
401 call messages_info(1)
402 end if
403
405 end subroutine ks_potential_load_vhxc
406
408 subroutine ks_potential_init_interpolation(this, vksold, order)
409 class(ks_potential_t), intent(in) :: this
410 type(potential_interpolation_t), intent(inout) :: vksold
411 integer, optional, intent(in) :: order
412
414
415 call potential_interpolation_init(vksold, this%np, this%nspin, this%needs_vtau, order=order)
416
419
421 subroutine ks_potential_run_zero_iter(this, vksold)
422 class(ks_potential_t), intent(in) :: this
423 type(potential_interpolation_t), intent(inout) :: vksold
424
426
427 call potential_interpolation_run_zero_iter(vksold, this%np, this%nspin, this%vhxc, vtau=this%vtau)
428
430 end subroutine ks_potential_run_zero_iter
431
433 subroutine ks_potential_interpolation_new(this, vksold, current_time, dt)
434 class(ks_potential_t), intent(in) :: this
435 type(potential_interpolation_t), intent(inout) :: vksold
436 real(real64), intent(in) :: current_time
437 real(real64), intent(in) :: dt
438
440
441 call potential_interpolation_new(vksold, this%np, this%nspin, current_time, dt, this%vhxc, vtau=this%vtau)
442
444 end subroutine ks_potential_interpolation_new
449 subroutine ks_potential_get_interpolated_potentials(this, vksold, history, storage)
450 class(ks_potential_t), intent(inout) :: this
451 type(potential_interpolation_t), intent(inout) :: vksold
452 integer, intent(in) :: history
453 type(xc_copied_potentials_t), optional, intent(inout) :: storage
454
455 if (this%theory_level == independent_particles) return
456
458
459 if (.not. present(storage)) then
460 if (this%needs_vtau) then
461 call potential_interpolation_get(vksold, this%np, this%nspin, history, this%vhxc, vtau = this%vtau)
463 else
464 call potential_interpolation_get(vksold, this%np, this%nspin, history, this%vhxc)
465 end if
466 else
467 call ks_potential_storage_allocate(this, storage)
468 if (this%needs_vtau) then
469 call potential_interpolation_get(vksold, this%np, this%nspin, history, storage%vhxc, vtau = storage%vtau)
470 else
471 call potential_interpolation_get(vksold, this%np, this%nspin, history, storage%vhxc)
472 end if
473 end if
474
477
479 subroutine ks_potential_set_interpolated_potentials(this, vksold, history)
480 class(ks_potential_t), intent(in) :: this
481 type(potential_interpolation_t), intent(inout) :: vksold
482 integer, intent(in) :: history
483
484 if (this%theory_level == independent_particles) return
485
487
488 if (this%needs_vtau) then
489 call potential_interpolation_set(vksold, this%np, this%nspin, history, this%vhxc, vtau = this%vtau)
490 else
491 call potential_interpolation_set(vksold, this%np, this%nspin, history, this%vhxc)
492 end if
493
496
497
499 subroutine ks_potential_interpolate_potentials(this, vksold, order, current_time, dt, interpolation_time)
500 class(ks_potential_t), intent(inout) :: this
501 type(potential_interpolation_t), intent(inout) :: vksold
502 integer, intent(in) :: order
503 real(real64), intent(in) :: current_time
504 real(real64), intent(in) :: dt
505 real(real64), intent(in) :: interpolation_time
506
507 if (this%theory_level == independent_particles) return
508
510
511 if (this%needs_vtau) then
512 call potential_interpolation_interpolate(vksold, order, current_time, dt, interpolation_time, &
513 this%vhxc, vtau = this%vtau)
515 else
516 call potential_interpolation_interpolate(vksold, order, current_time, dt, interpolation_time, &
517 this%vhxc)
518 end if
519
522
523 ! ---------------------------------------------------------
524 subroutine vtau_set_vout(field, this)
525 type(mixfield_t), intent(inout) :: field
526 type(ks_potential_t), intent(in) :: this
527
529
530 call mixfield_set_vout(field, this%vtau)
531
532 pop_sub(vtau_set_vout)
533 end subroutine vtau_set_vout
534
535 ! ---------------------------------------------------------
536 subroutine vtau_set_vin(field, this)
537 type(mixfield_t), intent(inout) :: field
538 type(ks_potential_t), intent(in) :: this
539
540 push_sub(vtau_set_vin)
541
542 call mixfield_set_vin(field, this%vtau)
543
545 end subroutine vtau_set_vin
546
547 ! ---------------------------------------------------------
548 subroutine vtau_get_vnew(field, this)
549 type(mixfield_t), intent(in) :: field
550 type(ks_potential_t), intent(inout) :: this
551
552 push_sub(vtau_get_vnew)
553
554 call mixfield_get_vnew(field, this%vtau)
555
557
558 pop_sub(vtau_get_vnew)
559 end subroutine vtau_get_vnew
560
562 subroutine ks_potential_output_potentials(this, namespace, how, dir, space, mesh, pos, atoms, grp)
563 class(ks_potential_t), intent(in) :: this
564 type(namespace_t), intent(in) :: namespace
565 integer(int64), intent(in) :: how
566 character(len=*), intent(in) :: dir
567 class(space_t), intent(in) :: space
568 class(mesh_t), intent(in) :: mesh
569 real(real64), optional, intent(in) :: pos(:,:)
570 type(atom_t), optional, intent(in) :: atoms(:)
571 type(mpi_grp_t), optional, intent(in) :: grp
572
573 character(len=MAX_PATH_LEN) :: fname
574 integer :: is, err
575
577
578 call dio_function_output(how, dir, 'vh', namespace, &
579 space, mesh, this%vhartree, units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
580
581 do is = 1, this%nspin
582 fname = get_filename_with_spin('vxc', this%nspin, is)
583 call dio_function_output(how, dir, fname, namespace, space, &
584 mesh, this%vxc(:, is), units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
585
586 if (this%needs_vtau) then
587 fname = get_filename_with_spin('vtau', this%nspin, is)
588 call dio_function_output(how, dir, fname, namespace, space, &
589 mesh, this%vtau(:, is), units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
590 end if
591 end do
592
595
597 subroutine ks_potential_storage_allocate(this, copy)
598 class(ks_potential_t), intent(in) :: this
599 type(xc_copied_potentials_t), intent(inout) :: copy
600
602
603 if (.not. allocated(copy%vhxc)) then
604 safe_allocate(copy%vhxc(1:this%np, 1:this%nspin))
605 end if
606 if (this%needs_vtau) then
607 if (.not. allocated(copy%vtau)) then
608 safe_allocate(copy%vtau(1:this%np, 1:this%nspin))
609 end if
610 end if
611
613 end subroutine ks_potential_storage_allocate
614
616 subroutine ks_potential_store_copy(this, copy)
617 class(ks_potential_t), intent(in) :: this
618 type(xc_copied_potentials_t), intent(inout) :: copy
620 if (this%theory_level == independent_particles) return
621
623
624 call ks_potential_storage_allocate(this, copy)
625 call lalg_copy(this%np, this%nspin, this%vhxc, copy%vhxc)
626 if (this%needs_vtau) then
627 call lalg_copy(this%np, this%nspin, this%vtau, copy%vtau)
628 end if
629
632
634 subroutine ks_potential_restore_copy(this, copy)
635 class(ks_potential_t), intent(inout) :: this
636 type(xc_copied_potentials_t), intent(in) :: copy
637
638 if (this%theory_level == independent_particles) return
639
641
642 assert(allocated(copy%vhxc))
643 call lalg_copy(this%np, this%nspin, copy%vhxc, this%vhxc)
644 if (this%needs_vtau) then
645 assert(allocated(copy%vtau))
646 call this%set_vtau(copy%vtau)
647 end if
648
650 end subroutine ks_potential_restore_copy
651
653 subroutine xc_copied_potentials_copy_vhxc_to_buffer(this, np, nspin, pnp, buffer)
654 class(xc_copied_potentials_t), intent(in) :: this
655 integer(int64), intent(in) :: np
656 integer, intent(in) :: nspin
657 integer(int64), intent(in) :: pnp
658 type(accel_mem_t), intent(inout) :: buffer
659
660 integer :: ispin
661
663
664 do ispin = 1, nspin
665 call accel_write_buffer(buffer, np, this%vhxc(:, ispin), offset = (ispin - 1)*pnp)
666 end do
667
670
672 subroutine xc_copied_potentials_end(this)
673 type(xc_copied_potentials_t), intent(inout) :: this
674
675 safe_deallocate_a(this%vhxc)
676 safe_deallocate_a(this%vtau)
677 end subroutine xc_copied_potentials_end
678
685 real(real64) function ks_potential_check_convergence(this, copy, mesh, rho, qtot) result(diff)
686 class(ks_potential_t), intent(in) :: this
687 type(xc_copied_potentials_t), intent(in) :: copy
688 class(mesh_t), intent(in) :: mesh
689 real(real64), intent(in) :: rho(:,:)
690 real(real64), intent(in) :: qtot
691
692 real(real64), allocatable :: diff_pot(:)
693 integer :: ip, is
694
695 push_sub(ks_potential_check_convergence)
696
697 safe_allocate(diff_pot(1:this%np))
698
699 !$omp parallel
700 !$omp do
701 do ip = 1, this%np
702 diff_pot(ip) = abs(copy%vhxc(ip, 1) - this%vhxc(ip, 1))*rho(ip, 1)
703 end do
704 do is = 2, this%nspin
705 !$omp do
706 do ip = 1, this%np
707 diff_pot(ip) = diff_pot(ip) + abs(copy%vhxc(ip, is) - this%vhxc(ip, is))*rho(ip, is)
708 end do
709 end do
710 !$omp end parallel
711 diff = dmf_integrate(mesh, diff_pot) / qtot
712
713 safe_deallocate_a(diff_pot)
714
715 pop_sub(ks_potential_check_convergence)
717
719 subroutine ks_potential_perform_interpolation(this, vksold, times, current_time)
720 class(ks_potential_t), intent(inout) :: this
721 type(potential_interpolation_t), intent(inout) :: vksold
722 real(real64), intent(in) :: times(:)
723 real(real64), intent(in) :: current_time
724
726
727 assert(size(times) == 3)
728
729 call interpolate(times, vksold%v_old(:, :, 1:3), current_time, vksold%v_old(:, :, 0))
730 if (this%needs_vtau) then
731 call interpolate(times, vksold%vtau_old(:, :, 1:3), current_time, vksold%vtau_old(:, :, 0))
732 end if
733
736
738 subroutine ks_potential_mix_potentials(this, vold, dt)
739 class(ks_potential_t), intent(in) :: this
740 type(xc_copied_potentials_t), intent(inout) :: vold
741 real(real64), intent(in) :: dt
742
743 integer :: ispin, ip
744
746
747 assert(not_in_openmp())
749 !$omp parallel private(ip)
750 do ispin = 1, this%nspin
751 !$omp do simd
752 do ip = 1, this%np
753 vold%vhxc(ip, ispin) = m_half*dt*(this%vhxc(ip, ispin) - vold%vhxc(ip, ispin))
754 end do
755
756 if (this%needs_vtau) then
757 !$omp do simd
758 do ip = 1, this%np
759 vold%vtau(ip, ispin) = m_half*dt*(this%vtau(ip, ispin) - vold%vtau(ip, ispin))
760 end do
761 end if
762 end do
763 !$omp end parallel
764
766 end subroutine ks_potential_mix_potentials
768
769#include "undef.F90"
770#include "complex.F90"
771#include "ks_potential_inc.F90"
772
773#include "undef.F90"
774#include "real.F90"
775#include "ks_potential_inc.F90"
776
777end module ks_potential_oct_m
778
779!! Local Variables:
780!! mode: f90
781!! coding: utf-8
782!! 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:918
pure logical function, public accel_is_enabled()
Definition: accel.F90:418
integer, parameter, public accel_mem_read_only
Definition: accel.F90:195
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:190
integer, parameter, public independent_particles
Theory level.
Definition: global.F90:236
logical pure function, public not_in_openmp()
Definition: global.F90:489
real(real64), parameter, public m_one
Definition: global.F90:191
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 ks_potential_load_vhxc(this, restart, space, mesh, ierr)
Loads the vhxc potentials.
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_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_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 ks_potential_dump_vhxc(this, restart, space, mesh, ierr)
Dumps the vhxc 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:600
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