(git:d5c4d39)
Loading...
Searching...
No Matches
realspace_grid_cube.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Generate Gaussian cube files
10! **************************************************************************************************
12 USE cp_files, ONLY: close_file,&
15 USE kinds, ONLY: dp
16 USE message_passing, ONLY: &
21 USE pw_types, ONLY: pw_r3d_rs_type
22#include "../base/base_uses.f90"
23
24 IMPLICIT NONE
25
26 PRIVATE
27
30
31 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'realspace_grid_cube'
32 INTEGER, PARAMETER, PRIVATE :: cube_entry_len = 13, &
33 cube_num_entries_line = 6
34 INTEGER, PARAMETER, PRIVATE :: cube_line_len = cube_entry_len*cube_num_entries_line
35 CHARACTER(len=*), PARAMETER :: cube_value_format = '(E13.5E3)', &
36 cube_values_format = '(6E13.5E3)'
37 LOGICAL, PARAMETER, PRIVATE :: debug_this_module = .false.
38
39CONTAINS
40
41! **************************************************************************************************
42!> \brief Read cube values from a character buffer.
43!> \param values value buffer from one cube line or one z-slice
44!> \param buffer parsed values
45! **************************************************************************************************
46 SUBROUTINE cube_read_values(values, buffer)
47 CHARACTER(LEN=*), INTENT(IN) :: values
48 REAL(kind=dp), DIMENSION(:), INTENT(OUT) :: buffer
49
50 CHARACTER(LEN=cube_entry_len) :: value
51 INTEGER :: i, pos, readstat
52
53 READ (values, *, iostat=readstat) buffer
54 IF (readstat == 0) RETURN
55
56 pos = 1
57 DO i = 1, SIZE(buffer)
58 IF (pos + cube_entry_len - 1 > len(values)) cpabort("Unexpected end of cube data.")
59 value = values(pos:pos + cube_entry_len - 1)
60 READ (value, '(E13.5)', iostat=readstat) buffer(i)
61 IF (readstat /= 0) cpabort("Bad value while reading cube data.")
62 pos = pos + cube_entry_len
63 IF (modulo(i, cube_num_entries_line) == 0) THEN
64 IF (pos <= len(values)) THEN
65 IF (values(pos:pos) == new_line('C')) pos = pos + 1
66 END IF
67 END IF
68 END DO
69
70 END SUBROUTINE cube_read_values
71
72! **************************************************************************************************
73!> \brief ...
74!> \param pw ...
75!> \param unit_nr ...
76!> \param title ...
77!> \param particles_r ...
78!> \param particles_z ...
79!> \param particles_zeff ...
80!> \param stride ...
81!> \param max_file_size_mb ...
82!> \param zero_tails ...
83!> \param silent ...
84!> \param mpi_io ...
85! **************************************************************************************************
86 SUBROUTINE pw_to_cube(pw, unit_nr, title, particles_r, particles_z, particles_zeff, &
87 stride, max_file_size_mb, zero_tails, silent, mpi_io)
88 TYPE(pw_r3d_rs_type), INTENT(IN) :: pw
89 INTEGER, INTENT(IN) :: unit_nr
90 CHARACTER(*), INTENT(IN), OPTIONAL :: title
91 REAL(kind=dp), DIMENSION(:, :), INTENT(IN), &
92 OPTIONAL :: particles_r
93 INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: particles_z
94 REAL(kind=dp), DIMENSION(:), INTENT(IN), OPTIONAL :: particles_zeff
95 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: stride
96 REAL(kind=dp), INTENT(IN), OPTIONAL :: max_file_size_mb
97 LOGICAL, INTENT(IN), OPTIONAL :: zero_tails, silent, mpi_io
98
99 CHARACTER(len=*), PARAMETER :: routinen = 'pw_to_cube'
100 INTEGER, PARAMETER :: entry_len = 13, num_entries_line = 6
101
102 INTEGER :: checksum, dest, handle, i, i1, i2, i3, iat, ip, l1, l2, l3, msglen, my_rank, &
103 my_stride(3), np, num_linebreak, num_pe, rank(2), size_of_z, source, tag, u1, u2, u3
104 LOGICAL :: be_silent, my_zero_tails, parallel_write
105 REAL(kind=dp) :: compression_factor, my_max_file_size_mb
106 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: buf
107 TYPE(mp_comm_type) :: gid
108 TYPE(mp_file_type) :: mp_unit
109
110 CALL timeset(routinen, handle)
111
112 my_zero_tails = .false.
113 be_silent = .false.
114 parallel_write = .false.
115 my_max_file_size_mb = 0.0_dp
116 IF (PRESENT(zero_tails)) my_zero_tails = zero_tails
117 IF (PRESENT(silent)) be_silent = silent
118 IF (PRESENT(mpi_io)) parallel_write = mpi_io
119 IF (PRESENT(max_file_size_mb)) my_max_file_size_mb = max_file_size_mb
120 cpassert(my_max_file_size_mb >= 0)
121
122 my_stride = 1
123 IF (PRESENT(stride)) THEN
124 IF (SIZE(stride) /= 1 .AND. SIZE(stride) /= 3) &
125 CALL cp_abort(__location__, "STRIDE keyword can accept only 1 "// &
126 "(the same for X,Y,Z) or 3 values. Correct your input file.")
127 IF (SIZE(stride) == 1) THEN
128 DO i = 1, 3
129 my_stride(i) = stride(1)
130 END DO
131 ELSE
132 my_stride = stride(1:3)
133 END IF
134 END IF
135
136 IF (my_max_file_size_mb > 0) THEN
137 ! A single grid point takes up 13 bytes, which is 1.3e-05 MB.
138 compression_factor = 1.3e-05_dp*product(real(pw%pw_grid%npts, dp))/max_file_size_mb
139 my_stride(:) = int(compression_factor**(1.0/3.0)) + 1
140 END IF
141
142 cpassert(my_stride(1) > 0)
143 cpassert(my_stride(2) > 0)
144 cpassert(my_stride(3) > 0)
145
146 IF (.NOT. parallel_write) THEN
147 IF (unit_nr > 0) THEN
148 ! this format seems to work for e.g. molekel and gOpenmol
149 ! latest version of VMD can read non orthorhombic cells
150 WRITE (unit_nr, '(a11)') "-Quickstep-"
151 IF (PRESENT(title)) THEN
152 WRITE (unit_nr, *) trim(title)
153 ELSE
154 WRITE (unit_nr, *) "No Title"
155 END IF
156
157 cpassert(PRESENT(particles_z) .EQV. PRESENT(particles_r))
158 np = 0
159 IF (PRESENT(particles_z)) THEN
160 cpassert(SIZE(particles_z) == SIZE(particles_r, dim=2))
161 ! cube files can only be written for 99999 particles due to a format limitation (I5)
162 ! so we limit the number of particles written.
163 np = min(99999, SIZE(particles_z))
164 END IF
165
166 WRITE (unit_nr, '(I5,3f12.6)') np, 0.0_dp, 0._dp, 0._dp !start of cube
167
168 WRITE (unit_nr, '(I5,3f12.6)') (pw%pw_grid%npts(1) + my_stride(1) - 1)/my_stride(1), &
169 pw%pw_grid%dh(1, 1)*real(my_stride(1), dp), pw%pw_grid%dh(2, 1)*real(my_stride(1), dp), &
170 pw%pw_grid%dh(3, 1)*real(my_stride(1), dp)
171 WRITE (unit_nr, '(I5,3f12.6)') (pw%pw_grid%npts(2) + my_stride(2) - 1)/my_stride(2), &
172 pw%pw_grid%dh(1, 2)*real(my_stride(2), dp), pw%pw_grid%dh(2, 2)*real(my_stride(2), dp), &
173 pw%pw_grid%dh(3, 2)*real(my_stride(2), dp)
174 WRITE (unit_nr, '(I5,3f12.6)') (pw%pw_grid%npts(3) + my_stride(3) - 1)/my_stride(3), &
175 pw%pw_grid%dh(1, 3)*real(my_stride(3), dp), pw%pw_grid%dh(2, 3)*real(my_stride(3), dp), &
176 pw%pw_grid%dh(3, 3)*real(my_stride(3), dp)
177
178 IF (PRESENT(particles_z)) THEN
179 IF (PRESENT(particles_zeff)) THEN
180 DO iat = 1, np
181 WRITE (unit_nr, '(I5,4f12.6)') particles_z(iat), particles_zeff(iat), particles_r(:, iat)
182 END DO
183 ELSE
184 DO iat = 1, np
185 WRITE (unit_nr, '(I5,4f12.6)') particles_z(iat), 0._dp, particles_r(:, iat)
186 END DO
187 END IF
188 END IF
189 END IF
190
191 ! shortcut
192 l1 = pw%pw_grid%bounds(1, 1)
193 l2 = pw%pw_grid%bounds(1, 2)
194 l3 = pw%pw_grid%bounds(1, 3)
195 u1 = pw%pw_grid%bounds(2, 1)
196 u2 = pw%pw_grid%bounds(2, 2)
197 u3 = pw%pw_grid%bounds(2, 3)
198
199 ALLOCATE (buf(l3:u3))
200
201 my_rank = pw%pw_grid%para%group%mepos
202 gid = pw%pw_grid%para%group
203 num_pe = pw%pw_grid%para%group%num_pe
204 tag = 1
205
206 rank(1) = unit_nr
207 rank(2) = my_rank
208 checksum = 0
209 IF (unit_nr > 0) checksum = 1
210
211 CALL gid%sum(checksum)
212 cpassert(checksum == 1)
213
214 CALL gid%maxloc(rank)
215 cpassert(rank(1) > 0)
216
217 dest = rank(2)
218 DO i1 = l1, u1, my_stride(1)
219 DO i2 = l2, u2, my_stride(2)
220
221 ! cycling through the CPUs, check if the current ray (I1,I2) is local to that CPU
222 IF (pw%pw_grid%para%mode /= pw_mode_local) THEN
223 DO ip = 0, num_pe - 1
224 IF (pw%pw_grid%para%bo(1, 1, ip, 1) <= i1 - l1 + 1 .AND. pw%pw_grid%para%bo(2, 1, ip, 1) >= i1 - l1 + 1 .AND. &
225 pw%pw_grid%para%bo(1, 2, ip, 1) <= i2 - l2 + 1 .AND. pw%pw_grid%para%bo(2, 2, ip, 1) >= i2 - l2 + 1) THEN
226 source = ip
227 END IF
228 END DO
229 ELSE
230 source = dest
231 END IF
232
233 IF (source == dest) THEN
234 IF (my_rank == source) THEN
235 buf(:) = pw%array(i1, i2, :)
236 END IF
237 ELSE
238 IF (my_rank == source) THEN
239 buf(:) = pw%array(i1, i2, :)
240 CALL gid%send(buf, dest, tag)
241 END IF
242 IF (my_rank == dest) THEN
243 CALL gid%recv(buf, source, tag)
244 END IF
245 END IF
246
247 IF (unit_nr > 0) THEN
248 IF (my_zero_tails) THEN
249 DO i3 = l3, u3
250 IF (buf(i3) < 1.e-7_dp) buf(i3) = 0.0_dp
251 END DO
252 END IF
253 WRITE (unit_nr, cube_values_format) (buf(i3), i3=l3, u3, my_stride(3))
254 END IF
255
256 ! this double loop generates so many messages that it can overload
257 ! the message passing system, e.g. on XT3
258 ! we therefore put a barrier here that limits the amount of message
259 ! that flies around at any given time.
260 ! if ever this routine becomes a bottleneck, we should go for a
261 ! more complicated rewrite
262 CALL gid%sync()
263
264 END DO
265 END DO
266
267 DEALLOCATE (buf)
268 ELSE
269 size_of_z = ceiling(real(pw%pw_grid%bounds(2, 3) - pw%pw_grid%bounds(1, 3) + 1, dp)/real(my_stride(3), dp))
270 num_linebreak = size_of_z/num_entries_line
271 IF (modulo(size_of_z, num_entries_line) /= 0) &
272 num_linebreak = num_linebreak + 1
273 msglen = (size_of_z*entry_len + num_linebreak)*mpi_character_size
274 CALL mp_unit%set_handle(unit_nr)
275 CALL pw_to_cube_parallel(pw, mp_unit, title, particles_r, particles_z, particles_zeff, &
276 my_stride, my_zero_tails, msglen)
277 END IF
278
279 CALL timestop(handle)
280
281 END SUBROUTINE pw_to_cube
282
283! **************************************************************************************************
284!> \brief Computes the external density on the grid
285!> hacked from external_read_density
286!> \param grid pw to read from cube file
287!> \param filename name of cube file
288!> \param scaling scale values before storing
289!> \param parallel_read ...
290!> \param silent ...
291!> \par History
292!> Created [M.Watkins] (01.2014)
293!> Use blocking, collective MPI read for parallel simulations [Nico Holmberg] (05.2017)
294! **************************************************************************************************
295 SUBROUTINE cube_to_pw(grid, filename, scaling, parallel_read, silent)
296
297 TYPE(pw_r3d_rs_type), INTENT(IN) :: grid
298 CHARACTER(len=*), INTENT(in) :: filename
299 REAL(kind=dp), INTENT(in) :: scaling
300 LOGICAL, INTENT(in) :: parallel_read
301 LOGICAL, INTENT(in), OPTIONAL :: silent
302
303 CHARACTER(len=*), PARAMETER :: routinen = 'cube_to_pw'
304 INTEGER, PARAMETER :: entry_len = 13, num_entries_line = 6
305
306 CHARACTER(LEN=cube_line_len) :: value_line
307 INTEGER :: extunit, handle, i, j, k, last_z, &
308 msglen, my_rank, nat, ndum, &
309 num_linebreak, num_pe, output_unit, &
310 size_of_z, tag
311 INTEGER, DIMENSION(3) :: lbounds, lbounds_local, npoints, &
312 npoints_local, ubounds, ubounds_local
313 LOGICAL :: be_silent
314 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: buffer
315 REAL(kind=dp), DIMENSION(3) :: dr, rdum
316 TYPE(mp_comm_type) :: gid
317
318 output_unit = cp_logger_get_default_io_unit()
319
320 CALL timeset(routinen, handle)
321
322 be_silent = .false.
323 IF (PRESENT(silent)) THEN
324 be_silent = silent
325 END IF
326 !get rs grids and parallel environment
327 gid = grid%pw_grid%para%group
328 my_rank = grid%pw_grid%para%group%mepos
329 num_pe = grid%pw_grid%para%group%num_pe
330 tag = 1
331
332 lbounds_local = grid%pw_grid%bounds_local(1, :)
333 ubounds_local = grid%pw_grid%bounds_local(2, :)
334 size_of_z = ubounds_local(3) - lbounds_local(3) + 1
335
336 IF (.NOT. parallel_read) THEN
337 npoints = grid%pw_grid%npts
338 lbounds = grid%pw_grid%bounds(1, :)
339 ubounds = grid%pw_grid%bounds(2, :)
340
341 DO i = 1, 3
342 dr(i) = grid%pw_grid%dh(i, i)
343 END DO
344
345 npoints_local = grid%pw_grid%npts_local
346 !pw grids at most pencils - all processors have a full set of z data for x,y
347 ALLOCATE (buffer(lbounds(3):ubounds(3)))
348
349 IF (my_rank == 0) THEN
350 IF (output_unit > 0 .AND. .NOT. be_silent) THEN
351 WRITE (output_unit, fmt="(/,T2,A,/,/,T2,A,/)") "Reading the cube file: ", trim(filename)
352 END IF
353
354 CALL open_file(file_name=filename, &
355 file_status="OLD", &
356 file_form="FORMATTED", &
357 file_action="READ", &
358 unit_number=extunit)
359
360 !skip header comments
361 DO i = 1, 2
362 READ (extunit, *)
363 END DO
364 READ (extunit, *) nat, rdum
365 DO i = 1, 3
366 READ (extunit, *) ndum, rdum
367 IF ((ndum /= npoints(i) .OR. (abs(rdum(i) - dr(i)) > 1e-4)) .AND. &
368 output_unit > 0) THEN
369 WRITE (output_unit, *) "Restart from density | ERROR! | CUBE FILE NOT COINCIDENT WITH INTERNAL GRID ", i
370 WRITE (output_unit, *) "Restart from density | ", ndum, " DIFFERS FROM ", npoints(i)
371 WRITE (output_unit, *) "Restart from density | ", rdum, " DIFFERS FROM ", dr(i)
372 END IF
373 END DO
374 !ignore atomic position data - read from coord or topology instead
375 DO i = 1, nat
376 READ (extunit, *)
377 END DO
378 END IF
379
380 !master sends all data to everyone
381 DO i = lbounds(1), ubounds(1)
382 DO j = lbounds(2), ubounds(2)
383 IF (my_rank == 0) THEN
384 DO k = lbounds(3), ubounds(3), cube_num_entries_line
385 last_z = min(k + cube_num_entries_line - 1, ubounds(3))
386 READ (extunit, '(A)') value_line
387 CALL cube_read_values(value_line, buffer(k:last_z))
388 END DO
389 END IF
390 CALL gid%bcast(buffer(lbounds(3):ubounds(3)), 0)
391
392 !only use data that is local to me - i.e. in slice of pencil I own
393 IF ((lbounds_local(1) <= i) .AND. (i <= ubounds_local(1)) .AND. (lbounds_local(2) <= j) &
394 .AND. (j <= ubounds_local(2))) THEN
395 !allow scaling of external potential values by factor 'scaling' (SCALING_FACTOR in input file)
396 grid%array(i, j, lbounds(3):ubounds(3)) = buffer(lbounds(3):ubounds(3))*scaling
397 END IF
398
399 END DO
400 END DO
401
402 IF (my_rank == 0) CALL close_file(unit_number=extunit)
403
404 CALL gid%sync()
405 ELSE
406 ! Parallel routine needs as input the byte size of each grid z-slice
407 ! This is a hack to prevent compilation errors with gcc -Wall (up to versions 6.3)
408 ! related to allocatable-length string declaration CHARACTER(LEN=:), ALLOCATABLE, DIMENSION(:) :: string
409 ! Each data line of a Gaussian cube contains max 6 entries with last line potentially containing less if nz % 6 /= 0
410 ! Thus, this size is simply the number of entries multiplied by the entry size + the number of line breaks
411 num_linebreak = size_of_z/num_entries_line
412 IF (modulo(size_of_z, num_entries_line) /= 0) &
413 num_linebreak = num_linebreak + 1
414 msglen = (size_of_z*entry_len + num_linebreak)*mpi_character_size
415 CALL cube_to_pw_parallel(grid, filename, scaling, msglen, silent=silent)
416 END IF
417
418 CALL timestop(handle)
419
420 END SUBROUTINE cube_to_pw
421
422! **************************************************************************************************
423!> \brief Reads a realspace potential/density from a cube file using collective MPI I/O and
424!> stores it in grid.
425!> \param grid pw to read from cube file
426!> \param filename name of cube file
427!> \param scaling scale values before storing
428!> \param msglen the size of each grid slice along z-axis in bytes
429!> \param silent ...
430!> \par History
431!> Created [Nico Holmberg] (05.2017)
432! **************************************************************************************************
433 SUBROUTINE cube_to_pw_parallel(grid, filename, scaling, msglen, silent)
434
435 TYPE(pw_r3d_rs_type), INTENT(IN) :: grid
436 CHARACTER(len=*), INTENT(in) :: filename
437 REAL(kind=dp), INTENT(in) :: scaling
438 INTEGER, INTENT(in) :: msglen
439 LOGICAL, INTENT(in), OPTIONAL :: silent
440
441 CHARACTER(LEN=cube_line_len) :: value_line
442 INTEGER, DIMENSION(3) :: lbounds, lbounds_local, npoints, &
443 npoints_local, ubounds, ubounds_local
444 INTEGER, ALLOCATABLE, DIMENSION(:), TARGET :: blocklengths
445 INTEGER(kind=file_offset), ALLOCATABLE, &
446 DIMENSION(:), TARGET :: displacements
447 INTEGER(kind=file_offset) :: bof
448 INTEGER :: extunit_handle, i, islice, j, k, last_z, &
449 my_rank, nat, ndum, nslices, num_pe, &
450 offset_global, output_unit, size_of_z, &
451 tag
452 CHARACTER(LEN=msglen), ALLOCATABLE, DIMENSION(:) :: readbuffer
453 LOGICAL :: be_silent, should_read(2)
454 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: buffer
455 REAL(kind=dp), DIMENSION(3) :: dr, rdum
456 TYPE(mp_comm_type) :: gid
457 TYPE(mp_file_descriptor_type) :: mp_file_desc
458 TYPE(mp_file_type) :: extunit
459
460 output_unit = cp_logger_get_default_io_unit()
461
462 be_silent = .false.
463 IF (PRESENT(silent)) THEN
464 be_silent = silent
465 END IF
466
467 !get rs grids and parallel envnment
468 gid = grid%pw_grid%para%group
469 my_rank = grid%pw_grid%para%group%mepos
470 num_pe = grid%pw_grid%para%group%num_pe
471 tag = 1
472
473 DO i = 1, 3
474 dr(i) = grid%pw_grid%dh(i, i)
475 END DO
476
477 npoints = grid%pw_grid%npts
478 lbounds = grid%pw_grid%bounds(1, :)
479 ubounds = grid%pw_grid%bounds(2, :)
480
481 npoints_local = grid%pw_grid%npts_local
482 lbounds_local = grid%pw_grid%bounds_local(1, :)
483 ubounds_local = grid%pw_grid%bounds_local(2, :)
484 size_of_z = ubounds_local(3) - lbounds_local(3) + 1
485 nslices = (ubounds_local(1) - lbounds_local(1) + 1)*(ubounds_local(2) - lbounds_local(2) + 1)
486 islice = 1
487
488 ! Read header information and determine byte offset of cube data on master process
489 IF (my_rank == 0) THEN
490 IF (output_unit > 0 .AND. .NOT. be_silent) THEN
491 WRITE (output_unit, fmt="(/,T2,A,/,/,T2,A,/)") "Reading the cube file: ", trim(filename)
492 END IF
493
494 CALL open_file(file_name=filename, &
495 file_status="OLD", &
496 file_form="FORMATTED", &
497 file_action="READ", &
498 file_access="STREAM", &
499 unit_number=extunit_handle)
500
501 !skip header comments
502 DO i = 1, 2
503 READ (extunit_handle, *)
504 END DO
505 READ (extunit_handle, *) nat, rdum
506 DO i = 1, 3
507 READ (extunit_handle, *) ndum, rdum
508 IF ((ndum /= npoints(i) .OR. (abs(rdum(i) - dr(i)) > 1e-4)) .AND. &
509 output_unit > 0) THEN
510 WRITE (output_unit, *) "Restart from density | ERROR! | CUBE FILE NOT COINCIDENT WITH INTERNAL GRID ", i
511 WRITE (output_unit, *) "Restart from density | ", ndum, " DIFFERS FROM ", npoints(i)
512 WRITE (output_unit, *) "Restart from density | ", rdum, " DIFFERS FROM ", dr(i)
513 END IF
514 END DO
515 !ignore atomic position data - read from coord or topology instead
516 DO i = 1, nat
517 READ (extunit_handle, *)
518 END DO
519 ! Get byte offset
520 INQUIRE (extunit_handle, pos=offset_global)
521 CALL close_file(unit_number=extunit_handle)
522 END IF
523 ! Sync offset and start parallel read
524 CALL gid%bcast(offset_global, grid%pw_grid%para%group%source)
525 ! INQUIRE(POS=...) returns a 1-based stream position, whereas MPI-IO
526 ! file offsets are 0-based.
527 bof = offset_global - 1
528 CALL extunit%open(groupid=gid, filepath=filename, amode_status=file_amode_rdonly)
529 ! Determine byte offsets for each grid z-slice which are local to a process
530 ALLOCATE (displacements(nslices))
531 displacements = 0
532 DO i = lbounds(1), ubounds(1)
533 should_read(:) = .true.
534 IF (i < lbounds_local(1)) THEN
535 should_read(1) = .false.
536 ELSE IF (i > ubounds_local(1)) THEN
537 EXIT
538 END IF
539 DO j = lbounds(2), ubounds(2)
540 should_read(2) = .true.
541 IF (j < lbounds_local(2) .OR. j > ubounds_local(2)) THEN
542 should_read(2) = .false.
543 END IF
544 IF (all(should_read .EQV. .true.)) THEN
545 IF (islice > nslices) cpabort("Index out of bounds.")
546 displacements(islice) = bof
547 islice = islice + 1
548 END IF
549 ! Update global byte offset
550 bof = bof + msglen
551 END DO
552 END DO
553 ! Size of each z-slice is msglen
554 ALLOCATE (blocklengths(nslices))
555 blocklengths(:) = msglen
556 ! Create indexed MPI type using calculated byte offsets as displacements and use it as a file view
557 mp_file_desc = mp_file_type_hindexed_make_chv(nslices, blocklengths, displacements)
558 bof = 0
559 CALL mp_file_type_set_view_chv(extunit, bof, mp_file_desc)
560 ! Collective read of cube
561 ALLOCATE (readbuffer(nslices))
562 readbuffer(:) = ''
563 CALL extunit%read_all(msglen, nslices, readbuffer, mp_file_desc)
564 CALL mp_file_type_free(mp_file_desc)
565 CALL extunit%close()
566 ! Convert cube values string -> real
567 i = lbounds_local(1)
568 j = lbounds_local(2)
569 ALLOCATE (buffer(lbounds(3):ubounds(3)))
570 buffer = 0.0_dp
571 DO islice = 1, nslices
572 CALL cube_read_values(readbuffer(islice), buffer(lbounds(3):ubounds(3)))
573 ! Optionally scale cube file values
574 grid%array(i, j, lbounds(3):ubounds(3)) = scaling*buffer(lbounds(3):ubounds(3))
575 j = j + 1
576 IF (j > ubounds_local(2)) THEN
577 j = lbounds_local(2)
578 i = i + 1
579 END IF
580 END DO
581 DEALLOCATE (readbuffer)
582 DEALLOCATE (blocklengths, displacements)
583 IF (debug_this_module) THEN
584 ! Check that cube was correctly read using intrinsic read on master who sends data to everyone
585 buffer = 0.0_dp
586 IF (my_rank == 0) THEN
587 IF (output_unit > 0 .AND. .NOT. be_silent) THEN
588 WRITE (output_unit, fmt="(/,T2,A,/,/,T2,A)") "Reading the cube file: ", filename
589 END IF
590
591 CALL open_file(file_name=filename, &
592 file_status="OLD", &
593 file_form="FORMATTED", &
594 file_action="READ", &
595 unit_number=extunit_handle)
596
597 !skip header comments
598 DO i = 1, 2
599 READ (extunit_handle, *)
600 END DO
601 READ (extunit_handle, *) nat, rdum
602 DO i = 1, 3
603 READ (extunit_handle, *) ndum, rdum
604 IF ((ndum /= npoints(i) .OR. (abs(rdum(i) - dr(i)) > 1e-4)) .AND. &
605 output_unit > 0) THEN
606 WRITE (output_unit, *) "Restart from density | ERROR! | CUBE FILE NOT COINCIDENT WITH INTERNAL GRID ", i
607 WRITE (output_unit, *) "Restart from density | ", ndum, " DIFFERS FROM ", npoints(i)
608 WRITE (output_unit, *) "Restart from density | ", rdum, " DIFFERS FROM ", dr(i)
609 END IF
610 END DO
611 !ignore atomic position data - read from coord or topology instead
612 DO i = 1, nat
613 READ (extunit_handle, *)
614 END DO
615 END IF
616
617 !master sends all data to everyone
618 DO i = lbounds(1), ubounds(1)
619 DO j = lbounds(2), ubounds(2)
620 IF (my_rank == 0) THEN
621 DO k = lbounds(3), ubounds(3), cube_num_entries_line
622 last_z = min(k + cube_num_entries_line - 1, ubounds(3))
623 READ (extunit_handle, '(A)') value_line
624 CALL cube_read_values(value_line, buffer(k:last_z))
625 END DO
626 END IF
627 CALL gid%bcast(buffer(lbounds(3):ubounds(3)), 0)
628
629 !only use data that is local to me - i.e. in slice of pencil I own
630 IF ((lbounds_local(1) <= i) .AND. (i <= ubounds_local(1)) .AND. (lbounds_local(2) <= j) &
631 .AND. (j <= ubounds_local(2))) THEN
632 !allow scaling of external potential values by factor 'scaling' (SCALING_FACTOR in input file)
633 IF (any(grid%array(i, j, lbounds(3):ubounds(3)) /= buffer(lbounds(3):ubounds(3))*scaling)) &
634 CALL cp_abort(__location__, &
635 "Error in parallel read of input cube file.")
636 END IF
637
638 END DO
639 END DO
640
641 IF (my_rank == 0) CALL close_file(unit_number=extunit_handle)
642
643 CALL gid%sync()
644 END IF
645 DEALLOCATE (buffer)
646
647 END SUBROUTINE cube_to_pw_parallel
648
649! **************************************************************************************************
650!> \brief Writes a realspace potential to a cube file using collective MPI I/O.
651!> \param grid the pw to output to the cube file
652!> \param unit_nr the handle associated with the cube file
653!> \param title title of the cube file
654!> \param particles_r Cartersian coordinates of the system
655!> \param particles_z atomic masses of atoms in the system
656!> \param particles_zeff effective atomic charges of atoms in the system
657!> \param stride every stride(i)th value of the potential is outputted (i=x,y,z)
658!> \param zero_tails flag that determines if small values of the potential should be zeroed
659!> \param msglen the size of each grid slice along z-axis in bytes
660!> \par History
661!> Created [Nico Holmberg] (11.2017)
662! **************************************************************************************************
663 SUBROUTINE pw_to_cube_parallel(grid, unit_nr, title, particles_r, particles_z, particles_zeff, &
664 stride, zero_tails, msglen)
665
666 TYPE(pw_r3d_rs_type), INTENT(IN) :: grid
667 TYPE(mp_file_type), INTENT(IN) :: unit_nr
668 CHARACTER(*), INTENT(IN), OPTIONAL :: title
669 REAL(kind=dp), DIMENSION(:, :), INTENT(IN), &
670 OPTIONAL :: particles_r
671 INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: particles_z
672 REAL(kind=dp), DIMENSION(:), INTENT(IN), OPTIONAL :: particles_zeff
673 INTEGER, INTENT(IN) :: stride(3)
674 LOGICAL, INTENT(IN) :: zero_tails
675 INTEGER, INTENT(IN) :: msglen
676
677 INTEGER, PARAMETER :: entry_len = 13, header_len = 41, &
678 header_len_z = 53, num_entries_line = 6
679
680 CHARACTER(LEN=entry_len) :: value
681 CHARACTER(LEN=header_len) :: header
682 CHARACTER(LEN=header_len_z) :: header_z
683 INTEGER, DIMENSION(3) :: lbounds, lbounds_local, ubounds, &
684 ubounds_local
685 INTEGER, ALLOCATABLE, DIMENSION(:), TARGET :: blocklengths
686 INTEGER(kind=file_offset), ALLOCATABLE, &
687 DIMENSION(:), TARGET :: displacements
688 INTEGER(kind=file_offset) :: bof
689 INTEGER :: counter, i, islice, j, k, last_z, &
690 my_rank, np, nslices, size_of_z
691 CHARACTER(LEN=msglen), ALLOCATABLE, DIMENSION(:) :: writebuffer
692 CHARACTER(LEN=msglen) :: tmp
693 LOGICAL :: should_write(2)
694 TYPE(mp_comm_type) :: gid
695 TYPE(mp_file_descriptor_type) :: mp_desc
696
697 !get rs grids and parallel envnment
698 gid = grid%pw_grid%para%group
699 my_rank = grid%pw_grid%para%group%mepos
700
701 ! Shortcut
702 lbounds = grid%pw_grid%bounds(1, :)
703 ubounds = grid%pw_grid%bounds(2, :)
704 lbounds_local = grid%pw_grid%bounds_local(1, :)
705 ubounds_local = grid%pw_grid%bounds_local(2, :)
706 ! Determine the total number of z-slices and the number of values per slice
707 size_of_z = ceiling(real(ubounds_local(3) - lbounds_local(3) + 1, dp)/real(stride(3), dp))
708 islice = 1
709 DO i = lbounds(1), ubounds(1), stride(1)
710 should_write(:) = .true.
711 IF (i < lbounds_local(1)) THEN
712 should_write(1) = .false.
713 ELSE IF (i > ubounds_local(1)) THEN
714 EXIT
715 END IF
716 DO j = lbounds(2), ubounds(2), stride(2)
717 should_write(2) = .true.
718 IF (j < lbounds_local(2) .OR. j > ubounds_local(2)) THEN
719 should_write(2) = .false.
720 END IF
721 IF (all(should_write .EQV. .true.)) THEN
722 islice = islice + 1
723 END IF
724 END DO
725 END DO
726 nslices = islice - 1
727 DO k = lbounds(3), ubounds(3), stride(3)
728 IF (k + stride(3) > ubounds(3)) last_z = k
729 END DO
730 islice = 1
731 ! Determine initial byte offset (0 or EOF if data is appended)
732 CALL unit_nr%get_position(bof)
733 ! Write header information on master process and update byte offset accordingly
734 IF (my_rank == 0) THEN
735 ! this format seems to work for e.g. molekel and gOpenmol
736 ! latest version of VMD can read non orthorhombic cells
737 CALL unit_nr%write_at(bof, "-Quickstep-"//new_line("C"))
738 bof = bof + len("-Quickstep-"//new_line("C"))*mpi_character_size
739 IF (PRESENT(title)) THEN
740 CALL unit_nr%write_at(bof, trim(title)//new_line("C"))
741 bof = bof + len(trim(title)//new_line("C"))*mpi_character_size
742 ELSE
743 CALL unit_nr%write_at(bof, "No Title"//new_line("C"))
744 bof = bof + len("No Title"//new_line("C"))*mpi_character_size
745 END IF
746
747 cpassert(PRESENT(particles_z) .EQV. PRESENT(particles_r))
748 np = 0
749 IF (PRESENT(particles_z)) THEN
750 cpassert(SIZE(particles_z) == SIZE(particles_r, dim=2))
751 ! cube files can only be written for 99999 particles due to a format limitation (I5)
752 ! so we limit the number of particles written.
753 np = min(99999, SIZE(particles_z))
754 END IF
755
756 WRITE (header, '(I5,3f12.6)') np, 0.0_dp, 0._dp, 0._dp !start of cube
757 CALL unit_nr%write_at(bof, header//new_line("C"))
758 bof = bof + len(header//new_line("C"))*mpi_character_size
759
760 WRITE (header, '(I5,3f12.6)') (grid%pw_grid%npts(1) + stride(1) - 1)/stride(1), &
761 grid%pw_grid%dh(1, 1)*real(stride(1), dp), grid%pw_grid%dh(2, 1)*real(stride(1), dp), &
762 grid%pw_grid%dh(3, 1)*real(stride(1), dp)
763 CALL unit_nr%write_at(bof, header//new_line("C"))
764 bof = bof + len(header//new_line("C"))*mpi_character_size
765
766 WRITE (header, '(I5,3f12.6)') (grid%pw_grid%npts(2) + stride(2) - 1)/stride(2), &
767 grid%pw_grid%dh(1, 2)*real(stride(2), dp), grid%pw_grid%dh(2, 2)*real(stride(2), dp), &
768 grid%pw_grid%dh(3, 2)*real(stride(2), dp)
769 CALL unit_nr%write_at(bof, header//new_line("C"))
770 bof = bof + len(header//new_line("C"))*mpi_character_size
771
772 WRITE (header, '(I5,3f12.6)') (grid%pw_grid%npts(3) + stride(3) - 1)/stride(3), &
773 grid%pw_grid%dh(1, 3)*real(stride(3), dp), grid%pw_grid%dh(2, 3)*real(stride(3), dp), &
774 grid%pw_grid%dh(3, 3)*real(stride(3), dp)
775 CALL unit_nr%write_at(bof, header//new_line("C"))
776 bof = bof + len(header//new_line("C"))*mpi_character_size
777
778 IF (PRESENT(particles_z)) THEN
779 IF (PRESENT(particles_zeff)) THEN
780 DO i = 1, np
781 WRITE (header_z, '(I5,4f12.6)') particles_z(i), particles_zeff(i), particles_r(:, i)
782 CALL unit_nr%write_at(bof, header_z//new_line("C"))
783 bof = bof + len(header_z//new_line("C"))*mpi_character_size
784 END DO
785 ELSE
786 DO i = 1, np
787 WRITE (header_z, '(I5,4f12.6)') particles_z(i), 0._dp, particles_r(:, i)
788 CALL unit_nr%write_at(bof, header_z//new_line("C"))
789 bof = bof + len(header_z//new_line("C"))*mpi_character_size
790 END DO
791 END IF
792 END IF
793 END IF
794 ! Sync offset
795 CALL gid%bcast(bof, grid%pw_grid%para%group%source)
796 ! Determine byte offsets for each grid z-slice which are local to a process
797 ! and convert z-slices to cube format compatible strings
798 ALLOCATE (displacements(nslices))
799 displacements = 0
800 ALLOCATE (writebuffer(nslices))
801 writebuffer(:) = ''
802 DO i = lbounds(1), ubounds(1), stride(1)
803 should_write(:) = .true.
804 IF (i < lbounds_local(1)) THEN
805 should_write(1) = .false.
806 ELSE IF (i > ubounds_local(1)) THEN
807 EXIT
808 END IF
809 DO j = lbounds(2), ubounds(2), stride(2)
810 should_write(2) = .true.
811 IF (j < lbounds_local(2) .OR. j > ubounds_local(2)) THEN
812 should_write(2) = .false.
813 END IF
814 IF (all(should_write .EQV. .true.)) THEN
815 IF (islice > nslices) cpabort("Index out of bounds.")
816 displacements(islice) = bof
817 tmp = ''
818 counter = 0
819 DO k = lbounds(3), ubounds(3), stride(3)
820 IF (zero_tails .AND. grid%array(i, j, k) < 1.e-7_dp) THEN
821 WRITE (value, cube_value_format) 0.0_dp
822 ELSE
823 WRITE (value, cube_value_format) grid%array(i, j, k)
824 END IF
825 tmp = trim(tmp)//trim(value)
826 counter = counter + 1
827 IF (modulo(counter, num_entries_line) == 0 .OR. k == last_z) &
828 tmp = trim(tmp)//new_line('C')
829 END DO
830 writebuffer(islice) = tmp
831 islice = islice + 1
832 END IF
833 ! Update global byte offset
834 bof = bof + msglen
835 END DO
836 END DO
837 ! Create indexed MPI type using calculated byte offsets as displacements
838 ! Size of each z-slice is msglen
839 ALLOCATE (blocklengths(nslices))
840 blocklengths(:) = msglen
841 mp_desc = mp_file_type_hindexed_make_chv(nslices, blocklengths, displacements)
842 ! Use the created type as a file view
843 ! NB. The vector 'displacements' contains the absolute offsets of each z-slice i.e.
844 ! they are given relative to the beginning of the file. The global offset to
845 ! set_view must therefore be set to 0
846 bof = 0
847 CALL mp_file_type_set_view_chv(unit_nr, bof, mp_desc)
848 ! Collective write of cube
849 CALL unit_nr%write_all(msglen, nslices, writebuffer, mp_desc)
850 ! Clean up
851 CALL mp_file_type_free(mp_desc)
852 DEALLOCATE (writebuffer)
853 DEALLOCATE (blocklengths, displacements)
854
855 END SUBROUTINE pw_to_cube_parallel
856
857! **************************************************************************************************
858!> \brief Prints a simple grid file: X Y Z value
859!> \param pw ...
860!> \param unit_nr ...
861!> \param stride ...
862!> \param pw2 ...
863!> \par History
864!> Created [Vladimir Rybkin] (08.2018)
865!> \author Vladimir Rybkin
866! **************************************************************************************************
867 SUBROUTINE pw_to_simple_volumetric(pw, unit_nr, stride, pw2)
868 TYPE(pw_r3d_rs_type), INTENT(IN) :: pw
869 INTEGER, INTENT(IN) :: unit_nr
870 INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: stride
871 TYPE(pw_r3d_rs_type), INTENT(IN), OPTIONAL :: pw2
872
873 CHARACTER(len=*), PARAMETER :: routinen = 'pw_to_simple_volumetric'
874
875 INTEGER :: checksum, dest, handle, i, i1, i2, i3, &
876 ip, l1, l2, l3, my_rank, my_stride(3), &
877 ngrids, npoints, num_pe, rank(2), &
878 source, tag, u1, u2, u3
879 LOGICAL :: double
880 REAL(kind=dp) :: x, y, z
881 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: buf, buf2
882 TYPE(mp_comm_type) :: gid
883
884 CALL timeset(routinen, handle)
885
886 ! Check if we write two grids
887 double = .false.
888 IF (PRESENT(pw2)) double = .true.
889
890 my_stride = 1
891 IF (PRESENT(stride)) THEN
892 IF (SIZE(stride) /= 1 .AND. SIZE(stride) /= 3) &
893 CALL cp_abort(__location__, "STRIDE keyword can accept only 1 "// &
894 "(the same for X,Y,Z) or 3 values. Correct your input file.")
895 IF (SIZE(stride) == 1) THEN
896 DO i = 1, 3
897 my_stride(i) = stride(1)
898 END DO
899 ELSE
900 my_stride = stride(1:3)
901 END IF
902 cpassert(my_stride(1) > 0)
903 cpassert(my_stride(2) > 0)
904 cpassert(my_stride(3) > 0)
905 END IF
906
907 ! shortcut
908 l1 = pw%pw_grid%bounds(1, 1)
909 l2 = pw%pw_grid%bounds(1, 2)
910 l3 = pw%pw_grid%bounds(1, 3)
911 u1 = pw%pw_grid%bounds(2, 1)
912 u2 = pw%pw_grid%bounds(2, 2)
913 u3 = pw%pw_grid%bounds(2, 3)
914
915 ! Write the header: number of points and number of spins
916 ngrids = 1
917 IF (double) ngrids = 2
918 npoints = ((pw%pw_grid%npts(1) + my_stride(1) - 1)/my_stride(1))* &
919 ((pw%pw_grid%npts(2) + my_stride(2) - 1)/my_stride(1))* &
920 ((pw%pw_grid%npts(3) + my_stride(3) - 1)/my_stride(1))
921 IF (unit_nr > 1) WRITE (unit_nr, '(I7,I5)') npoints, ngrids
922
923 ALLOCATE (buf(l3:u3))
924 IF (double) ALLOCATE (buf2(l3:u3))
925
926 my_rank = pw%pw_grid%para%group%mepos
927 gid = pw%pw_grid%para%group
928 num_pe = pw%pw_grid%para%group%num_pe
929 tag = 1
930
931 rank(1) = unit_nr
932 rank(2) = my_rank
933 checksum = 0
934 IF (unit_nr > 0) checksum = 1
935
936 CALL gid%sum(checksum)
937 cpassert(checksum == 1)
938
939 CALL gid%maxloc(rank)
940 cpassert(rank(1) > 0)
941
942 dest = rank(2)
943 DO i1 = l1, u1, my_stride(1)
944 DO i2 = l2, u2, my_stride(2)
945
946 ! cycling through the CPUs, check if the current ray (I1,I2) is local to that CPU
947 IF (pw%pw_grid%para%mode /= pw_mode_local) THEN
948 DO ip = 0, num_pe - 1
949 IF (pw%pw_grid%para%bo(1, 1, ip, 1) <= i1 - l1 + 1 .AND. pw%pw_grid%para%bo(2, 1, ip, 1) >= i1 - l1 + 1 .AND. &
950 pw%pw_grid%para%bo(1, 2, ip, 1) <= i2 - l2 + 1 .AND. pw%pw_grid%para%bo(2, 2, ip, 1) >= i2 - l2 + 1) THEN
951 source = ip
952 END IF
953 END DO
954 ELSE
955 source = dest
956 END IF
957
958 IF (source == dest) THEN
959 IF (my_rank == source) THEN
960 buf(:) = pw%array(i1, i2, :)
961 IF (double) buf2(:) = pw2%array(i1, i2, :)
962 END IF
963 ELSE
964 IF (my_rank == source) THEN
965 buf(:) = pw%array(i1, i2, :)
966 CALL gid%send(buf, dest, tag)
967 IF (double) THEN
968 buf2(:) = pw2%array(i1, i2, :)
969 CALL gid%send(buf2, dest, tag)
970 END IF
971 END IF
972 IF (my_rank == dest) THEN
973 CALL gid%recv(buf, source, tag)
974 IF (double) CALL gid%recv(buf2, source, tag)
975 END IF
976 END IF
977
978 IF (.NOT. double) THEN
979 DO i3 = l3, u3, my_stride(3)
980 x = pw%pw_grid%dh(1, 1)*i1 + &
981 pw%pw_grid%dh(2, 1)*i2 + &
982 pw%pw_grid%dh(3, 1)*i3
983
984 y = pw%pw_grid%dh(1, 2)*i1 + &
985 pw%pw_grid%dh(2, 2)*i2 + &
986 pw%pw_grid%dh(3, 2)*i3
987
988 z = pw%pw_grid%dh(1, 3)*i1 + &
989 pw%pw_grid%dh(2, 3)*i2 + &
990 pw%pw_grid%dh(3, 3)*i3
991
992 IF (unit_nr > 0) THEN
993 WRITE (unit_nr, '(6E13.5E3, 6E13.5E3, 6E13.5E3, 6E13.5E3)') x, y, z, buf(i3)
994 END IF
995 END DO
996
997 ELSE
998
999 DO i3 = l3, u3, my_stride(3)
1000 x = pw%pw_grid%dh(1, 1)*i1 + &
1001 pw%pw_grid%dh(2, 1)*i2 + &
1002 pw%pw_grid%dh(3, 1)*i3
1003
1004 y = pw%pw_grid%dh(1, 2)*i1 + &
1005 pw%pw_grid%dh(2, 2)*i2 + &
1006 pw%pw_grid%dh(3, 2)*i3
1007
1008 z = pw%pw_grid%dh(1, 3)*i1 + &
1009 pw%pw_grid%dh(2, 3)*i2 + &
1010 pw%pw_grid%dh(3, 3)*i3
1011
1012 IF (unit_nr > 0) THEN
1013 WRITE (unit_nr, '(6E13.5E3, 6E13.5E3, 6E13.5E3, 6E13.5E3)') x, y, z, buf(i3), buf2(i3)
1014 END IF
1015 END DO
1016
1017 END IF ! Double
1018
1019 ! this double loop generates so many messages that it can overload
1020 ! the message passing system, e.g. on XT3
1021 ! we therefore put a barrier here that limits the amount of message
1022 ! that flies around at any given time.
1023 ! if ever this routine becomes a bottleneck, we should go for a
1024 ! more complicated rewrite
1025 CALL gid%sync()
1026
1027 END DO
1028 END DO
1029
1030 DEALLOCATE (buf)
1031 IF (double) DEALLOCATE (buf2)
1032
1033 CALL timestop(handle)
1034
1035 END SUBROUTINE pw_to_simple_volumetric
1036
1037END MODULE realspace_grid_cube
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Utility routines to open and close files. Tracking of preconnections.
Definition cp_files.F:16
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
Definition cp_files.F:311
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Definition cp_files.F:122
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
type(mp_file_descriptor_type) function, public mp_file_type_hindexed_make_chv(count, lengths, displs)
Creates an indexed MPI type for arrays of strings using bytes for spacing (hindexed type)
subroutine, public mp_file_type_free(type_descriptor)
Releases the type used for MPI I/O.
integer, parameter, public mpi_character_size
integer, parameter, public file_offset
integer, parameter, public file_amode_rdonly
subroutine, public mp_file_type_set_view_chv(fh, offset, type_descriptor)
Uses a previously created indexed MPI character type to tell the MPI processes how to partition (set_...
integer, parameter, public pw_mode_local
Generate Gaussian cube files.
subroutine, public cube_to_pw(grid, filename, scaling, parallel_read, silent)
Computes the external density on the grid hacked from external_read_density.
character(len= *), parameter, public cube_values_format
subroutine, public pw_to_simple_volumetric(pw, unit_nr, stride, pw2)
Prints a simple grid file: X Y Z value.
character(len= *), parameter, public cube_value_format
subroutine, public pw_to_cube(pw, unit_nr, title, particles_r, particles_z, particles_zeff, stride, max_file_size_mb, zero_tails, silent, mpi_io)
...
subroutine, public cube_read_values(values, buffer)
Read cube values from a character buffer.
void writebuffer(int *psockfd, char *data, int *plen)
Writes to a socket.
Definition sockets.c:201
void readbuffer(int *psockfd, char *data, int *plen)
Reads from a socket.
Definition sockets.c:219