poisson_cutoff.F90

Go to the documentation of this file.

!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program 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 General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
!! 02110-1301, USA.
!!
 
#include "global.h"
 
module poisson_cutoff_oct_m
  use global_oct_m
  use, intrinsic :: iso_fortran_env
  use loct_math_oct_m
 
  implicit none
 
  private
  public ::                       &
    poisson_cutoff_3d_0d,         &
    poisson_cutoff_3d_1d,         &
    poisson_cutoff_3d_1d_finite,  &
    poisson_cutoff_3d_2d,         &
    poisson_cutoff_2d_0d,         &
    poisson_cutoff_2d_1d,         &
    poisson_cutoff_1d_0d,         &
    poisson_cutoff_intcoslog
 
 
  interface poisson_cutoff_3d_1d_finite
    real(real64) function c_poisson_cutoff_3d_1d_finite(gx, gperp, xsize, rsize)
      use, intrinsic :: iso_fortran_env
      implicit none
      real(real64), intent(in) :: gx, gperp, rsize, xsize
    end function c_poisson_cutoff_3d_1d_finite
  end interface poisson_cutoff_3d_1d_finite
 
  interface poisson_cutoff_2d_0d
    real(real64) function c_poisson_cutoff_2d_0d(x, y)
      use, intrinsic :: iso_fortran_env
      implicit none
      real(real64), intent(in) :: x, y
    end function c_poisson_cutoff_2d_0d
  end interface poisson_cutoff_2d_0d
 
  interface poisson_cutoff_2d_1d
    real(real64) function c_poisson_cutoff_2d_1d(gy, gx, r_c)
      use, intrinsic :: iso_fortran_env
      implicit none
      real(real64), intent(in) :: gy, gx, r_c
    end function c_poisson_cutoff_2d_1d
  end interface poisson_cutoff_2d_1d
 
  interface poisson_cutoff_1d_0d
    real(real64) function c_poisson_cutoff_1d_0d(g, a, r_c)
      use, intrinsic :: iso_fortran_env
      implicit none
      real(real64), intent(in) :: g, a, r_c
    end function c_poisson_cutoff_1d_0d
  end interface poisson_cutoff_1d_0d
 
  interface poisson_cutoff_intcoslog
    real(real64) function intcoslog(mu, gy, gx)
      use, intrinsic :: iso_fortran_env
      implicit none
      real(real64), intent(in) :: mu, gy, gx
    end function intcoslog
  end interface poisson_cutoff_intcoslog
 
 
contains
 
 
  ! ---------------------------------------------------------
  real(real64) function poisson_cutoff_3D_0D(x, r) result(cutoff)
    real(real64), intent(in) ::  x, r
 
    !no PUSH_SUB, called too often
 
    cutoff = m_one - cos(x*r)
 
  end function poisson_cutoff_3d_0d
  ! ---------------------------------------------------------
 
 
  ! ---------------------------------------------------------
  real(real64) function poisson_cutoff_3D_1D(x, p, rmax) result(cutoff)
    real(real64), intent(in) ::  x, p, rmax
 
    integer :: j
    real(real64) :: dr, r, sum
 
    integer, parameter :: nr = 1000_real64
 
    !no PUSH_SUB, called too often
 
    if (abs(x) < m_epsilon) then
      ! Simpson rule for the G_x = 0 contribution -log(r)
      dr = rmax/real(nr, real64)
      sum = m_zero
      do j = 1, nr - 1, 2
        r = j*dr
        sum = sum - m_four*r*loct_bessel_j0(p*r)*log(r)
        r = (j+1)*dr
        sum = sum - m_two*r*loct_bessel_j0(p*r)*log(r)
      end do
      sum = sum - rmax*loct_bessel_j0(p*rmax)*log(rmax)
      cutoff = (p**2)*m_third*sum*dr
    else
      cutoff = m_one + p*rmax*loct_bessel_j1(p*rmax)*loct_bessel_k0(x*rmax) &
        - x*rmax*loct_bessel_j0(p*rmax)*loct_bessel_k1(x*rmax)
    end if
 
  end function poisson_cutoff_3d_1d
 
  ! ---------------------------------------------------------
  real(real64) function poisson_cutoff_3D_2D(p, z, r) result(cutoff)
    real(real64), intent(in) ::  p, z, r
 
    !no PUSH_SUB, called too often
 
    if (abs(p) < m_epsilon) then
      cutoff = m_one - cos(z*r) - z*r*sin(z*r)
    else
      cutoff = m_one + exp(-p*r)*(z*sin(z*r)/p-cos(z*r))
    end if
 
  end function poisson_cutoff_3d_2d
  ! ---------------------------------------------------------
 
end module poisson_cutoff_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: