anchor.diff:           ! A NEGATIVE EXPONENTIAL DlSTRIBUTION WlTH DENSITY PROPORTIONAL
            ! TO EXP(-random_exponential), USING INVERSION.
   
  -         REAL  :: fn_val
  -         REAL  :: Sigma
  +         REAL(ReKi)  :: fn_val
  +         REAL(ReKi)  :: Sigma
   
            !     Local variable
  -         REAL  :: r
  +         REAL(ReKi)  :: r
   
            DO
              CALL RANDOM_NUMBER(r)
  
positive.code:         FUNCTION digamma(z) RESULT(phy)
  
           !Calculate the value of digamma function of z
           REAL, INTENT(IN) :: z
           REAL             :: phy
  
           phy = log(z) - 1/2/z - 1/12/z**2 + 1/120/z**4 - 1/252/z**6 + 1/240/z**8 - 5/660/z**10;
  
  
negative.code: 
           IMPLICIT NONE
  
           REAL(ReKi) :: mean
           REAL(ReKi) :: var
           REAL(ReKi) :: fn_val
  
           !Local variables                                                                                                                                                                     
           REAL(ReKi) :: k
           REAL(ReKi) :: Lambda
           REAL(ReKi) :: u
  
           k = wbpar (mean, var)
           lambda = mean / NWTC_gamma(1+1/k)
  