Skip to content

Commit

Permalink
Fixes for the Intel compiler related to #1664
Browse files Browse the repository at this point in the history
  • Loading branch information
@mkrack
mkrack committed Oct 5, 2021
1 parent 62b873b commit 2d2fde0
Show file tree
Hide file tree
Showing 3 changed files with 106 additions and 63 deletions.
96 changes: 64 additions & 32 deletions src/ec_orth_solver.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -75,11 +75,11 @@ MODULE ec_orth_solver

PRIVATE

! *** Global parameters ***
! Global parameters

CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'ec_orth_solver'

! *** Public subroutines ***
! Public subroutines

PUBLIC :: ec_response_ao

Expand Down Expand Up @@ -111,7 +111,7 @@ SUBROUTINE preconditioner(qs_env, matrix_ks, matrix_p, matrix_rhs, &
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: matrix_ks, matrix_p, matrix_rhs
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(OUT), &
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: matrix_cg_z
REAL(KIND=dp), INTENT(IN) :: eps_filter
INTEGER, INTENT(IN) :: iounit
Expand All @@ -129,6 +129,12 @@ SUBROUTINE preconditioner(qs_env, matrix_ks, matrix_p, matrix_rhs, &

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(qs_env))
CPASSERT(ASSOCIATED(matrix_ks))
CPASSERT(ASSOCIATED(matrix_p))
CPASSERT(ASSOCIATED(matrix_rhs))
CPASSERT(ASSOCIATED(matrix_cg_z))

NULLIFY (dft_control, linres_control)

t1 = m_walltime()
Expand Down Expand Up @@ -209,7 +215,7 @@ SUBROUTINE preconditioner(qs_env, matrix_ks, matrix_p, matrix_rhs, &
WRITE (iounit, "(/,T10,A,T25,A,T42,A,T62,A,/,T10,A)") &
"Iteration", "Stepsize", "Convergence", "Time", &
REPEAT("-", 58)
ENDIF
END IF

alpha(:) = 0.0_dp
max_iter = 200
Expand Down Expand Up @@ -283,7 +289,7 @@ SUBROUTINE preconditioner(qs_env, matrix_ks, matrix_p, matrix_rhs, &
! Convergence criteria
IF (norm_res .LT. linres_control%eps) THEN
converged = .TRUE.
ENDIF
END IF

t2 = m_walltime()
IF (i .EQ. 1 .OR. MOD(i, 1) .EQ. 0 .OR. converged) THEN
Expand All @@ -294,25 +300,25 @@ SUBROUTINE preconditioner(qs_env, matrix_ks, matrix_p, matrix_rhs, &
!WRITE (iounit, "(T10,I5,T25,1E8.2,T42,1E14.8,T58,F8.2)") &
! i, MAXVAL(alpha), norm_res, t2 - t1
CALL m_flush(iounit)
ENDIF
ENDIF
END IF
END IF
IF (converged) THEN
IF (iounit > 0) THEN
WRITE (iounit, "(/,T10,A,I4,A,/)") "The precon solver converged in ", i, " iterations."
CALL m_flush(iounit)
ENDIF
END IF
EXIT iteration
ENDIF
END IF

! Max number of iteration reached
IF (i == max_iter) THEN
IF (iounit > 0) THEN
WRITE (iounit, "(/,T10,A/)") &
"The precon solver didnt converge! Maximum number of iterations reached."
CALL m_flush(iounit)
ENDIF
END IF
converged = .FALSE.
ENDIF
END IF

END DO iteration

Expand Down Expand Up @@ -352,7 +358,7 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: matrix_hz
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(OUT), &
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: matrix_pz, matrix_wz
INTEGER, INTENT(IN) :: iounit
LOGICAL, INTENT(OUT) :: should_stop
Expand All @@ -378,6 +384,11 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(qs_env))
CPASSERT(ASSOCIATED(matrix_hz))
CPASSERT(ASSOCIATED(matrix_pz))
CPASSERT(ASSOCIATED(matrix_wz))

NULLIFY (dft_control, ksmat, matrix_s, linres_control, rho)

t1 = m_walltime()
Expand Down Expand Up @@ -583,7 +594,7 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul
WRITE (iounit, "(/,T3,A,T16,A,T25,A,T38,A,T52,A,/,T3,A)") &
"Iteration", "Method", "Stepsize", "Convergence", "Time", &
REPEAT("-", 80)
ENDIF
END IF

alpha(:) = 0.0_dp
restart = .FALSE.
Expand All @@ -597,7 +608,7 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul
! default for eps 10E-6 in MO_linres
IF (norm_res .LT. linres_control%eps) THEN
linres_control%converged = .TRUE.
ENDIF
END IF

t2 = m_walltime()
IF (i .EQ. 1 .OR. MOD(i, 1) .EQ. 0 .OR. linres_control%converged &
Expand All @@ -606,31 +617,31 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul
WRITE (iounit, "(T5,I5,T18,A3,T28,L1,T38,1E8.2,T48,F16.10,T68,F8.2)") &
i, linres_control%flag, restart, MAXVAL(alpha), norm_res, t2 - t1
CALL m_flush(iounit)
ENDIF
ENDIF
END IF
END IF
IF (linres_control%converged) THEN
IF (iounit > 0) THEN
WRITE (iounit, "(/,T2,A,I4,A,/)") "The linear solver converged in ", i, " iterations."
CALL m_flush(iounit)
ENDIF
END IF
EXIT iteration
ELSE IF (should_stop) THEN
IF (iounit > 0) THEN
WRITE (iounit, "(/,T2,A,I4,A,/)") "The linear solver did NOT converge! External stop"
CALL m_flush(iounit)
END IF
EXIT iteration
ENDIF
END IF

! Max number of iteration reached
IF (i == linres_control%max_iter) THEN
IF (iounit > 0) THEN
WRITE (iounit, "(/,T2,A/)") &
"The linear solver didnt converge! Maximum number of iterations reached."
CALL m_flush(iounit)
ENDIF
END IF
linres_control%converged = .FALSE.
ENDIF
END IF

! Hessian Ax = [F,B] + [G(B),P]
CALL build_hessian_op(qs_env=qs_env, &
Expand Down Expand Up @@ -730,9 +741,9 @@ SUBROUTINE ec_response_ao(qs_env, matrix_hz, matrix_pz, matrix_wz, iounit, shoul
! r_j+1 = r_j - alpha * Ap_j
DO ispin = 1, nspins
CALL dbcsr_add(matrix_res(ispin)%matrix, matrix_Ax(ispin)%matrix, 1.0_dp, -alpha(ispin))
ENDDO
END DO
restart = .FALSE.
ENDIF
END IF

! Preconditioner
linres_control%flag = ""
Expand Down Expand Up @@ -843,7 +854,7 @@ SUBROUTINE ec_wz_matrix(qs_env, matrix_z, matrix_wz, eps_filter)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: matrix_z
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(OUT), &
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: matrix_wz
REAL(KIND=dp), INTENT(IN) :: eps_filter

Expand All @@ -859,6 +870,10 @@ SUBROUTINE ec_wz_matrix(qs_env, matrix_z, matrix_wz, eps_filter)

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(qs_env))
CPASSERT(ASSOCIATED(matrix_z))
CPASSERT(ASSOCIATED(matrix_wz))

CALL get_qs_env(qs_env=qs_env, &
dft_control=dft_control, &
matrix_ks=matrix_ks, &
Expand Down Expand Up @@ -902,7 +917,7 @@ SUBROUTINE ec_wz_matrix(qs_env, matrix_z, matrix_wz, eps_filter)
! Whz = ZFP + PFZ
CALL dbcsr_copy(matrix_wz(ispin)%matrix, matrix_tmp, keep_sparsity=.TRUE.)

ENDDO
END DO

! Release matrices
CALL dbcsr_release(matrix_tmp)
Expand Down Expand Up @@ -933,7 +948,7 @@ SUBROUTINE hessian_op1(matrix_ks, matrix_p, matrix_cg, matrix_b, matrix_Ax, eps_

TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: matrix_ks, matrix_p, matrix_cg
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(OUT), &
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: matrix_b, matrix_Ax
REAL(KIND=dp), INTENT(IN) :: eps_filter

Expand All @@ -943,6 +958,12 @@ SUBROUTINE hessian_op1(matrix_ks, matrix_p, matrix_cg, matrix_b, matrix_Ax, eps_

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(matrix_ks))
CPASSERT(ASSOCIATED(matrix_p))
CPASSERT(ASSOCIATED(matrix_cg))
CPASSERT(ASSOCIATED(matrix_b))
CPASSERT(ASSOCIATED(matrix_Ax))

! Build intermediate density matrix
! B = [cg, P] = cg*P - P*cg = cg*P + (cg*P)^T
CALL commutator(matrix_cg, matrix_p, matrix_b, eps_filter, .TRUE.)
Expand All @@ -953,7 +974,7 @@ SUBROUTINE hessian_op1(matrix_ks, matrix_p, matrix_cg, matrix_b, matrix_Ax, eps_

CALL timestop(handle)

END SUBROUTINE
END SUBROUTINE hessian_op1

! **************************************************************************************************
!> \brief calculate lin transformation of Hessian matrix on a trial vector (matrix) matrix_cg
Expand All @@ -980,7 +1001,7 @@ SUBROUTINE build_hessian_op(qs_env, matrix_ks, matrix_p, matrix_s_sqrt_inv, &
TYPE(dbcsr_type), INTENT(IN) :: matrix_s_sqrt_inv
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: matrix_cg
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(OUT), &
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: matrix_Ax
REAL(KIND=dp), INTENT(IN) :: eps_filter

Expand All @@ -996,6 +1017,12 @@ SUBROUTINE build_hessian_op(qs_env, matrix_ks, matrix_p, matrix_s_sqrt_inv, &

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(qs_env))
CPASSERT(ASSOCIATED(matrix_ks))
CPASSERT(ASSOCIATED(matrix_p))
CPASSERT(ASSOCIATED(matrix_cg))
CPASSERT(ASSOCIATED(matrix_Ax))

CALL get_qs_env(qs_env=qs_env, &
dft_control=dft_control, &
para_env=para_env, &
Expand All @@ -1022,7 +1049,7 @@ SUBROUTINE build_hessian_op(qs_env, matrix_ks, matrix_p, matrix_s_sqrt_inv, &
chksum = 0.0_dp
DO ispin = 1, nspins
chksum = chksum + dbcsr_checksum(matrix_b(ispin)%matrix)
ENDDO
END DO

! skip the kernel if the DM is very small
IF (chksum .GT. 1.0E-14_dp) THEN
Expand Down Expand Up @@ -1065,7 +1092,7 @@ SUBROUTINE build_hessian_op(qs_env, matrix_ks, matrix_p, matrix_s_sqrt_inv, &

CALL timestop(handle)

END SUBROUTINE
END SUBROUTINE build_hessian_op

! **************************************************************************************************
!> \brief Calculate lin transformation of Hessian matrix on a trial vector (matrix) matrix_cg
Expand Down Expand Up @@ -1265,7 +1292,8 @@ END SUBROUTINE hessian_op2
! **************************************************************************************************
SUBROUTINE commutator(a, b, res, eps_filter, anticomm, alpha, beta)

TYPE(dbcsr_p_type), DIMENSION(:) :: a, b, res
TYPE(dbcsr_p_type), DIMENSION(:), INTENT(IN), &
POINTER :: a, b, res
REAL(KIND=dp) :: eps_filter
LOGICAL :: anticomm
REAL(KIND=dp), OPTIONAL :: alpha, beta
Expand All @@ -1278,6 +1306,10 @@ SUBROUTINE commutator(a, b, res, eps_filter, anticomm, alpha, beta)

CALL timeset(routineN, handle)

CPASSERT(ASSOCIATED(a))
CPASSERT(ASSOCIATED(b))
CPASSERT(ASSOCIATED(res))

CALL dbcsr_create(work, template=a(1)%matrix, matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_create(work2, template=a(1)%matrix, matrix_type=dbcsr_type_no_symmetry)

Expand Down Expand Up @@ -1382,7 +1414,7 @@ SUBROUTINE projector(qs_env, matrix_p, matrix_io, eps_filter)

CALL timestop(handle)

END SUBROUTINE
END SUBROUTINE projector

! **************************************************************************************************
!> \brief performs a tranformation of a matrix back to/into orthonormal basis
Expand Down Expand Up @@ -1427,6 +1459,6 @@ SUBROUTINE transform_m_orth(matrix, matrix_trafo, eps_filter)
CALL dbcsr_release(matrix_work)
CALL timestop(handle)

END SUBROUTINE
END SUBROUTINE transform_m_orth

END MODULE ec_orth_solver

0 comments on commit 2d2fde0

Please sign in to comment.