gaussian_torus Subroutine

private subroutine gaussian_torus(params, spp)

@brief Subroutine that generates a Gaussian radial distribution of particles in a circular cross-section torus of major and minor radi $R_0$ and $r_0$, respectively. @details We generate this exponentially decaying radial distribution $f(r)$ following the same approach as in \ref korc_spatial_distribution.disk, but this time, the radial distribution is given by:

$$f(r) = \left\{ \begin{array}{ll} \frac{1}{4\pi^2 \sigma^2 R_0}\frac{\exp{\left( -\frac{r^2}{2\sigma^2} \right)}}{1 - \exp{\left( -\frac{r_0^2}{2\sigma^2} \right)}} &\ r<r_0 \\ 0 & r>r_0 \end{array} \right.$$

The radial position of the particles $r$ is obtained using the Inverse Trasnform Sampling method, finding $r$ numerically through the Newton-Raphson method. First, we calculate the particles' radial distribution in a disk centered at $(R,Z) = (0,0)$. Then, we transfor to a new set of coordinates where the disk is centered at $(R,Z) = (R_0,Z_0)$. Finally, we generate the toroidal distribution by givin each particle a toroidal angle $\zeta$ which follows a uniform distribution in the interval $[0,2\pi]$.

@param[in] params Core KORC simulation parameters. @param[in,out] spp An instance of the derived type SPECIES containing all the parameters and simulation variables of the different species in the simulation. @param fl Variable used in the Newton-Raphson method for finding the radial position of each particle. @param fr Variable used in the Newton-Raphson method for finding the radial position of each particle. @param fm Variable used in the Newton-Raphson method for finding the radial position of each particle. @param rl Variable used in the Newton-Raphson method for finding the radial position of each particle. @param rr Variable used in the Newton-Raphson method for finding the radial position of each particle. @param rm Variable used in the Newton-Raphson method for finding the radial position of each particle. @param relerr Tolerance used to determine when to stop iterating the Newton-Raphson method for finding $r$. @param r Radial position of the particles $r$. @param theta Uniform deviates in the range $[0,2\pi]$ representing the uniform poloidal angle $\theta$ distribution of the particles. @param zeta Uniform deviates in the range $[0,2\pi]$ representing the uniform toroidal angle $\zeta$ distribution of the particles. @param pp Particle iterator.

Arguments

subroutine gaussian_torus(params,spp)
  TYPE(KORC_PARAMS), INTENT(IN) 			:: params
  TYPE(SPECIES), INTENT(INOUT) 			:: spp
  REAL(rp), DIMENSION(:), ALLOCATABLE 	:: theta
  REAL(rp), DIMENSION(:), ALLOCATABLE 	:: zeta
  REAL(rp), DIMENSION(:), ALLOCATABLE 	:: r ! temporary vars
  REAL(rp) 				:: sigma

  ALLOCATE( theta(spp%ppp) )
  ALLOCATE( zeta(spp%ppp) )
  ALLOCATE( r(spp%ppp) )

  ! Initial condition of uniformly distributed particles on a disk in the xz-plane
  ! A unique velocity direction
  call init_u_random(10986546_8)

  call init_random_seed()
  call RANDOM_NUMBER(theta)
  theta = 2.0_rp*C_PI*theta

  call init_random_seed()
  call RANDOM_NUMBER(zeta)
  zeta = 2.0_rp*C_PI*zeta

  ! Uniform distribution on a disk at a fixed azimuthal theta
  call init_random_seed()
  call RANDOM_NUMBER(r)

  sigma = 1.0_rp/SQRT(2.0_rp*(spp%falloff_rate/params%cpp%length))
  sigma = sigma/params%cpp%length

  r = sigma*SQRT(-2.0_rp*LOG(1.0_rp - (1.0_rp - &
       EXP(-0.5_rp*spp%r_outter**2/sigma**2))*r))
  spp%vars%X(:,1) = ( spp%Ro + r*COS(theta) )*SIN(zeta)
  spp%vars%X(:,2) = ( spp%Ro + r*COS(theta) )*COS(zeta)
  spp%vars%X(:,3) = spp%Zo + r*SIN(theta)

  DEALLOCATE(theta)
  DEALLOCATE(zeta)
  DEALLOCATE(r)
end subroutine gaussian_torus