spx_nnps_pairs.f90 Source File


Contents

Source Code


Source Code

!> NNPS 粒子对
!> @note (20230722) time used with build: 18.2%, time used with query: 52.7%
module spx_nnps_pairs

    use spx_kinds, only: rk
    use spx_math, only: dim
    use spx_smoothed_kernel_function, only: skf_obj
    use nnps_module, only: nnps_grid2d
    use spx_configuration, only: spc_obj
    use spx_contact_control, only: contact_control
    implicit none

    private
    public :: nnps_grid, nnps_pairs_type, pairs

    !> 交互粒子对
    !> @note 需要与 items 联系起来
    type nnps_pairs_type
        integer :: len                  !! 粒子对数量
        integer, private :: cap         !! 向量容量
        integer, pointer :: items(:)    !! 粒子对的2个粒子序号的数组 [(i, j), ...]
        real(rk), pointer :: rdx(:)     !! r, dx(:)
        real(rk), allocatable :: wdwdx(:, :)    !! w, dwdx(:)
        integer, allocatable :: contact_type(:) !! 接触控制类型
    contains
        procedure :: init, calc
        procedure, private :: extend
    end type nnps_pairs_type

    type(nnps_pairs_type) :: pairs  !! 交互粒子对
    type(nnps_grid2d) :: nnps_grid  !! NNPS 网格数据结构

contains

    !> 初始化内存空间
    subroutine init(self, n)
        class(nnps_pairs_type), intent(inout) :: self
        integer, intent(in) :: n    !! 粒子对容量

        allocate (self%wdwdx(dim + 1, n))
        allocate (self%contact_type(n))
        self%cap = n

    end subroutine init

    !> 拓展内存空间, @todo 优化效率
    subroutine extend(self, len)
        class(nnps_pairs_type), intent(inout) :: self
        integer, intent(in) :: len  !! 粒子对容量

        self%len = len
        if (self%cap < len) then
            self%cap = len*2
            deallocate (self%wdwdx, self%contact_type)
            allocate (self%wdwdx(1 + dim, self%cap), self%contact_type(self%cap))
        end if

    end subroutine extend

    !> 计算 dW/dX @todo 接触控制
    subroutine calc(self, itype)
        class(nnps_pairs_type), intent(inout) :: self
        integer, intent(in) :: itype(:)     !! 粒子类型
        integer :: i, ix

        call self%extend(size(self%items)/2)

        !$omp parallel do private(i, ix)
        do i = 1, self%len
            ix = 3*i
            call skf_obj%kernel(self%rdx(ix - 2), self%rdx(ix - 1:ix), &
                                self%wdwdx(1, i), self%wdwdx(2:3, i))
            ix = 2*i
            self%contact_type(i) = contact_control(itype, self%items(ix - 1), self%items(ix))
        end do

    end subroutine calc

end module spx_nnps_pairs