1!--------------------------------------------------------------------------------------------------! 2! Copyright (C) by the DBCSR developers group - All rights reserved ! 3! This file is part of the DBCSR library. ! 4! ! 5! For information on the license, see the LICENSE file. ! 6! For further information please visit https://dbcsr.cp2k.org ! 7! SPDX-License-Identifier: GPL-2.0+ ! 8!--------------------------------------------------------------------------------------------------! 9 10MODULE dbcsr_io 11 !! DBCSR input/output 12 USE dbcsr_array_types, ONLY: array_data 13 USE dbcsr_data_methods, ONLY: dbcsr_data_clear_pointer, & 14 dbcsr_data_get_size_referenced, & 15 dbcsr_data_init, & 16 dbcsr_data_new, & 17 dbcsr_get_data 18 USE dbcsr_dist_methods, ONLY: dbcsr_distribution_mp 19 USE dbcsr_kinds, ONLY: default_string_length, & 20 dp, & 21 real_4, & 22 real_4_size, & 23 real_8, & 24 real_8_size, & 25 sp 26 USE dbcsr_machine, ONLY: default_output_unit 27 USE dbcsr_methods, ONLY: & 28 dbcsr_data_area, dbcsr_distribution, dbcsr_get_data_size, dbcsr_get_data_type, & 29 dbcsr_get_matrix_type, dbcsr_get_num_blocks, dbcsr_name, dbcsr_nblkcols_total, & 30 dbcsr_nblkrows_total, dbcsr_valid_index 31 USE dbcsr_mp_methods, ONLY: dbcsr_mp_group, & 32 dbcsr_mp_pgrid 33 USE dbcsr_mp_operations, ONLY: dbcsr_mp_type_from_anytype 34 USE dbcsr_mpiwrap, ONLY: & 35 file_amode_create, file_amode_rdonly, file_amode_wronly, file_offset, mp_allgather, & 36 mp_environ, mp_file_close, mp_file_get_size, mp_file_open, mp_file_read_at_all, & 37 mp_file_write_at, mp_file_write_at_all, mp_sum, mp_type_descriptor_type, mp_type_size, & 38 mpi_character_size, mpi_integer_size 39 USE dbcsr_transformations, ONLY: dbcsr_datablock_redistribute 40 USE dbcsr_types, ONLY: dbcsr_data_obj, & 41 dbcsr_distribution_obj, & 42 dbcsr_mp_obj, & 43 dbcsr_type, & 44 dbcsr_type_complex_4, & 45 dbcsr_type_complex_8, & 46 dbcsr_type_real_4, & 47 dbcsr_type_real_8 48 USE dbcsr_work_operations, ONLY: dbcsr_create 49#include "base/dbcsr_base_uses.f90" 50 51!$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num, omp_get_num_threads 52 53 IMPLICIT NONE 54 55 PRIVATE 56 57 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dbcsr_io' 58 59 ! Main 60 PUBLIC :: dbcsr_print 61 PUBLIC :: dbcsr_print_block_sum 62 ! Low-level printing 63! PUBLIC :: dbcsr_printmat, dbcsr_print2dmat 64 ! Utility printing 65 PUBLIC :: dbcsr_binary_write 66 PUBLIC :: dbcsr_binary_read 67 68 LOGICAL, PARAMETER :: bcsr_debug = .TRUE. 69 LOGICAL, PARAMETER :: bcsr_info = .FALSE. 70 LOGICAL, PARAMETER :: bcsr_verbose = .FALSE. 71 72 INTERFACE dbcsr_printmat 73 MODULE PROCEDURE printmat_s, printmat_d, printmat_c, printmat_z 74 END INTERFACE 75 76CONTAINS 77 78 SUBROUTINE dbcsr_print(matrix, nodata, matlab_format, variable_name, unit_nr) 79 !! Prints a BCSR matrix (block-style, not full) 80 81 TYPE(dbcsr_type), INTENT(IN) :: matrix 82 !! matrix 83 LOGICAL, INTENT(IN), OPTIONAL :: nodata, matlab_format 84 !! don't print actual data 85 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: variable_name 86 INTEGER, INTENT(IN), OPTIONAL :: unit_nr 87 88 CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_print', routineP = moduleN//':'//routineN 89 90 COMPLEX(KIND=real_4), DIMENSION(:), POINTER :: c_sp 91 COMPLEX(KIND=real_8), DIMENSION(:), POINTER :: c_dp 92 INTEGER :: ablk_p, bc, blk, blk_p, br, ebr, fblk, & 93 handle, ibr, iunit, lblk, m, mn, n, & 94 sblk 95 INTEGER, DIMENSION(:), POINTER :: col_blk_offset, col_blk_size, & 96 local_cols, local_rows, & 97 row_blk_offset, row_blk_size 98 LOGICAL :: my_matlab_format, tr, yesprint 99 REAL(KIND=dp) :: blk_cs 100 REAL(KIND=real_4), DIMENSION(:), POINTER :: r_sp 101 REAL(KIND=real_8), DIMENSION(:), POINTER :: r_dp 102 103! --------------------------------------------------------------------------- 104 105 CALL timeset(routineN, handle) 106 IF (.NOT. dbcsr_valid_index(matrix)) & 107 DBCSR_WARN("Can not print invalid matrix.") 108 109 iunit = default_output_unit 110 IF (PRESENT(unit_nr)) iunit = unit_nr 111 112 my_matlab_format = .FALSE. 113 IF (PRESENT(matlab_format)) my_matlab_format = matlab_format 114 yesprint = .TRUE. 115 IF (PRESENT(nodata)) yesprint = .NOT. nodata 116 WRITE (iunit, *) routineP//' Contents of matrix named ', matrix%name 117 WRITE (iunit, *) routineP//' Flags ', matrix%symmetry, & 118 matrix%negate_real, matrix%negate_imaginary, "type", & 119 dbcsr_get_data_type(matrix), "serial", matrix%serial_number 120 WRITE (iunit, '(1X,A,3(1X,I9,1X,A))') routineP, matrix%nblks, "blocks", & 121 matrix%nze, "nzes,", dbcsr_get_data_size(matrix), "data els", & 122 dbcsr_data_get_size_referenced(matrix%data_area), "used" 123 WRITE (iunit, '(1X,A,I5,A,I5)') routineP//" Full size", & 124 matrix%nfullrows_total, "x", matrix%nfullcols_total 125 WRITE (iunit, '(1X,A,I5,A,I5)') routineP//" Blocked size", & 126 matrix%nblkrows_total, "x", matrix%nblkcols_total 127 SELECT CASE (matrix%data_type) 128 CASE (dbcsr_type_real_8) 129 CALL dbcsr_get_data(matrix%data_area, r_dp) 130 CASE (dbcsr_type_real_4) 131 CALL dbcsr_get_data(matrix%data_area, r_sp) 132 CASE (dbcsr_type_complex_8) 133 CALL dbcsr_get_data(matrix%data_area, c_dp) 134 CASE (dbcsr_type_complex_4) 135 CALL dbcsr_get_data(matrix%data_area, c_sp) 136 END SELECT 137 row_blk_size => array_data(matrix%row_blk_size) 138 col_blk_size => array_data(matrix%col_blk_size) 139 row_blk_offset => array_data(matrix%row_blk_offset) 140 col_blk_offset => array_data(matrix%col_blk_offset) 141 142 IF (matrix%nblks .GT. 0) THEN 143 IF (matrix%list_indexing) THEN 144 IF (SIZE(matrix%coo_l) .NE. 3*matrix%nblks) & 145 DBCSR_ABORT("Wrong list") 146 ebr = 1 147 sblk = 3 148 ELSE 149 ebr = matrix%nblkrows_total 150 sblk = 1 151 ENDIF 152 DO ibr = 1, ebr 153 IF (matrix%list_indexing) THEN 154 fblk = 1 155 lblk = SIZE(matrix%coo_l) 156 ELSE 157 br = ibr 158 fblk = matrix%row_p(br) + 1 159 lblk = matrix%row_p(br + 1) 160 m = row_blk_size(br) 161 ENDIF 162 DO blk = fblk, lblk, sblk 163 IF (matrix%list_indexing) THEN 164 br = matrix%coo_l(blk) 165 bc = matrix%coo_l(blk + 1) 166 IF (matrix%local_indexing) THEN 167 local_rows => array_data(matrix%local_rows) 168 local_cols => array_data(matrix%local_cols) 169 br = local_rows(br) 170 bc = local_cols(bc) 171 ENDIF 172 m = row_blk_size(br) 173 ablk_p = matrix%coo_l(blk + 2) 174 ELSE 175 bc = matrix%col_i(blk) 176 ablk_p = matrix%blk_p(blk) 177 ENDIF 178 n = col_blk_size(bc) 179 mn = m*n 180 blk_p = ABS(ablk_p) 181 tr = ablk_p .LT. 0 182 block_exists: IF (blk_p .NE. 0) THEN 183 IF (mn .GT. 0) THEN 184 SELECT CASE (matrix%data_type) 185 CASE (dbcsr_type_real_8) 186 blk_cs = REAL(DOT_PRODUCT(r_dp(blk_p:blk_p + mn - 1), & 187 r_dp(blk_p:blk_p + mn - 1)), KIND=dp) 188 !CALL & 189 ! dbcsr_printmat(r_dp(blk_p:blk_p+mn-1),m,n, tr=tr) 190 CASE (dbcsr_type_real_4) 191 blk_cs = REAL(DOT_PRODUCT(r_sp(blk_p:blk_p + mn - 1), & 192 r_sp(blk_p:blk_p + mn - 1)), KIND=dp) 193 !CALL & 194 ! dbcsr_printmat(r_sp(blk_p:blk_p+mn-1),m,n, tr=tr) 195 CASE (dbcsr_type_complex_8) 196 blk_cs = REAL(DOT_PRODUCT(c_dp(blk_p:blk_p + mn - 1), & 197 c_dp(blk_p:blk_p + mn - 1)), KIND=dp) 198 !CALL & 199 ! dbcsr_printmat(c_dp(blk_p:blk_p+mn-1),m,n, tr=tr) 200 CASE (dbcsr_type_complex_4) 201 blk_cs = REAL(DOT_PRODUCT(c_sp(blk_p:blk_p + mn - 1), & 202 c_sp(blk_p:blk_p + mn - 1)), KIND=dp) 203 !CALL & 204 ! dbcsr_printmat(c_sp(blk_p:blk_p+mn-1),m,n, tr=tr) 205 END SELECT 206 ELSE 207 blk_cs = 0.0_dp 208 ENDIF 209 !WRITE(iunit,*)routineP//' chksum for (',br,',',bc,') at',& 210 ! blk_p,'l',mn,'= ', blk_cs,'size',m,n 211 IF (.NOT. my_matlab_format) WRITE (iunit, '(A,I6,",",I6,A,I7,A,I6,I6,"=",I7,A,E12.3)') & 212 !" Checksum for (",br,bc,") at ",blk_p," size ",m,n,mn,& 213 " Checksum for (", br, bc, ") at ", ablk_p, " size ", m, n, mn, & 214 " checksum=", blk_cs 215 IF (yesprint .AND. blk_p .NE. 0) THEN 216 IF (mn .GT. 0) THEN 217 SELECT CASE (matrix%data_type) 218 CASE (dbcsr_type_real_8) 219 !WRITE(iunit,'(10(1X,F7.2))')r_dp(blk_p:blk_p+mn-1) 220 IF (my_matlab_format) THEN 221 CALL dbcsr_printmat_matlab_d(r_dp(blk_p:blk_p + mn - 1), m, n, & 222 row_blk_offset(br), col_blk_offset(bc), iunit, tr=tr, & 223 variable_name=variable_name) 224 ELSE 225 CALL dbcsr_printmat(r_dp(blk_p:blk_p + mn - 1), m, n, iunit=iunit, tr=tr) 226 ENDIF 227 CASE (dbcsr_type_real_4) 228 IF (my_matlab_format) THEN 229 CALL dbcsr_printmat_matlab_s(r_sp(blk_p:blk_p + mn - 1), m, n, & 230 row_blk_offset(br), col_blk_offset(bc), iunit, tr=tr, & 231 variable_name=variable_name) 232 ELSE 233 CALL dbcsr_printmat(r_sp(blk_p:blk_p + mn - 1), m, n, iunit=iunit, tr=tr) 234 ENDIF 235 CASE (dbcsr_type_complex_8) 236 IF (my_matlab_format) THEN 237 CALL dbcsr_printmat_matlab_z(c_dp(blk_p:blk_p + mn - 1), m, n, & 238 row_blk_offset(br), col_blk_offset(bc), iunit, tr=tr, & 239 variable_name=variable_name) 240 ELSE 241 CALL dbcsr_printmat(c_dp(blk_p:blk_p + mn - 1), m, n, iunit=iunit, tr=tr) 242 ENDIF 243 CASE (dbcsr_type_complex_4) 244 IF (my_matlab_format) THEN 245 CALL dbcsr_printmat_matlab_c(c_sp(blk_p:blk_p + mn - 1), m, n, & 246 row_blk_offset(br), col_blk_offset(bc), iunit, tr=tr, & 247 variable_name=variable_name) 248 ELSE 249 CALL dbcsr_printmat(c_sp(blk_p:blk_p + mn - 1), m, n, iunit=iunit, tr=tr) 250 ENDIF 251 END SELECT 252 ENDIF 253 ENDIF 254 ENDIF block_exists 255 ENDDO 256 ENDDO 257 ENDIF 258 CALL timestop(handle) 259 END SUBROUTINE dbcsr_print 260 261 SUBROUTINE dbcsr_printmat_matlab_d(matrix, rows, cols, r_offset, c_offset, iunit, tr, variable_name) 262 !! Prints the elements of a matrix. 263 264 REAL(KIND=real_8), DIMENSION(:), INTENT(IN) :: matrix 265 INTEGER, INTENT(IN) :: rows, cols, r_offset, c_offset, iunit 266 !! the logical (possibly detransposed) matrix size, not the stored size 267 !! the logical (possibly detransposed) matrix size, not the stored size 268 LOGICAL, INTENT(IN), OPTIONAL :: tr 269 !! specifies whether the elements are stored transposed 270 CHARACTER(len=*), INTENT(in), OPTIONAL :: variable_name 271 272 INTEGER :: c, c_off, m, n, r, r_off 273 LOGICAL :: t 274 275! --------------------------------------------------------------------------- 276 277 m = rows 278 n = cols 279 r_off = r_offset 280 c_off = c_offset 281 t = .FALSE. 282 IF (PRESENT(tr)) THEN 283 IF (tr) THEN 284 t = .TRUE. 285 m = cols 286 n = rows 287 r_off = c_offset 288 c_off = r_offset 289 ENDIF 290 ENDIF 291 292 DO c = 1, cols 293 DO r = 1, rows 294 IF (.NOT. t) THEN 295 IF (PRESENT(variable_name)) THEN 296 WRITE (iunit, '(A,I4,A,I4,A,E23.16,A)') & 297 variable_name//'(', r + r_offset - 1, ',', c + c_offset - 1, ')=', matrix(r + (c - 1)*rows), ';' 298 ELSE 299 WRITE (iunit, '(A,I4,A,I4,A,E23.16,A)') 'a(', r + r_offset - 1, ',', & 300 c + c_offset - 1, ')=', matrix(r + (c - 1)*rows), ';' 301 ENDIF 302 ELSE 303 IF (PRESENT(variable_name)) THEN 304 WRITE (iunit, '(A,I4,A,I4,A,E23.16,A)') & 305 variable_name//'(', r + r_offset - 1, ',', c + c_offset - 1, ')=', matrix((r - 1)*cols + c), ';' 306 ELSE 307 WRITE (iunit, '(A,I4,A,I4,A,E23.16,A)') 'a(', r + r_offset - 1, ',', & 308 c + c_offset - 1, ')=', matrix((r - 1)*cols + c), ';' 309 ENDIF 310 ENDIF 311 ENDDO 312 ENDDO 313 END SUBROUTINE dbcsr_printmat_matlab_d 314 315 SUBROUTINE dbcsr_printmat_matlab_s(matrix, rows, cols, r_offset, c_offset, iunit, tr, variable_name) 316 REAL(KIND=real_4), DIMENSION(:), INTENT(IN) :: matrix 317 INTEGER, INTENT(IN) :: rows, cols, r_offset, c_offset, iunit 318 LOGICAL, INTENT(IN), OPTIONAL :: tr 319 CHARACTER(len=*), INTENT(in), OPTIONAL :: variable_name 320 321 INTEGER :: c, c_off, m, n, r, r_off 322 LOGICAL :: t 323 324! --------------------------------------------------------------------------- 325 326 m = rows 327 n = cols 328 r_off = r_offset 329 c_off = c_offset 330 t = .FALSE. 331 IF (PRESENT(tr)) THEN 332 IF (tr) THEN 333 t = .TRUE. 334 m = cols 335 n = rows 336 r_off = c_offset 337 c_off = r_offset 338 ENDIF 339 ENDIF 340 341 DO c = 1, cols 342 DO r = 1, rows 343 IF (.NOT. t) THEN 344 IF (PRESENT(variable_name)) THEN 345 WRITE (iunit, '(A,I4,A,I4,A,E15.7,A)') & 346 variable_name//'(', r + r_offset - 1, ',', c + c_offset - 1, ')=', matrix(r + (c - 1)*rows), ';' 347 ELSE 348 WRITE (iunit, '(A,I4,A,I4,A,E15.7,A)') 'a(', r + r_offset - 1, ',', & 349 c + c_offset - 1, ')=', matrix(r + (c - 1)*rows), ';' 350 ENDIF 351 ELSE 352 IF (PRESENT(variable_name)) THEN 353 WRITE (iunit, '(A,I4,A,I4,A,E15.7,A)') & 354 variable_name//'(', r + r_offset - 1, ',', c + c_offset - 1, ')=', matrix((r - 1)*cols + c), ';' 355 ELSE 356 WRITE (iunit, '(A,I4,A,I4,A,E15.7,A)') 'a(', r + r_offset - 1, ',', & 357 c + c_offset - 1, ')=', matrix((r - 1)*cols + c), ';' 358 ENDIF 359 ENDIF 360 ENDDO 361 ENDDO 362 END SUBROUTINE dbcsr_printmat_matlab_s 363 364 SUBROUTINE dbcsr_printmat_matlab_z(matrix, rows, cols, r_offset, c_offset, iunit, tr, variable_name) 365 COMPLEX(KIND=real_8), DIMENSION(:), INTENT(IN) :: matrix 366 INTEGER, INTENT(IN) :: rows, cols, r_offset, c_offset, iunit 367 LOGICAL, INTENT(IN), OPTIONAL :: tr 368 CHARACTER(len=*), INTENT(in), OPTIONAL :: variable_name 369 370 INTEGER :: c, c_off, m, n, r, r_off 371 LOGICAL :: t 372 373! --------------------------------------------------------------------------- 374 375 m = rows 376 n = cols 377 r_off = r_offset 378 c_off = c_offset 379 t = .FALSE. 380 IF (PRESENT(tr)) THEN 381 IF (tr) THEN 382 t = .TRUE. 383 m = cols 384 n = rows 385 r_off = c_offset 386 c_off = r_offset 387 ENDIF 388 ENDIF 389 390 DO c = 1, cols 391 DO r = 1, rows 392 IF (.NOT. t) THEN 393 IF (PRESENT(variable_name)) THEN 394 WRITE (iunit, '(A,I3,A,I3,A,E23.16,A,E23.16,A)') variable_name//'(', r + r_offset - 1, ',', & 395 c + c_offset - 1, ')=', & 396 REAL(matrix(r + (c - 1)*rows)), '+', AIMAG(matrix(r + (c - 1)*rows)), 'i;' 397 ELSE 398 WRITE (iunit, '(A,I3,A,I3,A,E23.16,A,E23.16,A)') 'a(', r + r_offset - 1, ',', c + c_offset - 1, ')=', & 399 REAL(matrix(r + (c - 1)*rows)), '+', AIMAG(matrix(r + (c - 1)*rows)), 'i;' 400 ENDIF 401 ELSE 402 IF (PRESENT(variable_name)) THEN 403 WRITE (iunit, '(A,I3,A,I3,A,E23.16,A,E23.16,A)') variable_name//'(', r + r_offset - 1, ',', & 404 c + c_offset - 1, ')=', & 405 REAL(matrix((r - 1)*cols + c)), '+', AIMAG(matrix((r - 1)*cols + c)), 'i;' 406 ELSE 407 WRITE (iunit, '(A,I3,A,I3,A,E23.16,A,E23.16,A)') 'a(', r + r_offset - 1, ',', c + c_offset - 1, ')=', & 408 REAL(matrix((r - 1)*cols + c)), '+', AIMAG(matrix((r - 1)*cols + c)), 'i;' 409 ENDIF 410 ENDIF 411 ENDDO 412 ENDDO 413 END SUBROUTINE dbcsr_printmat_matlab_z 414 415 SUBROUTINE dbcsr_printmat_matlab_c(matrix, rows, cols, r_offset, c_offset, iunit, tr, variable_name) 416 COMPLEX(KIND=real_4), DIMENSION(:), INTENT(IN) :: matrix 417 INTEGER, INTENT(IN) :: rows, cols, r_offset, c_offset, iunit 418 LOGICAL, INTENT(IN), OPTIONAL :: tr 419 CHARACTER(len=*), INTENT(in), OPTIONAL :: variable_name 420 421 INTEGER :: c, c_off, m, n, r, r_off 422 LOGICAL :: t 423 424! --------------------------------------------------------------------------- 425 426 m = rows 427 n = cols 428 r_off = r_offset 429 c_off = c_offset 430 t = .FALSE. 431 IF (PRESENT(tr)) THEN 432 IF (tr) THEN 433 t = .TRUE. 434 m = cols 435 n = rows 436 r_off = c_offset 437 c_off = r_offset 438 ENDIF 439 ENDIF 440 441 DO c = 1, cols 442 DO r = 1, rows 443 IF (.NOT. t) THEN 444 IF (PRESENT(variable_name)) THEN 445 WRITE (iunit, '(A,I3,A,I3,A,E15.7,A,E15.7,A)') variable_name//'(', r + r_offset - 1, ',', & 446 c + c_offset - 1, ')=', & 447 REAL(matrix(r + (c - 1)*rows)), '+', AIMAG(matrix(r + (c - 1)*rows)), 'i;' 448 ELSE 449 WRITE (iunit, '(A,I3,A,I3,A,E15.7,A,E15.7,A)') 'a(', r + r_offset - 1, ',', c + c_offset - 1, ')=', & 450 REAL(matrix(r + (c - 1)*rows)), '+', AIMAG(matrix(r + (c - 1)*rows)), 'i;' 451 ENDIF 452 ELSE 453 IF (PRESENT(variable_name)) THEN 454 WRITE (iunit, '(A,I3,A,I3,A,E15.7,A,E15.7,A)') variable_name//'(', r + r_offset - 1, ',', & 455 c + c_offset - 1, ')=', & 456 REAL(matrix((r - 1)*cols + c)), '+', AIMAG(matrix((r - 1)*cols + c)), 'i;' 457 ELSE 458 WRITE (iunit, '(A,I3,A,I3,A,E15.7,A,E15.7,A)') 'a(', r + r_offset - 1, ',', c + c_offset - 1, ')=', & 459 REAL(matrix((r - 1)*cols + c)), '+', AIMAG(matrix((r - 1)*cols + c)), 'i;' 460 ENDIF 461 ENDIF 462 ENDDO 463 ENDDO 464 END SUBROUTINE dbcsr_printmat_matlab_c 465 466 SUBROUTINE printmat_s(matrix, rows, cols, iunit, title, tr) 467 !! Prints the elements of a matrix. 468 469 REAL(KIND=real_4), DIMENSION(:), INTENT(IN) :: matrix 470 INTEGER, INTENT(IN) :: rows, cols, iunit 471 !! the logical (possibly detransposed) matrix size, not the stored size 472 !! the logical (possibly detransposed) matrix size, not the stored size 473 CHARACTER(*), INTENT(IN), OPTIONAL :: title 474 LOGICAL, INTENT(IN), OPTIONAL :: tr 475 !! specifies whether the elements are stored transposed 476 477 CHARACTER(30) :: f 478 INTEGER :: m, n, r 479 LOGICAL :: t 480 REAL(KIND=dp) :: bit_bucket 481 482! --------------------------------------------------------------------------- 483 484 m = rows 485 n = cols 486 t = .FALSE. 487 IF (PRESENT(title)) WRITE (iunit, *) title 488 IF (PRESENT(tr)) THEN 489 IF (tr) THEN 490 t = .TRUE. 491 m = cols 492 n = rows 493 ENDIF 494 ENDIF 495 DO r = LBOUND(matrix, 1), UBOUND(matrix, 1) 496 bit_bucket = matrix(r) 497 ENDDO 498 bit_bucket = 0.0_dp 499 DO r = LBOUND(matrix, 1), UBOUND(matrix, 1) 500 bit_bucket = bit_bucket + matrix(r) 501 ENDDO 502 IF (m .GT. 10000) m = 0 503 IF (n .GT. 10000) n = 0 504 IF (m*n .LT. 1 .OR. m*n .GT. SIZE(matrix)) RETURN 505 WRITE (f, FMT="('(',I4,'(F9.4))')") cols 506 DO r = 1, rows 507 IF (.NOT. t) THEN 508 WRITE (iunit, FMT=f) matrix(r:r + (cols - 1)*rows:rows) 509 ELSE 510 WRITE (iunit, FMT=f) matrix((r - 1)*cols + 1:r*cols) 511 ENDIF 512 ENDDO 513 END SUBROUTINE printmat_s 514 515 SUBROUTINE printmat_d(matrix, rows, cols, iunit, title, tr) 516 REAL(KIND=real_8), DIMENSION(:), INTENT(IN) :: matrix 517 INTEGER, INTENT(IN) :: rows, cols, iunit 518 CHARACTER(*), INTENT(IN), OPTIONAL :: title 519 LOGICAL, INTENT(IN), OPTIONAL :: tr 520 521 IF (PRESENT(title)) THEN 522 IF (PRESENT(tr)) THEN 523 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title, tr) 524 ELSE 525 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title) 526 ENDIF 527 ELSE 528 IF (PRESENT(tr)) THEN 529 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, tr=tr) 530 ELSE 531 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit) 532 ENDIF 533 ENDIF 534 END SUBROUTINE printmat_d 535 536 SUBROUTINE printmat_c(matrix, rows, cols, iunit, title, tr) 537 COMPLEX(KIND=real_4), DIMENSION(:), INTENT(IN) :: matrix 538 INTEGER, INTENT(IN) :: rows, cols, iunit 539 CHARACTER(*), INTENT(IN), OPTIONAL :: title 540 LOGICAL, INTENT(IN), OPTIONAL :: tr 541 542 IF (PRESENT(title)) THEN 543 IF (PRESENT(tr)) THEN 544 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title, tr) 545 ELSE 546 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title) 547 ENDIF 548 ELSE 549 IF (PRESENT(tr)) THEN 550 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, tr=tr) 551 ELSE 552 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit) 553 ENDIF 554 ENDIF 555 END SUBROUTINE printmat_c 556 557 SUBROUTINE printmat_z(matrix, rows, cols, iunit, title, tr) 558 COMPLEX(KIND=real_8), DIMENSION(:), INTENT(IN) :: matrix 559 INTEGER, INTENT(IN) :: rows, cols, iunit 560 CHARACTER(*), INTENT(IN), OPTIONAL :: title 561 LOGICAL, INTENT(IN), OPTIONAL :: tr 562 563 IF (PRESENT(title)) THEN 564 IF (PRESENT(tr)) THEN 565 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title, tr) 566 ELSE 567 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, title) 568 ENDIF 569 ELSE 570 IF (PRESENT(tr)) THEN 571 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit, tr=tr) 572 ELSE 573 CALL printmat_s(REAL(matrix, KIND=sp), rows, cols, iunit) 574 ENDIF 575 ENDIF 576 END SUBROUTINE printmat_z 577 578 SUBROUTINE dbcsr_binary_write(matrix, filepath) 579 !! Writes a DBCSR matrix in a file 580 !! file's header consists of 3 sub-headers: 581 !! sub-header1 contains: 582 !! 1 string: (of length version_len) the current version of this routine, 583 !! 1 string: (of length default_string_length) matrix_name, 584 !! 1 character: matrix_type, 585 !! 4 integers: numnodes, data_type, nblkrows_total, nblkcols_total, 586 !! 2 vectors: row_blk_size (length = nblkrows_total), 587 !! col_blk_size (length = nblkcols_total), 588 !! sub-header2 contains: 589 !! 2 integers: nblks, data_area_size, 590 !! sub-header3 contains: 591 !! 3 vectors: row_p (length = nblkrows_total+1), 592 !! col_i (length = nblks), 593 !! blk_p (length = nblks); 594 !! and the file's body contains the block data 595 596 IMPLICIT NONE 597 598 TYPE(dbcsr_type), INTENT(IN) :: matrix 599 !! DBCSR matrix 600 CHARACTER(len=*), INTENT(IN) :: filepath 601 !! path to the file 602 603 CHARACTER(LEN=*), PARAMETER :: routineN = 'dbcsr_binary_write' 604 605 INTEGER :: nblkrows_total, nblkcols_total, & 606 nblks, size_of_pgrid, & 607 thefile, i, sendbuf, data_area_size, & 608 data_type, type_size, & 609 mp_group, mynode, numnodes, & 610 ginfo_size, linfo_size, handle 611 INTEGER, DIMENSION(:), POINTER :: row_p, col_i, blk_p, & 612 row_blk_size, col_blk_size 613 INTEGER, DIMENSION(:, :), POINTER :: pgrid 614 TYPE(mp_type_descriptor_type) :: mp_type 615 TYPE(dbcsr_mp_obj) :: mp_env 616 TYPE(dbcsr_distribution_obj) :: distribution 617 TYPE(dbcsr_data_obj) :: data_area 618 COMPLEX(sp), DIMENSION(:), POINTER :: c_sp 619 COMPLEX(dp), DIMENSION(:), POINTER :: c_dp 620 REAL(sp), DIMENSION(:), POINTER :: r_sp 621 REAL(dp), DIMENSION(:), POINTER :: r_dp 622 CHARACTER :: matrix_type 623 CHARACTER(LEN=80) :: matrix_name_v_1_0 624 CHARACTER(LEN=default_string_length) :: matrix_name 625 INTEGER, PARAMETER :: version_len = 10 626 CHARACTER(LEN=version_len), PARAMETER :: version = "DBCSRv_1.0" 627 INTEGER, ALLOCATABLE, DIMENSION(:) :: linfo_sizes, da_sizes 628 INTEGER(kind=file_offset), ALLOCATABLE, DIMENSION(:) :: bdata_disps, bdata_offsets, & 629 subh2_disps, subh2_offsets, & 630 subh3_disps, subh3_offsets 631 INTEGER(kind=file_offset), PARAMETER :: BOF = 0 632 INTEGER, PARAMETER :: char_count = version_len + default_string_length + 1 !version, matrix_name, matrix_type 633 634 CALL timeset(routineN, handle) 635 636 IF (default_string_length /= 80) & 637 CALL dbcsr_warn(__LOCATION__, "Changing the default string length affects "// & 638 "the format of the written matrix. Version needs to be adjusted") 639 640 nblkrows_total = dbcsr_nblkrows_total(matrix) 641 nblkcols_total = dbcsr_nblkcols_total(matrix) 642 distribution = dbcsr_distribution(matrix) 643 matrix_name = dbcsr_name(matrix) 644 data_area = dbcsr_data_area(matrix) 645 matrix_type = dbcsr_get_matrix_type(matrix) 646 data_type = dbcsr_get_data_type(matrix) 647 mp_env = dbcsr_distribution_mp(distribution) 648 mp_group = dbcsr_mp_group(mp_env) 649 nblks = dbcsr_get_num_blocks(matrix) 650 row_p => matrix%row_p 651 col_i => matrix%col_i 652 blk_p => matrix%blk_p 653 row_blk_size => array_data(matrix%row_blk_size) 654 col_blk_size => array_data(matrix%col_blk_size) 655 pgrid => dbcsr_mp_pgrid(mp_env) 656 size_of_pgrid = SIZE(pgrid) 657 658 CALL mp_environ(numnodes, mynode, mp_group) 659 660 ALLOCATE (linfo_sizes(numnodes), da_sizes(numnodes), & 661 subh2_disps(numnodes), subh2_offsets(numnodes), & 662 subh3_disps(numnodes), subh3_offsets(numnodes), & 663 bdata_disps(numnodes), bdata_offsets(numnodes)) 664 subh2_disps(:) = (/((i - 1)*2, i=1, numnodes)/) 665 subh3_disps = BOF 666 bdata_disps = BOF 667 linfo_sizes = BOF 668 subh2_offsets = BOF 669 subh3_offsets = BOF 670 bdata_offsets = BOF 671 da_sizes = BOF 672 673 ginfo_size = char_count + 4 + nblkrows_total + nblkcols_total 674 linfo_size = 1 + nblkrows_total + 2*nblks 675 676 sendbuf = linfo_size 677 CALL mp_allgather(sendbuf, linfo_sizes, mp_group) 678 CALL cumsum_l(INT(linfo_sizes, kind=file_offset), subh3_disps) 679 subh3_disps(:) = CSHIFT(subh3_disps, shift=-1) + ginfo_size + 2*numnodes 680 subh3_disps(1) = ginfo_size + 2*numnodes 681 682 data_area_size = dbcsr_data_get_size_referenced(matrix%data_area) 683 sendbuf = data_area_size 684 CALL mp_allgather(sendbuf, da_sizes, mp_group) 685 CALL cumsum_l(INT(da_sizes, kind=file_offset), bdata_disps) 686 bdata_disps(:) = CSHIFT(bdata_disps, shift=-1) + SUM(INT(linfo_sizes, KIND=file_offset)) + & 687 ginfo_size + numnodes*2 688 bdata_disps(1) = SUM(INT(linfo_sizes, KIND=file_offset)) + ginfo_size + numnodes*2 689 690 CALL mp_file_open(mp_group, thefile, filepath, file_amode_create + file_amode_wronly) 691 692 IF (mynode .EQ. 0) THEN 693 CALL mp_file_write_at(thefile, BOF, version) 694 matrix_name_v_1_0 = matrix_name 695 CALL mp_file_write_at(thefile, BOF + version_len*mpi_character_size, matrix_name_v_1_0) 696 CALL mp_file_write_at(thefile, BOF + (version_len + default_string_length)*mpi_character_size, matrix_type) 697 CALL mp_file_write_at(thefile, BOF + char_count*mpi_character_size, & 698 (/size_of_pgrid, data_type, & 699 nblkrows_total, nblkcols_total, & 700 row_blk_size, col_blk_size/)) 701 END IF 702! write sub-header2 703 subh2_disps(:) = subh2_disps(:) + ginfo_size 704 subh2_offsets(:) = BOF + (subh2_disps - char_count)*mpi_integer_size + & 705 char_count*mpi_character_size 706 CALL mp_file_write_at_all(thefile, subh2_offsets(mynode + 1), (/nblks, data_area_size/)) 707! write sub-header3 708 subh3_offsets(:) = BOF + (subh3_disps - char_count)*mpi_integer_size + & 709 char_count*mpi_character_size 710 CALL mp_file_write_at_all(thefile, subh3_offsets(mynode + 1), (/row_p, col_i, blk_p/)) 711! write block data 712 mp_type = dbcsr_mp_type_from_anytype(data_area) 713 CALL mp_type_size(mp_type, type_size) 714 bdata_offsets(:) = BOF + (/((bdata_disps(i) - bdata_disps(1))*type_size, i=1, numnodes)/) + & 715 (bdata_disps(1) - char_count)*mpi_integer_size + & 716 char_count*mpi_character_size 717 SELECT CASE (data_type) 718 CASE (dbcsr_type_real_4) 719 r_sp => data_area%d%r_sp 720 CALL mp_file_write_at_all(thefile, bdata_offsets(mynode + 1), r_sp, msglen=data_area_size) 721 CASE (dbcsr_type_real_8) 722 r_dp => data_area%d%r_dp 723 CALL mp_file_write_at_all(thefile, bdata_offsets(mynode + 1), r_dp, msglen=data_area_size) 724 CASE (dbcsr_type_complex_4) 725 c_sp => data_area%d%c_sp 726 CALL mp_file_write_at_all(thefile, bdata_offsets(mynode + 1), c_sp, msglen=data_area_size) 727 CASE (dbcsr_type_complex_8) 728 c_dp => data_area%d%c_dp 729 CALL mp_file_write_at_all(thefile, bdata_offsets(mynode + 1), c_dp, msglen=data_area_size) 730 END SELECT 731 CALL mp_file_close(thefile) 732 733 DEALLOCATE (linfo_sizes, da_sizes) 734 DEALLOCATE (subh2_disps, subh2_offsets, subh3_disps, subh3_offsets) 735 DEALLOCATE (bdata_disps, bdata_offsets) 736 737 CALL timestop(handle) 738 739 CONTAINS 740 SUBROUTINE cumsum_l(arr, cumsum) 741 INTEGER(kind=file_offset), DIMENSION(:), & 742 INTENT(IN) :: arr 743 INTEGER(kind=file_offset), DIMENSION(SIZE(arr)), & 744 INTENT(OUT) :: cumsum 745 746 INTEGER :: i 747 748 cumsum(1) = arr(1) 749 DO i = 2, SIZE(arr) 750 cumsum(i) = cumsum(i - 1) + arr(i) 751 END DO 752 END SUBROUTINE cumsum_l 753 END SUBROUTINE dbcsr_binary_write 754 755 SUBROUTINE dbcsr_binary_read(filepath, distribution, groupid, matrix_new) 756 !! Reads a DBCSR matrix from a file 757 758 IMPLICIT NONE 759 760 CHARACTER(len=*), INTENT(IN) :: filepath 761 !! path to the file 762 TYPE(dbcsr_distribution_obj), INTENT(IN) :: distribution 763 !! row and column distribution 764 INTEGER, INTENT(IN), OPTIONAL :: groupid 765 !! message passing environment identifier 766 TYPE(dbcsr_type), INTENT(INOUT) :: matrix_new 767 !! DBCSR matrix 768 769 CHARACTER(LEN=*), PARAMETER :: routineN = 'dbcsr_binary_read' 770 771 INTEGER :: nblkrows_total, nblkcols_total, & 772 nblks, darea_size, data_type, type_size, & 773 globalinfo_size, & 774 size_of_pgrid, & 775 thefile, i, j, & 776 nblocks, & 777 share_size, order, cur_blks, & 778 job_count, start_index, end_index, & 779 localinfo_length, blockdata_length, & 780 group_id, worker_id, group_list_size, handle, linfo_length 781 CHARACTER :: matrix_type 782 CHARACTER(LEN=default_string_length) :: matrix_name 783 INTEGER, PARAMETER :: version_len = 10 784 CHARACTER(LEN=version_len) :: version 785 CHARACTER(LEN=80) :: matrix_name_v_1_0 786 CHARACTER(LEN=version_len), PARAMETER :: version_v_1_0 = "DBCSRv_1.0" 787 788 INTEGER, DIMENSION(:), POINTER :: row_p, col_i, blk_p, & 789 proc_nblks, proc_darea_sizes 790 INTEGER, DIMENSION(4) :: values 791 INTEGER, ALLOCATABLE, DIMENSION(:) :: linfo_lens, bdata_lens 792 INTEGER, ALLOCATABLE, DIMENSION(:), TARGET :: ginfo_vec, linfo_vec, & 793 rowp, coli, blkp 794 INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: val_data 795 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: row_blk_size, col_blk_size 796 TYPE(dbcsr_mp_obj) :: mp_env 797 TYPE(dbcsr_data_obj) :: dblk 798 REAL(sp) :: rsp_dummy(1) 799 REAL(dp) :: rdp_dummy(1) 800 COMPLEX(sp) :: csp_dummy(1) 801 COMPLEX(dp) :: cdp_dummy(1) 802 REAL(sp), ALLOCATABLE, DIMENSION(:), TARGET :: rsp 803 REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rdp 804 COMPLEX(sp), ALLOCATABLE, DIMENSION(:), TARGET :: csp 805 COMPLEX(dp), ALLOCATABLE, DIMENSION(:), TARGET :: cdp 806 INTEGER(kind=file_offset), ALLOCATABLE, DIMENSION(:) :: subh2_offsets, & 807 subh3_disps, subh3_offsets, & 808 bdata_disps, bdata_offsets 809 INTEGER(kind=file_offset), PARAMETER :: BOF = 0 810 INTEGER(kind=file_offset) :: offset, subh2_start, subh3_start, bdata_start, file_size, & 811 localinfo_offset, blockdata_offset, sum_nblks, subh3_length, data_area_size 812 INTEGER, PARAMETER :: char_count = 1 + version_len + default_string_length 813 814 CALL timeset(routineN, handle) 815 816 IF (PRESENT(groupid)) THEN 817 group_id = groupid 818 ELSE 819 mp_env = dbcsr_distribution_mp(distribution) 820 group_id = dbcsr_mp_group(mp_env) 821 END IF 822 823 CALL mp_environ(group_list_size, worker_id, group_id) 824 CALL mp_file_open(group_id, thefile, filepath, file_amode_rdonly) 825 826! read version, matrix name and matrix type 827 CALL mp_file_read_at_all(thefile, BOF, version) 828 829 IF (version /= version_v_1_0) & 830 DBCSR_WARN("Trying to read an unknown version of the matrix data file. Good luck!") 831 832 CALL mp_file_read_at_all(thefile, BOF + version_len*mpi_character_size, matrix_name_v_1_0) 833 matrix_name = matrix_name_v_1_0 834 835 CALL mp_file_read_at_all(thefile, BOF + (version_len + default_string_length)*mpi_character_size, matrix_type) 836! read 4 integer values form sub-header1 837 CALL mp_file_read_at_all(thefile, BOF + char_count*mpi_character_size, values) 838 size_of_pgrid = values(1) 839 data_type = values(2) 840 nblkrows_total = values(3) 841 nblkcols_total = values(4) 842! read 2 vectors, row_blk_size and col_blk_size, from sub-header1 843 globalinfo_size = nblkrows_total + nblkcols_total 844 ALLOCATE (ginfo_vec(globalinfo_size)) 845 CALL mp_file_read_at_all(thefile, BOF + char_count*mpi_character_size + 4*mpi_integer_size, ginfo_vec) 846 row_blk_size => ginfo_vec(1:nblkrows_total) 847 col_blk_size => ginfo_vec(nblkrows_total + 1:globalinfo_size) 848 849! compute the offsets where sub-header2 and sub-header3 start 850 subh2_start = (4 + globalinfo_size)*mpi_integer_size + char_count*mpi_character_size 851 subh3_start = subh2_start + 2*size_of_pgrid*mpi_integer_size 852 853! compute the offsets in sub-header2 and read 2 integers nblocks, data_area_size 854 ! number of data chunks from sub-header 2 and 3 to be read by every node rounded up 855 ! to the next integer to make it even for all the nodes in the specified mpi group 856 share_size = CEILING(REAL(size_of_pgrid, KIND=dp)/group_list_size) 857 858 ALLOCATE (subh2_offsets(share_size)) 859 subh2_offsets = BOF 860 DO i = 1, share_size 861 offset = subh2_start + mpi_integer_size*2*(worker_id + (i - 1)*group_list_size) 862 IF (offset .GE. subh3_start) EXIT 863 subh2_offsets(i) = offset 864 END DO 865 866 ALLOCATE (val_data(3, share_size)) 867 val_data(:, :) = 0 868 DO i = 1, share_size 869 CALL mp_file_read_at_all(thefile, subh2_offsets(i), values, msglen=2) 870 nblocks = values(1) 871 data_area_size = values(2) 872 IF (subh2_offsets(i) .EQ. 0) EXIT 873 val_data(1, i) = nblocks 874 IF (data_area_size >= HUGE(val_data(2, i))) & 875 DBCSR_ABORT("Data area too large, fix code.") 876 val_data(2, i) = INT(data_area_size) 877 val_data(3, i) = worker_id + (i - 1)*group_list_size + 1 ! order 878 ! order = indices of an array of length size_of_pgrid to be accessed by the current node 879 END DO 880 nblks = SUM(val_data(1, :)) 881 darea_size = SUM(val_data(2, :)) 882 proc_nblks => val_data(1, :) ! to be passed to dbcsr_datablock_redistribute 883 proc_darea_sizes => val_data(2, :) ! to be passed to dbcsr_datablock_redistribute 884 885! compute the offsets in sub-header3 and read 3 vectors row_p, col_i, blk_p 886 ! actual number of chunks to be read by the current node 887 job_count = COUNT(val_data(3, :) .NE. 0) 888 CALL mp_file_get_size(thefile, file_size) 889 890 ALLOCATE (linfo_lens(size_of_pgrid)) 891 ALLOCATE (subh3_disps(size_of_pgrid)) 892 ALLOCATE (subh3_offsets(size_of_pgrid)) 893 linfo_lens = 0; subh3_disps = 0 894 DO i = 1, size_of_pgrid 895 DO j = 1, share_size 896 order = val_data(3, j) 897 IF (i .EQ. order) linfo_lens(order) = & 898 1 + nblkrows_total + 2*val_data(1, j) 899 END DO 900 END DO 901 CALL mp_sum(linfo_lens, group_id) 902 CALL cumsum_l(INT(linfo_lens, kind=file_offset), subh3_disps) 903 subh3_disps(:) = CSHIFT(subh3_disps, shift=-1) 904 subh3_disps(1) = BOF 905 subh3_offsets(:) = subh3_start + subh3_disps*mpi_integer_size 906 907 sum_nblks = INT(nblks, kind=file_offset) 908 CALL mp_sum(sum_nblks, group_id) 909 subh3_length = size_of_pgrid*INT(1 + nblkrows_total, KIND=file_offset) + 2*sum_nblks 910 911 linfo_length = nblkrows_total + 1 + 2*MAXVAL(val_data(1, :)) 912 913 ALLOCATE (linfo_vec(linfo_length)) 914 ALLOCATE (rowp((nblkrows_total + 1)*job_count)) 915 ALLOCATE (coli(nblks)) 916 ALLOCATE (blkp(nblks)) 917 DO i = 1, share_size 918 order = val_data(3, i) 919 cur_blks = val_data(1, i) 920 IF (order .EQ. 0) THEN 921 localinfo_offset = file_size 922 localinfo_length = 0 923 ELSE 924 localinfo_offset = subh3_offsets(order) 925 localinfo_length = linfo_lens(order) 926 END IF 927 CALL mp_file_read_at_all(thefile, localinfo_offset, linfo_vec, msglen=localinfo_length) 928 IF (localinfo_length .EQ. 0) EXIT 929 930 rowp((i - 1)*(nblkrows_total + 1) + 1:i*(nblkrows_total + 1)) = linfo_vec(1:nblkrows_total + 1) 931 start_index = SUM(val_data(1, 1:i - 1)) + 1 932 end_index = SUM(val_data(1, 1:i)) 933 coli(start_index:end_index) = & 934 linfo_vec(nblkrows_total + 2:cur_blks + nblkrows_total + 1) 935 blkp(start_index:end_index) = & 936 linfo_vec(cur_blks + nblkrows_total + 2:2*cur_blks + nblkrows_total + 1) 937 END DO 938 row_p => rowp 939 col_i => coli 940 blk_p => blkp 941 942! compute the offsets and read block data 943 ALLOCATE (bdata_lens(size_of_pgrid)) 944 ALLOCATE (bdata_disps(size_of_pgrid)) 945 ALLOCATE (bdata_offsets(size_of_pgrid)) 946 bdata_lens = 0 947 DO i = 1, size_of_pgrid 948 DO j = 1, share_size 949 order = val_data(3, j) 950 IF (i .EQ. order) bdata_lens(order) = val_data(2, j) 951 END DO 952 END DO 953 CALL mp_sum(bdata_lens, group_id) 954 CALL cumsum_l(INT(bdata_lens, kind=file_offset), bdata_disps) 955 bdata_disps(:) = CSHIFT(bdata_disps, shift=-1) 956 bdata_disps(1) = BOF 957 958 bdata_start = subh3_start + subh3_length*mpi_integer_size 959 SELECT CASE (data_type) 960 CASE (dbcsr_type_real_4) 961 type_size = real_4_size 962 CASE (dbcsr_type_real_8) 963 type_size = real_8_size 964 CASE (dbcsr_type_complex_4) 965 type_size = 2*real_4_size 966 CASE (dbcsr_type_complex_8) 967 type_size = 2*real_8_size 968 END SELECT 969 bdata_offsets(:) = bdata_start + bdata_disps*type_size 970 971 SELECT CASE (data_type) 972 CASE (dbcsr_type_real_4) 973 ALLOCATE (rsp(darea_size)) 974 DO i = 1, share_size 975 order = val_data(3, i) 976 ! use dummy one-sized data array as buffer in place of empty array 977 ! when nothing is supposed to be read (order = 0) 978 IF (order .EQ. 0) THEN 979 blockdata_offset = file_size 980 CALL mp_file_read_at_all(thefile, blockdata_offset, rsp_dummy) 981 ELSE 982 start_index = SUM(val_data(2, 1:i - 1)) + 1 983 end_index = SUM(val_data(2, 1:i)) 984 blockdata_length = bdata_lens(order) 985 blockdata_offset = bdata_offsets(order) 986 CALL mp_file_read_at_all(thefile, blockdata_offset, rsp(start_index:end_index), & 987 msglen=blockdata_length) 988 END IF 989 END DO 990 CASE (dbcsr_type_real_8) 991 ALLOCATE (rdp(darea_size)) 992 DO i = 1, share_size 993 order = val_data(3, i) 994 IF (order .EQ. 0) THEN 995 blockdata_offset = file_size 996 CALL mp_file_read_at_all(thefile, blockdata_offset, rdp_dummy) 997 ELSE 998 start_index = SUM(val_data(2, 1:i - 1)) + 1 999 end_index = SUM(val_data(2, 1:i)) 1000 blockdata_length = bdata_lens(order) 1001 blockdata_offset = bdata_offsets(order) 1002 CALL mp_file_read_at_all(thefile, blockdata_offset, rdp(start_index:end_index), & 1003 msglen=blockdata_length) 1004 END IF 1005 END DO 1006 CASE (dbcsr_type_complex_4) 1007 ALLOCATE (csp(darea_size)) 1008 DO i = 1, share_size 1009 order = val_data(3, i) 1010 IF (order .EQ. 0) THEN 1011 blockdata_offset = file_size 1012 CALL mp_file_read_at_all(thefile, blockdata_offset, csp_dummy) 1013 ELSE 1014 start_index = SUM(val_data(2, 1:i - 1)) + 1 1015 end_index = SUM(val_data(2, 1:i)) 1016 blockdata_length = bdata_lens(order) 1017 blockdata_offset = bdata_offsets(order) 1018 CALL mp_file_read_at_all(thefile, blockdata_offset, csp(start_index:end_index), & 1019 msglen=blockdata_length) 1020 END IF 1021 END DO 1022 CASE (dbcsr_type_complex_8) 1023 ALLOCATE (cdp(darea_size)) 1024 DO i = 1, share_size 1025 order = val_data(3, i) 1026 IF (order .EQ. 0) THEN 1027 blockdata_offset = file_size 1028 CALL mp_file_read_at_all(thefile, blockdata_offset, cdp_dummy) 1029 ELSE 1030 start_index = SUM(val_data(2, 1:i - 1)) + 1 1031 end_index = SUM(val_data(2, 1:i)) 1032 blockdata_length = bdata_lens(order) 1033 blockdata_offset = bdata_offsets(order) 1034 CALL mp_file_read_at_all(thefile, blockdata_offset, cdp(start_index:end_index), & 1035 msglen=blockdata_length) 1036 END IF 1037 END DO 1038 END SELECT 1039 CALL dbcsr_data_init(dblk) 1040 CALL dbcsr_data_new(dblk, data_type) 1041 IF (ALLOCATED(rdp)) dblk%d%r_dp => rdp 1042 IF (ALLOCATED(rsp)) dblk%d%r_sp => rsp 1043 IF (ALLOCATED(cdp)) dblk%d%c_dp => cdp 1044 IF (ALLOCATED(csp)) dblk%d%c_sp => csp 1045 1046 CALL mp_file_close(thefile) 1047 1048 CALL dbcsr_create(matrix_new, matrix_name, distribution, matrix_type, & 1049 row_blk_size, col_blk_size, nze=darea_size, & 1050 data_type=data_type) 1051 CALL dbcsr_datablock_redistribute(dblk, row_p, col_i, blk_p, proc_nblks, proc_darea_sizes, matrix_new) 1052 1053 DEALLOCATE (subh2_offsets, subh3_offsets, bdata_offsets) 1054 DEALLOCATE (subh3_disps, bdata_disps) 1055 DEALLOCATE (linfo_lens, bdata_lens) 1056 DEALLOCATE (val_data, ginfo_vec, linfo_vec) 1057 DEALLOCATE (rowp, coli, blkp) 1058 IF (ALLOCATED(rdp)) DEALLOCATE (rdp) 1059 IF (ALLOCATED(rsp)) DEALLOCATE (rsp) 1060 IF (ALLOCATED(cdp)) DEALLOCATE (cdp) 1061 IF (ALLOCATED(csp)) DEALLOCATE (csp) 1062 CALL dbcsr_data_clear_pointer(dblk) 1063 DEALLOCATE (dblk%d) 1064 1065 CALL timestop(handle) 1066 CONTAINS 1067 SUBROUTINE cumsum_l(arr, cumsum) 1068 INTEGER(kind=file_offset), DIMENSION(:), & 1069 INTENT(IN) :: arr 1070 INTEGER(kind=file_offset), DIMENSION(SIZE(arr)), & 1071 INTENT(OUT) :: cumsum 1072 1073 INTEGER :: i 1074 1075 cumsum(1) = arr(1) 1076 DO i = 2, SIZE(arr) 1077 cumsum(i) = cumsum(i - 1) + arr(i) 1078 END DO 1079 END SUBROUTINE cumsum_l 1080 1081 END SUBROUTINE dbcsr_binary_read 1082 1083 SUBROUTINE dbcsr_print_block_sum(matrix, unit_nr) 1084 !! Prints the sum of the elements for each block 1085 1086 TYPE(dbcsr_type), INTENT(IN) :: matrix 1087 !! matrix 1088 INTEGER, INTENT(IN), OPTIONAL :: unit_nr 1089 1090 CHARACTER(len=*), PARAMETER :: routineN = 'dbcsr_print_block_sum' 1091 1092 COMPLEX(KIND=real_4) :: blk_sum_c_sp 1093 COMPLEX(KIND=real_4), DIMENSION(:), POINTER :: c_sp 1094 COMPLEX(KIND=real_8) :: blk_sum_c_dp 1095 COMPLEX(KIND=real_8), DIMENSION(:), POINTER :: c_dp 1096 INTEGER :: bc, blk, blk_p, br, handle, iunit, m, & 1097 mn, n 1098 INTEGER, DIMENSION(:), POINTER :: col_blk_offset, col_blk_size, & 1099 row_blk_offset, row_blk_size 1100 REAL(KIND=real_4) :: blk_sum_r_sp 1101 REAL(KIND=real_4), DIMENSION(:), POINTER :: r_sp 1102 REAL(KIND=real_8) :: blk_sum_r_dp 1103 REAL(KIND=real_8), DIMENSION(:), POINTER :: r_dp 1104 1105! --------------------------------------------------------------------------- 1106 1107 CALL timeset(routineN, handle) 1108 IF (.NOT. dbcsr_valid_index(matrix)) & 1109 DBCSR_WARN("Can not print invalid matrix.") 1110 1111 iunit = default_output_unit 1112 IF (PRESENT(unit_nr)) iunit = unit_nr 1113 1114 IF (iunit > 0) THEN 1115 1116 SELECT CASE (matrix%data_type) 1117 CASE (dbcsr_type_real_8) 1118 CALL dbcsr_get_data(matrix%data_area, r_dp) 1119 CASE (dbcsr_type_real_4) 1120 CALL dbcsr_get_data(matrix%data_area, r_sp) 1121 CASE (dbcsr_type_complex_8) 1122 CALL dbcsr_get_data(matrix%data_area, c_dp) 1123 CASE (dbcsr_type_complex_4) 1124 CALL dbcsr_get_data(matrix%data_area, c_sp) 1125 END SELECT 1126 row_blk_size => array_data(matrix%row_blk_size) 1127 col_blk_size => array_data(matrix%col_blk_size) 1128 row_blk_offset => array_data(matrix%row_blk_offset) 1129 col_blk_offset => array_data(matrix%col_blk_offset) 1130 1131 IF (matrix%nblks .GT. 0) THEN 1132 DO br = 1, matrix%nblkrows_total 1133 m = row_blk_size(br) 1134 DO blk = matrix%row_p(br) + 1, matrix%row_p(br + 1) 1135 bc = matrix%col_i(blk) 1136 n = col_blk_size(bc) 1137 mn = m*n 1138 blk_p = ABS(matrix%blk_p(blk)) 1139 block_exists: IF (blk_p .NE. 0) THEN 1140 IF (mn .GT. 0) THEN 1141 SELECT CASE (matrix%data_type) 1142 CASE (dbcsr_type_real_8) 1143 blk_sum_r_dp = SUM(r_dp(blk_p:blk_p + mn - 1)) 1144 WRITE (iunit, '(I6,I6,ES18.9)') & 1145 br, bc, blk_sum_r_dp 1146 CASE (dbcsr_type_real_4) 1147 blk_sum_r_sp = SUM(r_sp(blk_p:blk_p + mn - 1)) 1148 WRITE (iunit, '(I6,I6,ES18.9)') & 1149 br, bc, blk_sum_r_sp 1150 CASE (dbcsr_type_complex_8) 1151 blk_sum_c_dp = SUM(c_dp(blk_p:blk_p + mn - 1)) 1152 WRITE (iunit, '(I6,I6,ES18.9," I*",ES18.9)') & 1153 br, bc, REAL(blk_sum_c_dp), AIMAG(blk_sum_c_dp) 1154 CASE (dbcsr_type_complex_4) 1155 blk_sum_c_sp = SUM(c_sp(blk_p:blk_p + mn - 1)) 1156 WRITE (iunit, '(I6,I6,ES18.9," I*",ES18.9)') & 1157 br, bc, REAL(blk_sum_c_sp), AIMAG(blk_sum_c_sp) 1158 END SELECT 1159 ELSE 1160 blk_sum_r_dp = 0.0_dp 1161 WRITE (iunit, '(I6,I6,ES18.9)') & 1162 br, bc, blk_sum_r_dp 1163 ENDIF 1164 ENDIF block_exists 1165 ENDDO 1166 ENDDO 1167 ENDIF 1168 1169 ENDIF ! unit > 0 1170 1171 CALL timestop(handle) 1172 1173 END SUBROUTINE dbcsr_print_block_sum 1174 1175END MODULE dbcsr_io 1176