diff --git a/source/KKRhost/bzkint0.f90 b/source/KKRhost/bzkint0.f90
index db35587ee33237ce888c9d66d06f699e5ccd08bc..eca5ea687a0e23abe48a13d469aad886de6ee888 100644
--- a/source/KKRhost/bzkint0.f90
+++ b/source/KKRhost/bzkint0.f90
@@ -25,7 +25,7 @@ contains
     nsheld, nembd, iemxd)
 
     use :: mod_constants, only: nsymaxd
-    use :: mod_runoptions, only: print_tau_structure, use_Chebychev_solver, use_full_BZ
+    use :: mod_runoptions, only: print_tau_structure, use_Chebychev_solver, use_full_BZ, decouple_spin_cheby, force_BZ_symm
     use :: mod_datatypes, only: dp
     use :: mod_gfshells, only: gfshells
     use :: mod_crtstar, only: crtstar
@@ -77,10 +77,15 @@ contains
     if (print_tau_structure) iprint = 2
 
     ! --> test: full BZ integration
-    if (use_full_BZ .or. use_Chebychev_solver) then
+    if (use_full_BZ .or. (use_Chebychev_solver .and. .not.decouple_spin_cheby) &
+        .or.(use_Chebychev_solver .and. .not.force_BZ_symm)) then
       nsymat = 1
       lirr = .false.
       write (1337, '(8X,2A,/)') 'Run option <use_full_BZ> or <use_Chebychev_solver>: ', ' overriding NSYMAT, generate full BZ k-mesh'
+    elseif (use_Chebychev_solver .and. force_BZ_symm) then
+      write (1337, '(8X,2A,/)') 'Run option <use_Chebychev_solver> and <force_BZ_symm>: ', ' use symmetries for k-mesh (can be deactivated with <use_full_BZ>)'
+    elseif (use_Chebychev_solver .and. decouple_spin_cheby) then
+      write (1337, '(8X,2A,/)') 'Run option <use_Chebychev_solver> and <decouple_spin_cheby>: ', ' use symmetries for k-mesh (can be deactivated with <use_full_BZ>)'
     end if
 
     ! --> generate BZ k-mesh
diff --git a/source/KKRhost/kkrmat01.F90 b/source/KKRhost/kkrmat01.F90
index 37f75774bd416b25cb4e3d5c51d3dcdd88be02df..ae581dc5ebb7e29d12a8d9770785aced61fa51ac 100644
--- a/source/KKRhost/kkrmat01.F90
+++ b/source/KKRhost/kkrmat01.F90
@@ -52,7 +52,7 @@ contains
 #endif
     use :: mod_types, only: t_inc
     use :: mod_runoptions, only: print_program_flow, use_Chebychev_solver, use_qdos, use_virtual_atoms, &
-      write_green_imp, write_rhoq_input, decouple_spins_cheby
+      write_green_imp, write_rhoq_input, decouple_spin_cheby
     use :: mod_rhoqtools, only: rhoq_find_kmask, rhoq_saveg, rhoq_write_tau0, rhoq_read_mu0_scoef
     use :: global_variables, only: nembd1, nembd2, nsheld, nclsd, naclsd, lmmaxd, nprincd, nrd, nrefd, lmgf0d, krel, ndim_slabinv, alm, almgf0 
     use :: mod_constants, only: czero, cone, nsymaxd, ci,pi
@@ -343,7 +343,7 @@ contains
         ! ----------------------------------------------------------------------
         ! LLY Lloyd
         ! --------------------------------------------------------------------
-        if (.not. use_Chebychev_solver .or. decouple_spins_cheby) then
+        if (.not. use_Chebychev_solver .or. decouple_spin_cheby) then
           ! $omp single
           do i2 = 1, alm
             do i1 = 1, alm
@@ -635,7 +635,16 @@ contains
         if ((((k_end-k_start)/50)==0 .or. mod(kpt-k_start,(k_end-k_start)/50)==0) .and. t_inc%i_write>0) write (1337, fmt=110, advance='no')
         ! $omp end critical
 #else
-        if (((nofks/50)==0 .or. mod(kpt,nofks/50)==0) .and. t_inc%i_write>0) write (1337, fmt=110, advance='no')
+        !if (((nofks/50)==0 .or. mod(kpt,nofks/50)==0) .and. t_inc%i_write>0) write (1337, fmt=110, advance='no')
+        if (t_inc%i_write>0) then
+          if (nofks/50==0) then
+            write (1337, fmt=110, advance='no')
+          else
+            if (mod(kpt,nofks/50)==0) then
+              write (1337, fmt=110, advance='no')
+            end if
+          end if
+        end if
 #endif
       end if                       ! mythread==0
 
diff --git a/source/KKRhost/main1b.F90 b/source/KKRhost/main1b.F90
index d03d3c1818365fb88696d9e69fb9aa1633b74848..657266646168633aa8bcf2f46018c4d04dedb70a 100644
--- a/source/KKRhost/main1b.F90
+++ b/source/KKRhost/main1b.F90
@@ -37,24 +37,24 @@ contains
     use :: mod_datatypes, only: dp, sp
     use :: mod_runoptions, only: calc_exchange_couplings, formatted_file, set_gmat_to_zero, use_Chebychev_solver, &
       use_qdos, use_readcpa, write_deci_tmat, write_gmat_plain, write_green_host, write_green_imp, write_kkrimp_input, &
-      write_pkkr_input, write_pkkr_operators, write_rhoq_input, write_gmat_ascii, decouple_spins_cheby
-    use :: mod_constants, only: czero, cone, pi
-    use :: mod_operators_for_fscode, only: operators_for_fscode
-    use :: mod_getscratch, only: opendafile
-    use :: mod_kloopz1, only: kloopz1_qdos
-    use :: mod_greenimp, only: greenimp
-    use :: mod_changerep, only: changerep
-    use :: mod_tmatimp_newsolver, only: tmatimp_newsolver
-    use :: mod_setfactl, only: setfactl
-    use :: mod_calctref13, only: calctref13
-    use :: mod_rotatespinframe, only: rotatematrix
-    use :: mod_types, only: t_tgmat, t_inc, t_lloyd, t_cpa, init_t_cpa, t_imp
-    use :: mod_timing, only: timing_start, timing_pause, timing_stop, timings_1b, print_time_and_date
-    use :: mod_wunfiles, only: get_params_1b, t_params !, read_angles
-    use :: mod_tbxccpljij, only: tbxccpljij
-    use :: mod_tbxccpljijdij, only: tbxccpljijdij
-    use :: mod_rhoqtools, only: rhoq_save_refpot
-    use :: mod_cinit, only: cinit
+      write_pkkr_input, write_pkkr_operators, write_rhoq_input, write_gmat_ascii, decouple_spin_cheby
+    use mod_constants, only: czero, cone, pi
+    use mod_operators_for_fscode, only: operators_for_fscode
+    use mod_getscratch, only: opendafile
+    use mod_kloopz1, only: kloopz1_qdos
+    use mod_greenimp, only: greenimp
+    use mod_changerep, only: changerep
+    use mod_tmatimp_newsolver, only: tmatimp_newsolver
+    use mod_setfactl, only: setfactl
+    use mod_calctref13, only: calctref13
+    use mod_rotatespinframe, only: rotatematrix
+    use mod_types, only: t_tgmat, t_inc, t_lloyd, t_cpa, init_t_cpa, t_imp
+    use mod_timing, only: timing_start, timing_pause, timing_stop, timings_1b, print_time_and_date
+    use mod_wunfiles, only: get_params_1b, t_params
+    use mod_tbxccpljij, only: tbxccpljij
+    use mod_tbxccpljijdij, only: tbxccpljijdij
+    use mod_rhoqtools, only: rhoq_save_refpot
+    use mod_cinit, only: cinit
     ! array dimensions
     use :: global_variables, only: maxmshd, iemxd, natypd, naezd, kpoibz, lmmaxd, lmgf0d, lmaxd, nrefd, nsheld, wlength, nofgij, &
       naclsd, nspind, nclsd, nembd, krel, korbit, natomimpd, nrd, nembd1, nspindd, nprincd, irmind, nspotd, irmd, lpotd, &
@@ -522,9 +522,9 @@ contains
             end if
             do i = 1, lm1
               trefll(i, i, i1) = wn1(i, i)
-              if (use_Chebychev_solver .and. .not.decouple_spins_cheby) trefll(lm1+i, lm1+i, i1) = wn1(i, i)
+              if (use_Chebychev_solver .and. .not.decouple_spin_cheby) trefll(lm1+i, lm1+i, i1) = wn1(i, i)
               dtrefll(i, i, i1) = wn2(i, i)                              ! LLY
-              if (use_Chebychev_solver .and. .not.decouple_spins_cheby) dtrefll(lm1+i, lm1+i, i1) = wn2(i, i) ! LLY
+              if (use_Chebychev_solver .and. .not.decouple_spin_cheby) dtrefll(lm1+i, lm1+i, i1) = wn2(i, i) ! LLY
             end do
 
             if (write_rhoq_input) then
@@ -567,7 +567,7 @@ contains
             tmat(:, :) = t_tgmat%tmat(:, :, irec)
           end if
 
-          if (use_Chebychev_solver .and. .not. decouple_spins_cheby) then
+          if (use_Chebychev_solver .and. .not. decouple_spin_cheby) then
             ! read in theta and phi for noncolinear
             theta = theta_at(i1)
             phi = phi_at(i1)
@@ -587,7 +587,7 @@ contains
               tmat(:, :) = t_lloyd%dtmat(:, :, irec)
             end if
 
-            if (use_Chebychev_solver .and. .not. decouple_spins_cheby) call rotatematrix(tmat, theta, phi, lmgf0d, 0) ! LLY
+            if (use_Chebychev_solver .and. .not. decouple_spin_cheby) call rotatematrix(tmat, theta, phi, lmgf0d, 0) ! LLY
 
             dtmatll(1:lmmaxd, 1:lmmaxd, i1) = tmat(1:lmmaxd, 1:lmmaxd) ! LLY
             if (t_lloyd%dtmat_to_file) then
@@ -752,7 +752,7 @@ contains
 
         if (lly/=0) then           ! LLY
 
-          if (use_Chebychev_solver .and. .not.decouple_spins_cheby) then
+          if (use_Chebychev_solver .and. .not.decouple_spin_cheby) then
             cdos_lly(ie, ispin) = tralpha(ie, ispin) - lly_grtr(ie, ispin)/volbz(1) + 2.0_dp*lly_g0tr(ie) ! LLY
           else
             if (lly/=2) then       ! LLY Lloyd
diff --git a/source/KKRhost/main1c.F90 b/source/KKRhost/main1c.F90
index 1a4c8f23e3cc3001fc36819f019b6cac1b9e4158..6452237519e26c3d68a71e97e8fd57815029e825 100644
--- a/source/KKRhost/main1c.F90
+++ b/source/KKRhost/main1c.F90
@@ -41,7 +41,7 @@ contains
     use mod_datatypes, only: dp
     use mod_runoptions, only: calc_gmat_lm_full, fix_nonco_angles, relax_SpinAngle_Dirac, use_Chebychev_solver, &
       use_decimation, use_qdos, write_DOS, write_complex_qdos, write_density_ascii, write_rho2ns, write_DOS_lm, &
-      decouple_spins_cheby, disable_print_serialnumber
+      decouple_spin_cheby, disable_print_serialnumber
     use mod_constants, only: czero, pi
     use mod_profiling, only: memocc
     use mod_mympi, only: myrank, master
@@ -811,8 +811,8 @@ contains
       ! Write out lm charges and moments
       ! -------------------------------------------------------------------
       withorbmom = krel+korbit
-      if (decouple_spins_cheby) withorbmom = nspin-1 !withorbmom+1
-      call wrmoms(withorbmom, natyp, nspinpot, charge, muorb, lmax, lmaxd1)
+      if (decouple_spin_cheby) withorbmom = nspin-1
+      if (.not. use_qdos) call wrmoms(withorbmom, natyp, nspinpot, charge, muorb, lmax, lmaxd1)
 
       ! ----------------------------------------------------------------------
       ! ITERMDIR
diff --git a/source/KKRhost/mixnocospin.f90 b/source/KKRhost/mixnocospin.f90
index 601e1957dabbc48df487b738a17510f8dac0e46c..bf5c23c48462dfe09df8ba4234ae863fdf177ecd 100644
--- a/source/KKRhost/mixnocospin.f90
+++ b/source/KKRhost/mixnocospin.f90
@@ -23,7 +23,7 @@ contains
   subroutine spinmix_noco(iter, natyp, theta, phi, fixdir, angles_new, totmoment, iounit)
     use mod_datatypes, only: dp
     use mod_constants, only: pi
-    use mod_runoptions, only: use_broyden_spinmix, write_angles_alliter, disable_print_serialnumber, fix_nonco_angles, decouple_spins_cheby
+    use mod_runoptions, only: use_broyden_spinmix, write_angles_alliter, disable_print_serialnumber, fix_nonco_angles, decouple_spin_cheby
     use global_variables, only: qbound_spin
     use mod_version_info, only: version_print_header
     use mod_wunfiles, only: t_params
@@ -45,7 +45,7 @@ contains
 
 
     ! MdSD,PR: write information on new angles to output file
-    if (.not.decouple_spins_cheby) then
+    if (.not.decouple_spin_cheby) then
 
       write (1337,*)
       write (1337, '("      I1    In/Out THETA[deg]       In/Out PHI[deg]        FIXDIR[boolean]   RMS(angles)[deg]")')
@@ -77,7 +77,7 @@ contains
         end do
       end if
 
-    end if ! .not.decouple_spins_cheby
+    end if ! .not.decouple_spin_cheby
   
     ! rewrite new theta and phi to nonco_angle_out.dat, nonco_angle.dat is the input
     if (.not. fix_nonco_angles) then
diff --git a/source/KKRhost/renorm_lly.f90 b/source/KKRhost/renorm_lly.f90
index 8543348d645fe74fe3ecbbbc2d726de52a85be51..b2fb75d9dfd37470a5d5fa06c66b4a36ba452244 100644
--- a/source/KKRhost/renorm_lly.f90
+++ b/source/KKRhost/renorm_lly.f90
@@ -29,7 +29,7 @@ contains
   subroutine renorm_lly(cdos_lly, ielast, nspin, natyp, cden, lmaxp1, conc, iestart, ieend, wez, ircut, ipan, ez, zat, rho2ns, r2nef, denef, denefat, espv)
 
     use :: mod_datatypes, only: dp
-    use :: mod_runoptions, only: use_Chebychev_solver, decouple_spins_cheby
+    use :: mod_runoptions, only: use_Chebychev_solver, decouple_spin_cheby
     use :: mod_constants, only: czero, pi
     use :: global_variables, only: ipand, natypd, lmaxd, npotd, iemxd, irmd, lmpotd, krel, nspind 
     implicit none
@@ -91,7 +91,7 @@ contains
       end do
     end do
     ! Now the locally-summed charge/energy is in cdos_loc, charge/energy/atom in chadd
-    if (.not. use_Chebychev_solver .or. decouple_spins_cheby) then
+    if (.not. use_Chebychev_solver .or. decouple_spin_cheby) then
       do ie = iestart, ieend
         do ispin = 1, nspin
           ! Renormalization factor per energy:
@@ -130,7 +130,7 @@ contains
     !    ENDDO
     !   ENDDO
 
-    if (nspin==1 .or. (use_Chebychev_solver .and. .not. decouple_spins_cheby) ) cren(:, 2) = cren(:, 1)
+    if (nspin==1 .or. (use_Chebychev_solver .and. .not. decouple_spin_cheby) ) cren(:, 2) = cren(:, 1)
 
 
     ! Now apply renormalization to energy-integrated density
diff --git a/source/KKRhost/rhooutnew.f90 b/source/KKRhost/rhooutnew.f90
index 3218e7895d2188b2e2a9f39a992d5ed0b6bab661..1ae6f9340d701225d10fdf2deac91a00ad5f2393 100644
--- a/source/KKRhost/rhooutnew.f90
+++ b/source/KKRhost/rhooutnew.f90
@@ -22,7 +22,7 @@ contains
      rll, ull, rllleft, sllleft, 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_spins_cheby
+    use :: mod_runoptions, only: calc_gmat_lm_full, use_ldau, decouple_spin_cheby
     use :: mod_profiling, only: memocc
     use :: global_variables, only: lmmaxd, ncleb, ntotd, nfund, korbit
     use :: mod_datatypes, only: dp
@@ -103,7 +103,7 @@ contains
     pnsi = czero
 
     ! set LMSHIFT value which is need to construct CDEN
-    if (decouple_spins_cheby) then
+    if (decouple_spin_cheby) then
       lmshift1(:) = 0 
       lmshift2(:) = 0 
     else
diff --git a/source/KKRhost/rhovalnew.F90 b/source/KKRhost/rhovalnew.F90
index 5441f6e64bf490ec210c8e788d7e0427f5e8e5b2..6ebde86ade65de2120cd2106f7b85380e56c6c65 100644
--- a/source/KKRhost/rhovalnew.F90
+++ b/source/KKRhost/rhovalnew.F90
@@ -41,7 +41,7 @@ contains
     use :: mod_save_wavefun, only: t_wavefunctions, read_wavefunc
     use :: mod_runoptions, only: calc_exchange_couplings, calc_gmat_lm_full, disable_tmat_sratrick, fix_nonco_angles, &
                                  use_qdos, write_complex_qdos, write_pkkr_operators, write_DOS_lm, set_cheby_nospeedup, &
-                                 decouple_spins_cheby, disable_print_serialnumber, set_gmat_to_zero
+                                 decouple_spin_cheby, disable_print_serialnumber, set_gmat_to_zero
     use :: mod_version_info, only: version_print_header
     use :: global_variables, only: lmmaxd, iemxd, ncleb, lmxspd, irmd, ntotd, nrmaxd, lmpotd, nspotd, nfund, korbit, mmaxd, nspind
     use :: mod_constants, only: czero, cvlight, cone, pi, ci
@@ -224,7 +224,7 @@ contains
     call memocc(i_stat, product(shape(vins))*kind(vins), 'VINS', 'RHOVALNEW')
     vins = 0.0_dp
     vins(1:irmdnew, 1:lmpotd, 1) = vinsnew(1:irmdnew, 1:lmpotd, ipot)
-    if (.not.decouple_spins_cheby)  vins(1:irmdnew, 1:lmpotd, nspin) = vinsnew(1:irmdnew, 1:lmpotd, ipot+nspin-1)
+    if (.not.decouple_spin_cheby)  vins(1:irmdnew, 1:lmpotd, nspin) = vinsnew(1:irmdnew, 1:lmpotd, ipot+nspin-1)
 
     ! set up the non-spherical ll' matrix for potential VLL'
     allocate (vnspll0(lmmaxd,lmmaxd,irmdnew), stat=i_stat)
@@ -477,7 +477,7 @@ contains
 
       ! recalculate wavefuntions, also include left solution
       ! contruct the spin-orbit coupling hamiltonian and add to potential
-      if ( .not. decouple_spins_cheby) then
+      if ( .not. decouple_spin_cheby) then
         call spinorbit_ham(lmax, lmmax0d, vins, rnew, eryd, zat, cvlight, socscale, nspin, lmpotd, theta, phi, ipan_intervall, &
           rpan_intervall, npan_tot, ncheb, irmdnew, nrmaxd, vnspll0, vnspll1(:,:,:,ith), '1')
       else
@@ -508,7 +508,7 @@ contains
         jlk2(:, :, ith) = czero
         gmatprefactor = czero
         jlk_index = 0
-        if (decouple_spins_cheby) then
+        if (decouple_spin_cheby) then
           use_fullgmat = 0
         else
           use_fullgmat = 1
@@ -554,7 +554,7 @@ contains
       !------------------------------------------------------------------------------
       if ((t_wavefunctions%nwfsavemax>0 .and. (.not. (rllleft_was_read_in .and. sllleft_was_read_in))) .or. (t_wavefunctions%nwfsavemax==0)) then
         ! read/recalc wavefunctions left contruct the TRANSPOSE spin-orbit coupling hamiltonian and add to potential
-        if ( .not. decouple_spins_cheby) then
+        if ( .not. decouple_spin_cheby) then
           call spinorbit_ham(lmax, lmmax0d, vins, rnew, eryd, zat, cvlight, socscale, nspin, lmpotd, theta, phi, ipan_intervall, rpan_intervall, npan_tot, ncheb, irmdnew, nrmaxd, &
             vnspll0, vnspll1(:,:,:,ith), 'rll')
         else
@@ -637,7 +637,7 @@ contains
         ! $omp end critical
 #endif
 
-        if ( .not. decouple_spins_cheby) then
+        if ( .not. decouple_spin_cheby) then
           ! rotate gmat from global frame to local frame
           call rotatematrix(gmat0, theta, phi, lmmax0d, 1)
         end if
@@ -706,7 +706,7 @@ contains
       !------------------------------------------------------------------------------
       ! Get orbital moment
       ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-      if (.not. decouple_spins_cheby) then
+      if (.not. decouple_spin_cheby) then
         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)
@@ -745,7 +745,7 @@ contains
             muorb(lmaxd1+1, 1:3) = muorb(lmaxd1+1, 1:3) + muorb(lm1, 1:3)
           end do
         end do ! IORB
-      end if ! .not. decouple_spins_cheby
+      end if ! .not. decouple_spin_cheby
 
     end do                         ! IE loop
 #ifdef CPP_OMP
@@ -1029,9 +1029,9 @@ contains
       ! MdSD: Should this also be corrected if the angles change?
       den_out(0:lmaxd1, 1:ielast, 1:nspin/(nspin-korbit)) = den(0:lmaxd1, 1:ielast, 1, 1:nspin/(nspin-korbit))
       ! calculate new THETA and PHI for non-colinear
-      ! if (.not. fix_nonco_angles .and. .not.decouple_spins_cheby) then
+      ! if (.not. fix_nonco_angles .and. .not.decouple_spin_cheby) then
       ! MdSD: now the new directions are always calculated, which can be useful for information purposes
-      if (.not.decouple_spins_cheby) then
+      if (.not.decouple_spin_cheby) then
         rho2ns_temp(1, 1) = rho2int(1)
         rho2ns_temp(2, 2) = rho2int(2)
         rho2ns_temp(1, 2) = rho2int(3)
diff --git a/source/KKRhost/rinput13.F90 b/source/KKRhost/rinput13.F90
index e9a1c912cf62d41b76209b9d801a8b3ff4633ee2..ce5db65a9bba6575536b0b7fc81a899e1476dd73 100644
--- a/source/KKRhost/rinput13.F90
+++ b/source/KKRhost/rinput13.F90
@@ -44,7 +44,7 @@ contains
       dirac_scale_SpeefOfLight, disable_charge_neutrality, disable_print_serialnumber, modify_soc_Dirac, relax_SpinAngle_Dirac, search_Efermi, &
     set_kmesh_large, stop_1b, stop_1c, use_BdG, use_Chebychev_solver, use_cond_LB, use_decimation, use_lloyd, use_qdos, &
     use_rigid_Efermi, use_semicore, use_virtual_atoms, write_green_host, write_green_imp, write_kkrimp_input, &
-    write_pkkr_input, write_pkkr_operators, use_ldau, set_cheby_nospeedup, decouple_spins_cheby, write_tb_coupling, set_cheby_nosoc
+    write_pkkr_input, write_pkkr_operators, use_ldau, set_cheby_nospeedup, decouple_spin_cheby, write_tb_coupling, set_cheby_nosoc
     use mod_constants, only: cvlight, ryd
     use mod_wunfiles, only: t_params
     use memoryhandling, only: allocate_semi_inf_host, allocate_magnetization, allocate_cell, allocate_cpa, allocate_soc, allocate_ldau
@@ -886,9 +886,9 @@ contains
 
     if (use_Chebychev_solver) korbit = 1
 
-    if (decouple_spins_cheby) then
-      write (*, '(A)') 'Warning: detected test option <decouple_spins_cheby>: use spin-decoupled radial equations with new solver'
-      write (1337, *)  'Warning: detected test option <decouple_spins_cheby>: reset KORBIT to zero but use NEWSOSOL for spin-decoupled matrices with explicit spin-loop'
+    if (use_Chebychev_solver .and. decouple_spin_cheby) then
+      write (*, '(A)') 'Warning: detected test option <decouple_spin_cheby>: use spin-decoupled radial equations with new solver'
+      write (1337, *)  'Warning: detected test option <decouple_spin_cheby>: reset KORBIT to zero but use NEWSOSOL for spin-decoupled matrices with explicit spin-loop'
       korbit = 0
     end if
 
@@ -1167,7 +1167,7 @@ contains
     ! End of allocation of SOC arrays
     !--------------------------------------------------------------------------------
     if (use_Chebychev_solver) then      ! Spin-orbit
-      if (use_Chebychev_solver .and. (nspin/=2) .and. .not.decouple_spins_cheby) stop ' set NSPIN = 2 for SOC solver in inputcard'
+      if (use_Chebychev_solver .and. (nspin/=2) .and. .not.decouple_spin_cheby) stop ' set NSPIN = 2 for SOC solver in inputcard'
       npan_log = 30
       npan_eq = 30
       ncheb = 10
@@ -1195,7 +1195,7 @@ contains
     end if
 
     call ioinput('<SOCSCL>        ', uio, 1, 7, ier)
-    if (ier==0 .and. .not.(set_cheby_nosoc .or. decouple_spins_cheby)) then
+    if (ier==0 .and. .not.(set_cheby_nosoc .or. decouple_spin_cheby)) then
       write (111, '(A10)') '<SOCSCL>  '
       do i = 1, natyp
         call ioinput('<SOCSCL>        ', uio, i, 7, ier)
@@ -1209,8 +1209,8 @@ contains
       ! !Bernd - old way
       ! write(111,FMT='(A10,50E10.2)') '<SOCSCL>= ',(SOCSCALE(I1),I1=1,NATYP)
       ! !Bernd - old way
-    elseif (set_cheby_nosoc .or. decouple_spins_cheby) then
-      write(*,*) 'Skipped reading <SOCSCL> because <set_cheby_nosoc>= T or <decouple_spins_cheby>= T. Automatically use <SOCSCL>=0.'
+    elseif (set_cheby_nosoc .or. decouple_spin_cheby) then
+      write(*,*) 'Skipped reading <SOCSCL> because <set_cheby_nosoc>= T or <decouple_spin_cheby>= T. Automatically use <SOCSCL>=0.'
       socscale(:) = 0.0_dp
       write (111, fmt='(A18,50E10.2)') '<SOCSCL>= ', (socscale(i1), i1=1, natyp)
     else
diff --git a/source/KKRhost/testdim.f90 b/source/KKRhost/testdim.f90
index 0d6082c6b72b94f50dde60a57bea80aa3309677a..6f737affb5e98c7311fe626e38d17f58924d959d 100644
--- a/source/KKRhost/testdim.f90
+++ b/source/KKRhost/testdim.f90
@@ -31,7 +31,7 @@ contains
   subroutine testdim(nspin,naez,nemb,natyp,ins,insref,nref,irns,nlayer,krel,nspind, &
     nprincd,knosph,irnsd,korbit,invmod)
 
-    use :: mod_runoptions, only: calc_complex_bandstructure, use_Chebychev_solver, use_cont, use_virtual_atoms, decouple_spins_cheby
+    use :: mod_runoptions, only: calc_complex_bandstructure, use_Chebychev_solver, use_cont, use_virtual_atoms, decouple_spin_cheby
     implicit none
 
     integer, intent (in) :: ins    !! 0 (MT), 1(ASA), 2(Full Potential)
@@ -103,7 +103,7 @@ contains
     end if
 
     if (use_Chebychev_solver .and. korbit==0) then
-      if (decouple_spins_cheby) then
+      if (decouple_spin_cheby) then
         write (6, *) 'Using NEWSOSOL for decoupled spin channels.'
       else
         write (6, *) 'Using option NEWSOSOL, change KORBIT in the inputcard from', korbit, 'to 1'
diff --git a/source/KKRhost/tmat_newsolver.F90 b/source/KKRhost/tmat_newsolver.F90
index 3b266c8b95feb74383458ab65287a7914a1efdf4..5927d95230fa7d49376ea75e8b5a1a9aacce7b2b 100644
--- a/source/KKRhost/tmat_newsolver.F90
+++ b/source/KKRhost/tmat_newsolver.F90
@@ -45,7 +45,7 @@ contains
 #endif
     use :: mod_datatypes, only: dp
     use :: mod_runoptions, only: calc_exchange_couplings, disable_tmat_sratrick, formatted_file, stop_1b, &
-      write_BdG_tests, write_pkkr_operators, write_rhoq_input, set_cheby_nospeedup, decouple_spins_cheby, calc_wronskian
+      write_BdG_tests, write_pkkr_operators, write_rhoq_input, set_cheby_nospeedup, decouple_spin_cheby, calc_wronskian
     use :: mod_constants, only: czero, cone, cvlight
     use :: global_variables, only: ntotd, ncleb, nrmaxd, mmaxd, nspind, nspotd, iemxd, lmmaxd, korbit
     use :: mod_wunfiles, only: t_params
@@ -190,7 +190,7 @@ contains
       vnspll0, vnspll1, hlk, jlk, hlk2, jlk2, tmatsph, ull, rll, sll, rllleft, sllleft)
 
     vins(1:irmdnew, 1:lmpot, 1) = vinsnew(1:irmdnew, 1:lmpot, ipot)
-    if (.not.decouple_spins_cheby) vins(1:irmdnew, 1:lmpot, nspin) = vinsnew(1:irmdnew, 1:lmpot, ipot+nspin-1)
+    if (.not.decouple_spin_cheby) vins(1:irmdnew, 1:lmpot, nspin) = vinsnew(1:irmdnew, 1:lmpot, ipot+nspin-1)
 
     KBdG = 0
 #ifdef CPP_BdG
@@ -345,7 +345,7 @@ contains
         !$omp end critical
 #endif
 
-        if ( .not. decouple_spins_cheby) then
+        if ( .not. decouple_spin_cheby) then
           ! Contruct the spin-orbit coupling hamiltonian and add to potential
           call spinorbit_ham(lmax, lmmax0d, vins, rnew, eryd, zat, cvlight, socscale, nspin, lmpot, theta, phi, ipan_intervall, rpan_intervall, npan_tot, ncheb, irmdnew, nrmaxd, &
             vnspll0(:,:,:), vnspll1(:,:,:,ith), '1')
@@ -411,7 +411,7 @@ contains
         hlk2(:, :, ith) = czero
         jlk2(:, :, ith) = czero
         gmatprefactor = czero
-        if (decouple_spins_cheby) then
+        if (decouple_spin_cheby) then
           use_fullgmat = 0
         else
           use_fullgmat = 1
diff --git a/source/KKRhost/wunfiles.F90 b/source/KKRhost/wunfiles.F90
index 13e58ac6b1e50cd99560d29ad835cb56641e9433..81d7179242e3c4953b8bfc98132a7309c4a593f5 100644
--- a/source/KKRhost/wunfiles.F90
+++ b/source/KKRhost/wunfiles.F90
@@ -366,7 +366,7 @@ contains
     use :: mod_runoptions, only: impurity_operator_only, relax_SpinAngle_Dirac, use_Chebychev_solver, use_qdos, &
       write_cpa_projection_file, write_deci_tmat, write_gmat_file, write_green_host, write_green_imp, write_gref_file, &
       write_kkrimp_input, write_lloyd_cdos_file, write_lloyd_dgref_file, write_lloyd_dtmat_file, write_lloyd_file, &
-      write_lloyd_g0tr_file, write_lloyd_tralpha_file, write_pkkr_input, write_pkkr_operators, write_tmat_file, decouple_spins_cheby
+      write_lloyd_g0tr_file, write_lloyd_tralpha_file, write_pkkr_input, write_pkkr_operators, write_tmat_file, decouple_spin_cheby
 
     implicit none
     ! ..
@@ -691,7 +691,7 @@ contains
     t_inc%naclsmax = naclsmax
     t_inc%nshell0 = nshell(0)
     if (use_Chebychev_solver) t_inc%newsosol = .true.
-    if (decouple_spins_cheby) t_inc%nosoc = .true.
+    if (decouple_spin_cheby) t_inc%nosoc = .true.
     if (write_deci_tmat) t_inc%deci_out = .true.
     !--------------------------------------------------------------------------------
     ! t_inc t_inc t_inc t_inc t_inc t_inc t_inc t_inc t_inc t_inc
diff --git a/source/common/radial_solver_Chebychev/sll_global_solutions.F90 b/source/common/radial_solver_Chebychev/sll_global_solutions.F90
index 1831d121d6067b7295aae3eeda594916d6a039b7..53f9e6bc8692868d489351dbee09f74aacd87e40 100644
--- a/source/common/radial_solver_Chebychev/sll_global_solutions.F90
+++ b/source/common/radial_solver_Chebychev/sll_global_solutions.F90
@@ -306,6 +306,7 @@ contains
     end do
     end if
  
+    ! test writeout
     !write(6,*) 'dllpmax',dllpmax,'iter_beta',iter_beta
 
     end do ! niter_beta
diff --git a/source/common/runoptions.F90 b/source/common/runoptions.F90
index 536ad1e866a60a6816978fab17999b84cd2010c1..c254f4d60ad0dc7dd735255fee5e68dd23bd6287 100644
--- a/source/common/runoptions.F90
+++ b/source/common/runoptions.F90
@@ -16,103 +16,109 @@ module mod_runoptions
   public
   save
 
-  logical :: calc_DOS_Efermi = .false.                 !!calculate DOS at Fermi energy only (former: 'DOS-EF')
-  logical :: calc_GF_Efermi = .false.                  !!calculation of cluster Green function at E Fermi (former: 'GF-EF')
-  logical :: set_cheby_nospeedup = .false.             !!always calculate irregular solution in Chebychev solver (even if not needed) (former: 'norllsll')
-  logical :: decouple_spins_cheby = .false.            !!decouple spin matrices in Chebychev solver neglecting SOC and for collinear calculations only
-  logical :: set_cheby_nosoc     = .false.             !!set SOC strength to 0 for all atoms (former: 'NOSOC')
-  logical :: use_cheby_quadprec  = .false.             !!use quadruple precision in Chebychev solver for irregular solution (former: 'quadprec')
-  logical :: calc_complex_bandstructure = .false.      !!complex band structure (former: 'COMPLEX')
-  logical :: calc_exchange_couplings = .false.         !!calculate magnetic exchange coupling parameters (former: 'XCPL')
-  logical :: calc_exchange_couplings_energy = .false.  !!write energy-resolved Jij-files also if npol/=0 (former: 'Jijenerg')
-  logical :: calc_gmat_lm_full = .false.               !!calculate all lm-lm components of systems greens function and store to file `gflle` (former: 'lmlm-dos')
-  logical :: dirac_scale_SpeefOfLight = .false.        !!scale the speed of light for Dirac solver (former: 'CSCALE')
-  logical :: disable_charge_neutrality = .false.       !!no charge neutrailty required: leaving Fermi level unaffected (former: 'no-neutr')
-  logical :: disable_print_serialnumber = .false.      !!deactivate writing of serial number and version information to files (for backwards compatibility) (former: 'noserial')
-  logical :: disable_reference_system = .false.        !!deactivate the tight-binding reference system (former: 'lrefsysf')
-  logical :: disable_tmat_sratrick = .false.           !!deactivate SRATRICK in solver for t-matirx (former: 'nosph')
-  logical :: fix_nonco_angles = .false.                !!fix direction of non-collinear magnetic moments (Chebychev solver) (former: 'FIXMOM')
-  logical :: formatted_file = .false.                  !!write files ascii-format. only effective with some other write-options (former: 'fileverb')
-  logical :: impurity_operator_only = .false.          !!only for `write_pkkr_operators`: disable costly recalculation of host operators (former: 'IMP_ONLY')
-  logical :: modify_soc_Dirac = .false.                !!modify SOC for Dirac solver (former: 'SOC')
-  logical :: no_madelung = .false.                     !!do not add some energy terms (coulomb, XC, eff. pot.) to total energy (former: 'NoMadel')
-  logical :: print_Gij = .false.                       !!print cluster G_ij matrices to outfile (former: 'Gmatij')
-  logical :: print_gmat = .false.                      !!print Gmat to outfile (former: 'Gmat')
-  logical :: print_ickeck = .false.                    !!enable test-output of ICHECK matrix from gfmask (former: 'ICHECK')
-  logical :: print_kmesh = .false.                     !!output of k-mesh (former: 'k-net')
-  logical :: print_kpoints = .false.                   !!print k-points to outfile (former: 'BZKP')
-  logical :: print_program_flow = .false.              !!monitor the program flow in some parts of the code (former: 'flow')
-  logical :: print_radial_mesh = .false.               !!write mesh information to output (former: 'RMESH')
-  logical :: print_refpot = .false.                    !!test output of refpot (former: 'REFPOT')
-  logical :: print_tau_structure = .false.             !!write extensive information about k-mesh symmetrization and structure of site-diagonal tau matrices to output (former: 'TAUSTRUC')
-  logical :: print_tmat = .false.                      !!print t-matrix to outfile (former: 'tmat')
-  logical :: relax_SpinAngle_Dirac = .false.           !!relax the spin angle in a SCF calculation [only DIRAC mode] (former: 'ITERMDIR')
-  logical :: search_Efermi = .false.                   !!modify convergence parameters to scan for fermi energy only (to reach charge neutrality). (former: 'SEARCHEF')
-  logical :: set_gmat_to_zero = .false.                !!set GMAT=0 in evaluation of density (former: 'GMAT=0')
-  logical :: set_empty_system = .false.                !!set potential and nuclear charge to zero (former: 'zeropot')
-  logical :: set_kmesh_large = .false.                 !!set equal k-mesh (largest) for all energy points (former: 'fix mesh')
-  logical :: set_kmesh_small = .false.                 !!set equal k-mesh (smallest) for all energy points (former: 'fix4mesh')
-  logical :: set_tmat_noinversion = .false.            !!do not perform inversion to get msst = Delta t^-1, but msst = Delta t. (former: 'testgmat')
-  logical :: simulate_asa = .false.                    !!set non-spherical potential to zero in full-potential calculation with Chebychev solver (former: 'simulasa')
-  logical :: slow_mixing_Efermi = .false.              !!renormalize Fermi-energy shift by mixing factor during mixing (former: 'slow-neu')
-  logical :: stop_1a = .false.                         !!stop after main1a (former: 'STOP1A')
-  logical :: stop_1b = .false.                         !!stop after main1b (former: 'STOP1B')
-  logical :: stop_1c = .false.                         !!stop after main1c (former: 'STOP1C')
-  logical :: symmetrize_gmat = .false.                 !!use symmetrization [G(k) + G(-k)]/2 in k-point loop (former: 'symG(k)')
-  logical :: symmetrize_potential_cubic= .false.       !!keep only symmetric part of potential (L=1,11,21,25,43,47). (former: 'potcubic')
-  logical :: symmetrize_potential_madelung = .false.   !!symmetrize potential in consistency to madelung potential (former: 'potsymm')
-  logical :: torque_operator_onlyMT = .false.          !!for torque operator: include only the part within the muffin tin (former: 'ONLYMT')
-  logical :: torque_operator_onlySph = .false.         !!for torque operator: include only the spherically symmetric part (former: 'ONLYSPH')
-  logical :: use_BdG = .false.                         !!use Bogoliubov-de-Gennes Formalism (former: 'useBdG')
-  logical :: use_Chebychev_solver = .false.            !!use the Chebychev solver (former: 'NEWSOSOL')
-  logical :: use_cond_LB = .false.                     !!perform calculation of conductance in Landauer-Büttiker formalism (former: 'CONDUCT')
-  logical :: use_cont = .false.                        !!no usage of embedding points. NEMB is set to 0. (former: 'CONT')
-  logical :: use_deci_onebulk= .false.                 !!in case of decimation: use same bulk on right and left. Speeds up calculations. (former: 'ONEBULK')
-  logical :: use_decimation = .false.                  !!use Decimation technique for semi-infinite systems (former: 'DECIMATE')
-  logical :: use_ewald_2d = .false.                    !!use 2D ewald sum instead of 3D sum (Attention: does not work always!) (former: 'ewald2d')
-  logical :: use_full_BZ = .false.                     !!use full Brillouin zone, i.e. switch off symmetries for k-space integration (former: 'fullBZ')
-  logical :: use_ldau = .false.                        !!use LDA+U as exchange-correlation potential (former: 'LDA+U')
-  logical :: use_lloyd = .false.                       !!use Lloyds formula to correct finite angular momentum cutoff (former: 'LLOYD')
-  logical :: use_qdos = .false.                        !!writes out qdos files for band structure calculations. (former: 'qdos')
-  logical :: use_readcpa = .false.                     !!read cpa t-matrix from file (former: 'readcpa')
-  logical :: use_rigid_Efermi = .false.                !!keep the Fermi energy fixed during self-consistency (former: 'rigid-ef')
-  logical :: use_semicore = .false.                    !!use semicore contour (former: 'SEMICORE')
-  logical :: use_spherical_potential_only = .false.    !!keeping only spherical component of potential (former: 'Vspher')
-  logical :: use_virtual_atoms = .false.               !!add virtual atoms (former: 'VIRATOMS')
-  logical :: write_BdG_tests = .false.                 !!test options for Bogouliubov-deGennes (former: 'BdG_dev')
-  logical :: write_DOS = .false.                       !!write out DOS files in any case (also if npol!=0) (former: 'DOS')
-  logical :: write_DOS_lm = .false.                    !!write out DOS files with decomposition into l and m components (former: 'lmdos')
-  logical :: write_gmat_plain = .false.                !!write out Green function as plain text file (former: 'GPLAIN')
-  logical :: write_green_host = .false.                !!write green function of the host to file `green_host` (former: 'WRTGREEN')
-  logical :: write_green_imp = .false.                 !!write out impurity Green function to GMATLL_GES (former: 'GREENIMP')
-  logical :: write_complex_qdos = .false.              !!write complex qdos to file (former: 'compqdos')
-  logical :: write_cpa_projection_file = .false.       !!write CPA projectors to file (former: 'projfile')
-  logical :: write_deci_pot = .false.                  !!write decimation-potential file (former: 'deci-pot')
-  logical :: write_deci_tmat = .false.                 !!write t-matrix to file 'decifile' (former: 'deci-out')
-  logical :: write_density_ascii = .false.             !!write density rho2ns to file densitydn.ascii (former: 'den-asci')
-  logical :: write_energy_mesh = .false.               !!write out the energy mesh to file `emesh.scf` (former: 'EMESH')
-  logical :: write_generalized_potential = .false.     !!write potential in general format. Usually prepares for running the VORONOI program. (former: 'GENPOT')
-  logical :: write_gmat_file = .false.                 !!write GMAT to file (former: 'gmatfile')
-  logical :: write_gref_file = .false.                 !!write GREF to file (former: 'greffile')
-  logical :: write_gmat_ascii = .false.                !!write GMAT to formatted file `gmat.ascii` (former: 'gmatasci')
-  logical :: write_kkrimp_input = .false.              !!write out files for KKRimp-code (former: 'KKRFLEX')
-  logical :: write_kkrsusc_input = .false.             !!write out files for KKRsusc-code (former: 'KKRSUSC')
-  logical :: write_kpts_file = .false.                 !!write and read k-mesh to/from file `kpoints` (former: 'kptsfile')
-  logical :: write_lloyd_cdos_file = .false.           !!write Lloyd array to file  (former: 'wrtcdos')
-  logical :: write_lloyd_dgref_file= .false.           !!write Lloyd array to file  (former: 'wrtdgref')
-  logical :: write_lloyd_dtmat_file= .false.           !!write Lloyd array to file  (former: 'wrtdtmat')
-  logical :: write_lloyd_file = .false.                !!write several Lloyd-arrays to files (former: 'llyfiles')
-  logical :: write_lloyd_g0tr_file = .false.           !!write Lloyd array to file  (former: 'wrtgotr')
-  logical :: write_lloyd_tralpha_file= .false.         !!write Lloyd array to file  (former: 'wrttral')
-  logical :: write_madelung_file = .false.             !!write madelung summation to file 'abvmad.unformatted' instead of keeping it in memory (former: 'madelfil')
-  logical :: write_pkkr_input = .false.                !!write out files for Pkkprime-code (former: 'FERMIOUT')
-  logical :: write_pkkr_operators = .false.            !!for Fermi-surface output: calculate various operators in KKR basis. (former: 'OPERATOR')
-  logical :: write_potential_tests = .false.           !!write potential at different steps in main2 to different files (former: 'vintrasp' and 'vpotout')
-  logical :: write_rho2ns = .false.                    !!write array rho2ns into file out_rhoval (from main1c) and out_rhotot (from main2) (former: 'RHOVALTW' and 'RHOVALW')
-  logical :: write_rhoq_input = .false.                !!write out files needed for rhoq module (Quasiparticle interference) (former: 'rhoqtest')
-  logical :: write_tmat_file = .false.                 !!write t-matix to file (former: 'tmatfile')
-  logical :: write_tb_coupling = .false.               !!write couplings in tight-binging reference system to file `couplings.dat` (former: 'godfrin')
-  logical :: calc_wronskian = .false.                  !!calculate the wronskian relations of first and second kind for the wavefunctions (see PhD Bauer pp 48)
+  logical :: calc_DOS_Efermi = .false.                 !! calculate DOS at Fermi energy only (former: 'DOS-EF')
+  logical :: calc_GF_Efermi = .false.                  !! calculation of cluster Green function at E Fermi (former: 'GF-EF')
+  logical :: set_cheby_nosoc = .false.                 !! set SOC strength to 0 for all atoms (former: 'NOSOC')
+  logical :: set_cheby_nospeedup = .false.             !! always calculate irregular solution in Chebychev solver (even if not needed) (former: 'norllsll')
+  logical :: decouple_spin_cheby = .false.             !! decouple spin matrices in Chebychev solver neglecting SOC and for collinear calculations only
+  logical :: use_cheby_quadprec  = .false.             !! use quadruple precision in Chebychev solver for irregular solution (former: 'quadprec')
+  logical :: calc_complex_bandstructure = .false.      !! complex band structure (former: 'COMPLEX')
+  logical :: calc_exchange_couplings = .false.         !! calculate magnetic exchange coupling parameters (former: 'XCPL')
+  logical :: calc_exchange_couplings_energy = .false.  !! write energy-resolved Jij-files also if npol/=0 (former: 'Jijenerg')
+  logical :: calc_gmat_lm_full = .false.               !! calculate all lm-lm components of systems greens function and store to file `gflle` (former: 'lmlm-dos')
+  logical :: dirac_scale_SpeefOfLight = .false.        !! scale the speed of light for Dirac solver (former: 'CSCALE')
+  logical :: disable_charge_neutrality = .false.       !! no charge neutrailty required: leaving Fermi level unaffected (former: 'no-neutr')
+  logical :: disable_print_serialnumber = .false.      !! deactivate writing of serial number and version information to files (for backwards compatibility) (former: 'noserial')
+  logical :: disable_reference_system = .false.        !! deactivate the tight-binding reference system (former: 'lrefsysf')
+  logical :: disable_tmat_sratrick = .false.           !! deactivate SRATRICK in solver for t-matirx (former: 'nosph')
+  logical :: fix_nonco_angles = .false.                !! fix direction of non-collinear magnetic moments (Chebychev solver) (former: 'FIXMOM')
+  logical :: formatted_file = .false.                  !! write files ascii-format. only effective with some other write-options (former: 'fileverb')
+  logical :: write_npy = .false.                       !! write files in npy format (python readable) (former: 'writenpy')
+  logical :: impurity_operator_only = .false.          !! only for `write_pkkr_operators`: disable costly recalculation of host operators (former: 'IMP_ONLY')
+  logical :: modify_soc_Dirac = .false.                !! modify SOC for Dirac solver (former: 'SOC')
+  logical :: no_madelung = .false.                     !! do not add some energy terms (coulomb, XC, eff. pot.) to total energy (former: 'NoMadel')
+  logical :: print_Gij = .false.                       !! print cluster G_ij matrices to outfile (former: 'Gmatij')
+  logical :: print_gmat = .false.                      !! print Gmat to outfile (former: 'Gmat')
+  logical :: print_ickeck = .false.                    !! enable test-output of ICHECK matrix from gfmask (former: 'ICHECK')
+  logical :: print_kmesh = .false.                     !! output of k-mesh (former: 'k-net')
+  logical :: print_kpoints = .false.                   !! print k-points to outfile (former: 'BZKP')
+  logical :: print_program_flow = .false.              !! monitor the program flow in some parts of the code (former: 'flow')
+  logical :: print_radial_mesh = .false.               !! write mesh information to output (former: 'RMESH')
+  logical :: print_refpot = .false.                    !! test output of refpot (former: 'REFPOT')
+  logical :: print_tau_structure = .false.             !! write extensive information about k-mesh symmetrization and structure of site-diagonal tau matrices to output (former: 'TAUSTRUC')
+  logical :: print_tmat = .false.                      !! print t-matrix to outfile (former: 'tmat')
+  logical :: relax_SpinAngle_Dirac = .false.           !! relax the spin angle in a SCF calculation [only DIRAC mode] (former: 'ITERMDIR')
+  logical :: search_Efermi = .false.                   !! modify convergence parameters to scan for fermi energy only (to reach charge neutrality). (former: 'SEARCHEF')
+  logical :: set_gmat_to_zero = .false.                !! set GMAT=0 in evaluation of density (former: 'GMAT=0')
+  logical :: set_empty_system = .false.                !! set potential and nuclear charge to zero (former: 'zeropot')
+  logical :: set_kmesh_large = .false.                 !! set equal k-mesh (largest) for all energy points (former: 'fix mesh')
+  logical :: set_kmesh_small = .false.                 !! set equal k-mesh (smallest) for all energy points (former: 'fix4mesh')
+  logical :: set_tmat_noinversion = .false.            !! do not perform inversion to get msst = Delta t^-1, but msst = Delta t. (former: 'testgmat')
+  logical :: simulate_asa = .false.                    !! set non-spherical potential to zero in full-potential calculation with Chebychev solver (former: 'simulasa')
+  logical :: slow_mixing_Efermi = .false.              !! renormalize Fermi-energy shift by mixing factor during mixing (former: 'slow-neu')
+  logical :: stop_1a = .false.                         !! stop after main1a (former: 'STOP1A')
+  logical :: stop_1b = .false.                         !! stop after main1b (former: 'STOP1B')
+  logical :: stop_1c = .false.                         !! stop after main1c (former: 'STOP1C')
+  logical :: symmetrize_gmat = .false.                 !! use symmetrization [G(k) + G(-k)]/2 in k-point loop (former: 'symG(k)')
+  logical :: symmetrize_potential_cubic= .false.       !! keep only symmetric part of potential (L=1,11,21,25,43,47). (former: 'potcubic')
+  logical :: symmetrize_potential_madelung = .false.   !! symmetrize potential in consistency to madelung potential (former: 'potsymm')
+  logical :: torque_operator_onlyMT = .false.          !! for torque operator: include only the part within the muffin tin (former: 'ONLYMT')
+  logical :: torque_operator_onlySph = .false.         !! for torque operator: include only the spherically symmetric part (former: 'ONLYSPH')
+  logical :: use_BdG = .false.                         !! use Bogoliubov-de-Gennes Formalism (former: 'useBdG')
+  logical :: use_Chebychev_solver = .false.            !! use the Chebychev solver (former: 'NEWSOSOL')
+  logical :: use_cond_LB = .false.                     !! perform calculation of conductance in Landauer-Büttiker formalism (former: 'CONDUCT')
+  logical :: use_cont = .false.                        !! no usage of embedding points. NEMB is set to 0. (former: 'CONT')
+  logical :: use_deci_onebulk= .false.                 !! in case of decimation: use same bulk on right and left. Speeds up calculations. (former: 'ONEBULK')
+  logical :: use_decimation = .false.                  !! use Decimation technique for semi-infinite systems (former: 'DECIMATE')
+  logical :: use_ewald_2d = .false.                    !! use 2D ewald sum instead of 3D sum (Attention: does not work always!) (former: 'ewald2d')
+  logical :: force_BZ_symm = .false.                   !! force the use of BZ symmetries for k-integration even with SOC (only overwritten by use_full_BZ option)
+  logical :: use_full_BZ = .false.                     !! use full Brillouin zone, i.e. switch off symmetries for k-space integration (former: 'fullBZ')
+  logical :: use_ldau = .false.                        !! use LDA+U as exchange-correlation potential (former: 'LDA+U')
+  logical :: use_lloyd = .false.                       !! use Lloyds formula to correct finite angular momentum cutoff (former: 'LLOYD')
+  logical :: use_qdos = .false.                        !! writes out qdos files for band structure calculations. (former: 'qdos')
+  logical :: use_readcpa = .false.                     !! read cpa t-matrix from file (former: 'readcpa')
+  logical :: use_rigid_Efermi = .false.                !! keep the Fermi energy fixed during self-consistency (former: 'rigid-ef')
+  logical :: use_semicore = .false.                    !! use semicore contour (former: 'SEMICORE')
+  logical :: use_spherical_potential_only = .false.    !! keeping only spherical component of potential (former: 'Vspher')
+  logical :: use_virtual_atoms = .false.               !! add virtual atoms (former: 'VIRATOMS')
+  logical :: write_BdG_tests = .false.                 !! test options for Bogouliubov-deGennes (former: 'BdG_dev')
+  logical :: write_DOS = .false.                       !! write out DOS files in any case (also if npol!=0) (former: 'DOS')
+  logical :: write_DOS_lm = .false.                    !! write out DOS files with decomposition into l and m components (former: 'lmdos')
+  logical :: write_gmat_plain = .false.                !! write out Green function as plain text file (former: 'GPLAIN')
+  logical :: write_green_host = .false.                !! write green function of the host to file `green_host` (former: 'WRTGREEN')
+  logical :: write_green_imp = .false.                 !! write out impurity Green function to GMATLL_GES (former: 'GREENIMP')
+  logical :: write_complex_qdos = .false.              !! write complex qdos to file (former: 'compqdos')
+  logical :: write_complex_lmdos = .false.             !! write complex lmdos to file
+  logical :: write_cpa_projection_file = .false.       !! write CPA projectors to file (former: 'projfile')
+  logical :: write_deci_pot = .false.                  !! write decimation-potential file (former: 'deci-pot')
+  logical :: write_deci_tmat = .false.                 !! write t-matrix to file 'decifile' (former: 'deci-out')
+  logical :: write_density_ascii = .false.             !! write density rho2ns to file densitydn.ascii (former: 'den-asci')
+  logical :: write_energy_mesh = .false.               !! write out the energy mesh to file `emesh.scf` (former: 'EMESH')
+  logical :: write_generalized_potential = .false.     !! write potential in general format. Usually prepares for running the VORONOI program. (former: 'GENPOT')
+  logical :: write_gmat_file = .false.                 !! write GMAT to file (former: 'gmatfile')
+  logical :: write_gref_file = .false.                 !! write GREF to file (former: 'greffile')
+  logical :: write_gmat_ascii = .false.                !! write GMAT to formatted file `gmat.ascii` (former: 'gmatasci')
+  logical :: write_kkrimp_input = .false.              !! write out files for KKRimp-code (former: 'KKRFLEX')
+  logical :: write_kkrsusc_input = .false.             !! write out files for KKRsusc-code (former: 'KKRSUSC')
+  logical :: write_kpts_file = .false.                 !! write and read k-mesh to/from file `kpoints` (former: 'kptsfile')
+  logical :: write_lloyd_cdos_file = .false.           !! write Lloyd array to file  (former: 'wrtcdos')
+  logical :: write_lloyd_dgref_file= .false.           !! write Lloyd array to file  (former: 'wrtdgref')
+  logical :: write_lloyd_dtmat_file= .false.           !! write Lloyd array to file  (former: 'wrtdtmat')
+  logical :: write_lloyd_file = .false.                !! write several Lloyd-arrays to files (former: 'llyfiles')
+  logical :: write_lloyd_g0tr_file = .false.           !! write Lloyd array to file  (former: 'wrtgotr')
+  logical :: write_lloyd_tralpha_file= .false.         !! write Lloyd array to file  (former: 'wrttral')
+  logical :: write_madelung_file = .false.             !! write madelung summation to file 'abvmad.unformatted' instead of keeping it in memory (former: 'madelfil')
+  logical :: write_pkkr_input = .false.                !! write out files for Pkkprime-code (former: 'FERMIOUT')
+  logical :: write_pkkr_operators = .false.            !! for Fermi-surface output: calculate various operators in KKR basis. (former: 'OPERATOR')
+  logical :: write_potential_tests = .false.           !! write potential at different steps in main2 to different files (former: 'vintrasp' and 'vpotout')
+  logical :: write_rho2ns = .false.                    !! write array rho2ns into file out_rhoval (from main1c) and out_rhotot (from main2) (former: 'RHOVALTW' and 'RHOVALW')
+  logical :: write_rhoq_input = .false.                !! write out files needed for rhoq module (Quasiparticle interference) (former: 'rhoqtest')
+  logical :: write_tmat_file = .false.                 !! write t-matix to file (former: 'tmatfile')
+  logical :: write_tb_coupling = .false.               !! write couplings in tight-binging reference system to file `couplings.dat` (former: 'godfrin')
+  logical :: calc_wronskian = .false.                  !! calculate the wronskian relations of first and second kind for the wavefunctions (see PhD Bauer pp 48)
+  logical :: calc_onsite_only = .false.                !! calculate not the full Green function for the density but take the onsite part alone
+  logical :: use_gmat_unity = .false.                  !! set the structural GF to the unity matrix for test purposes
+  logical :: nosmallcomp = .false.                     !! set small component of the wavefunction to zero
   logical :: use_broyden_spinmix = .false.             !! use broyden spin mixing for noncollinear angles
   logical :: write_angles_alliter= .false.             !! write out noncollinear angles for all iterations
 
@@ -205,7 +211,7 @@ module mod_runoptions
     call set_runoption(symmetrize_potential_madelung , '<symmetrize_potential_madelung>' , '<potsymm>' )
     call set_runoption(use_cheby_quadprec            , '<use_cheby_quadprec>'            , '<quadprec>')
     call set_runoption(set_cheby_nosoc               , '<set_cheby_nosoc>'               , '<NOSOC>'   )
-    call set_runoption(decouple_spins_cheby          , '<decouple_spins_cheby>')
+    call set_runoption(decouple_spin_cheby          , '<decouple_spin_cheby>')
     call set_runoption(set_empty_system              , '<set_empty_system>'              , '<zeropot>' )
     call set_runoption(write_tmat_file               , '<write_tmat_file>'               , '<tmatfile>')
     call set_runoption(write_tb_coupling             , '<write_tb_coupling>'             , '<godfrin>' )
@@ -748,7 +754,7 @@ module mod_runoptions
     call mpi_bcast(symmetrize_potential_madelung , 1, mpi_logical, master, mpi_comm_world, ierr)
     call mpi_bcast(use_cheby_quadprec            , 1, mpi_logical, master, mpi_comm_world, ierr)
     call mpi_bcast(set_cheby_nosoc               , 1, mpi_logical, master, mpi_comm_world, ierr)
-    call mpi_bcast(decouple_spins_cheby          , 1, mpi_logical, master, mpi_comm_world, ierr)
+    call mpi_bcast(decouple_spin_cheby          , 1, mpi_logical, master, mpi_comm_world, ierr)
     call mpi_bcast(set_empty_system              , 1, mpi_logical, master, mpi_comm_world, ierr)
     call mpi_bcast(write_tmat_file               , 1, mpi_logical, master, mpi_comm_world, ierr)
     call mpi_bcast(write_tb_coupling             , 1, mpi_logical, master, mpi_comm_world, ierr)
@@ -801,7 +807,7 @@ module mod_runoptions
     write(iounit, '(A35,1x,1L,3x,A)') '<calc_GF_Efermi>=', calc_GF_Efermi, "calculation of cluster Green function at E Fermi (former: 'GF-EF')"
     write(iounit, '(A35,1x,1L,3x,A)') '<set_cheby_nospeedup>=', set_cheby_nospeedup, "always calculate irregular solution in Chebychev solver (even if not needed) (former: 'norllsll')"
     write(iounit, '(A35,1x,1L,3x,A)') '<set_cheby_nosoc>=', set_cheby_nosoc, "set SOC strength to 0 for all atoms (former: 'NOSOC')"
-    write(iounit, '(A35,1x,1L,3x,A)') '<decouple_spins_cheby>=', decouple_spins_cheby, "decouple spin matrices in Chebychev solver neglecting SOC and for collinear calculations only"
+    write(iounit, '(A35,1x,1L,3x,A)') '<decouple_spin_cheby>=', decouple_spin_cheby, "decouple spin matrices in Chebychev solver neglecting SOC and for collinear calculations only"
     write(iounit, '(A35,1x,1L,3x,A)') '<calc_complex_bandstructure>=', calc_complex_bandstructure, "complex band structure (former: 'COMPLEX')"
     write(iounit, '(A35,1x,1L,3x,A)') '<calc_exchange_couplings>=', calc_exchange_couplings, "calculate magnetic exchange coupling parameters (former: 'XCPL')"
     write(iounit, '(A35,1x,1L,3x,A)') '<calc_exchange_couplings_energy>=', calc_exchange_couplings_energy, "write energy-resolved Jij-files also if npol/=0 (former: 'Jijenerg')"
diff --git a/source/common/spinorbit_ham.f90 b/source/common/spinorbit_ham.f90
index 5e0d00a2132c955b82a6c6387c0c47fda3fd4275..4f2829981526c6442fb15c3153497747984173bd 100644
--- a/source/common/spinorbit_ham.f90
+++ b/source/common/spinorbit_ham.f90
@@ -33,7 +33,7 @@ contains
 
     use :: mod_datatypes, only: dp
     use :: mod_constants, only: czero
-    use :: mod_runoptions, only: set_cheby_nosoc, decouple_spins_cheby
+    use :: mod_runoptions, only: set_cheby_nosoc, decouple_spin_cheby
     use :: mod_cheb, only: getclambdacinv
     use :: mod_spin_orbit_compl, only: spin_orbit_compl
     use :: mod_rotatespinframe, only: rotatematrix
@@ -136,7 +136,7 @@ contains
     ! contruct prefactor of spin-orbit hamiltonian
     hsofac = 0e0_dp
     vnspll1 = (0e0_dp, 0e0_dp)
-    if (set_cheby_nosoc .or. decouple_spins_cheby .or. zat<1e-6_dp) then
+    if (set_cheby_nosoc .or. decouple_spin_cheby .or. zat<1e-6_dp) then
       vnspll1(1:2*lmmaxd, 1:2*lmmaxd, 1:irmdnew) = vnspll(1:2*lmmaxd, 1:2*lmmaxd, 1:irmdnew)
     else
       do ir = 1, irmdnew
diff --git a/source/voronoi/maindriver12.f b/source/voronoi/maindriver12.f
index 28afb9e9bdd359d960735fc671812be58b7105bb..afd40002c12e3b05634e8baf9c5cc5e189d27738 100644
--- a/source/voronoi/maindriver12.f
+++ b/source/voronoi/maindriver12.f
@@ -828,7 +828,11 @@ c
   
 c           ! Redefine
             RMTCL(ISHAPE) = XRN_ALL(1   ,ISHAPE)*ALATC
-            RWSCL(ISHAPE) = XRN_ALL(NMESH+NRAD,ISHAPE)*ALATC
+            IF (.NOT.OPT('SIMULASA')) THEN
+              RWSCL(ISHAPE) = XRN_ALL(NMESH+NRAD,ISHAPE)*ALATC
+            ELSE
+              RWSCL(ISHAPE) = XRN_ALL(MESHN_ALL(ISHAPE),ISHAPE)*ALATC
+            END IF
             WRITE(6,*) 'rmt  = ',RMTCL(ISHAPE)
             WRITE(6,*) 'rmax = ',RWSCL(ISHAPE)
             LCONSTRUCTED(ISHAPE) = .TRUE.
diff --git a/source/voronoi/shape12.f90 b/source/voronoi/shape12.f90
index c90116c2577ec2db65e8659dc845d0ae6aa4bbcb..3ba4264ee094c81068c5e983427919339153d83a 100644
--- a/source/voronoi/shape12.f90
+++ b/source/voronoi/shape12.f90
@@ -2088,28 +2088,33 @@
       VRTY=V(2,IVERT)
       VRTZ=V(3,IVERT)
       INEW=1                          ! Save all different vertices
-      DO 13 IVRT=1,NVRT
-      T=(VRTX-VRT(1,IVRT))**2+(VRTY-VRT(2,IVRT))**2&
-     & +(VRTZ-VRT(3,IVRT))**2
-      IF(T.LT.TOLVDIST) INEW=0
-   13 CONTINUE
-      IF(INEW.EQ.1)                  THEN
-      NVRT=NVRT+1
-      IF(NVRT.GT.NVRTD) STOP 'INCREASE NVRTD'
-      VRT(1,NVRT)=V(1,IVERT)
-      VRT(2,NVRT)=V(2,IVERT)
-      VRT(3,NVRT)=V(3,IVERT)
-                                     END   IF
+      DO IVRT=1,NVRT
+        T=(VRTX-VRT(1,IVRT))**2+(VRTY-VRT(2,IVRT))**2&
+     &   +(VRTZ-VRT(3,IVRT))**2
+        IF(T.LT.TOLVDIST) INEW=0
+      END DO
+   
+
+      IF(INEW.EQ.1) THEN
+        NVRT=NVRT+1
+        IF(NVRT.GT.NVRTD) STOP 'INCREASE NVRTD'
+        VRT(1,NVRT)=V(1,IVERT)
+        VRT(2,NVRT)=V(2,IVERT)
+        VRT(3,NVRT)=V(3,IVERT)
+      END IF
+
       IVERTP=IVERT+1                  
       IF(IVERT.EQ.NVERT) IVERTP=1
       VRTPX=V(1,IVERTP)
       VRTPY=V(2,IVERTP)
       VRTPZ=V(3,IVERTP)
+
       IVERTM=IVERT-1
       IF(IVERT.EQ.1) IVERTM=NVERT
       VRTMX=V(1,IVERTM)
       VRTMY=V(2,IVERTM)               ! Check IF the  consecutive
       VRTMZ=V(3,IVERTM)               ! vertices define a polygon
+
       A1=SQRT((VRTPX-VRTX)**2+(VRTPY-VRTY)**2+(VRTPZ-VRTZ)**2)
       A2=SQRT((VRTMX-VRTX)**2+(VRTMY-VRTY)**2+(VRTMZ-VRTZ)**2)
       DOWN=A1*A2
@@ -2120,14 +2125,19 @@
       ! write(6,*) VRTMX,VRTMY,VRTMZ
       ! write(6,*) VRTPX,VRTPY,VRTPZ
       ! write(6,*) 'fisum ',ivert,a1,a2,up,arg,ACOS(ARG)
-      IF(DOWN.GE.TOLVDIST)   THEN
-      ARG=UP/DOWN
-      IF(ABS(ARG).GE.1.D0) ARG=SIGN(1.D0,ARG)
-      FISUM=FISUM-ACOS(ARG)
-      
-                             ELSE
-      STOP 'IDENTICAL CONSECUTIVE VERTICES'
-                             END    IF
+      IF(DOWN.GE.TOLVDIST) THEN
+        ARG=UP/DOWN
+        IF(ABS(ARG).GE.1.D0) ARG=SIGN(1.D0,ARG)
+        FISUM=FISUM-ACOS(ARG)
+      ELSE
+        write(*,*) 'ivert:', IVERT, IVERTP, IVERTM, NVERT
+        write(*,*) 'a1', VRTPX, VRTX, VRTPY, VRTY, VRTPZ, VRTZ
+        write(*,*) 'diffs:', VRTPX-VRTX, VRTPY-VRTY, VRTPZ-VRTZ
+        write(*,*) 'a2', VRTMX, VRTX, VRTMY, VRTY, VRTMZ, VRTZ
+        write(*,*) 'diffs:', VRTMX-VRTX, VRTMY-VRTY, VRTMZ-VRTZ
+        write(*,*) DOWN, TOLVDIST
+        stop 'IDENTICAL CONSECUTIVE VERTICES'
+      END IF
 !
 !------> T R E A T M E N T   O F   E D G E S 
 !