elpa_impl_generalized_transform_template.F90 5.83 KB
 Pavel Kus committed Aug 01, 2018 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 ``````! using elpa internal Hermitian multiply is faster then scalapack multiply, but we need an extra ! temporary matrix. ! using cannon algorithm should be the fastest. After this is verified, the other options should be removed ! however, we need the extra temporary matrix as well. #undef FORWARD_ELPA_CANNON #undef FORWARD_SCALAPACK #undef FORWARD_ELPA_HERMITIAN #if defined(REALCASE) && defined(DOUBLE_PRECISION) #define FORWARD_ELPA_CANNON !#define FORWARD_ELPA_HERMITIAN #else !TODO first just for real double... #define FORWARD_ELPA_HERMITIAN #endif #define BACKWARD_ELPA_CANNON #undef BACKWARD_SCALAPACK `````` Andreas Marek committed Mar 05, 2018 20 `````` subroutine elpa_transform_generalized_& `````` Pavel Kus committed Oct 30, 2017 21 `````` &ELPA_IMPL_SUFFIX& `````` Pavel Kus committed Feb 05, 2018 22 `````` &(self, a, b, is_already_decomposed, error) `````` Pavel Kus committed Oct 30, 2017 23 24 25 26 27 28 29 30 31 `````` implicit none #include "general/precision_kinds.F90" class(elpa_impl_t) :: self #ifdef USE_ASSUMED_SIZE MATH_DATATYPE(kind=rck) :: a(self%local_nrows, *), b(self%local_nrows, *) #else MATH_DATATYPE(kind=rck) :: a(self%local_nrows, self%local_ncols), b(self%local_nrows, self%local_ncols) #endif integer :: error `````` Pavel Kus committed Feb 05, 2018 32 `````` logical :: is_already_decomposed `````` Pavel Kus committed Feb 05, 2018 33 `````` integer :: sc_desc(SC_DESC_LEN) `````` Pavel Kus committed Aug 01, 2018 34 35 `````` integer(kind=ik) :: my_p, my_prow, my_pcol, np_rows, np_cols, mpierr, mpi_comm_rows, mpi_comm_cols, mpi_comm_all integer(kind=ik) :: BuffLevelInt `````` Pavel Kus committed Oct 30, 2017 36 `````` `````` Pavel Kus committed Aug 01, 2018 37 ``````#if defined(FORWARD_ELPA_HERMITIAN) || defined(FORWARD_ELPA_CANNON) `````` Pavel Kus committed Oct 30, 2017 38 `````` MATH_DATATYPE(kind=rck) :: tmp(self%local_nrows, self%local_ncols) `````` Pavel Kus committed Feb 05, 2018 39 ``````#endif `````` Pavel Kus committed Oct 30, 2017 40 `````` `````` Pavel Kus committed Feb 05, 2018 41 `````` call self%timer_start("transform_generalized()") `````` Pavel Kus committed Feb 05, 2018 42 `````` `````` Pavel Kus committed Feb 05, 2018 43 44 45 `````` error = self%construct_scalapack_descriptor(sc_desc) if(error .NE. ELPA_OK) return `````` Pavel Kus committed Feb 05, 2018 46 47 48 49 50 51 52 53 54 55 56 57 `````` if (.not. is_already_decomposed) then ! B = U^T*U, B<-U call self%elpa_cholesky_& &ELPA_IMPL_SUFFIX& &(b, error) if(error .NE. ELPA_OK) return ! B <- inv(U) call self%elpa_invert_trm_& &ELPA_IMPL_SUFFIX& &(b, error) if(error .NE. ELPA_OK) return end if `````` Pavel Kus committed Feb 05, 2018 58 `````` `````` Pavel Kus committed Aug 01, 2018 59 ``````#ifdef FORWARD_ELPA_HERMITIAN `````` Pavel Kus committed Feb 05, 2018 60 61 62 63 64 65 `````` ! tmp <- inv(U^T) * A (we have to use temporary variable) call self%elpa_hermitian_multiply_& &ELPA_IMPL_SUFFIX& &('U','F', self%na, b, a, self%local_nrows, self%local_ncols, tmp, & self%local_nrows, self%local_ncols, error) if(error .NE. ELPA_OK) return `````` Pavel Kus committed Feb 05, 2018 66 `````` `````` Pavel Kus committed Feb 05, 2018 67 68 `````` ! A <- inv(U)^T * A a(1:self%local_nrows, 1:self%local_ncols) = tmp(1:self%local_nrows, 1:self%local_ncols) `````` Pavel Kus committed Aug 01, 2018 69 70 ``````#endif #ifdef FORWARD_SCALAPACK `````` Pavel Kus committed Feb 05, 2018 71 72 `````` ! A <- inv(U)^T * A (using scalapack, we can directly update A) call self%timer_start("scalapack multiply inv(U)^T * A") `````` Pavel Kus committed Feb 05, 2018 73 ``````#ifdef WITH_MPI `````` Pavel Kus committed Feb 05, 2018 74 75 76 77 `````` call p& &BLAS_CHAR& &trmm("L", "U", BLAS_TRANS_OR_CONJ, "N", self%na, self%na, & ONE, b, 1, 1, sc_desc, a, 1, 1, sc_desc) `````` Pavel Kus committed Feb 05, 2018 78 ``````#else `````` Pavel Kus committed Feb 05, 2018 79 80 81 `````` call BLAS_CHAR& &trmm("L", "U", BLAS_TRANS_OR_CONJ, "N", self%na, self%na, & ONE, b, self%na, a, self%na) `````` Pavel Kus committed Feb 05, 2018 82 ``````#endif `````` Andreas Marek committed Mar 05, 2018 83 `````` `````` Pavel Kus committed Feb 05, 2018 84 `````` call self%timer_stop("scalapack multiply inv(U)^T * A") `````` Pavel Kus committed Aug 01, 2018 85 ``````#endif /* FORWARD_SCALAPACK */ `````` Pavel Kus committed Feb 05, 2018 86 `````` `````` Pavel Kus committed Aug 01, 2018 87 ``````#if defined(FORWARD_ELPA_HERMITIAN) || defined(FORWARD_SCALAPACK) `````` Pavel Kus committed Feb 05, 2018 88 89 90 `````` ! A <- inv(U)^T * A * inv(U) ! For this multiplication we do not have internal function in ELPA, ! so we have to call scalapack anyway `````` Pavel Kus committed Feb 05, 2018 91 92 `````` call self%timer_start("scalapack multiply A * inv(U)") #ifdef WITH_MPI `````` Pavel Kus committed Oct 30, 2017 93 94 95 96 97 98 99 100 101 `````` call p& &BLAS_CHAR& &trmm("R", "U", "N", "N", self%na, self%na, & ONE, b, 1, 1, sc_desc, a, 1, 1, sc_desc) #else call BLAS_CHAR& &trmm("R", "U", "N", "N", self%na, self%na, & ONE, b, self%na, a, self%na) #endif `````` Pavel Kus committed Feb 05, 2018 102 `````` call self%timer_stop("scalapack multiply A * inv(U)") `````` Pavel Kus committed Aug 01, 2018 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 ``````#endif /*(FORWARD_ELPA_HERMITIAN) || defined(FORWARD_SCALAPACK)*/ #ifdef FORWARD_ELPA_CANNON !TODO set the value properly !TODO tunable parameter? BuffLevelInt = 1 call self%get("mpi_comm_rows",mpi_comm_rows,error) call self%get("mpi_comm_cols",mpi_comm_cols,error) call self%get("mpi_comm_parent", mpi_comm_all,error) call mpi_comm_rank(mpi_comm_rows,my_prow,mpierr) call mpi_comm_size(mpi_comm_rows,np_rows,mpierr) call mpi_comm_rank(mpi_comm_cols,my_pcol,mpierr) call mpi_comm_size(mpi_comm_cols,np_cols,mpierr) call mpi_comm_rank(mpi_comm_all,my_p,mpierr) call cannons_reduction(a, b, self%local_nrows, self%local_ncols, np_rows, np_cols, my_prow, my_pcol, & sc_desc, tmp, BuffLevelInt, mpi_comm_rows, mpi_comm_cols) a(1:self%local_nrows, 1:self%local_ncols) = tmp(1:self%local_nrows, 1:self%local_ncols) #endif /*FORWARD_ELPA_CANNON*/ write(*, *) my_prow, my_pcol, "A(2,3)", a(2,3) `````` Pavel Kus committed Oct 30, 2017 126 `````` `````` Pavel Kus committed Feb 05, 2018 127 `````` call self%timer_stop("transform_generalized()") `````` Pavel Kus committed Oct 30, 2017 128 129 `````` end subroutine `````` Pavel Kus committed Nov 07, 2017 130 131 132 `````` subroutine elpa_transform_back_generalized_& &ELPA_IMPL_SUFFIX& `````` Pavel Kus committed Feb 05, 2018 133 `````` &(self, b, q, error) `````` Pavel Kus committed Nov 07, 2017 134 135 136 137 138 139 140 141 142 `````` implicit none #include "general/precision_kinds.F90" class(elpa_impl_t) :: self #ifdef USE_ASSUMED_SIZE MATH_DATATYPE(kind=rck) :: b(self%local_nrows, *), q(self%local_nrows, *) #else MATH_DATATYPE(kind=rck) :: b(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols) #endif integer :: error `````` Pavel Kus committed Feb 05, 2018 143 `````` integer :: sc_desc(SC_DESC_LEN) `````` Pavel Kus committed Nov 07, 2017 144 `````` `````` Pavel Kus committed Feb 05, 2018 145 146 `````` call self%timer_start("transform_back_generalized()") `````` Pavel Kus committed Feb 05, 2018 147 148 149 `````` error = self%construct_scalapack_descriptor(sc_desc) if(error .NE. ELPA_OK) return `````` Pavel Kus committed Feb 05, 2018 150 `````` call self%timer_start("scalapack multiply inv(U) * Q") `````` Pavel Kus committed Nov 07, 2017 151 ``````#ifdef WITH_MPI `````` Pavel Kus committed Feb 05, 2018 152 `````` ! Q <- inv(U) * Q `````` Pavel Kus committed Nov 07, 2017 153 154 `````` call p& &BLAS_CHAR& `````` Pavel Kus committed Mar 09, 2018 155 `````` &trmm("L", "U", "N", "N", self%na, self%nev, & `````` Pavel Kus committed Nov 07, 2017 156 157 158 `````` ONE, b, 1, 1, sc_desc, q, 1, 1, sc_desc) #else call BLAS_CHAR& `````` Pavel Kus committed Mar 09, 2018 159 `````` &trmm("L", "U", "N", "N", self%na, self%nev, & `````` Pavel Kus committed Nov 07, 2017 160 161 `````` ONE, b, self%na, q, self%na) #endif `````` Pavel Kus committed Feb 05, 2018 162 163 164 `````` call self%timer_stop("scalapack multiply inv(U) * Q") call self%timer_stop("transform_back_generalized()") `````` Pavel Kus committed Nov 07, 2017 165 166 167 `````` end subroutine ``````