diff --git a/src/GW/RGW_phBSE.f90 b/src/GW/RGW_phBSE.f90
index 32ddaae..4e1a926 100644
--- a/src/GW/RGW_phBSE.f90
+++ b/src/GW/RGW_phBSE.f90
@@ -89,8 +89,8 @@ subroutine RGW_phBSE(dophBSE2,exchange_kernel,TDA_W,TDA,dBSE,dTDA,singlet,triple
   call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
   call RGW_excitation_density(nOrb,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
 
-  call RGW_phBSE_static_kernel_A(eta,nOrb,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta)
-  call RGW_phBSE_static_kernel_B(eta,nOrb,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta)
+               call RGW_phBSE_static_kernel_A(eta,nOrb,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta)
+  if(.not.TDA) call RGW_phBSE_static_kernel_B(eta,nOrb,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta)
 
 !-------------------
 ! Singlet manifold
@@ -116,8 +116,8 @@ subroutine RGW_phBSE(dophBSE2,exchange_kernel,TDA_W,TDA,dBSE,dTDA,singlet,triple
       write(*,*) 
 
       call RGW_phBSE_static_kernel(eta,nOrb,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
-      call RGW_phBSE2_static_kernel_A(eta,nOrb,nC,nO,nV,nR,nS,1d0,eW,W,KA_sta)
 
+                   call RGW_phBSE2_static_kernel_A(eta,nOrb,nC,nO,nV,nR,nS,1d0,eW,W,KA_sta)
       if(.not.TDA) call RGW_phBSE2_static_kernel_B(eta,nOrb,nC,nO,nV,nR,nS,1d0,eW,W,KB_sta)
 
       deallocate(W)
diff --git a/src/GW/UGW_phBSE_static_kernel_A.f90 b/src/GW/UGW_phBSE_static_kernel_A.f90
index 3939043..6b2f0b7 100644
--- a/src/GW/UGW_phBSE_static_kernel_A.f90
+++ b/src/GW/UGW_phBSE_static_kernel_A.f90
@@ -32,12 +32,14 @@ subroutine UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
 
   double precision,intent(out)  :: KA(nSt,nSt)
 
+! Initialization
+
+  KA(:,:) = 0d0
+
 !--------------------------------------------------!
 ! Build BSE matrix for spin-conserving transitions !
 !--------------------------------------------------!
 
-  KA(:,:) = 0d0
-
   if(ispin == 1) then
 
     ! aaaa block
diff --git a/src/GW/UGW_phBSE_static_kernel_B.f90 b/src/GW/UGW_phBSE_static_kernel_B.f90
index d0bb9e1..d6c6d21 100644
--- a/src/GW/UGW_phBSE_static_kernel_B.f90
+++ b/src/GW/UGW_phBSE_static_kernel_B.f90
@@ -32,12 +32,14 @@ subroutine UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
 
   double precision,intent(out)  :: KB(nSt,nSt)
 
+! Initialization
+
+  KB(:,:) = 0d0
+
 !--------------------------------------------------!
 ! Build BSE matrix for spin-conserving transitions !
 !--------------------------------------------------!
 
-  KB(:,:) = 0d0
-
   if(ispin == 1) then
 
     ! aaaa block
diff --git a/src/LR/phULR_A.f90 b/src/LR/phULR_A.f90
index 19fed4e..e3fe0d6 100644
--- a/src/LR/phULR_A.f90
+++ b/src/LR/phULR_A.f90
@@ -39,6 +39,10 @@ subroutine phULR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda,eHF,ERI_aaaa,E
   delta_dRPA = 0d0
   if(dRPA) delta_dRPA = 1d0
 
+! Initialization
+
+  Aph(:,:) = 0d0
+
 !----------------------------------------------
 ! Build A matrix for spin-conserved transitions
 !----------------------------------------------
@@ -127,8 +131,6 @@ subroutine phULR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda,eHF,ERI_aaaa,E
 
   if(ispin == 2) then 
 
-    Aph(:,:) = 0d0
-
     ! abab block
 
     ia = 0
diff --git a/src/LR/phULR_B.f90 b/src/LR/phULR_B.f90
index a9a231e..69c0e8c 100644
--- a/src/LR/phULR_B.f90
+++ b/src/LR/phULR_B.f90
@@ -38,6 +38,10 @@ subroutine phULR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda,ERI_aaaa,ERI_a
   delta_dRPA = 0d0
   if(dRPA) delta_dRPA = 1d0
 
+! Initialization
+
+  Bph(:,:) = 0d0
+
 !-----------------------------------------------
 ! Build B matrix for spin-conserving transitions
 !-----------------------------------------------
@@ -124,8 +128,6 @@ subroutine phULR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,lambda,ERI_aaaa,ERI_a
 
   if(ispin == 2) then
 
-    Bph(:,:) = 0d0
-
     ! abba block
 
     ia = 0