Skip to content

Commit

Permalink
Fixed indentation
Browse files Browse the repository at this point in the history
  • Loading branch information
Linnea Andersson committed Jul 17, 2021
1 parent f5eeefd commit 6304e52
Show file tree
Hide file tree
Showing 5 changed files with 72 additions and 72 deletions.
20 changes: 10 additions & 10 deletions src/prg_ewald_mod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -436,8 +436,8 @@ subroutine Ewald_Real_Space_Test(COULOMBV,I,RX,RY,RZ,LBox, &
Ra(2) = RY(I)
Ra(3) = RZ(I)

! !$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(nnI,J,Rb,Rab,dR,MAGR,MAGR2,TJ,DC,Z,NUMREP_ERFC,CA) &
! !$OMP REDUCTION(+:COULOMBV)
! !$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(nnI,J,Rb,Rab,dR,MAGR,MAGR2,TJ,DC,Z,NUMREP_ERFC,CA) &
! !$OMP REDUCTION(+:COULOMBV)
do nnI = 1,nrnnlist(I)
Rb(1) = nnRx(I,nnI)
Rb(2) = nnRy(I,nnI)
Expand Down Expand Up @@ -488,7 +488,7 @@ subroutine Ewald_Real_Space_Test(COULOMBV,I,RX,RY,RZ,LBox, &

endif
enddo
! !$OMP END PARALLEL DO
! !$OMP END PARALLEL DO
COULOMBV = KECONST*COULOMBV

end subroutine Ewald_Real_Space_Test
Expand Down Expand Up @@ -815,16 +815,16 @@ subroutine Ewald_k_Space_latte_single(COULOMBV,J,RXYZ,Box,DELTAQ,Nr_atoms,COULAC
M22 = M*RECIPVECS(2,2)

do N = NMIN,NMAX
K(1) = L11
K(2) = M22
K(1) = L11
K(2) = M22
K(3) = N*RECIPVECS(3,3)
K2 = K(1)*K(1) + K(2)*K(2) + K(3)*K(3)

PREFACTOR = EIGHT*pi*exp(-K2/(4.D0*CALPHA2))/(COULVOL*K2)
KEPREF = KECONST*PREFACTOR
JDOT = K(1)*RXYZ(1,J) + K(2)*RXYZ(2,J) + K(3)*RXYZ(3,J)
SINJDOT = sin(JDOT)
COSJDOT = cos(JDOT)
COSJDOT = cos(JDOT)
do I = 1,Nr_atoms
IDOT = K(1)*RXYZ(1,I) + K(2)*RXYZ(2,I) + K(3)*RXYZ(3,I)
COULOMBV(I) = COULOMBV(I) + KEPREF*DELTAQ(J)*(COSJDOT*cos(IDOT)+SINJDOT*sin(IDOT))
Expand All @@ -841,9 +841,9 @@ subroutine Ewald_k_Space_latte_single(COULOMBV,J,RXYZ,Box,DELTAQ,Nr_atoms,COULAC
end subroutine Ewald_k_Space_latte_single

subroutine Ewald_k_Space_Test(COULOMBV,RX,RY,RZ,LBox,DELTAQ,Nr_atoms,COULACC,Max_Nr_Neigh)
!
!
implicit none
!
!
integer, parameter :: PREC = 8
real(PREC), parameter :: ONE = 1.D0, TWO = 2.D0, ZERO = 0.D0, SIX = 6.D0, THREE = 3.D0, FOUR = 4.D0
real(PREC), parameter :: FOURTYEIGHT = 48.D0, ELEVEN = 11.D0, SIXTEEN = 16.D0, EIGHT = 8.D0
Expand All @@ -864,7 +864,7 @@ subroutine Ewald_k_Space_Test(COULOMBV,RX,RY,RZ,LBox,DELTAQ,Nr_atoms,COULACC,Max
real(PREC) :: PREVIR, COSSUM,SINSUM,DOT, KEPREF, COSSUM2, SINSUM2

integer :: I,J,L,M,N, ccnt, nnI, LMAX,MMAX,NMAX,NMIN,MMIN
!
!
COULVOL = LBox(1)*LBox(2)*LBox(3)
SQRTX = sqrt(-log(COULACC))

Expand Down Expand Up @@ -946,7 +946,7 @@ subroutine Ewald_k_Space_Test(COULOMBV,RX,RY,RZ,LBox,DELTAQ,Nr_atoms,COULACC,Max
SINSUM2 = SINSUM*SINSUM

! Add up energy and force contributions
!
!
KEPREF = KECONST*PREFACTOR
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(I)
do I = 1,Nr_atoms
Expand Down
92 changes: 46 additions & 46 deletions src/prg_implicit_fermi_mod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ subroutine prg_implicit_fermi_save_inverse(Inv_bml, h_bml, p_bml, nsteps, nocc,
call prg_conjgrad(y_bml, ai_bml, I_bml, aux_bml, d_bml, w_bml, 0.0001_dp, threshold)
do i = 1, nsteps
call bml_copy(I_bml, Inv_bml(i))
enddo
enddo
else
! Otherwise use previous inverse as starting guess
call bml_copy(Inv_bml(1),ai_bml)
Expand All @@ -123,7 +123,7 @@ subroutine prg_implicit_fermi_save_inverse(Inv_bml, h_bml, p_bml, nsteps, nocc,
call bml_scale_add_identity(y_bml, 2.0_dp, 1.0_dp, threshold)
! Find inverse ai = (2*(P2-P)+I)^-1
!call prg_conjgrad(y_bml, Inv_bml(i), I_bml, aux_bml, d_bml, w_bml, tol, threshold)
!call bml_copy(Inv_bml(i),ai_bml)
!call bml_copy(Inv_bml(i),ai_bml)
call prg_newtonschulz(y_bml, Inv_bml(i), d_bml, w_bml, aux_bml, I_bml, tol, threshold)
call bml_multiply(Inv_bml(i), p2_bml, p_bml, 1.0_dp, 0.0_dp, threshold)
!call bml_copy(ai_bml, Inv_bml(i)) ! Save inverses for use in perturbation response calculation
Expand All @@ -133,37 +133,37 @@ subroutine prg_implicit_fermi_save_inverse(Inv_bml, h_bml, p_bml, nsteps, nocc,
trdPdmu = beta*(trP0 - bml_sum_squares(p_bml)) ! sum p(i,j)**2
occErr = abs(trP0 - nocc)
write(*,*) 'occerr =', occErr

! If occupation error is too large, do bisection method
if (occerr > 1.0_dp) then
! if (newerr > occerr) then
! if (newerr > occerr) then
if (nocc-trP0 < 0.0_dp .and. prev .eq. -1) then
prev = -1
else if (nocc-trP0 > 0.0_dp .and. prev .eq. 1) then
else if (nocc-trP0 > 0.0_dp .and. prev .eq. 1) then
prev = 1
else if (nocc-trP0 > 0.0_dp .and. prev .eq. -1) then
prev = 1
alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! do while(newerr > occerr)
! alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! enddo
else
prev = -1
alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! do while(newerr > occerr)
! alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! enddo
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! do while(newerr > occerr)
! alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! enddo
else
prev = -1
alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! do while(newerr > occerr)
! alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
! enddo
endif
!newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
!do while(newerr > occerr .and. abs(occerr-newerr)/occerr<0.1 )
! alpha = alpha/2
! newerr = abs(trP0+prev*alpha*trdPdmu-nocc)
!enddo
! write(*,*) 'newerr =', newerr
! write(*,*) 'newerr =', newerr
mu = mu + prev*alpha
muadj = 1

Expand All @@ -180,9 +180,9 @@ subroutine prg_implicit_fermi_save_inverse(Inv_bml, h_bml, p_bml, nsteps, nocc,
end if
enddo

if (iter .ge. maxiter) then
write(*,*) 'Could not converge chemical potential in prg_impplicit_fermi_save_inverse'
end if
if (iter .ge. maxiter) then
write(*,*) 'Could not converge chemical potential in prg_impplicit_fermi_save_inverse'
end if
! Adjusting the occupation sometimes causes the perturbation calculation to not converge.
! For now we recompute the DM one extra time if mu was adjusted.
!if (muadj .eq. 1) then
Expand Down Expand Up @@ -317,9 +317,9 @@ subroutine prg_implicit_fermi(h_bml, p_bml, nsteps, k, nocc, &
trdPdmu = trdPdmu - bml_sum_squares(p_bml) ! sum p(i,j)**2
trdPdmu = beta * trdPdmu
occErr = abs(trP0 - nocc)
if (occErr .gt. occErrLimit) then
if (occErr .gt. occErrLimit) then
mu = mu + (nocc - trP0)/trdPdmu
end if
end if
write(*,*) "mu =", mu
enddo

Expand Down Expand Up @@ -503,25 +503,25 @@ subroutine prg_implicit_fermi_first_order_response(H0_bml, H1_bml, P0_bml, P1_bm
call bml_multiply(Inv_bml(i), C_bml, P1_bml, 1.0_dp, 0.0_dp, threshold)
enddo

! do i = 1, nsteps-1
! D = A^-1*P0
! call bml_multiply(Inv_bml(i), B0_bml, C0_bml, 1.0_dp, 0.0_dp, threshold)
! call bml_multiply(C0_bml, B0_bml, B_bml, 1.0_dp, 0.0_dp, threshold)
! B0 = A^-1*P0^2
! call bml_copy(B_bml,B0_bml)
! B = I + D -P0*D
! call bml_add(B_bml, C0_bml, -1.0_dp, 1.0_dp, threshold)
! call bml_scale_add_identity(B_bml, 1.0_dp, 1.0_dp, threshold)
! P1 = 2D*P1(I+D-P0*D)
! call bml_multiply(C0_bml, P1_bml, C_bml, 1.0_dp, 0.0_dp, threshold)
! call bml_multiply(C_bml, B_bml, P1_bml, 2.0_dp, 0.0_dp, threshold)
! enddo
! call bml_multiply(B0_bml, P1_bml, C_bml, 2.0_dp, 0.0_dp, threshold)
! call bml_copy(P1_bml, B_bml)
! call bml_add(B_bml, C_bml, 2.0_dp, -2.0_dp, threshold)
! Get next P1
! call bml_multiply(B_bml, P0_bml, C_bml, 1.0_dp, 1.0_dp, threshold)
! call bml_multiply(Inv_bml(i), C_bml, P1_bml, 1.0_dp, 0.0_dp, threshold)
! do i = 1, nsteps-1
! D = A^-1*P0
! call bml_multiply(Inv_bml(i), B0_bml, C0_bml, 1.0_dp, 0.0_dp, threshold)
! call bml_multiply(C0_bml, B0_bml, B_bml, 1.0_dp, 0.0_dp, threshold)
! B0 = A^-1*P0^2
! call bml_copy(B_bml,B0_bml)
! B = I + D -P0*D
! call bml_add(B_bml, C0_bml, -1.0_dp, 1.0_dp, threshold)
! call bml_scale_add_identity(B_bml, 1.0_dp, 1.0_dp, threshold)
! P1 = 2D*P1(I+D-P0*D)
! call bml_multiply(C0_bml, P1_bml, C_bml, 1.0_dp, 0.0_dp, threshold)
! call bml_multiply(C_bml, B_bml, P1_bml, 2.0_dp, 0.0_dp, threshold)
! enddo
! call bml_multiply(B0_bml, P1_bml, C_bml, 2.0_dp, 0.0_dp, threshold)
! call bml_copy(P1_bml, B_bml)
! call bml_add(B_bml, C_bml, 2.0_dp, -2.0_dp, threshold)
! Get next P1
! call bml_multiply(B_bml, P0_bml, C_bml, 1.0_dp, 1.0_dp, threshold)
! call bml_multiply(Inv_bml(i), C_bml, P1_bml, 1.0_dp, 0.0_dp, threshold)


! dPdmu = beta*P0(I-P0)
Expand Down Expand Up @@ -1035,7 +1035,7 @@ subroutine prg_conjgrad(A_bml, p_bml, p2_bml, tmp_bml, d_bml, w_bml, cg_tol, thr

do while (r_norm_new .gt. cg_tol)

! write(*,*) r_norm_new
! write(*,*) r_norm_new
k = k + 1
if (k .eq. 1) then
call bml_copy(tmp_bml, d_bml)
Expand Down Expand Up @@ -1173,8 +1173,8 @@ subroutine prg_test_density_matrix(ham_bml, p_bml, beta, mu, nocc, osteps, occEr
trdPdmu = beta * trdPdmu
occErr = abs(trP0 - nocc)
if (occErr .gt. occErrLimit) then
mu = mu + (nocc - trP0)/trdPdmu
end if
mu = mu + (nocc - trP0)/trdPdmu
end if
!write(*,*) "mu = ", mu
enddo

Expand Down
4 changes: 2 additions & 2 deletions src/prg_syrotation_mod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -190,7 +190,7 @@ subroutine prg_rotate(rot,r,verbose)
v2(2)=rot%v2(2)
v2(3)=rot%v2(3)

vQ(1)=rot%vQ(1) !Rotation center
vQ(1)=rot%vQ(1) !Rotation center
vQ(2)=rot%vQ(2)
vQ(3)=rot%vQ(3)

Expand Down Expand Up @@ -232,7 +232,7 @@ subroutine prg_rotate(rot,r,verbose)
v2=pq2-vQ
endif

vtr(1)=0.0_dp !Translation
vtr(1)=0.0_dp !Translation
vtr(2)=0.0_dp
vtr(3)=0.0_dp

Expand Down
14 changes: 7 additions & 7 deletions src/prg_xlbokernel_mod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -332,8 +332,8 @@ subroutine prg_kernel_multirank(KRes,KK0_bml,Res,FelTol,L,LMAX,HO_bml,mu,beta,RX
enddo
!$OMP END PARALLEL DO

! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc,TIMERATIO,Max_Nr_Neigh)
! Coulomb_Pot_dq_v = Coulomb_Pot_Real+Coulomb_Pot_k
! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc,TIMERATIO,Max_Nr_Neigh)
! Coulomb_Pot_dq_v = Coulomb_Pot_Real+Coulomb_Pot_k

row1 = 0.0_dp
do J = 1,HDIM
Expand Down Expand Up @@ -557,7 +557,7 @@ subroutine prg_full_kernel_latte(KK,DO_bml,mu0,RXYZ,Box,Hubbard_U, &
enddo

write(*,*) 'ewaldracc =', ewaldracc
write(*,*) 'ewaldkacc =', ewaldkacc
write(*,*) 'ewaldkacc =', ewaldkacc
write(*,*) 'implcit response time =', respacc
call Invert(JJ,KK,Nr_atoms)

Expand Down Expand Up @@ -648,8 +648,8 @@ subroutine prg_full_kernel(KK,DO_bml,mu0,RX,RY,RZ,LBox,Hubbard_U,Element_Type, &
Coulomb_Pot_Real(I) = Coulomb_Pot_Real_I
enddo
!$OMP END PARALLEL DO
! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc, &
! TIMERATIO,Max_Nr_Neigh)
! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc, &
! TIMERATIO,Max_Nr_Neigh)
Coulomb_Pot_dq_v = Coulomb_Pot_Real+Coulomb_Pot_k

diagonal = 0.0_dp
Expand Down Expand Up @@ -766,8 +766,8 @@ subroutine prg_kernel_matrix_multirank(KRes,KK0_bml,Res,FelTol,L,LMAX,HO_bml,mu,
enddo
!$OMP END PARALLEL DO

! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc,TIMERATIO,Max_Nr_Neigh)
! Coulomb_Pot_dq_v = Coulomb_Pot_Real+Coulomb_Pot_k
! call Ewald_k_Space(Coulomb_Pot_k,Coulomb_Force_k,RX,RY,RZ,LBox,dq_v,Nr_atoms,Coulomb_acc,TIMERATIO,Max_Nr_Neigh)
! Coulomb_Pot_dq_v = Coulomb_Pot_Real+Coulomb_Pot_k

row1 = 0.0_dp
do J = 1,HDIM
Expand Down
14 changes: 7 additions & 7 deletions tests/src/main.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ program main
implicit none

integer :: norb, mdim, verbose
type(bml_matrix_t) :: inv_bml(10)
type(bml_matrix_t) :: inv_bml(10)
type(bml_matrix_t) :: ham_bml
type(bml_matrix_t) :: rho_bml, rho1_bml
type(bml_matrix_t) :: rho_ortho_bml
Expand Down Expand Up @@ -215,7 +215,7 @@ program main
error stop
endif

case("prg_implicit_fermi_save_inverse")
case("prg_implicit_fermi_save_inverse")

mu = 0.2_dp
beta = 4.0_dp !nocc,osteps,occerrlimit
Expand All @@ -224,18 +224,18 @@ program main

do i = 1,norecs
call bml_identity_matrix(bml_type,bml_element_real,dp,norb,norb,inv_bml(i))
enddo
enddo

call bml_zero_matrix(bml_type,bml_element_real,dp,norb,norb,rho_bml)
call bml_zero_matrix(bml_type,bml_element_real,dp,norb,norb,ham_bml)

call bml_read_matrix(ham_bml,'hamiltonian_ortho.mtx')

call bml_zero_matrix(bml_type,bml_element_real,dp,norb,norb,rho1_bml)
mu = 0.2_dp

mu = 0.2_dp
call prg_implicit_fermi_save_inverse(inv_bml,ham_bml,rho_bml,norecs,nocc,mu,beta,1e-4_dp, threshold, 1e-5_dp, 1,occiter)
call prg_test_density_matrix(ham_bml,rho1_bml,beta,mu,nocc,1,1e-4_dp,threshold)
call prg_test_density_matrix(ham_bml,rho1_bml,beta,mu,nocc,1,1e-4_dp,threshold)
write(*,*) mu
call bml_scale(0.5_dp,rho_bml)
call bml_add(rho1_bml,rho_bml,1.0_dp,-1.0_dp,threshold)
Expand Down Expand Up @@ -1132,7 +1132,7 @@ program main
call bml_add_deprecated(-1.0_dp,rho_bml,1.0_dp,rho1_bml,0.0_dp)
error_calc = bml_fnorm(rho_bml)
write(*,*)error_calc


if(error_calc.gt.error_tol)then
write(*,*) "Error is too high", error_calc
Expand Down

0 comments on commit 6304e52

Please sign in to comment.