diff --git a/src/qs_cdft_methods.F b/src/qs_cdft_methods.F
index f2fd34d428..8c59c6b077 100644
--- a/src/qs_cdft_methods.F
+++ b/src/qs_cdft_methods.F
@@ -755,8 +755,8 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
          CALL pw_set(cdft_control%group(igroup)%hw_rho_total_constraint, 1.0_dp)
 
          IF (igroup == 1) THEN
-            CALL auxbas_pw_pool%create_pw(cdft_control%group(1)%hw_rho_total)
-            CALL pw_set(cdft_control%group(1)%hw_rho_total, 1.0_dp)
+            CALL auxbas_pw_pool%create_pw(cdft_control%hw_rho_total)
+            CALL pw_set(cdft_control%hw_rho_total, 1.0_dp)
 
             IF (hirshfeld_control%print_density) THEN
                DO iatom = 1, cdft_control%natoms
@@ -878,7 +878,7 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
          DEALLOCATE (pab)
 
          IF (igroup == 1) THEN
-            CALL transfer_rs2pw(rs_rho_all, cdft_control%group(1)%hw_rho_total)
+            CALL transfer_rs2pw(rs_rho_all, cdft_control%hw_rho_total)
          END IF
 
          CALL transfer_rs2pw(rs_rho_constr, cdft_control%group(igroup)%hw_rho_total_constraint)
@@ -887,7 +887,7 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
          ! Calculate weight function
          CALL hfun_scale(cdft_control%group(igroup)%weight%array, &
                          cdft_control%group(igroup)%hw_rho_total_constraint%array, &
-                         cdft_control%group(1)%hw_rho_total%array, divide=.TRUE., &
+                         cdft_control%hw_rho_total%array, divide=.TRUE., &
                          small=hirshfeld_control%eps_cutoff)
 
          ! Calculate charges
@@ -896,7 +896,7 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
                CALL transfer_rs2pw(rs_single_charge(i), cdft_control%group(igroup)%hw_rho_atomic_charge(i))
                CALL hfun_scale(cdft_control%charge(i)%array, &
                                cdft_control%group(igroup)%hw_rho_atomic_charge(i)%array, &
-                               cdft_control%group(1)%hw_rho_total%array, divide=.TRUE., &
+                               cdft_control%hw_rho_total%array, divide=.TRUE., &
                                small=hirshfeld_control%eps_cutoff)
             END DO
          END IF
@@ -916,7 +916,7 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
          CALL auxbas_pw_pool%give_back_pw(cdft_control%group(igroup)%hw_rho_total_constraint)
 
          IF (.NOT. cdft_control%in_memory .AND. igroup == 1) THEN
-            CALL auxbas_pw_pool%give_back_pw(cdft_control%group(1)%hw_rho_total)
+            CALL auxbas_pw_pool%give_back_pw(cdft_control%hw_rho_total)
          END IF
 
          IF (hirshfeld_control%print_density .AND. igroup == 1) THEN
@@ -1094,11 +1094,11 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
 
                      ra(:) = particle_set(iatom)%r
 
-                     IF (cdft_control%group(1)%hw_rho_total%array(i, j, k) > hirshfeld_control%eps_cutoff) THEN
+                     IF (cdft_control%hw_rho_total%array(i, j, k) > hirshfeld_control%eps_cutoff) THEN
 
                         exp_eval = (coefficients(iatom) - &
                                     cdft_control%group(igroup)%weight%array(i, j, k))/ &
-                                   cdft_control%group(1)%hw_rho_total%array(i, j, k)
+                                   cdft_control%hw_rho_total%array(i, j, k)
 
                         r2 = (/i*dr_pw(1), j*dr_pw(2), k*dr_pw(3)/) + origin
                         r_pbc = pbc(ra, r2, cell)
@@ -1124,7 +1124,7 @@ SUBROUTINE hirshfeld_constraint_low(qs_env, just_gradients)
                END DO
             END DO
          END DO
-         CALL auxbas_pw_pool%give_back_pw(cdft_control%group(1)%hw_rho_total)
+         CALL auxbas_pw_pool%give_back_pw(cdft_control%hw_rho_total)
       END IF
 
       CALL rs_grid_release(rs_rho_all)
diff --git a/src/qs_cdft_types.F b/src/qs_cdft_types.F
index a7d01b19f5..a5cf4601a6 100644
--- a/src/qs_cdft_types.F
+++ b/src/qs_cdft_types.F
@@ -214,7 +214,6 @@ MODULE qs_cdft_types
       TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER             :: hw_rho_atomic_dr => NULL()
       TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER             :: hw_rho_atomic_charge => NULL()
       TYPE(pw_r3d_rs_type)             :: hw_rho_total_constraint = pw_r3d_rs_type()
-      TYPE(pw_r3d_rs_type)          :: hw_rho_total = pw_r3d_rs_type()
    END TYPE cdft_group_type
 
    TYPE cdft_control_type
@@ -256,6 +255,7 @@ MODULE qs_cdft_types
       TYPE(pw_r3d_rs_type), POINTER, &
          DIMENSION(:, :)                   :: fragments => NULL()
       TYPE(qs_outer_scf_type)              :: constraint = qs_outer_scf_type()
+      TYPE(pw_r3d_rs_type)                 :: hw_rho_total = pw_r3d_rs_type()
    END TYPE cdft_control_type
 
    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_cdft_types'
diff --git a/src/qs_cdft_utils.F b/src/qs_cdft_utils.F
index 421355405f..6cf9ad4a37 100644
--- a/src/qs_cdft_utils.F
+++ b/src/qs_cdft_utils.F
@@ -1341,7 +1341,7 @@ SUBROUTINE cdft_print_hirshfeld_density(qs_env)
          unit_nr = cp_print_key_unit_nr(logger, cdft_constraint_section, "PROGRAM_RUN_INFO", mpi_io=mpi_io, &
                                         file_position="REWIND", middle_name=middle_name, extension=".cube")
 
-         CALL cp_pw_to_cube(cdft_control%group(igroup)%hw_rho_total, unit_nr, "CDFT Weight Function", mpi_io=mpi_io, &
+         CALL cp_pw_to_cube(cdft_control%hw_rho_total, unit_nr, "CDFT Weight Function", mpi_io=mpi_io, &
                   particles=particles, stride=section_get_ivals(cdft_constraint_section, "PROGRAM_RUN_INFO%WEIGHT_FUNCTION%STRIDE"))
 
          CALL cp_print_key_finished_output(unit_nr, logger, cdft_constraint_section, "PROGRAM_RUN_INFO", mpi_io=mpi_io)