453 INTEGER,
INTENT(IN),
OPTIONAL :: ncol, start_col
455 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_fm_init_random'
457 INTEGER :: handle, icol_global, icol_local, irow_local, my_ncol, my_start_col, ncol_global, &
458 ncol_local, nrow_global, nrow_local
459 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
460 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: buff
461 REAL(kind=
dp),
CONTIGUOUS,
DIMENSION(:, :), &
462 POINTER :: local_data
463 REAL(kind=
dp),
DIMENSION(3, 2),
SAVE :: &
464 seed = reshape((/1.0_dp, 2.0_dp, 3.0_dp, 4.0_dp, 5.0_dp, 6.0_dp/), (/3, 2/))
467 CALL timeset(routinen, handle)
470 CALL matrix%matrix_struct%para_env%bcast(seed, 0)
473 extended_precision=.true., seed=seed)
475 cpassert(.NOT. matrix%use_sp)
477 CALL cp_fm_get_info(matrix, nrow_global=nrow_global, ncol_global=ncol_global, &
478 nrow_local=nrow_local, ncol_local=ncol_local, &
479 local_data=local_data, &
480 row_indices=row_indices, col_indices=col_indices)
483 IF (
PRESENT(start_col)) my_start_col = start_col
484 my_ncol = matrix%matrix_struct%ncol_global
485 IF (
PRESENT(ncol)) my_ncol = ncol
487 IF (ncol_global < (my_start_col + my_ncol - 1)) &
488 cpabort(
"ncol_global>=(my_start_col+my_ncol-1)")
490 ALLOCATE (buff(nrow_global))
496 DO icol_local = 1, ncol_local
497 cpassert(col_indices(icol_local) > icol_global)
499 CALL rng%reset_to_next_substream()
500 icol_global = icol_global + 1
501 IF (icol_global == col_indices(icol_local))
EXIT
504 DO irow_local = 1, nrow_local
505 local_data(irow_local, icol_local) = buff(row_indices(irow_local))
518 CALL rng%get(ig=seed)
520 CALL timestop(handle)
767 start_col, n_rows, n_cols, alpha, beta, transpose)
769 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(in) :: new_values
770 INTEGER,
INTENT(in),
OPTIONAL :: start_row, start_col, n_rows, n_cols
771 REAL(kind=
dp),
INTENT(in),
OPTIONAL :: alpha, beta
772 LOGICAL,
INTENT(in),
OPTIONAL :: transpose
774 INTEGER :: i, i0, j, j0, ncol, ncol_block, &
775 ncol_global, ncol_local, nrow, &
776 nrow_block, nrow_global, nrow_local, &
778 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
780 REAL(kind=
dp) :: al, be
781 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: full_block
783 al = 1.0_dp; be = 0.0_dp; i0 = 1; j0 = 1; tr_a = .false.
787 cpassert(.NOT. fm%use_sp)
789 IF (
PRESENT(alpha)) al = alpha
790 IF (
PRESENT(beta)) be = beta
791 IF (
PRESENT(start_row)) i0 = start_row
792 IF (
PRESENT(start_col)) j0 = start_col
793 IF (
PRESENT(transpose)) tr_a = transpose
795 nrow =
SIZE(new_values, 2)
796 ncol =
SIZE(new_values, 1)
798 nrow =
SIZE(new_values, 1)
799 ncol =
SIZE(new_values, 2)
801 IF (
PRESENT(n_rows)) nrow = n_rows
802 IF (
PRESENT(n_cols)) ncol = n_cols
804 full_block => fm%local_data
807 nrow_global=nrow_global, ncol_global=ncol_global, &
808 nrow_block=nrow_block, ncol_block=ncol_block, &
809 nrow_local=nrow_local, ncol_local=ncol_local, &
810 row_indices=row_indices, col_indices=col_indices)
812 IF (al == 1.0 .AND. be == 0.0)
THEN
814 this_col = col_indices(j) - j0 + 1
815 IF (this_col .GE. 1 .AND. this_col .LE. ncol)
THEN
817 IF (i0 == 1 .AND. nrow_global == nrow)
THEN
819 full_block(i, j) = new_values(this_col, row_indices(i))
823 this_row = row_indices(i) - i0 + 1
824 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
825 full_block(i, j) = new_values(this_col, this_row)
830 IF (i0 == 1 .AND. nrow_global == nrow)
THEN
832 full_block(i, j) = new_values(row_indices(i), this_col)
836 this_row = row_indices(i) - i0 + 1
837 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
838 full_block(i, j) = new_values(this_row, this_col)
847 this_col = col_indices(j) - j0 + 1
848 IF (this_col .GE. 1 .AND. this_col .LE. ncol)
THEN
851 this_row = row_indices(i) - i0 + 1
852 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
853 full_block(i, j) = al*new_values(this_col, this_row) + &
859 this_row = row_indices(i) - i0 + 1
860 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
861 full_block(i, j) = al*new_values(this_row, this_col) + &
900 start_col, n_rows, n_cols, transpose)
902 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(out) :: target_m
903 INTEGER,
INTENT(in),
OPTIONAL :: start_row, start_col, n_rows, n_cols
904 LOGICAL,
INTENT(in),
OPTIONAL :: transpose
906 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_fm_get_submatrix'
908 INTEGER :: handle, i, i0, j, j0, ncol, ncol_global, &
909 ncol_local, nrow, nrow_global, &
910 nrow_local, this_col, this_row
911 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
913 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: full_block
916 CALL timeset(routinen, handle)
918 i0 = 1; j0 = 1; tr_a = .false.
920 cpassert(.NOT. fm%use_sp)
922 IF (
PRESENT(start_row)) i0 = start_row
923 IF (
PRESENT(start_col)) j0 = start_col
924 IF (
PRESENT(transpose)) tr_a = transpose
926 nrow =
SIZE(target_m, 2)
927 ncol =
SIZE(target_m, 1)
929 nrow =
SIZE(target_m, 1)
930 ncol =
SIZE(target_m, 2)
932 IF (
PRESENT(n_rows)) nrow = n_rows
933 IF (
PRESENT(n_cols)) ncol = n_cols
935 para_env => fm%matrix_struct%para_env
937 full_block => fm%local_data
938#if defined(__SCALAPACK)
940 IF (
SIZE(target_m, 1)*
SIZE(target_m, 2) /= 0)
THEN
941 CALL dcopy(
SIZE(target_m, 1)*
SIZE(target_m, 2), [0.0_dp], 0, target_m, 1)
946 nrow_global=nrow_global, ncol_global=ncol_global, &
947 nrow_local=nrow_local, ncol_local=ncol_local, &
948 row_indices=row_indices, col_indices=col_indices)
951 this_col = col_indices(j) - j0 + 1
952 IF (this_col .GE. 1 .AND. this_col .LE. ncol)
THEN
954 IF (i0 == 1 .AND. nrow_global == nrow)
THEN
956 target_m(this_col, row_indices(i)) = full_block(i, j)
960 this_row = row_indices(i) - i0 + 1
961 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
962 target_m(this_col, this_row) = full_block(i, j)
967 IF (i0 == 1 .AND. nrow_global == nrow)
THEN
969 target_m(row_indices(i), this_col) = full_block(i, j)
973 this_row = row_indices(i) - i0 + 1
974 IF (this_row >= 1 .AND. this_row <= nrow)
THEN
975 target_m(this_row, this_col) = full_block(i, j)
983 CALL para_env%sum(target_m)
985 CALL timestop(handle)
1065 REAL(kind=
dp),
INTENT(OUT) :: a_max
1066 INTEGER,
INTENT(OUT),
OPTIONAL :: ir_max, ic_max
1068 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_fm_maxabsval'
1070 INTEGER :: handle, i, ic_max_local, ir_max_local, &
1071 j, mepos, ncol_local, nrow_local, &
1073 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: ic_max_vec, ir_max_vec
1074 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
1076 REAL(
dp),
ALLOCATABLE,
DIMENSION(:) :: a_max_vec
1077 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: my_block
1078 REAL(kind=
sp),
DIMENSION(:, :),
POINTER :: my_block_sp
1080 CALL timeset(routinen, handle)
1082 my_block => matrix%local_data
1083 my_block_sp => matrix%local_data_sp
1085 CALL cp_fm_get_info(matrix, nrow_local=nrow_local, ncol_local=ncol_local, &
1086 row_indices=row_indices, col_indices=col_indices)
1088 IF (matrix%use_sp)
THEN
1089 a_max = real(maxval(abs(my_block_sp(1:nrow_local, 1:ncol_local))),
dp)
1091 a_max = maxval(abs(my_block(1:nrow_local, 1:ncol_local)))
1094 IF (
PRESENT(ir_max))
THEN
1095 num_pe = matrix%matrix_struct%para_env%num_pe
1096 mepos = matrix%matrix_struct%para_env%mepos
1097 ALLOCATE (ir_max_vec(0:num_pe - 1))
1098 ir_max_vec(0:num_pe - 1) = 0
1099 ALLOCATE (ic_max_vec(0:num_pe - 1))
1100 ic_max_vec(0:num_pe - 1) = 0
1101 ALLOCATE (a_max_vec(0:num_pe - 1))
1102 a_max_vec(0:num_pe - 1) = 0.0_dp
1105 IF ((ncol_local > 0) .AND. (nrow_local > 0))
THEN
1106 DO i = 1, ncol_local
1107 DO j = 1, nrow_local
1108 IF (matrix%use_sp)
THEN
1109 IF (abs(my_block_sp(j, i)) > my_max)
THEN
1110 my_max = real(my_block_sp(j, i),
dp)
1115 IF (abs(my_block(j, i)) > my_max)
THEN
1116 my_max = my_block(j, i)
1124 a_max_vec(mepos) = my_max
1125 ir_max_vec(mepos) = row_indices(ir_max_local)
1126 ic_max_vec(mepos) = col_indices(ic_max_local)
1130 CALL matrix%matrix_struct%para_env%sum(a_max_vec)
1131 CALL matrix%matrix_struct%para_env%sum(ir_max_vec)
1132 CALL matrix%matrix_struct%para_env%sum(ic_max_vec)
1135 DO i = 0, num_pe - 1
1136 IF (a_max_vec(i) > my_max)
THEN
1137 ir_max = ir_max_vec(i)
1138 ic_max = ic_max_vec(i)
1142 DEALLOCATE (ir_max_vec, ic_max_vec, a_max_vec)
1143 cpassert(ic_max > 0)
1144 cpassert(ir_max > 0)
1148 CALL matrix%matrix_struct%para_env%max(a_max)
1150 CALL timestop(handle)
1311 SUBROUTINE cp_fm_to_fm_matrix(source, destination)
1313 TYPE(
cp_fm_type),
INTENT(IN) :: source, destination
1315 CHARACTER(LEN=*),
PARAMETER :: routinen =
'cp_fm_to_fm_matrix'
1317 INTEGER :: handle, npcol, nprow
1319 CALL timeset(routinen, handle)
1321 nprow = source%matrix_struct%context%num_pe(1)
1322 npcol = source%matrix_struct%context%num_pe(2)
1326 destination%matrix_struct))
THEN
1327 IF (source%use_sp .AND. destination%use_sp)
THEN
1328 IF (
SIZE(source%local_data_sp, 1) /=
SIZE(destination%local_data_sp, 1) .OR. &
1329 SIZE(source%local_data_sp, 2) /=
SIZE(destination%local_data_sp, 2)) &
1330 CALL cp_abort(__location__, &
1331 "Cannot copy full matrix <"//trim(source%name)// &
1332 "> to full matrix <"//trim(destination%name)// &
1333 ">. The local_data blocks have different sizes.")
1334 CALL scopy(
SIZE(source%local_data_sp, 1)*
SIZE(source%local_data_sp, 2), &
1335 source%local_data_sp, 1, destination%local_data_sp, 1)
1336 ELSE IF (source%use_sp .AND. .NOT. destination%use_sp)
THEN
1337 IF (
SIZE(source%local_data_sp, 1) /=
SIZE(destination%local_data, 1) .OR. &
1338 SIZE(source%local_data_sp, 2) /=
SIZE(destination%local_data, 2)) &
1339 CALL cp_abort(__location__, &
1340 "Cannot copy full matrix <"//trim(source%name)// &
1341 "> to full matrix <"//trim(destination%name)// &
1342 ">. The local_data blocks have different sizes.")
1343 destination%local_data = real(source%local_data_sp,
dp)
1344 ELSE IF (.NOT. source%use_sp .AND. destination%use_sp)
THEN
1345 IF (
SIZE(source%local_data, 1) /=
SIZE(destination%local_data_sp, 1) .OR. &
1346 SIZE(source%local_data, 2) /=
SIZE(destination%local_data_sp, 2)) &
1347 CALL cp_abort(__location__, &
1348 "Cannot copy full matrix <"//trim(source%name)// &
1349 "> to full matrix <"//trim(destination%name)// &
1350 ">. The local_data blocks have different sizes.")
1351 destination%local_data_sp = real(source%local_data,
sp)
1353 IF (
SIZE(source%local_data, 1) /=
SIZE(destination%local_data, 1) .OR. &
1354 SIZE(source%local_data, 2) /=
SIZE(destination%local_data, 2)) &
1355 CALL cp_abort(__location__, &
1356 "Cannot copy full matrix <"//trim(source%name)// &
1357 "> to full matrix <"//trim(destination%name)// &
1358 ">. The local_data blocks have different sizes.")
1359 CALL dcopy(
SIZE(source%local_data, 1)*
SIZE(source%local_data, 2), &
1360 source%local_data, 1, destination%local_data, 1)
1363 cpabort(
"Data structures of source and target full matrix are not equivalent")
1366 CALL timestop(handle)
1568 TYPE(
cp_fm_type),
INTENT(IN) :: source, destination
1572 CHARACTER(LEN=*),
PARAMETER :: routinen =
'cp_fm_start_copy_general'
1574 INTEGER :: dest_p_i, dest_q_j, global_rank, global_size, handle, i, j, k, mpi_rank, &
1575 ncol_block_dest, ncol_block_src, ncol_local_recv, ncol_local_send, ncols, &
1576 nrow_block_dest, nrow_block_src, nrow_local_recv, nrow_local_send, nrows, p, q, &
1577 recv_rank, recv_size, send_rank, send_size
1578 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: all_ranks, dest2global, dest_p, dest_q, &
1579 recv_count, send_count, send_disp, &
1580 source2global, src_p, src_q
1581 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: dest_blacs2mpi
1582 INTEGER,
DIMENSION(2) :: dest_block, dest_block_tmp, dest_num_pe, &
1583 src_block, src_block_tmp, src_num_pe
1584 INTEGER,
DIMENSION(:),
POINTER :: recv_col_indices, recv_row_indices, &
1585 send_col_indices, send_row_indices
1589 CALL timeset(routinen, handle)
1593 nrow_local_send =
SIZE(source%local_data, 1)
1594 ncol_local_send =
SIZE(source%local_data, 2)
1595 ALLOCATE (info%send_buf(nrow_local_send*ncol_local_send))
1597 DO j = 1, ncol_local_send
1598 DO i = 1, nrow_local_send
1600 info%send_buf(k) = source%local_data(i, j)
1605 IF (source%use_sp) cpabort(
"only DP kind implemented")
1606 IF (destination%use_sp) cpabort(
"only DP kind implemented")
1608 NULLIFY (recv_dist, send_dist)
1609 NULLIFY (recv_col_indices, recv_row_indices, send_col_indices, send_row_indices)
1612 global_size = para_env%num_pe
1613 global_rank = para_env%mepos
1619 IF (
ASSOCIATED(destination%matrix_struct))
THEN
1620 recv_dist => destination%matrix_struct
1621 recv_rank = recv_dist%para_env%mepos
1626 IF (
ASSOCIATED(source%matrix_struct))
THEN
1627 send_dist => source%matrix_struct
1628 send_rank = send_dist%para_env%mepos
1634 ALLOCATE (all_ranks(0:global_size - 1))
1636 CALL para_env%allgather(send_rank, all_ranks)
1637 IF (
ASSOCIATED(recv_dist))
THEN
1638 ALLOCATE (source2global(0:count(all_ranks .NE.
mp_proc_null) - 1))
1639 DO i = 0, global_size - 1
1641 source2global(all_ranks(i)) = i
1646 CALL para_env%allgather(recv_rank, all_ranks)
1647 IF (
ASSOCIATED(send_dist))
THEN
1648 ALLOCATE (dest2global(0:count(all_ranks .NE.
mp_proc_null) - 1))
1649 DO i = 0, global_size - 1
1651 dest2global(all_ranks(i)) = i
1655 DEALLOCATE (all_ranks)
1662 IF (global_rank == 0)
THEN
1664 CALL para_env%irecv(src_block,
mp_any_source, recv_req(1), tag=src_tag)
1665 CALL para_env%irecv(dest_block,
mp_any_source, recv_req(2), tag=dest_tag)
1666 CALL para_env%irecv(src_num_pe,
mp_any_source, recv_req(3), tag=src_tag)
1667 CALL para_env%irecv(dest_num_pe,
mp_any_source, recv_req(4), tag=dest_tag)
1670 IF (
ASSOCIATED(send_dist))
THEN
1671 IF ((send_rank .EQ. 0))
THEN
1673 src_block_tmp = (/send_dist%nrow_block, send_dist%ncol_block/)
1674 CALL para_env%isend(src_block_tmp, 0, send_req(1), tag=src_tag)
1675 CALL para_env%isend(send_dist%context%num_pe, 0, send_req(2), tag=src_tag)
1679 IF (
ASSOCIATED(recv_dist))
THEN
1680 IF ((recv_rank .EQ. 0))
THEN
1681 dest_block_tmp = (/recv_dist%nrow_block, recv_dist%ncol_block/)
1682 CALL para_env%isend(dest_block_tmp, 0, send_req(3), tag=dest_tag)
1683 CALL para_env%isend(recv_dist%context%num_pe, 0, send_req(4), tag=dest_tag)
1687 IF (global_rank == 0)
THEN
1690 ALLOCATE (info%src_blacs2mpi(0:src_num_pe(1) - 1, 0:src_num_pe(2) - 1), &
1691 dest_blacs2mpi(0:dest_num_pe(1) - 1, 0:dest_num_pe(2) - 1) &
1693 CALL para_env%irecv(info%src_blacs2mpi,
mp_any_source, recv_req(5), tag=src_tag)
1694 CALL para_env%irecv(dest_blacs2mpi,
mp_any_source, recv_req(6), tag=dest_tag)
1697 IF (
ASSOCIATED(send_dist))
THEN
1698 IF ((send_rank .EQ. 0))
THEN
1699 CALL para_env%isend(send_dist%context%blacs2mpi(:, :), 0, send_req(5), tag=src_tag)
1703 IF (
ASSOCIATED(recv_dist))
THEN
1704 IF ((recv_rank .EQ. 0))
THEN
1705 CALL para_env%isend(recv_dist%context%blacs2mpi(:, :), 0, send_req(6), tag=dest_tag)
1709 IF (global_rank == 0)
THEN
1714 CALL para_env%bcast(src_block, 0)
1715 CALL para_env%bcast(dest_block, 0)
1716 CALL para_env%bcast(src_num_pe, 0)
1717 CALL para_env%bcast(dest_num_pe, 0)
1718 info%src_num_pe(1:2) = src_num_pe(1:2)
1719 info%nblock_src(1:2) = src_block(1:2)
1720 IF (global_rank .NE. 0)
THEN
1721 ALLOCATE (info%src_blacs2mpi(0:src_num_pe(1) - 1, 0:src_num_pe(2) - 1), &
1722 dest_blacs2mpi(0:dest_num_pe(1) - 1, 0:dest_num_pe(2) - 1) &
1725 CALL para_env%bcast(info%src_blacs2mpi, 0)
1726 CALL para_env%bcast(dest_blacs2mpi, 0)
1728 recv_size = dest_num_pe(1)*dest_num_pe(2)
1729 send_size = src_num_pe(1)*src_num_pe(2)
1730 info%send_size = send_size
1748 IF (
ASSOCIATED(recv_dist))
THEN
1750 col_indices=recv_col_indices &
1752 info%recv_col_indices => recv_col_indices
1753 info%recv_row_indices => recv_row_indices
1754 nrow_block_src = src_block(1)
1755 ncol_block_src = src_block(2)
1756 ALLOCATE (recv_count(0:send_size - 1), info%recv_disp(0:send_size - 1), info%recv_request(0:send_size - 1))
1759 nrow_local_recv = recv_dist%nrow_locals(recv_dist%context%mepos(1))
1760 ncol_local_recv = recv_dist%ncol_locals(recv_dist%context%mepos(2))
1761 info%nlocal_recv(1) = nrow_local_recv
1762 info%nlocal_recv(2) = ncol_local_recv
1764 ALLOCATE (src_p(nrow_local_recv), src_q(ncol_local_recv))
1765 DO i = 1, nrow_local_recv
1768 src_p(i) = mod(((recv_row_indices(i) - 1)/nrow_block_src), src_num_pe(1))
1770 DO j = 1, ncol_local_recv
1772 src_q(j) = mod(((recv_col_indices(j) - 1)/ncol_block_src), src_num_pe(2))
1776 DO q = 0, src_num_pe(2) - 1
1777 ncols = count(src_q .EQ. q)
1778 DO p = 0, src_num_pe(1) - 1
1779 nrows = count(src_p .EQ. p)
1781 recv_count(info%src_blacs2mpi(p, q)) = nrows*ncols
1784 DEALLOCATE (src_p, src_q)
1788 ALLOCATE (info%recv_buf(sum(recv_count(:))))
1789 info%recv_disp(0) = 0
1790 DO i = 1, send_size - 1
1791 info%recv_disp(i) = info%recv_disp(i - 1) + recv_count(i - 1)
1795 DO k = 0, send_size - 1
1796 IF (recv_count(k) .GT. 0)
THEN
1797 CALL para_env%irecv(info%recv_buf(info%recv_disp(k) + 1:info%recv_disp(k) + recv_count(k)), &
1798 source2global(k), info%recv_request(k))
1801 DEALLOCATE (source2global)
1805 IF (
ASSOCIATED(send_dist))
THEN
1807 col_indices=send_col_indices &
1809 nrow_block_dest = dest_block(1)
1810 ncol_block_dest = dest_block(2)
1811 ALLOCATE (send_count(0:recv_size - 1), send_disp(0:recv_size - 1), info%send_request(0:recv_size - 1))
1814 nrow_local_send = send_dist%nrow_locals(send_dist%context%mepos(1))
1815 ncol_local_send = send_dist%ncol_locals(send_dist%context%mepos(2))
1819 ALLOCATE (dest_p(nrow_local_send), dest_q(ncol_local_send))
1821 DO i = 1, nrow_local_send
1823 dest_p(i) = mod(((send_row_indices(i) - 1)/nrow_block_dest), dest_num_pe(1))
1825 DO j = 1, ncol_local_send
1826 dest_q(j) = mod(((send_col_indices(j) - 1)/ncol_block_dest), dest_num_pe(2))
1830 DO q = 0, dest_num_pe(2) - 1
1831 ncols = count(dest_q .EQ. q)
1832 DO p = 0, dest_num_pe(1) - 1
1833 nrows = count(dest_p .EQ. p)
1834 send_count(dest_blacs2mpi(p, q)) = nrows*ncols
1837 DEALLOCATE (dest_p, dest_q)
1840 ALLOCATE (info%send_buf(sum(send_count(:))))
1842 DO k = 1, recv_size - 1
1843 send_disp(k) = send_disp(k - 1) + send_count(k - 1)
1848 DO j = 1, ncol_local_send
1850 dest_q_j = mod(((send_col_indices(j) - 1)/ncol_block_dest), dest_num_pe(2))
1851 DO i = 1, nrow_local_send
1852 dest_p_i = mod(((send_row_indices(i) - 1)/nrow_block_dest), dest_num_pe(1))
1853 mpi_rank = dest_blacs2mpi(dest_p_i, dest_q_j)
1854 send_count(mpi_rank) = send_count(mpi_rank) + 1
1855 info%send_buf(send_disp(mpi_rank) + send_count(mpi_rank)) = source%local_data(i, j)
1860 DO k = 0, recv_size - 1
1861 IF (send_count(k) .GT. 0)
THEN
1862 CALL para_env%isend(info%send_buf(send_disp(k) + 1:send_disp(k) + send_count(k)), &
1863 dest2global(k), info%send_request(k))
1866 DEALLOCATE (send_count, send_disp, dest2global)
1868 DEALLOCATE (dest_blacs2mpi)
1872 CALL timestop(handle)
2148 INTEGER,
INTENT(IN) :: unit
2150 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_fm_write_unformatted'
2152 INTEGER :: handle, j, max_block, &
2153 ncol_global, nrow_global
2155#if defined(__SCALAPACK)
2156 INTEGER :: i, i_block, icol_local, &
2158 iprow, irow_local, &
2161 INTEGER,
DIMENSION(9) :: desc
2162 REAL(kind=
dp),
DIMENSION(:),
POINTER :: vecbuf
2163 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: newdat
2165 INTEGER,
EXTERNAL :: numroc
2168 CALL timeset(routinen, handle)
2169 CALL cp_fm_get_info(fm, nrow_global=nrow_global, ncol_global=ncol_global, ncol_block=max_block, &
2172#if defined(__SCALAPACK)
2173 num_pe = para_env%num_pe
2174 mepos = para_env%mepos
2178 CALL ictxt_loc%gridinit(para_env,
'R', 1, num_pe)
2179 CALL descinit(desc, nrow_global, ncol_global, rb, max_block, 0, 0, ictxt_loc%get_handle(), nrow_global, info)
2181 associate(nprow => ictxt_loc%num_pe(1), npcol => ictxt_loc%num_pe(2), &
2182 myprow => ictxt_loc%mepos(1), mypcol => ictxt_loc%mepos(2))
2183 in = numroc(ncol_global, max_block, mypcol, 0, npcol)
2185 ALLOCATE (newdat(nrow_global, max(1, in)))
2188 CALL pdgemr2d(nrow_global, ncol_global, fm%local_data, 1, 1, &
2189 fm%matrix_struct%descriptor, &
2190 newdat, 1, 1, desc, ictxt_loc%get_handle())
2192 ALLOCATE (vecbuf(nrow_global*max_block))
2193 vecbuf = huge(1.0_dp)
2195 DO i = 1, ncol_global, max(max_block, 1)
2196 i_block = min(max_block, ncol_global - i + 1)
2197 CALL infog2l(1, i, desc, nprow, npcol, myprow, mypcol, &
2198 irow_local, icol_local, iprow, ipcol)
2199 IF (ipcol == mypcol)
THEN
2201 vecbuf((j - 1)*nrow_global + 1:nrow_global*j) = newdat(:, icol_local + j - 1)
2205 IF (ipcol == 0)
THEN
2208 IF (ipcol == mypcol)
THEN
2209 CALL para_env%send(vecbuf(:), 0, tag)
2211 IF (mypcol == 0)
THEN
2212 CALL para_env%recv(vecbuf(:), ipcol, tag)
2218 WRITE (unit) vecbuf((j - 1)*nrow_global + 1:nrow_global*j)
2226 CALL ictxt_loc%gridexit()
2233 DO j = 1, ncol_global
2234 WRITE (unit) fm%local_data(:, j)
2239 CALL timestop(handle)
2252 INTEGER,
INTENT(IN) :: unit
2253 CHARACTER(LEN=*),
INTENT(IN),
OPTIONAL ::
header, value_format
2255 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_fm_write_formatted'
2257 CHARACTER(LEN=21) :: my_value_format
2258 INTEGER :: handle, i, j, max_block, &
2259 ncol_global, nrow_global, &
2261 TYPE(mp_para_env_type),
POINTER :: para_env
2262#if defined(__SCALAPACK)
2263 INTEGER :: i_block, icol_local, &
2265 iprow, irow_local, &
2266 mepos, num_pe, rb, tag, k, &
2268 INTEGER,
DIMENSION(9) :: desc
2269 REAL(kind=dp),
DIMENSION(:),
POINTER :: vecbuf
2270 REAL(kind=dp),
DIMENSION(:, :),
POINTER :: newdat
2271 TYPE(cp_blacs_type) :: ictxt_loc
2272 INTEGER,
EXTERNAL :: numroc
2275 CALL timeset(routinen, handle)
2276 CALL cp_fm_get_info(fm, nrow_global=nrow_global, ncol_global=ncol_global, ncol_block=max_block, &
2277 nrow_block=nrow_block, para_env=para_env)
2279 IF (
PRESENT(value_format))
THEN
2280 cpassert(len_trim(adjustl(value_format)) < 11)
2281 my_value_format =
"(I10, I10, "//trim(adjustl(value_format))//
")"
2283 my_value_format =
"(I10, I10, ES24.12)"
2288 WRITE (unit,
"(A2, A8, A10, A24)")
"#",
"Row",
"Column", adjustl(
"Value")
2291#if defined(__SCALAPACK)
2292 num_pe = para_env%num_pe
2293 mepos = para_env%mepos
2297 CALL ictxt_loc%gridinit(para_env,
'R', 1, num_pe)
2298 CALL descinit(desc, nrow_global, ncol_global, rb, max_block, 0, 0, ictxt_loc%get_handle(), nrow_global, info)
2300 associate(nprow => ictxt_loc%num_pe(1), npcol => ictxt_loc%num_pe(2), &
2301 myprow => ictxt_loc%mepos(1), mypcol => ictxt_loc%mepos(2))
2302 in = numroc(ncol_global, max_block, mypcol, 0, npcol)
2304 ALLOCATE (newdat(nrow_global, max(1, in)))
2307 CALL pdgemr2d(nrow_global, ncol_global, fm%local_data, 1, 1, &
2308 fm%matrix_struct%descriptor, &
2309 newdat, 1, 1, desc, ictxt_loc%get_handle())
2311 ALLOCATE (vecbuf(nrow_global*max_block))
2312 vecbuf = huge(1.0_dp)
2316 DO i = 1, ncol_global, max(max_block, 1)
2317 i_block = min(max_block, ncol_global - i + 1)
2318 CALL infog2l(1, i, desc, nprow, npcol, myprow, mypcol, &
2319 irow_local, icol_local, iprow, ipcol)
2320 IF (ipcol == mypcol)
THEN
2322 vecbuf((j - 1)*nrow_global + 1:nrow_global*j) = newdat(:, icol_local + j - 1)
2326 IF (ipcol == 0)
THEN
2329 IF (ipcol == mypcol)
THEN
2330 CALL para_env%send(vecbuf(:), 0, tag)
2332 IF (mypcol == 0)
THEN
2333 CALL para_env%recv(vecbuf(:), ipcol, tag)
2339 DO k = (j - 1)*nrow_global + 1, nrow_global*j
2340 WRITE (unit=unit, fmt=my_value_format) irow, icol, vecbuf(k)
2342 IF (irow > nrow_global)
THEN
2354 CALL ictxt_loc%gridexit()
2361 DO j = 1, ncol_global
2362 DO i = 1, nrow_global
2363 WRITE (unit=unit, fmt=my_value_format) i, j, fm%local_data(i, j)
2369 CALL timestop(handle)