23 #include "../base/base_uses.f90"
31 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'particle_types'
36 TYPE(atomic_kind_type),
POINTER :: atomic_kind => null()
37 REAL(KIND=
dp),
DIMENSION(3) :: f = 0.0_dp, &
41 INTEGER :: atom_index = 0, &
44 END TYPE particle_type
48 PUBLIC :: particle_type
69 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
70 INTEGER,
INTENT(IN) :: nparticle
72 IF (
ASSOCIATED(particle_set))
THEN
75 ALLOCATE (particle_set(nparticle))
87 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
89 IF (
ASSOCIATED(particle_set))
THEN
90 DEALLOCATE (particle_set)
91 NULLIFY (particle_set)
107 TYPE(particle_type),
INTENT(INOUT) :: particle_set(:)
109 CLASS(mp_comm_type),
INTENT(IN) :: int_group
110 REAL(kind=
dp),
INTENT(INOUT),
OPTIONAL :: pos(:, :), vel(:, :),
for(:, :)
111 LOGICAL,
INTENT(IN),
OPTIONAL :: add
113 CHARACTER(len=*),
PARAMETER :: routinen =
'update_particle_set'
115 INTEGER :: handle, iparticle, nparticle
116 LOGICAL :: my_add, update_for, update_pos, &
119 CALL timeset(routinen, handle)
121 nparticle =
SIZE(particle_set)
122 update_pos =
PRESENT(pos)
123 update_vel =
PRESENT(vel)
124 update_for =
PRESENT(
for)
126 IF (
PRESENT(add)) my_add = add
129 CALL int_group%sum(pos)
131 DO iparticle = 1, nparticle
132 particle_set(iparticle)%r(:) = particle_set(iparticle)%r(:) + pos(:, iparticle)
135 DO iparticle = 1, nparticle
136 particle_set(iparticle)%r(:) = pos(:, iparticle)
141 CALL int_group%sum(vel)
143 DO iparticle = 1, nparticle
144 particle_set(iparticle)%v(:) = particle_set(iparticle)%v(:) + vel(:, iparticle)
147 DO iparticle = 1, nparticle
148 particle_set(iparticle)%v(:) = vel(:, iparticle)
153 CALL int_group%sum(
for)
155 DO iparticle = 1, nparticle
156 particle_set(iparticle)%f(:) = particle_set(iparticle)%f(:) +
for(:, iparticle)
159 DO iparticle = 1, nparticle
160 particle_set(iparticle)%f(:) =
for(:, iparticle)
165 CALL timestop(handle)
182 INTEGER,
INTENT(IN) :: iatom
183 TYPE(particle_type),
DIMENSION(:),
INTENT(IN) :: particle_set
184 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: vector
185 REAL(kind=
dp),
DIMENSION(3) :: x
188 REAL(kind=
dp) :: fc, fs, mass
191 IF (particle_set(iatom)%shell_index == 0)
THEN
192 x(1:3) = vector(ic + 1:ic + 3)
194 is = 3*(
SIZE(particle_set) + particle_set(iatom)%shell_index - 1)
195 mass = particle_set(iatom)%atomic_kind%mass
196 fc = particle_set(iatom)%atomic_kind%shell%mass_core/mass
197 fs = particle_set(iatom)%atomic_kind%shell%mass_shell/mass
198 x(1:3) = fc*vector(ic + 1:ic + 3) + fs*vector(is + 1:is + 3)
219 INTEGER,
INTENT(IN) :: iatom
220 TYPE(particle_type),
DIMENSION(:),
INTENT(IN) :: particle_set
221 REAL(kind=
dp),
DIMENSION(3),
INTENT(INOUT) :: x
222 REAL(kind=
dp),
DIMENSION(:),
INTENT(INOUT) :: vector
225 REAL(kind=
dp) :: fc, fs, mass
228 IF (particle_set(iatom)%shell_index == 0)
THEN
229 vector(ic + 1:ic + 3) = vector(ic + 1:ic + 3) + x(1:3)
230 x(1:3) = vector(ic + 1:ic + 3)
232 is = 3*(
SIZE(particle_set) + particle_set(iatom)%shell_index - 1)
233 mass = particle_set(iatom)%atomic_kind%mass
234 fc = particle_set(iatom)%atomic_kind%shell%mass_core/mass
235 fs = particle_set(iatom)%atomic_kind%shell%mass_shell/mass
236 vector(ic + 1:ic + 3) = vector(ic + 1:ic + 3) + x(1:3)
237 vector(is + 1:is + 3) = vector(is + 1:is + 3) + x(1:3)
238 x(1:3) = fc*vector(ic + 1:ic + 3) + fs*vector(is + 1:is + 3)
for(int lxp=0;lxp<=lp;lxp++)
Define the atomic kind types and their sub types.
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
Define the data structure for the particle information.
pure subroutine, public update_particle_pos_or_vel(iatom, particle_set, x, vector)
Update the atomic position or velocity by x and return the updated atomic position or velocity in x e...
subroutine, public deallocate_particle_set(particle_set)
Deallocate a particle set.
pure real(kind=dp) function, dimension(3), public get_particle_pos_or_vel(iatom, particle_set, vector)
Return the atomic position or velocity of atom iatom in x from a packed vector even if core-shell par...
subroutine, public update_particle_set(particle_set, int_group, pos, vel, for, add)
...
subroutine, public allocate_particle_set(particle_set, nparticle)
Allocate a particle set.