initialize_system_geometry Subroutine

public subroutine initialize_system_geometry(full_geom_file, monomer_file, sys_geom, error)

Read full geometry and monomer template, initialize system_geometry_t

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: full_geom_file
character(len=*), intent(in) :: monomer_file
type(system_geometry_t), intent(out) :: sys_geom
type(error_t), intent(out) :: error

Calls

proc~~initialize_system_geometry~~CallsGraph proc~initialize_system_geometry initialize_system_geometry proc~element_symbol_to_number element_symbol_to_number proc~initialize_system_geometry->proc~element_symbol_to_number proc~error_add_context error_t%error_add_context proc~initialize_system_geometry->proc~error_add_context proc~error_has_error error_t%error_has_error proc~initialize_system_geometry->proc~error_has_error proc~error_set error_t%error_set proc~initialize_system_geometry->proc~error_set proc~geometry_destroy geometry_type%geometry_destroy proc~initialize_system_geometry->proc~geometry_destroy proc~read_xyz_file read_xyz_file proc~initialize_system_geometry->proc~read_xyz_file proc~to_bohr to_bohr proc~initialize_system_geometry->proc~to_bohr to_lower to_lower proc~element_symbol_to_number->to_lower to_upper to_upper proc~element_symbol_to_number->to_upper proc~read_xyz_file->proc~error_add_context proc~read_xyz_file->proc~error_has_error proc~read_xyz_file->proc~error_set proc~read_xyz_string read_xyz_string proc~read_xyz_file->proc~read_xyz_string proc~read_xyz_string->proc~error_set proc~int_to_string int_to_string proc~read_xyz_string->proc~int_to_string proc~split_lines split_lines proc~read_xyz_string->proc~split_lines

Variables

Type Visibility Attributes Name Initial
type(geometry_type), private :: full_geom
integer, private :: i
type(geometry_type), private :: monomer_geom

Source Code

   subroutine initialize_system_geometry(full_geom_file, monomer_file, sys_geom, error)
      !! Read full geometry and monomer template, initialize system_geometry_t
      character(len=*), intent(in) :: full_geom_file, monomer_file
      type(system_geometry_t), intent(out) :: sys_geom
      type(error_t), intent(out) :: error

      type(geometry_type) :: full_geom, monomer_geom
      integer :: i

      call read_xyz_file(full_geom_file, full_geom, error)
      if (error%has_error()) then
         call error%add_context("mqc_physical_fragment:initialize_system_geometry")
         return
      end if

      ! Read monomer template
      ! this will be changed once we have a proper input file parsing
      call read_xyz_file(monomer_file, monomer_geom, error)
      if (error%has_error()) then
         call error%add_context("mqc_physical_fragment:initialize_system_geometry")
         call full_geom%destroy()
         return
      end if

      ! Validate that full geometry is a multiple of monomer size
      sys_geom%atoms_per_monomer = monomer_geom%natoms
      sys_geom%total_atoms = full_geom%natoms

      if (mod(sys_geom%total_atoms, sys_geom%atoms_per_monomer) /= 0) then
         call error%set(ERROR_VALIDATION, "Full geometry atoms not a multiple of monomer atoms")
         call full_geom%destroy()
         call monomer_geom%destroy()
         return
      end if

      sys_geom%n_monomers = sys_geom%total_atoms/sys_geom%atoms_per_monomer

      ! TODO JORGE: this can be a sys_geom%allocate()
      allocate (sys_geom%element_numbers(sys_geom%total_atoms))
      allocate (sys_geom%coordinates(3, sys_geom%total_atoms))

      do i = 1, sys_geom%total_atoms
         sys_geom%element_numbers(i) = element_symbol_to_number(full_geom%elements(i))
      end do

      ! Store coordinates in Bohr (convert from Angstroms)
      ! TODO JORGE: need a way to handle units
      sys_geom%coordinates = to_bohr(full_geom%coords)

      call full_geom%destroy()
      call monomer_geom%destroy()

   end subroutine initialize_system_geometry