two_normals Function

private function two_normals(self) result(rands)

Return two normal random variates with mean 0 and variance 1. http://en.wikipedia.org/wiki/Marsaglia_polar_method

Type Bound

rng_t

Arguments

Type IntentOptional Attributes Name
class(rng_t), intent(inout) :: self

Return Value real(kind=dp), (2)


Calls

proc~~two_normals~~CallsGraph proc~two_normals rng_t%two_normals proc~unif_01 rng_t%unif_01 proc~two_normals->proc~unif_01 proc~next rng_t%next proc~unif_01->proc~next proc~rotl rotl proc~next->proc~rotl

Called by

proc~~two_normals~~CalledByGraph proc~two_normals rng_t%two_normals proc~normal rng_t%normal proc~normal->proc~two_normals proc~get_maxwellian_velocity get_maxwellian_velocity proc~get_maxwellian_velocity->proc~normal

Source Code

  function two_normals(self) result(rands)
    class(rng_t), intent(inout) :: self
    real(dp)                    :: rands(2), sum_sq

    do
       rands(1) = 2 * self%unif_01() - 1
       rands(2) = 2 * self%unif_01() - 1
       sum_sq = sum(rands**2)
       if (sum_sq < 1.0_dp .and. sum_sq > 0.0_dp) exit
    end do
    rands = rands * sqrt(-2 * log(sum_sq) / sum_sq)
  end function two_normals