subroutine initialize_fragmented_system(nfrag, geom, fragments, charge, multiplicity, &
allow_overlapping, use_angstrom, sys_geom, error)
!! Shared helper to initialize system_geometry_t for fragmented calculations
!! Handles fragment allocation, size checking, and overlap validation
use mqc_geometry, only: geometry_type
use mqc_config_parser, only: input_fragment_t
integer, intent(in) :: nfrag
type(geometry_type), intent(in) :: geom
type(input_fragment_t), intent(in) :: fragments(:)
integer, intent(in) :: charge, multiplicity
logical, intent(in) :: allow_overlapping
logical, intent(in) :: use_angstrom
type(system_geometry_t), intent(out) :: sys_geom
type(error_t), intent(out) :: error
integer :: i, j, atoms_in_first_frag, max_frag_size
logical :: all_same_size
! Set up basic system geometry
sys_geom%n_monomers = nfrag
sys_geom%total_atoms = geom%natoms
sys_geom%charge = charge
sys_geom%multiplicity = multiplicity
! Allocate fragment info arrays
allocate (sys_geom%fragment_sizes(nfrag))
allocate (sys_geom%fragment_charges(nfrag))
allocate (sys_geom%fragment_multiplicities(nfrag))
! Get fragment sizes
max_frag_size = 0
atoms_in_first_frag = size(fragments(1)%indices)
all_same_size = .true.
do i = 1, nfrag
sys_geom%fragment_sizes(i) = size(fragments(i)%indices)
sys_geom%fragment_charges(i) = fragments(i)%charge
sys_geom%fragment_multiplicities(i) = fragments(i)%multiplicity
max_frag_size = max(max_frag_size, sys_geom%fragment_sizes(i))
if (sys_geom%fragment_sizes(i) /= atoms_in_first_frag) then
all_same_size = .false.
end if
end do
! Allocate fragment_atoms array
allocate (sys_geom%fragment_atoms(max_frag_size, nfrag))
sys_geom%fragment_atoms = -1 ! Initialize with invalid index
! Store fragment atom indices (0-indexed from input file)
do i = 1, nfrag
do j = 1, sys_geom%fragment_sizes(i)
sys_geom%fragment_atoms(j, i) = fragments(i)%indices(j)
end do
end do
! Check for overlapping fragments if not allowed
if (.not. allow_overlapping) then
call check_fragment_overlap(fragments, nfrag, error)
if (error%has_error()) then
call error%add_context("mqc_config_adapter:geometry_to_system_fragmented")
return
end if
end if
! Set atoms_per_monomer: use common size if identical, else 0
if (all_same_size) then
sys_geom%atoms_per_monomer = atoms_in_first_frag
else
sys_geom%atoms_per_monomer = 0 ! Signal variable-sized fragments
end if
allocate (sys_geom%element_numbers(sys_geom%total_atoms))
allocate (sys_geom%coordinates(3, sys_geom%total_atoms))
! Convert element symbols to atomic numbers
do i = 1, sys_geom%total_atoms
sys_geom%element_numbers(i) = element_symbol_to_number(geom%elements(i))
end do
! Store coordinates (convert to Bohr if needed)
if (use_angstrom) then
sys_geom%coordinates = to_bohr(geom%coords)
else
sys_geom%coordinates = geom%coords
end if
end subroutine initialize_fragmented_system