1 |
chuckv |
4 |
module xyz_io |
2 |
|
|
use definitions, ONLY : DP,machdep_lnblnk,machdep_flush |
3 |
|
|
use parameter |
4 |
|
|
use simulation |
5 |
|
|
use model_module |
6 |
|
|
use status |
7 |
|
|
use file_units, ONLY : next_unit |
8 |
|
|
#ifdef MPI |
9 |
|
|
use mpi_module |
10 |
|
|
#endif |
11 |
|
|
implicit none |
12 |
|
|
PRIVATE |
13 |
|
|
|
14 |
|
|
|
15 |
|
|
!various common formats |
16 |
|
|
character(len=*), parameter :: xyz_line1 = "(i11)" |
17 |
|
|
character(len=*), parameter :: xyz_line2 = "(i11,3es12.3)" |
18 |
|
|
character(len=*), parameter :: xyz_coord = "(a3,6es12.3,i3)" |
19 |
|
|
|
20 |
|
|
|
21 |
|
|
|
22 |
|
|
! atype interface |
23 |
|
|
character(len=2) :: charident |
24 |
|
|
character(len=2) :: atid |
25 |
|
|
|
26 |
|
|
integer :: dump_unit |
27 |
|
|
integer :: eor_unit |
28 |
|
|
integer :: force_unit |
29 |
|
|
|
30 |
|
|
#ifdef MPI |
31 |
|
|
integer, dimension(mpi_status_size) :: istatus |
32 |
|
|
#endif |
33 |
|
|
|
34 |
|
|
public :: open_xyz_files |
35 |
|
|
public :: read_initial_config |
36 |
|
|
public :: dump_config |
37 |
|
|
public :: write_final_config |
38 |
|
|
public :: dump_forces |
39 |
|
|
#ifdef MPI |
40 |
|
|
public :: mpi_dump_q |
41 |
|
|
#endif |
42 |
|
|
|
43 |
|
|
|
44 |
|
|
contains |
45 |
|
|
|
46 |
|
|
!! open input and output files for simulation |
47 |
|
|
subroutine open_xyz_files(this_error) |
48 |
|
|
character(len=80) :: dump_file, msg |
49 |
|
|
character(len=80) :: force_file |
50 |
|
|
character(len=*), parameter :: routine_name = "OPEN_FILES" |
51 |
|
|
integer, optional :: this_error |
52 |
|
|
integer :: err |
53 |
|
|
#ifndef MPI |
54 |
|
|
integer, parameter :: node = 0 |
55 |
|
|
#endif |
56 |
|
|
|
57 |
|
|
if (present(this_error)) this_error = 0 |
58 |
|
|
|
59 |
|
|
if (node == 0) then |
60 |
|
|
|
61 |
|
|
dump_unit = next_unit() |
62 |
|
|
|
63 |
|
|
dump_file = file_prefix(1:machdep_lnblnk(file_prefix)) // '.dump' |
64 |
|
|
open(file=dump_file, unit=dump_unit, status='replace', iostat=err, & |
65 |
|
|
action='write') |
66 |
|
|
if (err > 0) then |
67 |
|
|
write(msg,*) "Error opening dump file", trim(dump_file) |
68 |
|
|
call info(routine_name, msg) |
69 |
|
|
if (present(this_error)) this_error = -1 |
70 |
|
|
end if |
71 |
|
|
write(msg,*) 'configurations will dump to: ', trim(dump_file) |
72 |
|
|
call info(routine_name, msg) |
73 |
|
|
|
74 |
|
|
|
75 |
|
|
|
76 |
|
|
!!$ force_file = file_prefix(1:machdep_lnblnk(file_prefix)) // '.force' |
77 |
|
|
!!$ open(file=force_file, unit=force_unit, status='replace', iostat= err, & |
78 |
|
|
!!$ action='write') |
79 |
|
|
!!$ |
80 |
|
|
!!$ if (err > 0) then |
81 |
|
|
!!$ write(msg,*) "Error opening force file", trim(force_file) |
82 |
|
|
!!$ call info(routine_name, msg) |
83 |
|
|
!!$ end if |
84 |
|
|
!!$ write(msg,*) 'Force dump will be in: ', trim(force_file) |
85 |
|
|
!!$ call info(routine_name, msg) |
86 |
|
|
|
87 |
|
|
end if |
88 |
|
|
end subroutine open_xyz_files |
89 |
|
|
|
90 |
|
|
|
91 |
|
|
|
92 |
|
|
subroutine read_initial_config(this_error) |
93 |
|
|
! Reads In the current configuration in XYZ file format. The next |
94 |
|
|
! frame in XYZ is simply appended to the end of the current frame. |
95 |
|
|
|
96 |
|
|
!standard io |
97 |
|
|
integer :: fu, i, a, atom, nats, id, junk,j |
98 |
|
|
integer :: natoms_tmp, id_tmp |
99 |
|
|
integer, intent(out), optional :: this_error |
100 |
|
|
character(len=2), allocatable, dimension(:) :: atid_tmp, atid_array |
101 |
|
|
character(len=120) :: msg |
102 |
|
|
integer :: read_error |
103 |
|
|
integer :: xyz_input_unit |
104 |
|
|
integer :: err |
105 |
|
|
character(len=80) :: input_file |
106 |
|
|
|
107 |
|
|
!mpi stuff |
108 |
|
|
integer :: iproc, ncount |
109 |
|
|
#ifndef MPI |
110 |
|
|
integer,parameter :: node = 0 |
111 |
|
|
#endif |
112 |
|
|
|
113 |
|
|
natoms_tmp = 0 |
114 |
|
|
|
115 |
|
|
if (present(this_error)) this_error = 0 |
116 |
|
|
|
117 |
|
|
call info("READ_INITIAL_CONFIG","Reading old trajectory file") |
118 |
|
|
#ifdef MPI |
119 |
|
|
if (node == 0) then |
120 |
|
|
#endif |
121 |
|
|
input_file = file_prefix(1:machdep_lnblnk(file_prefix)) // '.in' |
122 |
|
|
xyz_input_unit = next_unit() |
123 |
|
|
open(file=input_file, unit=xyz_input_unit, status='old',iostat=err, & |
124 |
|
|
action='read') |
125 |
|
|
if (err > 0 ) then |
126 |
|
|
write(msg,*) "Error opening config file", trim(input_file) |
127 |
|
|
call info("READ_INITIAL_CONFIG", msg) |
128 |
|
|
if (present(this_error)) this_error = -1 |
129 |
|
|
end if |
130 |
|
|
|
131 |
|
|
write(msg,*) 'reading initial config from: ', trim(input_file) |
132 |
|
|
call info("READ_INITIAL_CONFIG", msg) |
133 |
|
|
#ifdef MPI |
134 |
|
|
endif |
135 |
|
|
#endif |
136 |
|
|
|
137 |
|
|
|
138 |
|
|
#ifndef MPI |
139 |
|
|
|
140 |
|
|
read(unit=xyz_input_unit,fmt=*,iostat=read_error) natoms_tmp |
141 |
|
|
read(xyz_input_unit,*) |
142 |
|
|
|
143 |
|
|
if (read_error /= 0) then |
144 |
|
|
write(msg,*) "Error reading particle configuration: ", i, "Atom: ",atom |
145 |
|
|
call info("READ_INITIAL_CONFIG",msg) |
146 |
|
|
if (present(this_error)) this_error = -1 |
147 |
|
|
return |
148 |
|
|
else if (natoms /= natoms_tmp) then |
149 |
|
|
write(msg,*) 'Atom number mismatch read: ',natoms_tmp, & |
150 |
|
|
'Allocated: ', natoms |
151 |
|
|
call info('READ_INITIAL_CONFIG', msg) |
152 |
|
|
if (present(this_error)) this_error = -1 |
153 |
|
|
return |
154 |
|
|
end if |
155 |
|
|
|
156 |
|
|
|
157 |
|
|
|
158 |
|
|
! this block is only true for monatomic systems |
159 |
|
|
!! fix me for general case |
160 |
|
|
nmol = natoms |
161 |
|
|
do i = 1, nmol |
162 |
|
|
na(i) = 1 |
163 |
|
|
atom_index(i,1) = i |
164 |
|
|
enddo |
165 |
|
|
|
166 |
|
|
do i = 1, nmol |
167 |
|
|
do a = 1, na(i) |
168 |
|
|
atom = atom_index(i,a) |
169 |
|
|
read(unit=xyz_input_unit,fmt=*,iostat=read_error) atid, q(1,atom), q(2,atom), q(3,atom), & |
170 |
|
|
v(1,atom), v(2,atom), v(3,atom) |
171 |
|
|
if (read_error /= 0) then |
172 |
|
|
write(msg,*) "Error reading particle configuration: ", i, "Atom: ",atom |
173 |
|
|
call info("READ_INITIAL_CONFIG",msg) |
174 |
|
|
if (present(this_error)) this_error = -1 |
175 |
|
|
return |
176 |
|
|
end if |
177 |
|
|
id = atype_atoi(atid) |
178 |
|
|
ident(atom) = id |
179 |
|
|
! this line only works for monatomic systems: |
180 |
|
|
model(i) = atype_identifier(id) |
181 |
|
|
mass(atom) = atype_mass(id) |
182 |
|
|
end do |
183 |
|
|
end do |
184 |
|
|
|
185 |
|
|
#else |
186 |
|
|
if (node == 0) then |
187 |
|
|
allocate(atid_tmp(max_local)) |
188 |
|
|
else |
189 |
|
|
allocate(atid_array(max_local)) |
190 |
|
|
end if |
191 |
|
|
|
192 |
|
|
|
193 |
|
|
if (node == 0) then |
194 |
|
|
read(xyz_input_unit,*) natoms_tmp |
195 |
|
|
read(xyz_input_unit,*) |
196 |
|
|
end if |
197 |
|
|
|
198 |
|
|
call mpi_bcast(natoms_tmp,1,mpi_integer,0,mpi_comm_world,mpi_err) |
199 |
|
|
if (mpi_err /= 0) then |
200 |
|
|
call error("READ_INITIAL_CONFIG","Error bcast natoms_tmp") |
201 |
|
|
endif |
202 |
|
|
|
203 |
|
|
if (natoms_tmp /= natoms) then |
204 |
|
|
write(msg,*) 'Atom number mismatch read: ',natoms_tmp, & |
205 |
|
|
'Allocated: ', natoms |
206 |
|
|
call error('READ_INITIAL_CONFIG', msg) |
207 |
|
|
end if |
208 |
|
|
|
209 |
|
|
call mpi_bcast(natoms,1,mpi_integer,0,mpi_comm_world,mpi_err) |
210 |
|
|
if (mpi_err /= 0) then |
211 |
|
|
call error("READ_INITIAL_CONFIG","Error bcast natoms") |
212 |
|
|
endif |
213 |
|
|
|
214 |
|
|
! this block is only true for monatomic systems |
215 |
|
|
nmol = natoms |
216 |
|
|
|
217 |
|
|
do i = 1, maxmol_local |
218 |
|
|
na(i) = 1 |
219 |
|
|
atom_index(i,1) = i |
220 |
|
|
end do |
221 |
|
|
|
222 |
|
|
if (node == 0) then |
223 |
|
|
do iproc = 0, nprocs - 1 |
224 |
|
|
|
225 |
|
|
ncount = nint(float(iproc+1)/nprocs*nmol) - & |
226 |
|
|
nint(float(iproc)/nprocs*nmol) |
227 |
|
|
do i = 1,ncount |
228 |
|
|
do a = 1, na(i) |
229 |
|
|
atom = atom_index(i,a) |
230 |
|
|
read(xyz_input_unit,*) atid_tmp(i), q_tmp(1,atom), q_tmp(2,atom), & |
231 |
|
|
q_tmp(3,atom), & |
232 |
|
|
v_tmp(1,atom), v_tmp(2,atom), v_tmp(3,atom) |
233 |
|
|
end do |
234 |
|
|
end do |
235 |
|
|
|
236 |
|
|
if (iproc.eq.0) then |
237 |
|
|
do i = 1,ncount |
238 |
|
|
do a = 1, na(i) |
239 |
|
|
atom = atom_index(i,a) |
240 |
|
|
q(1,atom) = q_tmp(1,atom) |
241 |
|
|
q(2,atom) = q_tmp(2,atom) |
242 |
|
|
q(3,atom) = q_tmp(3,atom) |
243 |
|
|
v(1,atom) = v_tmp(1,atom) |
244 |
|
|
v(2,atom) = v_tmp(2,atom) |
245 |
|
|
v(3,atom) = v_tmp(3,atom) |
246 |
|
|
id = atype_atoi(atid_tmp(i)) |
247 |
|
|
ident(atom) = id |
248 |
|
|
! this line only works for monatomic systems: |
249 |
|
|
model(i) = atype_identifier(id) |
250 |
|
|
mass(atom) = atype_mass(id) |
251 |
|
|
end do |
252 |
|
|
enddo |
253 |
|
|
else |
254 |
|
|
|
255 |
|
|
call mpi_send(atid_tmp,2*ncount,mpi_character,iproc,0, & |
256 |
|
|
mpi_comm_world,mpi_err) |
257 |
|
|
if (mpi_err /= 0) then |
258 |
|
|
call error("READ_INITIAL_CONFIG","Error mpi_send atid_tmp") |
259 |
|
|
endif |
260 |
|
|
call mpi_send(q_tmp,3*ncount,mpi_double_precision, & |
261 |
|
|
iproc,0,mpi_comm_world,mpi_err) |
262 |
|
|
if (mpi_err /= 0) then |
263 |
|
|
call error("READ_INITIAL_CONFIG","Error mpi_send q_tmp") |
264 |
|
|
endif |
265 |
|
|
call mpi_send(v_tmp,3*ncount,mpi_double_precision, & |
266 |
|
|
iproc,0,mpi_comm_world,mpi_err) |
267 |
|
|
if (mpi_err /= 0) then |
268 |
|
|
call error("READ_INITIAL_CONFIG","Error mpi_send v_tmp") |
269 |
|
|
endif |
270 |
|
|
|
271 |
|
|
endif |
272 |
|
|
|
273 |
|
|
end do |
274 |
|
|
close(xyz_input_unit) |
275 |
|
|
else |
276 |
|
|
call mpi_recv(atid_array,2*max_local,mpi_character,0,0, & |
277 |
|
|
mpi_comm_world,istatus,mpi_err) |
278 |
|
|
if (mpi_err /= 0) then |
279 |
|
|
call error("READ_INITIAL_CONFIG","Error reading atid_array") |
280 |
|
|
endif |
281 |
|
|
call mpi_recv(q,3*max_local, & |
282 |
|
|
mpi_double_precision,0,0,mpi_comm_world,istatus,mpi_err) |
283 |
|
|
if (mpi_err /= 0) then |
284 |
|
|
call error("READ_INITIAL_CONFIG","Error reading q_array") |
285 |
|
|
endif |
286 |
|
|
call mpi_recv(v,3*max_local, & |
287 |
|
|
mpi_double_precision,0,0,mpi_comm_world,istatus,mpi_err) |
288 |
|
|
if (mpi_err /= 0) then |
289 |
|
|
call error("READ_INITIAL_CONFIG","Error reading v_array") |
290 |
|
|
endif |
291 |
|
|
|
292 |
|
|
do j = 1,nlocal |
293 |
|
|
! array error here..... |
294 |
|
|
id = atype_atoi(atid_array(j)) |
295 |
|
|
ident(j) = id |
296 |
|
|
! this line only works for monatomic systems: |
297 |
|
|
model(j) = atype_identifier(id) |
298 |
|
|
mass(j) = atype_mass(id) |
299 |
|
|
end do |
300 |
|
|
|
301 |
|
|
end if |
302 |
|
|
call mpi_barrier(mpi_comm_world,mpi_err) |
303 |
|
|
|
304 |
|
|
|
305 |
|
|
|
306 |
|
|
if (node == 0) then |
307 |
|
|
deallocate(atid_tmp) |
308 |
|
|
close(xyz_input_unit) |
309 |
|
|
else |
310 |
|
|
deallocate(atid_array) |
311 |
|
|
end if |
312 |
|
|
#endif |
313 |
|
|
end subroutine read_initial_config |
314 |
|
|
|
315 |
|
|
subroutine dump_config(time) |
316 |
|
|
include 'headers/atom.h' |
317 |
|
|
|
318 |
|
|
real(kind=DP), intent(in) :: time |
319 |
|
|
#ifdef MPI |
320 |
|
|
integer, dimension(maxmol_local) :: na_tmp |
321 |
|
|
integer, dimension(maxmol_local,3) :: atom_index_tmp |
322 |
|
|
integer, dimension(max_local) :: ident_tmp |
323 |
|
|
#endif |
324 |
|
|
integer ::i,j, iproc, a |
325 |
|
|
integer :: atom |
326 |
|
|
integer :: icnt, i_a_cnt |
327 |
|
|
character(len=10) :: chartmp |
328 |
|
|
|
329 |
|
|
#ifndef MPI |
330 |
|
|
write(dump_unit,*) natoms |
331 |
|
|
write(dump_unit,*) time |
332 |
|
|
|
333 |
|
|
do i = 1, natoms |
334 |
|
|
charident = atype_identifier(ident(i)) // ' ' |
335 |
|
|
write(unit=dump_unit,fmt=xyz_coord) charident, (q(j,i),j=1,3), & |
336 |
|
|
(v(j,i),j=1,3), ident(i) |
337 |
|
|
end do |
338 |
|
|
|
339 |
|
|
call machdep_flush(dump_unit) |
340 |
|
|
#else |
341 |
|
|
if (node /= 0) then |
342 |
|
|
call mpi_send(na,nmol_local,mpi_integer,0,0, & |
343 |
|
|
mpi_comm_world,mpi_err) |
344 |
|
|
if (mpi_err /= 0) then |
345 |
|
|
call error("DUMP_CONFIG","Error sending na") |
346 |
|
|
endif |
347 |
|
|
call mpi_send(atom_index,3*nmol_local,mpi_integer,0,0, & |
348 |
|
|
mpi_comm_world,mpi_err) |
349 |
|
|
if (mpi_err /= 0) then |
350 |
|
|
call error("DUMP_CONFIG","Error sending atom_index") |
351 |
|
|
endif |
352 |
|
|
call mpi_send(ident,nlocal,mpi_integer,0,0, & |
353 |
|
|
mpi_comm_world,mpi_err) |
354 |
|
|
if (mpi_err /= 0) then |
355 |
|
|
call error("DUMP_CONFIG","Error sending ident") |
356 |
|
|
endif |
357 |
|
|
call mpi_send(q,3*nlocal,mpi_double_precision,0,0, & |
358 |
|
|
mpi_comm_world,mpi_err) |
359 |
|
|
if (mpi_err /= 0) then |
360 |
|
|
call error("DUMP_CONFIG","Error sending q") |
361 |
|
|
endif |
362 |
|
|
call mpi_send(v,3*nlocal,mpi_double_precision,0,0, & |
363 |
|
|
mpi_comm_world,mpi_err) |
364 |
|
|
if (mpi_err /= 0) then |
365 |
|
|
call error("DUMP_CONFIG","Error sending v") |
366 |
|
|
endif |
367 |
|
|
|
368 |
|
|
else |
369 |
|
|
|
370 |
|
|
|
371 |
|
|
write(dump_unit,*) natoms |
372 |
|
|
write(dump_unit,*) time |
373 |
|
|
|
374 |
|
|
|
375 |
|
|
do i = 1, nlocal |
376 |
|
|
charident = atype_identifier(ident(i)) // ' ' |
377 |
|
|
write(unit=dump_unit,fmt=xyz_coord) charident, (q(j,i),j=1,3), & |
378 |
|
|
(v(j,i),j=1,3), ident(i) |
379 |
|
|
end do |
380 |
|
|
|
381 |
|
|
do iproc = 1,nprocs - 1 |
382 |
|
|
icnt = 0 |
383 |
|
|
call mpi_recv(na_tmp, maxmol_local, mpi_integer,iproc,0, & |
384 |
|
|
mpi_comm_world,istatus,mpi_err) |
385 |
|
|
|
386 |
|
|
call mpi_get_count(istatus,mpi_integer,icnt,mpi_err) |
387 |
|
|
|
388 |
|
|
call mpi_recv(atom_index_tmp,3*icnt,mpi_integer,iproc,0, & |
389 |
|
|
mpi_comm_world,istatus,mpi_err) |
390 |
|
|
if (mpi_err /= 0) then |
391 |
|
|
call error("DUMP_CONFIG","Error recieving atom_index") |
392 |
|
|
endif |
393 |
|
|
|
394 |
|
|
call mpi_recv(ident_tmp,max_local,mpi_integer,iproc,0, & |
395 |
|
|
mpi_comm_world,istatus,mpi_err) |
396 |
|
|
if (mpi_err /= 0) then |
397 |
|
|
call error("DUMP_CONFIG","Error recieving ident_tmp") |
398 |
|
|
endif |
399 |
|
|
|
400 |
|
|
call mpi_get_count(istatus,mpi_integer,i_a_cnt,mpi_err) |
401 |
|
|
|
402 |
|
|
call mpi_recv(q_tmp,3*i_a_cnt,mpi_double_precision,iproc,0, & |
403 |
|
|
mpi_comm_world,istatus,mpi_err) |
404 |
|
|
if (mpi_err /= 0) then |
405 |
|
|
call error("DUMP_CONFIG","Error recieving q_tmp") |
406 |
|
|
endif |
407 |
|
|
|
408 |
|
|
call mpi_recv(v_tmp,3*i_a_cnt,mpi_double_precision,iproc,0, & |
409 |
|
|
mpi_comm_world,istatus,mpi_err) |
410 |
|
|
if (mpi_err /= 0) then |
411 |
|
|
call error("DUMP_CONFIG","Error recieving v_tmp") |
412 |
|
|
endif |
413 |
|
|
|
414 |
|
|
|
415 |
|
|
do i = 1,icnt |
416 |
|
|
do a = 1, na_tmp(i) |
417 |
|
|
atom = atom_index_tmp(i,a) |
418 |
|
|
atid = atype_identifier(ident_tmp(atom)) // ' ' |
419 |
|
|
write(dump_unit,fmt=xyz_coord) atid, q_tmp(1,atom), & |
420 |
|
|
q_tmp(2,atom), q_tmp(3,atom), v_tmp(1,atom), & |
421 |
|
|
v_tmp(2,atom), v_tmp(3,atom), ident_tmp(atom) |
422 |
|
|
end do |
423 |
|
|
enddo |
424 |
|
|
|
425 |
|
|
end do |
426 |
|
|
endif |
427 |
|
|
|
428 |
|
|
call mpi_barrier(mpi_comm_world,mpi_err) |
429 |
|
|
|
430 |
|
|
if (node == 0) & |
431 |
|
|
call machdep_flush(dump_unit) |
432 |
|
|
#endif |
433 |
|
|
end subroutine dump_config |
434 |
|
|
|
435 |
|
|
|
436 |
|
|
subroutine write_final_config() |
437 |
|
|
! Reads In the current configuration in XYZ file format. The next |
438 |
|
|
! frame in XYZ is simply appended to the end of the current frame. |
439 |
|
|
|
440 |
|
|
include 'headers/atom.h' |
441 |
|
|
|
442 |
|
|
character(len = 80) :: eor_file |
443 |
|
|
character(len = 120) :: msg |
444 |
|
|
integer err |
445 |
|
|
|
446 |
|
|
#ifdef MPI |
447 |
|
|
integer, dimension(maxmol_local) :: na_tmp |
448 |
|
|
integer, dimension(maxmol_local,3) :: atom_index_tmp |
449 |
|
|
integer, dimension(max_local) :: ident_tmp |
450 |
|
|
#endif |
451 |
|
|
integer ::i,j, iproc, a |
452 |
|
|
integer :: atom |
453 |
|
|
integer :: icnt, i_a_cnt |
454 |
|
|
|
455 |
|
|
character(len=10) :: chartmp |
456 |
|
|
|
457 |
|
|
|
458 |
|
|
|
459 |
|
|
|
460 |
|
|
#ifdef MPI |
461 |
|
|
if (node == 0) then |
462 |
|
|
#endif |
463 |
|
|
eor_unit = next_unit() |
464 |
|
|
|
465 |
|
|
eor_file = file_prefix(1:machdep_lnblnk(file_prefix)) // '.eor' |
466 |
|
|
open(file=eor_file, unit=eor_unit, status='replace', iostat= err, & |
467 |
|
|
action='write') |
468 |
|
|
|
469 |
|
|
if (err > 0) then |
470 |
|
|
write(msg,*) "Error opening eor file", trim(eor_file) |
471 |
|
|
call info("DUMP_FINAL_CONFIG", msg) |
472 |
|
|
end if |
473 |
|
|
write(msg,*) 'final configuration will be in: ', trim(eor_file) |
474 |
|
|
call info("DUMP_FINAL_CONFIG", msg) |
475 |
|
|
|
476 |
|
|
#ifdef MPI |
477 |
|
|
endif |
478 |
|
|
#endif |
479 |
|
|
|
480 |
|
|
|
481 |
|
|
|
482 |
|
|
|
483 |
|
|
|
484 |
|
|
|
485 |
|
|
|
486 |
|
|
|
487 |
|
|
#ifndef MPI |
488 |
|
|
|
489 |
|
|
|
490 |
|
|
write(eor_unit,fmt=xyz_line1) natoms |
491 |
|
|
write(eor_unit,*) |
492 |
|
|
|
493 |
|
|
do i = 1, natoms |
494 |
|
|
do a = 1, na(i) |
495 |
|
|
atom = atom_index(i,a) |
496 |
|
|
atid = atype_identifier(ident(atom)) // ' ' |
497 |
|
|
write(eor_unit,fmt=xyz_coord) atid, q(1,atom), q(2,atom), q(3,atom), & |
498 |
|
|
v(1,atom), v(2,atom), v(3,atom), ident(i) |
499 |
|
|
end do |
500 |
|
|
end do |
501 |
|
|
|
502 |
|
|
close(eor_unit) |
503 |
|
|
#else |
504 |
|
|
|
505 |
|
|
if (node /= 0) then |
506 |
|
|
|
507 |
|
|
call mpi_send(na,nmol_local,mpi_integer,0,0, & |
508 |
|
|
mpi_comm_world,mpi_err) |
509 |
|
|
if (mpi_err /= 0) then |
510 |
|
|
call error("WRITE_FINALCONFIG","Error sending na") |
511 |
|
|
endif |
512 |
|
|
call mpi_send(atom_index,3*nmol_local,mpi_integer,0,0, & |
513 |
|
|
mpi_comm_world,mpi_err) |
514 |
|
|
if (mpi_err /= 0) then |
515 |
|
|
call error("WRITE_FINALCONFIG","Error sending atom_index") |
516 |
|
|
endif |
517 |
|
|
call mpi_send(ident,nlocal,mpi_integer,0,0, & |
518 |
|
|
mpi_comm_world,mpi_err) |
519 |
|
|
if (mpi_err /= 0) then |
520 |
|
|
call error("WRITE_FINALCONFIG","Error sending ident") |
521 |
|
|
endif |
522 |
|
|
call mpi_send(q,3*nlocal,mpi_double_precision,0,0, & |
523 |
|
|
mpi_comm_world,mpi_err) |
524 |
|
|
if (mpi_err /= 0) then |
525 |
|
|
call error("WRITE_FINALCONFIG","Error sending q") |
526 |
|
|
endif |
527 |
|
|
call mpi_send(v,3*nlocal,mpi_double_precision,0,0, & |
528 |
|
|
mpi_comm_world,mpi_err) |
529 |
|
|
if (mpi_err /= 0) then |
530 |
|
|
call error("WRITE_FINALCONFIG","Error sending v") |
531 |
|
|
endif |
532 |
|
|
|
533 |
|
|
else |
534 |
|
|
|
535 |
|
|
|
536 |
|
|
write(eor_unit,*) natoms |
537 |
|
|
write(eor_unit,*) |
538 |
|
|
|
539 |
|
|
|
540 |
|
|
do i = 1, nlocal |
541 |
|
|
charident = atype_identifier(ident(i)) // ' ' |
542 |
|
|
write(unit=eor_unit,fmt=xyz_coord) charident, (q(j,i),j=1,3), & |
543 |
|
|
(v(j,i),j=1,3), ident(i) |
544 |
|
|
end do |
545 |
|
|
|
546 |
|
|
do iproc = 1,nprocs - 1 |
547 |
|
|
icnt = 0 |
548 |
|
|
|
549 |
|
|
call mpi_recv(na_tmp, maxmol_local, mpi_integer,iproc,0, & |
550 |
|
|
mpi_comm_world,istatus,mpi_err) |
551 |
|
|
if (mpi_err /= 0) then |
552 |
|
|
call error("WRITE_FINALCONFIG","Error receiving na_tmp") |
553 |
|
|
endif |
554 |
|
|
|
555 |
|
|
call mpi_get_count(istatus,mpi_integer,icnt,mpi_err) |
556 |
|
|
|
557 |
|
|
call mpi_recv(atom_index_tmp,3*icnt,mpi_integer,iproc,0, & |
558 |
|
|
mpi_comm_world,istatus,mpi_err) |
559 |
|
|
if (mpi_err /= 0) then |
560 |
|
|
call error("WRITE_FINALCONFIG","Error receiving atom_index_tmp") |
561 |
|
|
endif |
562 |
|
|
|
563 |
|
|
call mpi_recv(ident_tmp,max_local,mpi_integer,iproc,0, & |
564 |
|
|
mpi_comm_world,istatus,mpi_err) |
565 |
|
|
if (mpi_err /= 0) then |
566 |
|
|
call error("WRITE_FINAL_CONFIG","Error receiving ident_tmp") |
567 |
|
|
endif |
568 |
|
|
|
569 |
|
|
call mpi_get_count(istatus,mpi_integer,i_a_cnt,mpi_err) |
570 |
|
|
call mpi_recv(q_tmp,3*i_a_cnt,mpi_double_precision,iproc,0, & |
571 |
|
|
mpi_comm_world,istatus,mpi_err) |
572 |
|
|
if (mpi_err /= 0) then |
573 |
|
|
call error("WRITE_FINAL_CONFIG","Error receiving q_tmp") |
574 |
|
|
endif |
575 |
|
|
call mpi_recv(v_tmp,3*i_a_cnt,mpi_double_precision,iproc,0, & |
576 |
|
|
mpi_comm_world,istatus,mpi_err) |
577 |
|
|
if (mpi_err /= 0) then |
578 |
|
|
call error("WRITE_FINAL_CONFIG","Error receiving v_tmp") |
579 |
|
|
endif |
580 |
|
|
|
581 |
|
|
|
582 |
|
|
do i = 1,icnt |
583 |
|
|
do a = 1, na_tmp(i) |
584 |
|
|
atom = atom_index_tmp(i,a) |
585 |
|
|
atid = atype_identifier(ident_tmp(atom)) // ' ' |
586 |
|
|
write(eor_unit,fmt=xyz_coord) atid, q_tmp(1,atom), & |
587 |
|
|
q_tmp(2,atom), q_tmp(3,atom), v_tmp(1,atom), & |
588 |
|
|
v_tmp(2,atom), v_tmp(3,atom), ident_tmp(atom) |
589 |
|
|
end do |
590 |
|
|
enddo |
591 |
|
|
|
592 |
|
|
end do |
593 |
|
|
endif |
594 |
|
|
|
595 |
|
|
call mpi_barrier(mpi_comm_world,mpi_err) |
596 |
|
|
|
597 |
|
|
if (node == 0) then |
598 |
|
|
call machdep_flush(eor_unit) |
599 |
|
|
close(eor_unit) |
600 |
|
|
end if |
601 |
|
|
#endif |
602 |
|
|
end subroutine write_final_config |
603 |
|
|
|
604 |
|
|
|
605 |
|
|
subroutine dump_forces(this_time) |
606 |
|
|
! Reads In the current configuration in XYZ file format. The next |
607 |
|
|
! frame in XYZ is simply appended to the end of the current frame. |
608 |
|
|
|
609 |
|
|
! include 'headers/fileio.h' |
610 |
|
|
! include 'headers/atom.h' |
611 |
|
|
|
612 |
|
|
real( kind=DP ), intent(in) :: this_time |
613 |
|
|
integer ::i,j, iproc, a |
614 |
|
|
#ifdef MPI |
615 |
|
|
integer, dimension(max_local) :: tag_tmp |
616 |
|
|
#endif |
617 |
|
|
integer :: icnt |
618 |
|
|
character(Len=*), parameter :: forcefmt = '(I3,3es18.8)' |
619 |
|
|
|
620 |
|
|
#ifdef MPI |
621 |
|
|
if (node /= 0) then |
622 |
|
|
|
623 |
|
|
call mpi_send(tag_local,nlocal,mpi_integer,0,0, & |
624 |
|
|
mpi_comm_world,mpi_err) |
625 |
|
|
|
626 |
|
|
call mpi_send(f,3*nlocal,mpi_double_precision,0,0, & |
627 |
|
|
mpi_comm_world,mpi_err) |
628 |
|
|
|
629 |
|
|
call mpi_send(rho,nlocal,mpi_double_precision,0,0, & |
630 |
|
|
mpi_comm_world,mpi_err) |
631 |
|
|
|
632 |
|
|
! call mpi_send(frho,nlocal,mpi_double_precision,0,0, & |
633 |
|
|
! mpi_comm_world,mpi_err) |
634 |
|
|
|
635 |
|
|
else |
636 |
|
|
|
637 |
|
|
write(force_unit,'(es12.4)') this_time |
638 |
|
|
|
639 |
|
|
|
640 |
|
|
do i = 1, nlocal |
641 |
|
|
write(unit=force_unit,fmt=*) tag_local(i),rho(i),f(1,i), & |
642 |
|
|
f(2,i),f(3,i) |
643 |
|
|
end do |
644 |
|
|
|
645 |
|
|
do iproc = 1,nprocs - 1 |
646 |
|
|
icnt = 0 |
647 |
|
|
f_tmp = 0.0d0 |
648 |
|
|
|
649 |
|
|
|
650 |
|
|
call mpi_recv(tag_tmp,max_local,mpi_integer,iproc,0, & |
651 |
|
|
mpi_comm_world,istatus,mpi_err) |
652 |
|
|
call mpi_get_count(istatus,mpi_integer,icnt,mpi_err) |
653 |
|
|
call mpi_recv(f_tmp,3*icnt,mpi_double_precision,iproc,0, & |
654 |
|
|
mpi_comm_world,istatus,mpi_err) |
655 |
|
|
call mpi_recv(rho_tmp,icnt,mpi_double_precision,iproc,0, & |
656 |
|
|
mpi_comm_world,istatus,mpi_err) |
657 |
|
|
! call mpi_recv(frho_tmp,icnt,mpi_double_precision,iproc,0, & |
658 |
|
|
! mpi_comm_world,istatus,mpi_err) |
659 |
|
|
|
660 |
|
|
|
661 |
|
|
do i = 1,icnt |
662 |
|
|
write(unit=force_unit,fmt=*) tag_tmp(i),rho_tmp(i), & |
663 |
|
|
f_tmp(1,i),f_tmp(2,i),f_tmp(3,i) |
664 |
|
|
enddo |
665 |
|
|
|
666 |
|
|
end do |
667 |
|
|
endif |
668 |
|
|
|
669 |
|
|
call mpi_barrier(mpi_comm_world,mpi_err) |
670 |
|
|
|
671 |
|
|
if (node == 0) then |
672 |
|
|
call machdep_flush(dump_unit) |
673 |
|
|
end if |
674 |
|
|
#else |
675 |
|
|
write(force_unit,'(es12.4)') this_time |
676 |
|
|
do i = 1, natoms |
677 |
|
|
write(unit=force_unit,fmt=*) i,rho(i),f(1,i), & |
678 |
|
|
f(2,i),f(3,i) |
679 |
|
|
end do |
680 |
|
|
|
681 |
|
|
#endif |
682 |
|
|
|
683 |
|
|
end subroutine dump_forces |
684 |
|
|
|
685 |
|
|
|
686 |
|
|
!! remove me |
687 |
|
|
#ifdef MPI |
688 |
|
|
subroutine mpi_dump_q() |
689 |
|
|
integer, parameter :: q_row_unit = 110 |
690 |
|
|
integer, parameter :: q_col_unit = 111 |
691 |
|
|
integer :: processor |
692 |
|
|
integer :: i, err |
693 |
|
|
integer :: nrow_count,ncol_count |
694 |
|
|
character(len=60) :: q_row_file,q_col_file,q_file |
695 |
|
|
character(len=5) :: char_node |
696 |
|
|
real(kind=8), dimension(3,max_row) :: q_row_tmp |
697 |
|
|
real(kind=8), dimension(3,max_col) :: q_col_tmp |
698 |
|
|
integer, dimension(max_row) :: tag_row_tmp,ident_row_tmp |
699 |
|
|
integer, dimension(max_col) :: tag_col_tmp,ident_col_tmp |
700 |
|
|
integer :: nrow_tmp, ncol_tmp, nlocal_tmp |
701 |
|
|
|
702 |
|
|
if (node /= 0 ) then |
703 |
|
|
call mpi_send(nlocal,1,mpi_integer,0,0,mpi_comm_world,mpi_err) |
704 |
|
|
call mpi_send(nrow,1,mpi_integer,0,0,mpi_comm_world,mpi_err) |
705 |
|
|
call mpi_send(ncol,1,mpi_integer,0,0,mpi_comm_world,mpi_err) |
706 |
|
|
call mpi_send(ident_row,nrow,mpi_integer,0,0, & |
707 |
|
|
mpi_comm_world,mpi_err) |
708 |
|
|
call mpi_send(ident_col,ncol,mpi_integer,0,0, & |
709 |
|
|
mpi_comm_world,mpi_err) |
710 |
|
|
call mpi_send(q_row,3*nrow,mpi_double_precision,0,0, & |
711 |
|
|
mpi_comm_world,mpi_err) |
712 |
|
|
call mpi_send(q_col,3*ncol,mpi_double_precision,0,0, & |
713 |
|
|
mpi_comm_world,mpi_err) |
714 |
|
|
call mpi_send(tag_row,nrow,mpi_integer,0,0, & |
715 |
|
|
mpi_comm_world,mpi_err) |
716 |
|
|
call mpi_send(tag_col,ncol,mpi_integer,0,0, & |
717 |
|
|
mpi_comm_world,mpi_err) |
718 |
|
|
else |
719 |
|
|
|
720 |
|
|
!! node 0 writes its arrays first |
721 |
|
|
write(char_node,"(a1,i1,a1)") "_",node,"_" |
722 |
|
|
char_node = trim(char_node) |
723 |
|
|
q_row_file = file_prefix(1:machdep_lnblnk(file_prefix)) // & |
724 |
|
|
char_node(1:machdep_lnblnk(char_node)) // & |
725 |
|
|
'.row' |
726 |
|
|
open(file=q_row_file, unit=q_row_unit, status='replace', iostat= err, & |
727 |
|
|
action='write') |
728 |
|
|
write(unit=q_row_unit,fmt="(i5,i5)") node,nrow |
729 |
|
|
do i = 1, nrow |
730 |
|
|
write(unit=q_row_unit,fmt="(i5,3es12.4,i4)") tag_row(i),q_row(1,i), & |
731 |
|
|
q_row(2,i), q_row(3,i), ident_row(i) |
732 |
|
|
end do |
733 |
|
|
|
734 |
|
|
close(unit=q_row_unit) |
735 |
|
|
|
736 |
|
|
q_col_file = file_prefix(1:machdep_lnblnk(file_prefix)) // & |
737 |
|
|
char_node(1:machdep_lnblnk(char_node)) // & |
738 |
|
|
'.col' |
739 |
|
|
open(file=q_col_file, unit=q_col_unit, status='replace', iostat= err, & |
740 |
|
|
action='write') |
741 |
|
|
write(unit=q_col_unit,fmt="(i5,i5)") node,ncol |
742 |
|
|
do i = 1, ncol |
743 |
|
|
write(unit=q_col_unit,fmt="(i5,3es12.4,i4)") tag_col(i),q_col(1,i), q_col(2,i), & |
744 |
|
|
q_col(3,i), ident_col(i) |
745 |
|
|
end do |
746 |
|
|
|
747 |
|
|
close(unit=q_col_unit) |
748 |
|
|
write(*,*) node, nlocal, nrow, ncol |
749 |
|
|
do processor = 1, nprocs - 1 |
750 |
|
|
!! recv data from other procs |
751 |
|
|
call mpi_recv(nlocal_tmp,1,mpi_integer,processor,0,mpi_comm_world,& |
752 |
|
|
istatus,mpi_err) |
753 |
|
|
call mpi_recv(nrow_tmp,1,mpi_integer,processor,0,mpi_comm_world,& |
754 |
|
|
istatus,mpi_err) |
755 |
|
|
call mpi_recv(ncol_tmp,1,mpi_integer,processor,0,mpi_comm_world,& |
756 |
|
|
istatus,mpi_err) |
757 |
|
|
write(*,*) processor, nlocal_tmp, nrow_tmp, ncol_tmp |
758 |
|
|
call mpi_recv(ident_row_tmp,nrow_tmp,mpi_integer,processor,0, & |
759 |
|
|
mpi_comm_world,istatus,mpi_err) |
760 |
|
|
write(*,*) "Past ident row" |
761 |
|
|
! call mpi_get_count(istatus,mpi_integer, nrow_count,mpi_err) |
762 |
|
|
call mpi_recv(ident_col_tmp,ncol_tmp,mpi_integer,processor,0, & |
763 |
|
|
mpi_comm_world,istatus,mpi_err) |
764 |
|
|
! call mpi_get_count(istatus,mpi_integer,ncol_count,mpi_err) |
765 |
|
|
write(*,*) "Past ident_col " |
766 |
|
|
call mpi_recv(q_row_tmp,3*nrow_tmp,mpi_double_precision,processor,0,& |
767 |
|
|
mpi_comm_world,istatus,mpi_err) |
768 |
|
|
write(*,*) "Past q_row" |
769 |
|
|
call mpi_recv(q_col_tmp,3*ncol_tmp,mpi_double_precision,processor,0,& |
770 |
|
|
mpi_comm_world,istatus,mpi_err) |
771 |
|
|
write(*,*) "Past q_col" |
772 |
|
|
call mpi_recv(tag_row_tmp,nrow_tmp,mpi_integer,processor,0, & |
773 |
|
|
mpi_comm_world,istatus,mpi_err) |
774 |
|
|
write(*,*) "Past tag row" |
775 |
|
|
call mpi_recv(tag_col_tmp,ncol_tmp,mpi_integer,processor,0, & |
776 |
|
|
mpi_comm_world,istatus,mpi_err) |
777 |
|
|
write(*,*) "Past tag col" |
778 |
|
|
write(char_node,"(a1,i1,a1)") "_",processor,"_" |
779 |
|
|
char_node = trim(char_node) |
780 |
|
|
q_row_file = file_prefix(1:machdep_lnblnk(file_prefix)) // & |
781 |
|
|
char_node(1:machdep_lnblnk(char_node)) // '.row' |
782 |
|
|
open(file=q_row_file, unit=q_row_unit, status='replace', iostat= err, & |
783 |
|
|
action='write') |
784 |
|
|
write(unit=q_row_unit,fmt="(i5,i5)") node,nrow |
785 |
|
|
do i = 1, nrow_tmp |
786 |
|
|
write(unit=q_row_unit,fmt="(i5,3es12.4,i4)") tag_row_tmp(i),q_row_tmp(1,i), & |
787 |
|
|
q_row_tmp(2,i),q_row_tmp(3,i), ident_row_tmp(i) |
788 |
|
|
end do |
789 |
|
|
|
790 |
|
|
close(unit=q_row_unit) |
791 |
|
|
|
792 |
|
|
q_col_file = file_prefix(1:machdep_lnblnk(file_prefix)) // & |
793 |
|
|
char_node(1:machdep_lnblnk(char_node)) // & |
794 |
|
|
'.col' |
795 |
|
|
open(file=q_col_file, unit=q_col_unit, status='replace', iostat= err, & |
796 |
|
|
action='write') |
797 |
|
|
write(unit=q_col_unit,fmt="(i5,i5)") node,ncol |
798 |
|
|
do i = 1, ncol_tmp |
799 |
|
|
write(unit=q_col_unit,fmt="(i5,3es12.4,i4)") tag_col_tmp(i),q_col_tmp(1,i), & |
800 |
|
|
q_col_tmp(2,i), q_col_tmp(3,i), ident_col_tmp(i) |
801 |
|
|
end do |
802 |
|
|
|
803 |
|
|
close(unit=q_col_unit) |
804 |
|
|
enddo |
805 |
|
|
|
806 |
|
|
endif |
807 |
|
|
call mpi_barrier(mpi_comm_world,mpi_err) |
808 |
|
|
end subroutine mpi_dump_q |
809 |
|
|
|
810 |
|
|
#endif |
811 |
|
|
|
812 |
|
|
end module xyz_io |