! ! MiMiC: A Framework for Multiscale Modeling in Computational Chemistry ! Copyright (C) 2015-2022 The MiMiC Contributors (see CONTRIBUTORS file for details). ! ! This file is part of MiMiC. ! ! MiMiC is free software: you can redistribute it and/or modify ! it under the terms of the GNU Lesser General Public License as ! published by the Free Software Foundation, either version 3 of ! the License, or (at your option) any later version. ! ! MiMiC is distributed in the hope that it will be useful, but ! WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU Lesser General Public License for more details. ! ! You should have received a copy of the GNU Lesser General Public License ! along with this program. If not, see <http://www.gnu.org/licenses/>. ! !> contains particle types module mimic_particles use mimic_precision, only: dp use mimic_properties implicit none private public :: nucleus_type public :: atom_type public :: site_type public :: compute_minimum_distance public :: compute_centroid public :: compute_center_of_mass public :: compute_center_of_charge !> base type for particles type, abstract :: particle_type !> particle id integer, public :: id !> particle charge real(dp), public :: charge !> particle mass real(dp), public :: mass !> pointer to coordinate of the particle real(dp), dimension(:), pointer, public :: coordinate => null() !> pointer to force acting on particle real(dp), dimension(:), pointer, public :: force => null() !> flag indicating that the particle is treated by some other code logical :: overlapped end type particle_type !> type representing a nucleus type, extends(particle_type) :: nucleus_type private !> id of CPMD species integer, public :: species_id !> id of CPMD atom integer, public :: atom_id contains private procedure :: init_nucleus generic, public :: init => init_nucleus end type nucleus_type !> type representing an atom type, extends(nucleus_type) :: atom_type private !> covalent radius real(dp), public :: radius !> multipoles type(multipoles_type), public :: multipoles contains private procedure :: init_atom generic, public :: init => init_atom end type atom_type !> non-atomic site (e.g. bond midpoints) type, extends(particle_type) :: site_type private !> coordinate of site real(dp), dimension(3) :: internal_coordinate !> force acting on site real(dp), dimension(3) :: internal_force !> multipoles associated with the site type(multipoles_type), public :: multipoles contains private procedure, public :: init => init_site ! procedure :: update_coordinate ! procedure :: update_force end type site_type contains subroutine init_nucleus(this, id, & species_id, atom_id, & charge, coordinate, & force, overlapped) class(nucleus_type), intent(inout) :: this integer, intent(in) :: id, species_id, atom_id real(dp), intent(in) :: charge real(dp), dimension(:), intent(in), target :: coordinate, force logical, intent(in) :: overlapped this%id = id this%species_id = species_id this%atom_id = atom_id this%charge = charge this%coordinate => coordinate this%force => force this%overlapped = overlapped end subroutine init_nucleus subroutine init_atom(this, id, & species_id, atom_id, & charge, radius, & coordinate, force, & multipoles, overlapped) class(atom_type), intent(inout) :: this integer, intent(in) :: id, species_id, atom_id real(dp), intent(in) :: charge, radius real(dp), dimension(:), intent(in), target :: coordinate, force type(multipoles_type), intent(in) :: multipoles logical, intent(in) :: overlapped call this%nucleus_type%init(id, species_id, atom_id, charge, coordinate, force, overlapped) this%radius = radius this%multipoles = multipoles end subroutine init_atom subroutine init_site(this, id, coordinate, force, multipoles, overlapped) class(site_type), intent(inout) :: this integer, intent(in) :: id real(dp), dimension(:), intent(in) :: coordinate real(dp), dimension(:), intent(in) :: force type(multipoles_type), intent(in), optional :: multipoles logical, intent(in) :: overlapped this%id = id allocate(this%coordinate(3)) this%coordinate = coordinate allocate(this%force(3)) this%force = force if (present(multipoles)) then this%multipoles = multipoles end if this%overlapped = overlapped end subroutine init_site !> computes minimum distance between two sets of particles pure function compute_minimum_distance(particles_one, particles_two) result(minimum_distance) class(particle_type), dimension(:), intent(in) :: particles_one class(particle_type), dimension(:), intent(in) :: particles_two integer :: i, j real(dp) :: current_distance real(dp) :: minimum_distance real(dp), dimension(3) :: coordinate_one real(dp), dimension(3) :: coordinate_two minimum_distance = huge(dp) do i = 1, size(particles_one) coordinate_one = particles_one(i)%coordinate do j = 1, size(particles_two) coordinate_two = particles_two(j)%coordinate current_distance = norm2(coordinate_two - coordinate_one) if (current_distance < minimum_distance) then minimum_distance = current_distance end if end do end do end function compute_minimum_distance !> compute centroid of a set of particles pure function compute_centroid(particles) result(centroid) class(particle_type), dimension(:), intent(in) :: particles integer :: i real(dp), dimension(3) :: centroid centroid = 0.0_dp do i = 1, size(particles) centroid = centroid + particles(i)%coordinate end do centroid = centroid / real(size(particles), dp) end function compute_centroid !> compute center-of-mass of a set of particles pure function compute_center_of_mass(particles) result(center_of_mass) class(particle_type), dimension(:), intent(in) :: particles integer :: i real(dp) :: mass real(dp) :: total_mass real(dp), dimension(3) :: coordinate real(dp), dimension(3) :: center_of_mass total_mass = 0.0_dp center_of_mass = 0.0_dp do i = 1, size(particles) mass = particles(i)%mass coordinate = particles(i)%coordinate total_mass = total_mass + mass center_of_mass = center_of_mass + mass * coordinate end do center_of_mass = center_of_mass / total_mass end function compute_center_of_mass !> compute center-of-charge of a set particles pure function compute_center_of_charge(particles) result(center_of_charge) class(particle_type), dimension(:), intent(in) :: particles integer :: i real(dp) :: charge real(dp) :: total_charge real(dp), dimension(3) :: coordinate real(dp), dimension(3) :: center_of_charge total_charge = 0.0_dp center_of_charge = 0.0_dp do i = 1, size(particles) charge = particles(i)%charge coordinate = particles(i)%coordinate total_charge = total_charge + charge center_of_charge = center_of_charge + charge * coordinate end do center_of_charge = center_of_charge / total_charge end function compute_center_of_charge end module mimic_particles