34 use,
intrinsic :: iso_fortran_env
69 logical :: needs_vtau = .false.
71 real(real64),
public,
allocatable :: vhartree(:)
72 real(real64),
public,
allocatable :: vxc(:,:)
73 real(real64),
public,
allocatable :: vhxc(:,:)
76 real(real64),
allocatable :: vtau(:,:)
77 real(real64),
allocatable :: lapl_vtau(:,:)
79 type(accel_mem_t) :: vtau_accel
80 type(accel_mem_t) :: lapl_vtau_accel
82 type(derivatives_t),
pointer :: der => null()
117 real(real64),
public,
allocatable :: vhxc(:,:)
118 real(real64),
allocatable :: vtau(:,:)
119 real(real64),
allocatable :: lapl_vtau(:,:)
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
138 this%theory_level = theory_level
139 this%needs_vtau = needs_vtau
146 safe_allocate(this%vhxc(1:np, 1:nspin))
147 this%vhxc(1:np, 1:nspin) =
m_zero
151 safe_allocate(this%vhartree(1:np_part))
154 safe_allocate(this%vxc(1:np, 1:nspin))
158 safe_allocate(this%vtau(1:np_part, 1:nspin))
160 safe_allocate(this%lapl_vtau(1:np, 1:nspin))
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)
198 real(real64),
contiguous,
intent(in) :: vtau(:,:)
202 assert(this%needs_vtau)
215 integer :: offset, ispin
219 assert(this%needs_vtau)
221 do ispin = 1, this%nspin
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)
240 real(real64),
contiguous,
intent(inout) :: pot(:,:)
241 integer,
optional,
intent(in) :: nspin
249 assert(
size(pot, dim=1) >= this%np)
250 assert(
size(pot, dim=2) >= nspin_)
264 class(
mesh_t),
intent(in) :: mesh
265 integer,
intent(out) :: ierr
267 integer :: iunit, err, err2(2), isp
268 character(len=MAX_PATH_LEN) :: filename
269 character(len=100) :: lines(2)
281 message(1) =
"Debug: Writing Vhxc restart."
286 iunit = restart%open(
'vhxc')
287 lines(1) =
'# #spin #nspin filename'
289 call restart%write(iunit, lines, 2, err)
290 if (err /= 0) ierr = ierr + 1
293 do isp = 1, this%nspin
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
299 call restart%write_mesh_function(filename, mesh, this%vhxc(:,isp), err)
300 if (err /= 0) err2(2) = err2(2) + 1
303 if (err2(1) /= 0) ierr = ierr + 2
304 if (err2(2) /= 0) ierr = ierr + 4
307 call restart%write(iunit, lines, 1, err)
308 if (err /= 0) ierr = ierr + 4
311 if (this%needs_vtau)
then
312 lines(1) =
'# #spin #nspin filename'
314 call restart%write(iunit, lines, 2, err)
315 if (err /= 0) ierr = ierr + 8
318 do isp = 1, this%nspin
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
324 call restart%write_mesh_function(filename, mesh, this%vtau(:,isp), err)
325 if (err /= 0) err2(1) = err2(1) + 1
328 if (err2(1) /= 0) ierr = ierr + 32
329 if (err2(2) /= 0) ierr = ierr + 64
332 call restart%write(iunit, lines, 1, err)
333 if (err /= 0) ierr = ierr + 128
336 call restart%close(iunit)
339 message(1) =
"Debug: Writing Vhxc restart done."
352 class(
mesh_t),
intent(in) :: mesh
353 integer,
intent(out) :: ierr
355 integer :: err, err2, isp
356 character(len=MAX_PATH_LEN) :: filename
369 message(1) =
"Debug: Reading Vhxc restart."
374 do isp = 1, this%nspin
377 call restart%read_mesh_function(filename, mesh, this%vhxc(:,isp), err)
378 if (err /= 0) err2 = err2 + 1
381 if (err2 /= 0) ierr = ierr + 1
385 if (this%needs_vtau)
then
386 do isp = 1, this%nspin
389 call restart%read_mesh_function(filename, mesh, this%vtau(:,isp), err)
390 if (err /= 0) err2 = err2 + 1
394 if (err2 /= 0) ierr = ierr + 2
398 message(1) =
"Debug: Reading Vhxc restart done."
409 integer,
optional,
intent(in) :: order
434 real(real64),
intent(in) :: current_time
435 real(real64),
intent(in) :: dt
450 integer,
intent(in) :: history
457 if (.not.
present(storage))
then
458 if (this%needs_vtau)
then
466 if (this%needs_vtau)
then
480 integer,
intent(in) :: history
486 if (this%needs_vtau)
then
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
509 if (this%needs_vtau)
then
511 this%vhxc, vtau = this%vtau)
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
571 character(len=MAX_PATH_LEN) :: fname
577 space, mesh, this%vhartree,
units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
579 do is = 1, this%nspin
582 mesh, this%vxc(:, is),
units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
584 if (this%needs_vtau)
then
587 mesh, this%vtau(:, is),
units_out%energy, err, pos=pos, atoms=atoms, grp = grp)
601 if (.not.
allocated(copy%vhxc))
then
602 safe_allocate(copy%vhxc(1:this%np, 1:this%nspin))
604 if (this%needs_vtau)
then
605 if (.not.
allocated(copy%vtau))
then
606 safe_allocate(copy%vtau(1:this%np, 1:this%nspin))
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)
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)
653 integer(int64),
intent(in) :: np
654 integer,
intent(in) :: nspin
655 integer(int64),
intent(in) :: pnp
673 safe_deallocate_a(this%vhxc)
674 safe_deallocate_a(this%vtau)
683 real(real64) function ks_potential_check_convergence(this, copy, mesh, rho, qtot) result(diff)
686 class(
mesh_t),
intent(in) :: mesh
687 real(real64),
intent(in) :: rho(:,:)
688 real(real64),
intent(in) :: qtot
690 real(real64),
allocatable :: diff_pot(:)
693 push_sub(ks_potential_check_convergence)
695 safe_allocate(diff_pot(1:this%np))
700 diff_pot(ip) = abs(copy%vhxc(ip, 1) - this%vhxc(ip, 1))*rho(ip, 1)
702 do is = 2, this%nspin
705 diff_pot(ip) = diff_pot(ip) + abs(copy%vhxc(ip, is) - this%vhxc(ip, is))*rho(ip, is)
711 safe_deallocate_a(diff_pot)
713 pop_sub(ks_potential_check_convergence)
719 type(potential_interpolation_t),
intent(inout) :: vksold
720 real(real64),
intent(in) :: times(:)
721 real(real64),
intent(in) :: current_time
725 assert(
size(times) == 3)
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))
739 real(real64),
intent(in) :: dt
745 assert(not_in_openmp())
748 do ispin = 1, this%nspin
751 vold%vhxc(ip, ispin) = m_half*dt*(this%vhxc(ip, ispin) - vold%vhxc(ip, ispin))
754 if (this%needs_vtau)
then
757 vold%vtau(ip, ispin) = m_half*dt*(this%vtau(ip, ispin) - vold%vtau(ip, ispin))
768#include "complex.F90"
769#include "ks_potential_inc.F90"
773#include "ks_potential_inc.F90"
constant times a vector plus a vector
Copies a vector x, to a vector y.
subroutine, public accel_release_buffer(this, async)
pure logical function, public accel_is_enabled()
integer, parameter, public accel_mem_read_only
type(debug_t), save, public debug
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
integer, parameter, public independent_particles
Theory level.
logical pure function, public not_in_openmp()
real(real64), parameter, public m_one
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.
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
this module contains the low-level part of the output system
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.
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
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
Describes mesh distribution to nodes.
Quantities used in mixing: Input, output and new potentials, and the residuals.
This is defined even when running serial.