Eigenvalue SoLvers for Petaflop-Applications (ELPA) 2024.05.001
Loading...
Searching...
No Matches
transform_columns_template.F90
Go to the documentation of this file.
1subroutine transform_columns_&
2&precision&
3&(obj, col1, col2, na, tmp, l_rqs, l_rqe, q, ldq, matrixcols, &
4 l_rows, mpi_comm_cols, p_col, l_col, qtrans)
5 use precision
7#ifdef WITH_OPENMP_TRADITIONAL
8 use elpa_omp
9#endif
10 use elpa_mpi
11 implicit none
12 class(elpa_abstract_impl_t), intent(inout) :: obj
13 integer(kind=ik), intent(in) :: na, l_rqs, l_rqe, ldq, matrixCols
14 integer(kind=ik), intent(in) :: l_rows, mpi_comm_cols
15 integer(kind=ik), intent(in) :: p_col(na), l_col(na)
16#ifdef USE_ASSUMED_SIZE
17 real(kind=real_datatype), intent(inout) :: q(ldq,*)
18#else
19 real(kind=real_datatype), intent(inout) :: q(ldq,matrixcols)
20#endif
21 real(kind=real_datatype), intent(in) :: qtrans(2,2)
22#ifdef WITH_MPI
23 integer(kind=MPI_KIND) :: mpierrMPI, my_pcolMPI
24 integer(kind=ik) :: mpierr
25#endif
26 integer(kind=ik) :: my_pcol
27 integer(kind=ik) :: col1, col2
28 real(kind=real_datatype) :: tmp(na)
29 integer(kind=ik) :: pc1, pc2, lc1, lc2
30
31 if (l_rows==0) return ! My processor column has no work to do
32
33#ifdef WITH_MPI
34 call obj%timer%start("mpi_communication")
35 call mpi_comm_rank(int(mpi_comm_cols,kind=mpi_kind) ,my_pcolmpi, mpierr)
36 !call mpi_comm_size(int(mpi_comm_cols,kind=MPI_KIND) ,np_colsMPI, mpierr)
37
38 my_pcol = int(my_pcolmpi,kind=c_int)
39 !np_cols = int(np_colsMPI,kind=c_int)
40
41 call obj%timer%stop("mpi_communication")
42#else
43#endif
44 pc1 = p_col(col1)
45 lc1 = l_col(col1)
46 pc2 = p_col(col2)
47 lc2 = l_col(col2)
48
49 if (pc1==my_pcol) then
50 if (pc2==my_pcol) then
51 ! both columns are local
52 tmp(1:l_rows) = q(l_rqs:l_rqe,lc1)*qtrans(1,1) + q(l_rqs:l_rqe,lc2)*qtrans(2,1)
53 q(l_rqs:l_rqe,lc2) = q(l_rqs:l_rqe,lc1)*qtrans(1,2) + q(l_rqs:l_rqe,lc2)*qtrans(2,2)
54 q(l_rqs:l_rqe,lc1) = tmp(1:l_rows)
55 else
56#ifdef WITH_MPI
57 call obj%timer%start("mpi_communication")
58 call mpi_sendrecv(q(l_rqs,lc1), int(l_rows,kind=mpi_kind), mpi_real_precision, pc2, 1_mpi_kind, &
59 tmp, int(l_rows,kind=mpi_kind), mpi_real_precision, pc2, 1_mpi_kind, &
60 int(mpi_comm_cols,kind=mpi_kind), mpi_status_ignore, mpierr)
61 call obj%timer%stop("mpi_communication")
62#else /* WITH_MPI */
63#endif /* WITH_MPI */
64 q(l_rqs:l_rqe,lc1) = q(l_rqs:l_rqe,lc1)*qtrans(1,1) + tmp(1:l_rows)*qtrans(2,1)
65 endif
66 else if (pc2==my_pcol) then
67#ifdef WITH_MPI
68 call obj%timer%start("mpi_communication")
69 call mpi_sendrecv(q(l_rqs,lc2), int(l_rows,kind=mpi_kind), mpi_real_precision, pc1, 1_mpi_kind, &
70 tmp, int(l_rows,kind=mpi_kind), mpi_real_precision, pc1, 1_mpi_kind, &
71 int(mpi_comm_cols,kind=mpi_kind), mpi_status_ignore, mpierr)
72 call obj%timer%stop("mpi_communication")
73#else /* WITH_MPI */
74 tmp(1:l_rows) = q(l_rqs:l_rqe,lc2)
75#endif /* WITH_MPI */
76
77 q(l_rqs:l_rqe,lc2) = tmp(1:l_rows)*qtrans(1,2) + q(l_rqs:l_rqe,lc2)*qtrans(2,2)
78 endif
79end subroutine transform_columns_&
80 &precision
81
Fortran module to provide an abstract definition of the implementation. Do not use directly....
Definition elpa_abstract_impl.F90:50
Definition elpa_abstract_impl.F90:73