diff --git a/src/GW/RGW_phACFDT.f90 b/src/GW/RGW_phACFDT.f90
index e2b2254..154fa6f 100644
--- a/src/GW/RGW_phACFDT.f90
+++ b/src/GW/RGW_phACFDT.f90
@@ -99,8 +99,8 @@ subroutine RGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,singlet,triplet,eta,nBas,
   call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
   call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
 
-  call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA)
-  call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB)
+               call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA)
+  if(.not.TDA) call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB)
 
 ! Singlet manifold
 
@@ -129,8 +129,8 @@ subroutine RGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,singlet,triplet,eta,nBas,
         call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
         call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
    
-        call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
-        call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
+                     call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
+        if(.not.TDA) call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
 
       end if
 
@@ -186,8 +186,8 @@ subroutine RGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,singlet,triplet,eta,nBas,
         call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
         call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
   
-        call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
-        call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
+                     call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
+        if(.not.TDA) call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
 
       end if  
 
diff --git a/src/GW/UGW_phACFDT.f90 b/src/GW/UGW_phACFDT.f90
index 9801c0c..9693d8d 100644
--- a/src/GW/UGW_phACFDT.f90
+++ b/src/GW/UGW_phACFDT.f90
@@ -112,14 +112,14 @@ subroutine UGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,spin_conserved,spin_flip,
   allocate(OmRPA(nS_sc),XpY_RPA(nS_sc,nS_sc),XmY_RPA(nS_sc,nS_sc),rho_RPA(nBas,nBas,nS_sc,nspin))
   allocate(Aph(nS_sc,nS_sc),Bph(nS_sc,nS_sc),KA(nS_sc,nS_sc),KB(nS_sc,nS_sc))
 
-               call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,1d0,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
-  if(.not.TDA) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,1d0,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
+                 call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,1d0,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
+  if(.not.TDA_W) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,1d0,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
 
   call phULR(TDA_W,nS_aa,nS_bb,nS_sc,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
   call UGW_excitation_density(nBas,nC,nO,nR,nS_aa,nS_bb,nS_sc,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY_RPA,rho_RPA)
 
-  call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,OmRPA,rho_RPA,KA)
-  call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,OmRPA,rho_RPA,KB)
+               call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,OmRPA,rho_RPA,KA)
+  if(.not.TDA) call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0,OmRPA,rho_RPA,KB)
 
 ! Spin-conserved manifold
 
@@ -144,24 +144,24 @@ subroutine UGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,spin_conserved,spin_flip,
 
       if(doXBS) then
 
-                     call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
-        if(.not.TDA) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
+                       call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
+        if(.not.TDA_W) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
    
         call phULR(TDA_W,nS_aa,nS_bb,nS_sc,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
         call UGW_excitation_density(nBas,nC,nO,nR,nS_aa,nS_bb,nS_sc,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY_RPA,rho_RPA)
 
-        call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KA)
-        call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KB)
+                     call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KA)
+        if(.not.TDA) call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KB)
 
       end if
 
-                   call phULR_A(isp_W,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
-      if(.not.TDA) call phULR_B(isp_W,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
+                   call phULR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
+      if(.not.TDA) call phULR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
 
                    Aph(:,:) = Aph(:,:) + KA(:,:)
       if(.not.TDA) Bph(:,:) = Bph(:,:) + KB(:,:)
 
-      call phULR(TDA_W,nS_aa,nS_bb,nS_sc,Aph,Bph,EcAC(ispin),OmRPA,XpY_RPA,XmY_RPA)
+      call phULR(TDA,nS_aa,nS_bb,nS_sc,Aph,Bph,EcAC(ispin),Om,XpY,XmY)
 
       call phUACFDT_correlation_energy(ispin,exchange_kernel,nBas,nC,nO,nV,nR,nS,nS_aa,nS_bb,nS_sc, &
                                        ERI_aaaa,ERI_aabb,ERI_bbbb,XpY,XmY,Ec(iAC,ispin))
@@ -208,24 +208,24 @@ subroutine UGW_phACFDT(exchange_kernel,doXBS,TDA_W,TDA,spin_conserved,spin_flip,
 
       if(doXBS) then
 
-                     call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
-        if(.not.TDA) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
+                       call phULR_A(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
+        if(.not.TDA_W) call phULR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
 
         call phULR(TDA_W,nS_aa,nS_bb,nS_sc,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
         call UGW_excitation_density(nBas,nC,nO,nR,nS_aa,nS_bb,nS_sc,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY_RPA,rho_RPA)
 
-        call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KA)
-        call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KB)
+                     call UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KA)
+        if(.not.TDA) call UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,lambda,OmRPA,rho_RPA,KB)
 
       end if
 
-                   call phULR_A(isp_W,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sf,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
-      if(.not.TDA) call phULR_B(isp_W,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sf,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
+                   call phULR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sf,lambda,eW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
+      if(.not.TDA) call phULR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sf,lambda,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
 
                    Aph(:,:) = Aph(:,:) + KA(:,:)
       if(.not.TDA) Bph(:,:) = Bph(:,:) + KB(:,:)
 
-      call phULR(TDA_W,nS_aa,nS_bb,nS_sf,Aph,Bph,EcAC(ispin),OmRPA,XpY_RPA,XmY_RPA)
+      call phULR(TDA,nS_aa,nS_bb,nS_sf,Aph,Bph,EcAC(ispin),Om,XpY,XmY)
 
       call phUACFDT_correlation_energy(ispin,exchange_kernel,nBas,nC,nO,nV,nR,nS,nS_ab,nS_ba,nS_sf, &
                                        ERI_aaaa,ERI_aabb,ERI_bbbb,XpY,XmY,Ec(iAC,ispin))
diff --git a/src/GW/UGW_phBSE_static_kernel_A.f90 b/src/GW/UGW_phBSE_static_kernel_A.f90
index 71eb7cf..3939043 100644
--- a/src/GW/UGW_phBSE_static_kernel_A.f90
+++ b/src/GW/UGW_phBSE_static_kernel_A.f90
@@ -57,7 +57,7 @@ subroutine UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,j,kc,1)*rho(a,b,kc,1)*Om(kc)/eps
             end do
  
-            KA(ia,jb) = KA(ia,jb) + 2d0*lambda*chi
+            KA(ia,jb) = 2d0*lambda**2*chi
  
           end do
         end do
@@ -81,7 +81,7 @@ subroutine UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,j,kc,2)*rho(a,b,kc,2)*Om(kc)/eps 
             end do
  
-            KA(nSa+ia,nSa+jb) = KA(nSa+ia,nSa+jb) + 2d0*lambda*chi
+            KA(nSa+ia,nSa+jb) = 2d0*lambda**2*chi
  
           end do
         end do
@@ -113,7 +113,7 @@ subroutine UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,j,kc,1)*rho(a,b,kc,2)*Om(kc)/eps
             end do
 
-            KA(ia,jb) = KA(ia,jb) + 2d0*lambda*chi
+            KA(ia,jb) = 2d0*lambda**2*chi
 
           end  do
         end  do
@@ -137,7 +137,7 @@ subroutine UGW_phBSE_static_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,j,kc,2)*rho(a,b,kc,1)*Om(kc)/eps
             end do
 
-            KA(nSa+ia,nSa+jb) = KA(nSa+ia,nSa+jb) + 2d0*lambda*chi
+            KA(nSa+ia,nSa+jb) = 2d0*lambda**2*chi
 
           end  do
         end  do
diff --git a/src/GW/UGW_phBSE_static_kernel_B.f90 b/src/GW/UGW_phBSE_static_kernel_B.f90
index 9d63068..d0bb9e1 100644
--- a/src/GW/UGW_phBSE_static_kernel_B.f90
+++ b/src/GW/UGW_phBSE_static_kernel_B.f90
@@ -57,7 +57,7 @@ subroutine UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,b,kc,1)*rho(a,j,kc,1)*Om(kc)/eps   
             end do
  
-            KB(ia,jb) = KB(ia,jb) + 2d0*lambda*chi
+            KB(ia,jb) = 2d0*lambda**2*chi
  
           end do
         end do
@@ -82,7 +82,7 @@ subroutine UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,b,kc,2)*rho(a,j,kc,2)*Om(kc)/eps
             end do
  
-            KB(nSa+ia,nSa+jb) = KB(nSa+ia,nSa+jb) + 2d0*lambda*chi
+            KB(nSa+ia,nSa+jb) = 2d0*lambda**2*chi
  
           end do
         end do
@@ -115,7 +115,7 @@ subroutine UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,b,kc,1)*rho(a,j,kc,2)*Om(kc)/eps
             end do
 
-            KB(ia,nSa+jb) = KB(ia,nSa+jb) + 2d0*lambda*chi
+            KB(ia,nSa+jb) = 2d0*lambda**2*chi
 
           end  do
         end  do
@@ -139,7 +139,7 @@ subroutine UGW_phBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nSa,nSb,nSt,nS_s
               chi = chi + rho(i,b,kc,2)*rho(a,j,kc,1)*Om(kc)/eps
             end do
 
-            KB(nSa+ia,jb) =  KB(nSa+ia,jb) + 2d0*lambda*chi
+            KB(nSa+ia,jb) =  2d0*lambda**2*chi
 
           end  do
         end  do