Octopus
nlcc.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!! Copyright (C) 2021 Nicolas Tancogne-Dejean
3!!
4!! This program is free software; you can redistribute it and/or modify
5!! it under the terms of the GNU General Public License as published by
6!! the Free Software Foundation; either version 2, or (at your option)
7!! any later version.
8!!
9!! This program is distributed in the hope that it will be useful,
10!! but WITHOUT ANY WARRANTY; without even the implied warranty of
11!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12!! GNU General Public License for more details.
13!!
14!! You should have received a copy of the GNU General Public License
15!! along with this program; if not, write to the Free Software
16!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
17!! 02110-1301, USA.
18!!
19
20#include "global.h"
21
22module nlcc_oct_m
23 use atom_oct_m
24 use comm_oct_m
25 use debug_oct_m
28 use global_oct_m
31 use ions_oct_m
34 use mesh_oct_m
39 use ps_oct_m
42 use space_oct_m
47
48 implicit none
49
50 private
51 public :: &
52 nlcc_t
53
54 type, extends(density_interaction_t) :: nlcc_t
55 private
56
57 type(mesh_t), pointer :: mesh => null()
58 type(space_t), pointer :: space => null()
59
60 ! This information should be copied, and not obtained thanks to pointer
61 ! This is a temporary change here
62 type(distributed_t), pointer, public :: atoms_dist => null()
63 type(atom_t), pointer, public :: atom(:) => null()
64 real(real64), pointer, public :: pos(:,:) => null()
65 type(lattice_vectors_t), pointer :: latt => null()
66
67 contains
68 procedure :: init => nlcc_init
69 procedure :: calculate => nlcc_calculate
70 procedure :: calculate_energy => nlcc_calculate_energy
71 procedure :: end => nlcc_end
72 final :: nlcc_finalize
73 end type nlcc_t
74
75 interface nlcc_t
76 module procedure nlcc_constructor
77 end interface nlcc_t
78
79contains
80
81 ! ---------------------------------------------------------
82 function nlcc_constructor(partner) result(this)
83 class(interaction_partner_t), target, intent(inout) :: partner
84 class(nlcc_t), pointer :: this
85
86 push_sub(nlcc_constructor)
87
88 allocate(this)
89
90 this%label = "nlcc"
91
92 this%partner => partner
93
94 pop_sub(nlcc_constructor)
95 end function nlcc_constructor
96
97 ! ---------------------------------------------------------
98 subroutine nlcc_init(this, mesh, ions)
99 class(nlcc_t), intent(inout) :: this
100 class(mesh_t), target, intent(in) :: mesh
101 type(ions_t), target, intent(in) :: ions
102
103
104 push_sub(nlcc_init)
105
106 this%mesh => mesh
107
108 safe_allocate(this%density(1:mesh%np,1))
109
110 this%atoms_dist => ions%atoms_dist
111 this%atom => ions%atom
112 this%space => ions%space
113 this%pos => ions%pos
114 this%latt => ions%latt
116 pop_sub(nlcc_init)
117 end subroutine nlcc_init
118
119 ! ---------------------------------------------------------
120 subroutine nlcc_calculate(this)
121 class(nlcc_t), intent(inout) :: this
122
123 integer :: ia
124
125 push_sub(nlcc_calculate)
126
127 call profiling_in("NLCC_INTER")
128
129 this%density(:,:) = m_zero
130
131 do ia = this%atoms_dist%start, this%atoms_dist%end
132 if (this%atom(ia)%species%is_ps_with_nlcc()) then
133 call species_get_nlcc(this%atom(ia)%species, this%space, this%latt, this%pos(:, ia), this%mesh, &
134 this%density(:,1), accumulate=.true.)
135 end if
136 end do
137
138 ! reduce over atoms if required
139 if (this%atoms_dist%parallel) then
140 call comm_allreduce(this%atoms_dist%mpi_grp, this%density(:,1))
141 end if
142
143 call profiling_out("NLCC_INTER")
144
145 pop_sub(nlcc_calculate)
146 end subroutine nlcc_calculate
148 ! ---------------------------------------------------------
149 subroutine nlcc_end(this)
150 class(nlcc_t), intent(inout) :: this
152 push_sub(nlcc_end)
153
154 safe_deallocate_a(this%density)
155 nullify(this%mesh)
159 pop_sub(nlcc_end)
160 end subroutine nlcc_end
163 ! ---------------------------------------------------------
164 subroutine nlcc_finalize(this)
165 type(nlcc_t), intent(inout) :: this
166
167 push_sub(nlcc_finalize)
168
169 call nlcc_end(this)
170
171 pop_sub(nlcc_finalize)
172 end subroutine nlcc_finalize
173
174 ! ---------------------------------------------------------
175 subroutine nlcc_calculate_energy(this)
176 class(nlcc_t), intent(inout) :: this
177
178 push_sub(nlcc_calculate_energy)
179
180 !NLCC has no energy by itself. The energy is obtained through exchange and correlation
181
182 this%energy = m_zero
183
184 pop_sub(nlcc_calculate_energy)
185 end subroutine nlcc_calculate_energy
186
187end module nlcc_oct_m
188
189!! Local Variables:
190!! mode: f90
191!! coding: utf-8
192!! End:
real(real64), parameter, public m_zero
Definition: global.F90:187
This module defines the abstract interaction_t class, and some auxiliary classes for interactions.
subroutine, public interaction_end(this)
This module defines classes and functions for interaction partners.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
subroutine nlcc_end(this)
Definition: nlcc.F90:243
subroutine nlcc_calculate(this)
Definition: nlcc.F90:214
subroutine nlcc_finalize(this)
Definition: nlcc.F90:258
subroutine nlcc_init(this, mesh, ions)
Definition: nlcc.F90:192
class(nlcc_t) function, pointer nlcc_constructor(partner)
Definition: nlcc.F90:176
subroutine nlcc_calculate_energy(this)
Definition: nlcc.F90:269
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
Definition: profiling.F90:623
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
Definition: profiling.F90:552
Definition: ps.F90:114
This module defines the quantity_t class and the IDs for quantities, which can be exposed by a system...
Definition: quantity.F90:137
subroutine, public species_get_nlcc(species, space, latt, pos, mesh, rho_core, accumulate)
int true(void)