white_noise.F

Go to the documentation of this file.

#include "cppdefs.h"
      MODULE white_noise_mod
!
!git $Id$
!================================================== Hernan G. Arango ===
!  Copyright (c) 2002-2025 The ROMS Group                              !
!    Licensed under a MIT/X style license                              !
!    See License_ROMS.md                                               !
!=======================================================================
!                                                                      !
!  These routines generates "white noise" arrays with random numbers   !
!  over specified range.  These random numbers are scaled to insure,   !
!  approximately, zero mean expectation and unit variance.             !
!                                                                      !
!  Several random number generation schemes are allowed:               !
!                                                                      !
!  Rscheme = 0         Use F90 intrinsic random numbers generator,     !
!                        0 <= R < 1                                    !
!  Rscheme = 1         Random deviates with Gaussian distribution,     !
!                        -1 < R < 1                                    !
!                                                                      !
!  Routines:                                                           !
!                                                                      !
!  white_noise1d       Random numbers for 1D arrays                    !
!  white_noise2d       Random numbers for 2D arrays                    !
!  white_noise2d_bry   Random numbers for 2D boundary arrays           !
!  white_noise3d       Random numbers for 3D arrays                    !
!  white_noise3d_bry   Random numbers for 3D boundary arrays           !
!                                                                      !
!=======================================================================
!
      USE mod_kinds
 
      implicit none
 
      PRIVATE
 
      PUBLIC :: white_noise1d
      PUBLIC :: white_noise2d
      PUBLIC :: white_noise2d_bry
#ifdef SOLVE3D
      PUBLIC :: white_noise3d
      PUBLIC :: white_noise3d_bry
#endif
 
      CONTAINS
!
!***********************************************************************
      SUBROUTINE white_noise1d (ng, model, Rscheme,                     &
     &                          Imin, Imax, LBi, UBi,                   &
     &                          Rmin, Rmax, R)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_scalars
!
      USE nrutil, ONLY : gasdev
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, model, rscheme
      integer, intent(in) :: imin, imax, lbi, ubi
 
      real(r8), intent(out) :: rmin, rmax
 
#ifdef ASSUMED_SHAPE
      real(r8), intent(out) :: r(lbi:)
#else
      real(r8), intent(out) :: r(lbi:ubi)
#endif
!
!  Local variable declarations.
!
      integer :: i, ic
 
      real(r8), parameter :: fac = 2.0_r8 * 1.73205080756887720_r8
 
      real(r8), dimension((UBi-LBi+1)) :: random
!
!-----------------------------------------------------------------------
!  Generate an array of random numbers.
!-----------------------------------------------------------------------
!
!  Initialize output array.
!
      DO i=lbi,ubi
        r(i)=0.0_r8
      END DO
!
!  F90 intrinsic pseudorandom numbers with an uniform distribution
!  over the range 0 <= R < 1.
!
      IF (rscheme.eq.0) THEN
        CALL random_number (r(lbi:))
!
!  Scale (fac=2*SQRT(3)) the random number to insure the expectation
!  mean to be approximately zero, E(R)=0, and the expectation variance
!  to be approximately unity, E(R^2)=1 (Bennett, 2002; book page 72).
!
!  Recall that,
!
!      E(R) = sum[R f(R)]
!
!  where F(R) is the random variable probability distribution.
!
        DO i=imin,imax
          r(i)=fac*(r(i)-0.5_r8)
        END DO
        rmin=minval(r)
        rmax=maxval(r)
!
!  Random deviates with a Gaussian (normal) distribuiton over the
!  range from -1 to 1.
!
      ELSE IF (rscheme.eq.1) THEN
        IF (master) THEN
          CALL gasdev (random)
          rmin=minval(random)
          rmax=maxval(random)
        END IF
        ic=0
        DO i=imin,imax
          ic=ic+1
          r(i)=random(ic)
        END DO
      END IF
 
      RETURN
      END SUBROUTINE white_noise1d
 
!
!***********************************************************************
      SUBROUTINE white_noise2d (ng, model, gtype, Rscheme,              &
     &                          Imin, Imax, Jmin, Jmax,                 &
     &                          LBi, UBi, LBj, UBj,                     &
     &                          Rmin, Rmax, R)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_scalars
!
#ifdef DISTRIBUTE
      USE distribute_mod, ONLY : mp_scatter2d
#endif
      USE nrutil, ONLY : gasdev
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, model, gtype, rscheme
      integer, intent(in) :: imin, imax, jmin, jmax
      integer, intent(in) :: lbi, ubi, lbj, ubj
 
      real(r8), intent(out) :: rmin, rmax
 
      real(r8), intent(out) :: r(lbi:ubi,lbj:ubj)
!
!  Local variable declarations.
!
      integer :: npts, i, ic, j
 
      real(r8), parameter :: fac = 2.0_r8 * 1.73205080756887720_r8
 
      real(r8), allocatable :: random(:)
!
!-----------------------------------------------------------------------
!  Generate an array of random numbers.
!-----------------------------------------------------------------------
!
!  Initialize output array.
!
      DO j=lbj,ubj
        DO i=lbi,ubi
          r(i,j)=0.0_r8
        END DO
      END DO
!
!  F90 intrinsic pseudorandom numbers with an uniform distribution
!  over the range 0 <= R < 1.
!
      IF (rscheme.eq.0) THEN
        CALL random_number (r(lbi:,lbj:))
!
!  Scale (fac=2*SQRT(3)) the random number to insure the expectation
!  mean to be approximately zero, E(R)=0, and the expectation variance
!  to be approximately unity, E(R^2)=1 (Bennett, 2002; book page 72).
!
!  Recall that,
!
!      E(R) = sum[R f(R)]
!
!  where F(R) is the random variable probability distribution.
!
        DO j=jmin,jmax
          DO i=imin,imax
            r(i,j)=fac*(r(i,j)-0.5_r8)
          END DO
        END DO
        rmin=minval(r)
        rmax=maxval(r)
!
!  Random deviates with a Gaussian (normal) distribuiton with zero
!  mean and unit variance.
!
      ELSE IF (rscheme.eq.1) THEN
        npts=(lm(ng)+2)*(mm(ng)+2)
        IF (.not.allocated(random)) THEN
          allocate ( random(npts+2) )
        END IF
        IF (master) THEN
          CALL gasdev (random)
          rmin=random(1)
          rmax=random(1)
          DO ic=1,npts
            rmin=min(rmin,random(ic))
            rmax=max(rmax,random(ic))
          END DO
        END IF
#ifdef DISTRIBUTE
        npts=SIZE(random)
        CALL mp_scatter2d (ng, model, lbi, ubi, lbj, ubj,               &
     &                     nghostpoints, gtype, rmin, rmax,             &
# if defined READ_WATER && defined MASKING
     &                     (lm(ng)+2)*(mm(ng)+2),                       &
     &                     scalars(ng)%IJwater(:,gtype),                &
# endif
     &                     npts, random, r)
#else
        ic=0
        DO j=jmin,jmax
          DO i=imin,imax
            ic=ic+1
            r(i,j)=random(ic)
          END DO
        END DO
#endif
        IF (allocated(random)) THEN
          deallocate (random)
        END IF
      END IF
 
      RETURN
      END SUBROUTINE white_noise2d
 
!
!***********************************************************************
      SUBROUTINE white_noise2d_bry (ng, tile, model, boundary, Rscheme, &
     &                              Imin, Imax,                         &
     &                              LBij, UBij,                         &
     &                              Rmin, Rmax, R)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_scalars
!
#ifdef DISTRIBUTE
      USE distribute_mod, ONLY : mp_bcastf
      USE mp_exchange_mod, ONLY : mp_exchange2d_bry
#endif
      USE nrutil, ONLY : gasdev
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, tile, model, boundary, rscheme
      integer, intent(in) :: imin, imax
      integer, intent(in) :: lbij, ubij
 
      real(r8), intent(out) :: rmin, rmax
 
#ifdef ASSUMED_SHAPE
      real(r8), intent(out) :: r(lbij:)
#else
      real(r8), intent(out) :: r(lbij:ubij)
#endif
!
!  Local variable declarations.
!
      integer :: ioff, i, ic
 
      real(r8), parameter :: fac = 2.0_r8 * 1.73205080756887720_r8
 
      real(r8), dimension(UBij-LBij+1) :: random
!
!-----------------------------------------------------------------------
!  Generate an array of random numbers.
!-----------------------------------------------------------------------
!
!  Initialize output array.
!
      DO i=lbij,ubij
        r(i)=0.0_r8
      END DO
!
!  F90 intrinsic pseudorandom numbers with an uniform distribution
!  over the range 0 <= R < 1.
!
      IF (rscheme.eq.0) THEN
        CALL random_number (r(lbij:))
!
!  Scale (fac=2*SQRT(3)) the random number to insure the expectation
!  mean to be approximately zero, E(R)=0, and the expectation variance
!  to be approximately unity, E(R^2)=1 (Bennett, 2002; book page 72).
!
!  Recall that,
!
!      E(R) = sum[R f(R)]
!
!  where F(R) is the random variable probability distribution.
!
        DO i=imin,imax
          r(i)=fac*(r(i)-0.5_r8)
        END DO
        rmin=minval(r)
        rmax=maxval(r)
!
!  Random deviates with a Gaussian (normal) distribuiton with zero
!  mean and unit variance.
!
      ELSE IF (rscheme.eq.1) THEN
        IF (master) THEN
          CALL gasdev (random)
        END IF
#ifdef DISTRIBUTE
        CALL mp_bcastf (ng, model, random)
#endif
        rmin=minval(random)
        rmax=maxval(random)
 
        ioff=1-lbij
        DO i=imin,imax
          ic=i+ioff
          r(i)=random(ic)
        END DO
#ifdef DISTRIBUTE
        CALL mp_exchange2d_bry (ng, tile, model, 1, boundary,           &
     &                          lbij, ubij,                             &
     &                          nghostpoints,                           &
     &                          ewperiodic(ng), nsperiodic(ng),         &
     &                          r)
#endif
      END IF
 
      RETURN
      END SUBROUTINE white_noise2d_bry
 
 
#ifdef SOLVE3D
!
!***********************************************************************
      SUBROUTINE white_noise3d (ng, model, gtype, Rscheme,              &
     &                          Imin, Imax, Jmin, Jmax,                 &
     &                          LBi, UBi, LBj, UBj, LBk, UBk,           &
     &                          Rmin, Rmax, R)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_scalars
!
# ifdef DISTRIBUTE
      USE distribute_mod, ONLY : mp_scatter3d
# endif
      USE nrutil, ONLY : gasdev
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, model, gtype, rscheme
      integer, intent(in) :: imin, imax, jmin, jmax
      integer, intent(in) :: lbi, ubi, lbj, ubj, lbk, ubk
 
      real(r8), intent(out) :: rmin, rmax
 
# ifdef ASSUMED_SHAPE
      real(r8), intent(out) :: r(lbi:,lbj:,lbk:)
# else
      real(r8), intent(out) :: r(lbi:ubi,lbj:ubj,lbk:ubk)
# endif
!
!  Local variable declarations.
!
      integer :: npts, i, ic, j, k
 
      real(r8), parameter :: fac = 2.0_r8 * 1.73205080756887720_r8
 
      real(r8), allocatable :: random(:)
!
!-----------------------------------------------------------------------
!  Generate an array of random numbers.
!-----------------------------------------------------------------------
!
!  Initialize output array.
!
      DO k=lbk,ubk
        DO j=lbj,ubj
          DO i=lbi,ubi
            r(i,j,k)=0.0_r8
          END DO
        END DO
      END DO
!
!  F90 intrinsic pseudorandom numbers with an uniform distribution
!  over the range 0 <= R < 1.
!
      IF (rscheme.eq.0) THEN
        CALL random_number (r(lbi:,lbj:,lbk:))
!
!  Scale (fac=2*SQRT(3)) the random number to insure the expectation
!  mean to be approximately zero, E(R)=0, and the expectation variance
!  to be approximately unity, E(R^2)=1.
!
!  Recall that,
!
!      E(R) = sum[R f(R)]
!
!  where F(R) is the random variable probability distribution
!
        DO k=lbk,ubk
          DO j=jmin,jmax
            DO i=imin,imax
              r(i,j,k)=fac*(r(i,j,k)-0.5_r8)
            END DO
          END DO
        END DO
        rmin=minval(r)
        rmax=maxval(r)
!
!  Random deviates with a Gaussian (normal) distribution with zero
!  mean and unit variance.
!
      ELSE IF (rscheme.eq.1) THEN
        npts=(lm(ng)+2)*(mm(ng)+2)*(ubk-lbk+1)
        IF (.not.allocated(random)) THEN
          allocate ( random(npts+2) )
        END IF
        IF (master) THEN
          CALL gasdev (random)
          rmin=random(1)
          rmax=random(1)
          DO ic=1,npts
            rmin=min(rmin,random(ic))
            rmax=max(rmax,random(ic))
          END DO
        END IF
# ifdef DISTRIBUTE
        npts=SIZE(random)
        CALL mp_scatter3d (ng, model, lbi, ubi, lbj, ubj, lbk, ubk,     &
     &                     nghostpoints, gtype, rmin, rmax,             &
#  if defined READ_WATER && defined MASKING
     &                     (lm(ng)+2)*(mm(ng)+2),                       &
     &                     scalars(ng)%IJwater(:,gtype),                &
#  endif
     &                     npts, random, r)
# else
        ic=0
        DO k=lbk,ubk
          DO j=jmin,jmax
            DO i=imin,imax
              ic=ic+1
              r(i,j,k)=random(ic)
            END DO
          END DO
        END DO
# endif
        IF (allocated(random)) THEN
          deallocate (random)
        END IF
      END IF
 
      RETURN
      END SUBROUTINE white_noise3d
 
!
!***********************************************************************
      SUBROUTINE white_noise3d_bry (ng, tile, model, boundary, Rscheme, &
     &                              Imin, Imax,                         &
     &                              LBij, UBij, LBk, UBk,               &
     &                              Rmin, Rmax, R)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_scalars
!
# ifdef DISTRIBUTE
      USE distribute_mod, ONLY : mp_bcastf
      USE mp_exchange_mod, ONLY : mp_exchange3d_bry
# endif
      USE nrutil, ONLY : gasdev
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, tile, model, boundary, rscheme
      integer, intent(in) :: imin, imax
      integer, intent(in) :: lbij, ubij, lbk, ubk
 
      real(r8), intent(out) :: rmin, rmax
 
# ifdef ASSUMED_SHAPE
      real(r8), intent(out) :: r(lbij:,lbk:)
# else
      real(r8), intent(out) :: r(lbij:ubij,lbk:ubk)
# endif
!
!  Local variable declarations.
!
      integer :: ilen, ioff, i, ic, k, kc
 
      real(r8), parameter :: fac = 2.0_r8 * 1.73205080756887720_r8
 
      real(r8), dimension((UBij-LBij+1)*(UBk-LBk+1)) :: random
!
!-----------------------------------------------------------------------
!  Generate an array of random numbers.
!-----------------------------------------------------------------------
!
!  Initialize output array.
!
      DO k=lbk,ubk
        DO i=lbij,ubij
          r(i,k)=0.0_r8
        END DO
      END DO
!
!  F90 intrinsic pseudorandom numbers with an uniform distribution
!  over the range 0 <= R < 1.
!
      IF (rscheme.eq.0) THEN
        CALL random_number (r(lbij:,lbk:))
!
!  Scale (fac=2*SQRT(3)) the random number to insure the expectation
!  mean to be approximately zero, E(R)=0, and the expectation variance
!  to be approximately unity, E(R^2)=1 (Bennett, 2002; book page 72).
!
!  Recall that,
!
!      E(R) = sum[R f(R)]
!
!  where F(R) is the random variable probability distribution.
!
        DO k=lbk,ubk
          DO i=imin,imax
            r(i,k)=fac*(r(i,k)-0.5_r8)
          END DO
        END DO
        rmin=minval(r)
        rmax=maxval(r)
!
!  Random deviates with a Gaussian (normal) distribuiton with zero
!  mean and unit variance.
!
      ELSE IF (rscheme.eq.1) THEN
        IF (master) THEN
          CALL gasdev (random)
        END IF
# ifdef DISTRIBUTE
        CALL mp_bcastf (ng, model, random)
# endif
        rmin=minval(random)
        rmax=maxval(random)
 
        ilen=ubij-lbij+1
        ioff=1-lbij
 
        DO k=lbk,ubk
          kc=(k-lbk)*ilen
          DO i=imin,imax
            ic=i+ioff+kc
            r(i,k)=random(ic)
          END DO
        END DO
# ifdef DISTRIBUTE
        CALL mp_exchange3d_bry (ng, tile, model, 1, boundary,           &
     &                          lbij, ubij, lbk, ubk,                   &
     &                          nghostpoints,                           &
     &                          ewperiodic(ng), nsperiodic(ng),         &
     &                          r)
# endif
      END IF
 
      RETURN
      END SUBROUTINE white_noise3d_bry
#endif
 
      END MODULE white_noise_mod