diff --git a/source/KKRhost/rhooutnew.f90 b/source/KKRhost/rhooutnew.f90
index 1ae6f9340d701225d10fdf2deac91a00ad5f2393..7bd9eee288c31de05c7808154a01391ff3e1dbbf 100644
--- a/source/KKRhost/rhooutnew.f90
+++ b/source/KKRhost/rhooutnew.f90
@@ -19,7 +19,7 @@ contains
!> Deprecated: False
!> Calculation of the density for the new solver
subroutine rhooutnew(nsra, lmax, gmatll, ek, lmpot, df, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, &
- rll, ull, rllleft, sllleft, cden, cdenlm, cdenns, rho2nsc, corbital, gflle_part, rpan_intervall, ipan_intervall, nspin)
+ rll, sll, ull, rllleft, sllleft, ullleft, cden, cdenlm, cdenns, rho2nsc, corbital, gflle_part, rpan_intervall, ipan_intervall, nspin)
use :: mod_constants, only: cone,czero,pi
use :: mod_runoptions, only: calc_gmat_lm_full, use_ldau, decouple_spin_cheby
@@ -51,10 +51,12 @@ contains
real (kind=dp), dimension (0:ntotd), intent (in) :: rpan_intervall
real (kind=dp), dimension (ntotd*(ncheb+1), nfund), intent (in) :: thetasnew
complex (kind=dp), dimension (lmmaxd, lmmaxd), intent (in) :: gmatll !! GMATLL=diagonal elements of the G matrix (system) Note that SLL is not needed for calculation of density, only needed for calculation of Green function
- complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: ull
complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: rll
+ complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: sll
+ complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: ull
complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: rllleft
complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: sllleft
+ complex (kind=dp), dimension (nsra*lmmaxd, lmmaxd, irmdnew), intent (in) :: ullleft
! .. Output variables
complex (kind=dp), dimension (irmdnew, nspin*(1+korbit)), intent (out) :: cdenns
diff --git a/source/KKRhost/rhovalnew.F90 b/source/KKRhost/rhovalnew.F90
index 0517062ded50b85798f88965e0c4f239e5c807c5..d4678a6f879176e13145e7648fbba8b57c5159f9 100644
--- a/source/KKRhost/rhovalnew.F90
+++ b/source/KKRhost/rhovalnew.F90
@@ -486,6 +486,19 @@ contains
lmmaxd,irmdnew,ith,nth,rll_was_read_in,sll_was_read_in, &
rllleft_was_read_in,sllleft_was_read_in)
end if
+ ! MdSD: TEST
+ ! if (myrank == 0) then
+ ! write(*,'("In rhovalnew:")')
+ ! write(*,'("nwfsavemax=",i4)') t_wavefunctions%nwfsavemax
+ ! write(*,'("save_rll=",l4)') t_wavefunctions%save_rll
+ ! write(*,'("save_sll=",l4)') t_wavefunctions%save_sll
+ ! write(*,'("save_rllleft=",l4)') t_wavefunctions%save_rllleft
+ ! write(*,'("save_sllleft=",l4)') t_wavefunctions%save_sllleft
+ ! write(*,'("rll_was_read_in=",l4)') rll_was_read_in
+ ! write(*,'("sll_was_read_in=",l4)') sll_was_read_in
+ ! write(*,'("rllleft_was_read_in=",l4)') rllleft_was_read_in
+ ! write(*,'("sllleft_was_read_in=",l4)') sllleft_was_read_in
+ ! end if
! recalculate wavefuntions, also include left solution
! contruct the spin-orbit coupling hamiltonian and add to potential
@@ -497,7 +510,7 @@ contains
end if
! Add magnetic field
- if ( t_params%bfield%lbfield .or. t_params%bfield%lbfield_constr ) then
+ if (t_params%bfield%lbfield) then
imt1 = ipan_intervall(t_params%npan_log+t_params%npan_eq) + 1
call add_bfield(t_params%bfield,i1,lmax,nspin,irmdnew,imt1,iend,ncheb,theta,phi,t_params%ifunm1(:,t_params%ntcell(i1)),&
t_params%icleb,t_params%cleb(:,1),t_params%thetasnew(1:irmdnew,:,t_params%ntcell(i1)),'1',vnspll2(:,:,:), &
@@ -554,15 +567,19 @@ contains
! faster calculation of RLL.
! no irregular solutions SLL are needed in self-consistent iterations
! because the density depends only on RLL, RLLLEFT and SLLLEFT
- if (.not.set_cheby_nospeedup .and. .not. (calc_exchange_couplings .or. write_pkkr_operators)) then
- call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
- hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
-! call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), rll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
-! hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
- else
- call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rll(:,:,:,ith), sll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
- hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alphall)
- end if
+ call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
+ ! MdSD: right now it seems that sll is not used for anything
+ ! if (set_cheby_nospeedup .or. calc_exchange_couplings .or. write_pkkr_operators) then
+ ! call sll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), sll(:,:,:,ith), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ ! nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick)
+ ! MdSD: this is the old interface
+ ! call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rll(:,:,:,ith), sll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
+ ! hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alphall)
+ ! MdSD: if using the old rllsll this is needed for rhooutnew
+ ! ull(:,:,:,ith) = rll(:,:,:,ith)
+ ! end if
+ ! MdSD: if using the old rllsll check if this is needed
if (nsra==2) then
ull(lmmaxd+1:nvec*lmmaxd, :, :, ith) = ull(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
rll(lmmaxd+1:nvec*lmmaxd, :, :, ith) = rll(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
@@ -584,7 +601,7 @@ contains
end if
! Add magnetic field
- if ( t_params%bfield%lbfield .or. t_params%bfield%lbfield_constr ) then
+ if (t_params%bfield%lbfield) then
call add_bfield(t_params%bfield,i1,lmax,nspin,irmdnew,imt1,iend,ncheb,theta,phi,t_params%ifunm1(:,t_params%ntcell(i1)),&
t_params%icleb,t_params%cleb(:,1),t_params%thetasnew(1:irmdnew,:,t_params%ntcell(i1)),'transpose',vnspll2(:,:,:), &
vnspll1(:,:,:,ith),t_params%bfield%thetallmat(:,:,1:irmdnew,t_params%ntcell(i1)))
@@ -632,15 +649,14 @@ contains
tmattemp = czero
alphall = czero
! faster calculation of RLLLEFT and SLLLEFT.
- if (.not.set_cheby_nospeedup .and. .not. (calc_exchange_couplings .or. write_pkkr_operators)) then
- call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ullleft(:,:,:,ith), rllleft(:,:,:,ith), tmattemp, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
- jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
- call sll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), sllleft(:,:,:,ith), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
- hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick)
- else
- call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), tmattemp, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
- jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alphall)
- end if
+ call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ullleft(:,:,:,ith), rllleft(:,:,:,ith), tmattemp, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
+ call sll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), sllleft(:,:,:,ith), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick)
+ ! MdSD: this is the old interface
+ ! call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), tmattemp, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
+ ! jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alphall)
+ ! MdSD: if using the old rllsll check if this is needed
if (nsra==2) then
ullleft(lmmaxd+1:nvec*lmmaxd, :, :, ith) = ullleft(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
rllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith) = rllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
@@ -680,9 +696,10 @@ contains
end do
! calculate density
-
- call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, df, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, rll(:,:,:,ith), & ! SLL(:,:,:,ith), commented out since sll is not used in rhooutnew
- ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), rho2nsc_loop(:,:,:,ie), 0, gflle(:,:,ie,iq), rpan_intervall, ipan_intervall, nspin/(nspin-korbit))
+ ! MdSD: modified the interface to allow for all types of scattering solutions; some can be dummy arguments
+ call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, df, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, &
+ rll(:,:,:,ith), sll(:,:,:,ith), ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), ullleft(:,:,:,ith), &
+ cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), rho2nsc_loop(:,:,:,ie), 0, gflle(:,:,ie,iq), rpan_intervall, ipan_intervall, nspin/(nspin-korbit))
do jspin = 1, nspin/(nspin-korbit)*(1+korbit)
@@ -730,9 +747,10 @@ contains
! Get charge at the Fermi energy (IELAST)
!------------------------------------------------------------------------------
if (ie==ielast .and. ldorhoef) then
- call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, cone, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, rll(:,:,:,ith), & ! SLL(:,:,:,ith), ! commented out since sll is not used in rhooutnew
- ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), r2nefc_loop(:,:,:,ith), 0, gflle_part(:,:,ith), rpan_intervall, &
- ipan_intervall, nspin/(nspin-korbit))
+ ! MdSD: modified the interface to allow for all types of scattering solutions; some can be dummy arguments
+ call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, cone, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, &
+ rll(:,:,:,ith), sll(:,:,:,ith), ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), ullleft(:,:,:,ith), &
+ cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), r2nefc_loop(:,:,:,ith), 0, gflle_part(:,:,ith), rpan_intervall, ipan_intervall, nspin/(nspin-korbit))
end if
!------------------------------------------------------------------------------
@@ -744,8 +762,10 @@ contains
proj(2) = sin(phi)*sin(theta)
proj(3) = cos(theta)
do iorb = 1, 3
- call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, cone, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, rll(:,:,:,ith), & ! SLL(:,:,:,ith), ! commented out since sll is not used in rhooutnew
- ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), r2orbc(:,:,:,ith), iorb, gflle_part(:,:,ith), rpan_intervall, ipan_intervall, nspin)
+ ! MdSD: modified the interface to allow for all types of scattering solutions; some can be dummy arguments
+ call rhooutnew(nsra, lmax, gmatll(1,1,ie), ek, lmpotd, cone, npan_tot, ncheb, cleb, icleb, iend, irmdnew, thetasnew, ifunm, imt1, lmsp, &
+ rll(:,:,:,ith), sll(:,:,:,ith), ull(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), ullleft(:,:,:,ith), &
+ cden(:,:,:,ith), cdenlm(:,:,:,ith), cdenns(:,:,ith), r2orbc(:,:,:,ith), iorb, gflle_part(:,:,ith), rpan_intervall, ipan_intervall, nspin)
do jspin = 1, nspin*(1+korbit)
if (jspin<=2) then
do lm1 = 0, lmax
diff --git a/source/KKRhost/tmat_newsolver.F90 b/source/KKRhost/tmat_newsolver.F90
index 9d4a327850c90b979bfc528e1432b0178c8f9945..2cfac3cb47ff60fcf4df26c7a39ec413e542582b 100644
--- a/source/KKRhost/tmat_newsolver.F90
+++ b/source/KKRhost/tmat_newsolver.F90
@@ -56,6 +56,7 @@ contains
use :: mod_jijhelp, only: calc_dtmatjij
use :: mod_calcsph, only: calcsph
use :: mod_rll_global_solutions, only: rll_global_solutions
+ use :: mod_sll_global_solutions, only: sll_global_solutions ! MdSD: TEST
use :: mod_rllsllsourceterms, only: rllsllsourceterms
use :: mod_rllsll, only: rllsll
use :: mod_spinorbit_ham, only: spinorbit_ham
@@ -278,6 +279,15 @@ contains
if (stop_1b .and. .not. write_pkkr_operators) then
t_wavefunctions%nwfsavemax = 0
end if
+ ! MdSD: TEST
+ ! if (myrank == 0) then
+ ! write(*,'("In tmat_newsolver:")')
+ ! write(*,'("nwfsavemax=",i4)') t_wavefunctions%nwfsavemax
+ ! write(*,'("save_rll=",l4)') t_wavefunctions%save_rll
+ ! write(*,'("save_sll=",l4)') t_wavefunctions%save_sll
+ ! write(*,'("save_rllleft=",l4)') t_wavefunctions%save_rllleft
+ ! write(*,'("save_sllleft=",l4)') t_wavefunctions%save_sllleft
+ ! end if
end if
#ifdef CPP_OMP
@@ -356,7 +366,7 @@ contains
end if
! Add magnetic field
- if ( t_params%bfield%lbfield .or. t_params%bfield%lbfield_constr ) then
+ if (t_params%bfield%lbfield) then
! MdSD: constraining fields
if (t_inc%i_write>1) then
write (1337,'("tmat_newsolver: myrank=",i8," iatom=",i8)') myrank, i1
@@ -513,16 +523,17 @@ contains
! faster calculation of RLL.
! no irregular solutions are needed in self-consistent iterations
! because the t-matrix depends only on RLL
- if (.not. set_cheby_nospeedup .and. .not. (calc_exchange_couplings .or. write_pkkr_operators) .and. .not.calc_wronskian) then
-! call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rll(:,:,:,ith), tmat0(:,:), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
-! jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alpha0(:,:))
- call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, &
- hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
- else
- call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rll(:,:,:,ith), sll(:,:,:,ith), tmat0(:,:), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
- jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alpha0(:,:))
+ call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rll(:,:,:,ith), tmatll, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick, alphall)
+ ! MdSD: check if these are actually needed except in rhovalnew
+ if (calc_wronskian) then
+ call sll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), sll(:,:,:,ith), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', use_sratrick)
end if
-
+ ! MdSD: this is the old interface
+ ! call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rll(:,:,:,ith), sll(:,:,:,ith), tmat0(:,:), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
+ ! jlk_index, hlk(:,:,ith), jlk(:,:,ith), hlk2(:,:,ith), jlk2(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alpha0(:,:))
+ ! MdSD: if using the old rllsll check if this is needed
if (nsra==2) then
rll(lmmaxd+1:nvec*lmmaxd, :, :, ith) = rll(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
sll(lmmaxd+1:nvec*lmmaxd, :, :, ith) = sll(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
@@ -622,7 +633,7 @@ contains
if ( t_dtmatjij_at%calculate .or. (t_wavefunctions%isave_wavefun(i1,ie)>0 .and. &
(t_wavefunctions%save_rllleft .or. t_wavefunctions%save_sllleft)) .or. &
((write_rhoq_input .and. ie==2) .and. (i1==mu0)) .or. & ! rhoqtest
- calc_wronskian ) then
+ calc_exchange_couplings .or. write_pkkr_operators .or. calc_wronskian ) then ! MdSD: seems to make more sense to check here than below
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Calculate the left-hand side solution this needs to be done for the
! calculation of t-matrices for Jij tensor or if wavefunctions should be saved
@@ -634,7 +645,7 @@ contains
nrmaxd,vnspll0(:,:,:),vnspll2(:,:,:),'transpose')
! Add magnetic field
- if ( t_params%bfield%lbfield .or. t_params%bfield%lbfield_constr ) then
+ if (t_params%bfield%lbfield) then
call add_bfield(t_params%bfield,i1,lmax,nspin,irmdnew,imt1,iend,ncheb,theta,phi,t_params%ifunm1(:,t_params%ntcell(i1)),&
t_params%icleb,t_params%cleb(:,1),t_params%thetasnew(1:irmdnew,:,t_params%ntcell(i1)),'transpose',vnspll2(:,:,:), &
vnspll1(:,:,:,ith),t_params%bfield%thetallmat(:,:,1:irmdnew,t_params%ntcell(i1)))
@@ -683,15 +694,17 @@ contains
! notice that exchange the order of left and right hankel/bessel functions
tmat0 = czero
alpha0 = czero ! LLY
- ! faster calculation of RLL.
- ! no left solutions are needed in self-consistent iterations
- ! because the t-matrix depends only on RLL
- if (.not. set_cheby_nospeedup .and. .not. ( calc_exchange_couplings .or. write_pkkr_operators) .and. .not.calc_wronskian) then
- ! do nothing
- else
- call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), tmat0, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
- jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alpha0)
+ call rll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), ull(:,:,:,ith), rllleft(:,:,:,ith), tmat0, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick, alpha0)
+ ! MdSD: check if these are actually needed except in rhovalnew
+ if (calc_wronskian) then
+ call sll_global_solutions(rpan_intervall, rnew, vnspll(:,:,:,ith), sllleft(:,:,:,ith), ncheb, npan_tot, lmmaxd, nvec*lmmaxd, &
+ nsra*(1+korbit)*(lmax+1), irmdnew, nsra, jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', use_sratrick)
end if
+ ! MdSD: this is the old interface
+ ! call rllsll(rpan_intervall, rnew, vnspll(:,:,:,ith), rllleft(:,:,:,ith), sllleft(:,:,:,ith), tmat0, ncheb, npan_tot, lmmaxd, nvec*lmmaxd, nsra*(1+korbit)*(lmax+1), irmdnew, nsra, &
+ ! jlk_index, hlk2(:,:,ith), jlk2(:,:,ith), hlk(:,:,ith), jlk(:,:,ith), gmatprefactor, '1', '1', '0', use_sratrick, alpha0)
+ ! MdSD: if using the old rllsll check if this is needed
if (nsra==2) then
rllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith) = rllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight
sllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith) = sllleft(lmmaxd+1:nvec*lmmaxd, :, :, ith)/cvlight