source: palm/trunk/SOURCE/random_gauss.f90 @ 3768

Last change on this file since 3768 was 3655, checked in by knoop, 3 years ago

Bugfix: made "unit" and "found" intend INOUT in module interface subroutines + automatic copyright update

  • Property svn:keywords set to Id
File size: 3.4 KB
Line 
1!> @file random_gauss.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2019 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! -----------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: random_gauss.f90 3655 2019-01-07 16:51:22Z raasch $
27! Corrected "Former revisions" section
28!
29! 2696 2017-12-14 17:12:51Z kanani
30! Change in file header (GPL part)
31!
32! 2101 2017-01-05 16:42:31Z suehring
33!
34! 2000 2016-08-20 18:09:15Z knoop
35! Forced header and separation lines into 80 columns
36!
37! 1682 2015-10-07 23:56:08Z knoop
38! Code annotations made doxygen readable
39!
40! 1342 2014-03-26 17:04:47Z kanani
41! REAL constants defined as wp-kind
42!
43! 1320 2014-03-20 08:40:49Z raasch
44! ONLY-attribute added to USE-statements,
45! kind-parameters added to all INTEGER and REAL declaration statements,
46! kinds are defined in new module kinds,
47! old module precision_kind is removed,
48! revision history before 2012 removed,
49! comment fields (!:) to be used for variable explanations added to
50! all variable declaration statements
51!
52! 1036 2012-10-22 13:43:42Z raasch
53! code put under GPL (PALM 3.9)
54!
55! RCS Log replace by Id keyword, revision history cleaned up
56!
57! Revision 1.1  1998/03/25 20:09:47  raasch
58! Initial revision
59!
60!
61! Description:
62! ------------
63!> Generates a gaussian distributed random number (mean value 1, sigma = 1)
64!> This routine is taken from the "numerical recipies".
65!------------------------------------------------------------------------------!
66 FUNCTION random_gauss( idum, upper_limit )
67 
68
69    USE kinds
70
71    USE random_function_mod,                                                   &
72        ONLY:  random_function
73
74    IMPLICIT NONE
75
76    INTEGER(iwp) ::  idum          !<
77    INTEGER(iwp) ::  iset          !<
78
79    REAL(wp)     ::  fac           !<
80    REAL(wp)     ::  gset          !<
81    REAL(wp)     ::  random_gauss  !<
82    REAL(wp)     ::  rsq           !<
83    REAL(wp)     ::  upper_limit   !<
84    REAL(wp)     ::  v1            !<
85    REAL(wp)     ::  v2            !<
86
87    SAVE  iset, gset
88
89    DATA  iset /0/
90
91!
92!-- Random numbers are created as long as they do not fall below the given
93!-- upper limit
94    DO
95
96       IF ( iset == 0 )  THEN
97          rsq = 0.0_wp
98          DO  WHILE ( rsq >= 1.0_wp  .OR.  rsq == 0.0_wp )
99             v1  = 2.0_wp * random_function( idum ) - 1.0_wp
100             v2  = 2.0_wp * random_function( idum ) - 1.0_wp
101             rsq = v1**2 + v2**2
102          ENDDO
103          fac          = SQRT( -2.0_wp * LOG( rsq ) / rsq )
104          gset         = v1 * fac
105          random_gauss = v2 * fac + 1.0_wp
106          iset         = 1
107       ELSE
108          random_gauss = gset + 1.0_wp
109          iset         = 0
110       ENDIF
111
112       IF ( ABS( random_gauss - 1.0_wp ) <= upper_limit )  EXIT
113
114    ENDDO
115
116 END FUNCTION random_gauss
Note: See TracBrowser for help on using the repository browser.