mimic_particles.f90 Source File


Contents

Source Code


Source Code

!
!    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