From 2e8ced0eef799853a03b0c08062b14089498f587 Mon Sep 17 00:00:00 2001 From: Abdallah Ammar Date: Mon, 10 Apr 2023 16:11:52 +0200 Subject: [PATCH] overlap with dgemm --- src/bi_ortho_mos/overlap.irp.f | 92 ++++++++++++++++++++-------------- 1 file changed, 54 insertions(+), 38 deletions(-) diff --git a/src/bi_ortho_mos/overlap.irp.f b/src/bi_ortho_mos/overlap.irp.f index d7f45c94..ff5d5c84 100644 --- a/src/bi_ortho_mos/overlap.irp.f +++ b/src/bi_ortho_mos/overlap.irp.f @@ -12,32 +12,27 @@ double precision :: accu_d, accu_nd double precision, allocatable :: tmp(:,:) - ! TODO : re do the DEGEMM +! overlap_bi_ortho = 0.d0 +! do i = 1, mo_num +! do k = 1, mo_num +! do m = 1, ao_num +! do n = 1, ao_num +! overlap_bi_ortho(k,i) += ao_overlap(n,m) * mo_l_coef(n,k) * mo_r_coef(m,i) +! enddo +! enddo +! enddo +! enddo - overlap_bi_ortho = 0.d0 - do i = 1, mo_num - do k = 1, mo_num - do m = 1, ao_num - do n = 1, ao_num - overlap_bi_ortho(k,i) += ao_overlap(n,m) * mo_l_coef(n,k) * mo_r_coef(m,i) - enddo - enddo - enddo - enddo - -! allocate( tmp(mo_num,ao_num) ) -! -! ! tmp <-- L.T x S_ao -! call dgemm( "T", "N", mo_num, ao_num, ao_num, 1.d0 & -! , mo_l_coef, size(mo_l_coef, 1), ao_overlap, size(ao_overlap, 1) & -! , 0.d0, tmp, size(tmp, 1) ) -! -! ! S <-- tmp x R -! call dgemm( "N", "N", mo_num, mo_num, ao_num, 1.d0 & -! , tmp, size(tmp, 1), mo_r_coef, size(mo_r_coef, 1) & -! , 0.d0, overlap_bi_ortho, size(overlap_bi_ortho, 1) ) -! -! deallocate( tmp ) + allocate( tmp(mo_num,ao_num) ) + ! tmp <-- L.T x S_ao + call dgemm( "T", "N", mo_num, ao_num, ao_num, 1.d0 & + , mo_l_coef(1,1), size(mo_l_coef, 1), ao_overlap(1,1), size(ao_overlap, 1) & + , 0.d0, tmp(1,1), size(tmp, 1) ) + ! S <-- tmp x R + call dgemm( "N", "N", mo_num, mo_num, ao_num, 1.d0 & + , tmp(1,1), size(tmp, 1), mo_r_coef(1,1), size(mo_r_coef, 1) & + , 0.d0, overlap_bi_ortho(1,1), size(overlap_bi_ortho, 1) ) + deallocate(tmp) do i = 1, mo_num overlap_diag_bi_ortho(i) = overlap_bi_ortho(i,i) @@ -84,20 +79,41 @@ END_PROVIDER END_DOC implicit none - integer :: i, j, p, q + integer :: i, j, p, q + double precision, allocatable :: tmp(:,:) - overlap_mo_r = 0.d0 - overlap_mo_l = 0.d0 - do i = 1, mo_num - do j = 1, mo_num - do p = 1, ao_num - do q = 1, ao_num - overlap_mo_r(j,i) += mo_r_coef(q,i) * mo_r_coef(p,j) * ao_overlap(q,p) - overlap_mo_l(j,i) += mo_l_coef(q,i) * mo_l_coef(p,j) * ao_overlap(q,p) - enddo - enddo - enddo - enddo + !overlap_mo_r = 0.d0 + !overlap_mo_l = 0.d0 + !do i = 1, mo_num + ! do j = 1, mo_num + ! do p = 1, ao_num + ! do q = 1, ao_num + ! overlap_mo_r(j,i) += mo_r_coef(q,i) * mo_r_coef(p,j) * ao_overlap(q,p) + ! overlap_mo_l(j,i) += mo_l_coef(q,i) * mo_l_coef(p,j) * ao_overlap(q,p) + ! enddo + ! enddo + ! enddo + !enddo + + allocate( tmp(mo_num,ao_num) ) + + tmp = 0.d0 + call dgemm( "T", "N", mo_num, ao_num, ao_num, 1.d0 & + , mo_r_coef(1,1), size(mo_r_coef, 1), ao_overlap(1,1), size(ao_overlap, 1) & + , 0.d0, tmp(1,1), size(tmp, 1) ) + call dgemm( "N", "N", mo_num, mo_num, ao_num, 1.d0 & + , tmp(1,1), size(tmp, 1), mo_r_coef(1,1), size(mo_r_coef, 1) & + , 0.d0, overlap_mo_r(1,1), size(overlap_mo_r, 1) ) + + tmp = 0.d0 + call dgemm( "T", "N", mo_num, ao_num, ao_num, 1.d0 & + , mo_l_coef(1,1), size(mo_l_coef, 1), ao_overlap(1,1), size(ao_overlap, 1) & + , 0.d0, tmp(1,1), size(tmp, 1) ) + call dgemm( "N", "N", mo_num, mo_num, ao_num, 1.d0 & + , tmp(1,1), size(tmp, 1), mo_l_coef(1,1), size(mo_l_coef, 1) & + , 0.d0, overlap_mo_l(1,1), size(overlap_mo_l, 1) ) + + deallocate(tmp) END_PROVIDER