Octopus
atom.F90
Go to the documentation of this file.
1#include "global.h"
2
3module atom_oct_m
10 use unit_oct_m
12
13 implicit none
14
15 private
16 public :: &
17 atom_init, &
23
24 type, public :: atom_t
25 !private
26 character(len=LABEL_LEN) :: label = ""
27 class(species_t), pointer :: species => null()
28 integer, allocatable :: c(:)
29
30 !Components of the force
31 real(real64), allocatable :: f_ii(:)
32 real(real64), allocatable :: f_vdw(:)
33 real(real64), allocatable :: f_loc(:)
34 real(real64), allocatable :: f_nl(:)
35 real(real64), allocatable :: f_fields(:)
36 real(real64), allocatable :: f_u(:)
37 real(real64), allocatable :: f_scf(:)
38 real(real64), allocatable :: f_nlcc(:)
39 real(real64), allocatable :: f_photons(:)
40 contains
41 procedure :: copy => atom_copy
42 generic :: assignment(=) => copy
43 procedure :: get_label => atom_get_label
44 final :: atom_finalize
45 end type atom_t
46
47 interface atom_same_species
48 module procedure atom_same_species_aa
49 module procedure atom_same_species_as
50 end interface atom_same_species
51
52contains
53
54 ! ---------------------------------------------------------
55 subroutine atom_init(this, dim, label, species)
56 type(atom_t), intent(out) :: this
57 integer, intent(in) :: dim
58 character(len=*), intent(in) :: label
59 class(species_t), target, optional, intent(in) :: species
60
61 push_sub(atom_init)
62
63 this%label = trim(adjustl(label))
64 this%species => null()
65 if (present(species)) this%species => species
66
67 safe_allocate(this%c(1:dim))
68 this%c = 0
69
70 safe_allocate(this%f_ii(1:dim))
71 safe_allocate(this%f_vdw(1:dim))
72 safe_allocate(this%f_loc(1:dim))
73 safe_allocate(this%f_nl(1:dim))
74 safe_allocate(this%f_fields(1:dim))
75 safe_allocate(this%f_u(1:dim))
76 safe_allocate(this%f_scf(1:dim))
77 safe_allocate(this%f_nlcc(1:dim))
78 safe_allocate(this%f_photons(1:dim))
79 this%f_ii = m_zero
80 this%f_vdw = m_zero
81 this%f_loc = m_zero
82 this%f_nl = m_zero
83 this%f_fields = m_zero
84 this%f_u = m_zero
85 this%f_scf = m_zero
86 this%f_nlcc = m_zero
87 this%f_photons = m_zero
88
89 pop_sub(atom_init)
90 end subroutine atom_init
91
92 ! ---------------------------------------------------------
93 subroutine atom_copy(atom_out, atom_in)
94 class(atom_t), intent(out) :: atom_out
95 class(atom_t), intent(in) :: atom_in
97 push_sub(atom_copy)
98
99 atom_out%label = atom_in%label
100 atom_out%species => atom_in%species
101
102 safe_allocate_source_a(atom_out%c, atom_in%c)
103 safe_allocate_source_a(atom_out%f_ii, atom_in%f_ii)
104 safe_allocate_source_a(atom_out%f_vdw, atom_in%f_vdw)
105 safe_allocate_source_a(atom_out%f_loc, atom_in%f_loc)
106 safe_allocate_source_a(atom_out%f_nl, atom_in%f_nl)
107 safe_allocate_source_a(atom_out%f_fields, atom_in%f_fields)
108 safe_allocate_source_a(atom_out%f_u, atom_in%f_u)
109 safe_allocate_source_a(atom_out%f_scf, atom_in%f_scf)
110 safe_allocate_source_a(atom_out%f_nlcc, atom_in%f_nlcc)
111 safe_allocate_source_a(atom_out%f_photons, atom_in%f_photons)
112
113 pop_sub(atom_copy)
114 end subroutine atom_copy
115
116 ! ---------------------------------------------------------
117 impure elemental subroutine atom_finalize(this)
118 type(atom_t), intent(inout) :: this
122 this%label = ""
123 this%species => null()
125 safe_deallocate_a(this%c)
127 safe_deallocate_a(this%f_ii)
128 safe_deallocate_a(this%f_vdw)
129 safe_deallocate_a(this%f_loc)
130 safe_deallocate_a(this%f_nl)
131 safe_deallocate_a(this%f_fields)
132 safe_deallocate_a(this%f_u)
133 safe_deallocate_a(this%f_scf)
134 safe_deallocate_a(this%f_nlcc)
135 safe_deallocate_a(this%f_photons)
138 end subroutine atom_finalize
139
140 ! ---------------------------------------------------------
141
143 pure function atom_get_label(this) result(label)
144 class(atom_t), intent(in) :: this
145
146 character(len=len_trim(adjustl(this%label))) :: label
147
148 label=trim(adjustl(this%label))
149
150 end function atom_get_label
151
152 ! ---------------------------------------------------------
153 subroutine atom_set_species(this, species)
154 type(atom_t), intent(inout) :: this
155 class(species_t), target, intent(in) :: species
156
157 push_sub(atom_set_species)
158
159 this%species => species
160 pop_sub(atom_set_species)
161
162 end subroutine atom_set_species
163
164 ! ---------------------------------------------------------
165 subroutine atom_get_species(this, species)
166 type(atom_t), target, intent(in) :: this
167 class(species_t), pointer, intent(out) :: species
168
169 ! NO PUSH_SUB, called too often
170
171 species => null()
172 if (associated(this%species)) species => this%species
173
174 end subroutine atom_get_species
175
176 ! ---------------------------------------------------------
177 elemental function atom_same_species_aa(this, that) result(is)
178 type(atom_t), intent(in) :: this
179 type(atom_t), intent(in) :: that
180
181 logical :: is
182
183 is = this%label == that%label
184
185 end function atom_same_species_aa
187 ! ---------------------------------------------------------
188 elemental function atom_same_species_as(this, species) result(is)
189 type(atom_t), intent(in) :: this
190 class(species_t), intent(in) :: species
191
192 logical :: is
193
194 is = this%label == species%get_label()
195
196 end function atom_same_species_as
197
199 pure logical function all_species_are_jellium_slab(atom)
200 type(atom_t), intent(in) :: atom(:)
201
202 integer :: i
203
205 do i = 1, size(atom)
206 select type(spec => atom(i)%species)
207 type is(jellium_slab_t)
208 class default
210 end select
211 enddo
212
214
216 pure logical function any_species_is_jellium_sphere(atom)
217 type(atom_t), intent(in) :: atom(:)
218
219 integer :: i
220
222 do i = 1, size(atom)
223 select type(spec => atom(i)%species)
224 type is(jellium_sphere_t)
226 end select
227 enddo
228
230
231end module atom_oct_m
232
233!! Local Variables:
234!! mode: f90
235!! coding: utf-8
236!! End:
subroutine atom_copy(atom_out, atom_in)
Definition: atom.F90:187
pure character(len=len_trim(adjustl(this%label))) function atom_get_label(this)
Getter for label attribute.
Definition: atom.F90:237
pure logical function, public all_species_are_jellium_slab(atom)
Check if all species are jellium slab.
Definition: atom.F90:293
elemental logical function atom_same_species_aa(this, that)
Definition: atom.F90:271
impure elemental subroutine atom_finalize(this)
Definition: atom.F90:211
pure logical function, public any_species_is_jellium_sphere(atom)
Check if any species is a jellium sphere.
Definition: atom.F90:310
elemental logical function atom_same_species_as(this, species)
Definition: atom.F90:282
subroutine, public atom_init(this, dim, label, species)
Definition: atom.F90:149
subroutine, public atom_get_species(this, species)
Definition: atom.F90:259
subroutine, public atom_set_species(this, species)
Definition: atom.F90:247
real(real64), parameter, public m_zero
Definition: global.F90:188
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:132
This module defines the unit system, used for input and output.
An abstract class for species. Derived classes include jellium, all electron, and pseudopotential spe...
Definition: species.F90:143
int true(void)