Ignore:
Timestamp:
Mar 25, 2014 7:45:13 PM (10 years ago)
Author:
kanani
Message:

REAL constants defined as wp-kind

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/prandtl_fluxes.f90

    r1321 r1340  
    2020! Current revisions:
    2121! -----------------
    22 !
     22! REAL constants defined as wp-kind
    2323!
    2424! Former revisions:
     
    105105       DO  i = nxlg, nxrg
    106106          DO  j = nysg, nyng
    107              ts(j,i) = -shf(j,i) / ( us(j,i) + 1E-30 )
     107             ts(j,i) = -shf(j,i) / ( us(j,i) + 1E-30_wp )
    108108!
    109109!--          ts must be limited, because otherwise overflow may occur in case of
    110110!--          us=0 when computing rif further below
    111              IF ( ts(j,i) < -1.05E5 )  ts(j,i) = -1.0E5
    112              IF ( ts(j,i) >   1.0E5 )  ts(j,i) =  1.0E5
     111             IF ( ts(j,i) < -1.05E5_wp )  ts(j,i) = -1.0E5_wp
     112             IF ( ts(j,i) >   1.0E5_wp )  ts(j,i) =  1.0E5_wp
    113113          ENDDO
    114114       ENDDO
     
    131131             z_p = zu(k+1) - zw(k)
    132132
    133              IF ( rif(j,i) >= 0.0 )  THEN
     133             IF ( rif(j,i) >= 0.0_wp )  THEN
    134134!
    135135!--             Stable stratification
    136136                ts(j,i) = kappa * ( pt(k+1,j,i) - pt(k,j,i) ) / (           &
    137137                                  LOG( z_p / z0h(j,i) ) +                   &
    138                                   5.0 * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &
     138                                  5.0_wp * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &
    139139                                                                )
    140140             ELSE
    141141!
    142142!--             Unstable stratification
    143                 a = SQRT( 1.0 - 16.0 * rif(j,i) )
    144                 b = SQRT( 1.0 - 16.0 * rif(j,i) * z0h(j,i) / z_p )
     143                a = SQRT( 1.0_wp - 16.0_wp * rif(j,i) )
     144                b = SQRT( 1.0_wp - 16.0_wp * rif(j,i) * z0h(j,i) / z_p )
    145145
    146146                ts(j,i) = kappa * ( pt(k+1,j,i) - pt(k,j,i) ) /  (          &
    147147                          LOG( z_p / z0h(j,i) ) -                           &
    148                           2.0 * LOG( ( 1.0 + a ) / ( 1.0 + b ) ) )
     148                          2.0_wp * LOG( ( 1.0_wp + a ) / ( 1.0_wp + b ) ) )
    149149             ENDIF
    150150
     
    163163             z_p = zu(k+1) - zw(k)
    164164             rif(j,i) = z_p * kappa * g * ts(j,i) / &
    165                         ( pt(k+1,j,i) * ( us(j,i)**2 + 1E-30 ) )
     165                        ( pt(k+1,j,i) * ( us(j,i)**2 + 1E-30_wp ) )
    166166!
    167167!--          Limit the value range of the Richardson numbers.
     
    182182             z_p = zu(k+1) - zw(k)
    183183             rif(j,i) = z_p * kappa * g *                            &
    184                         ( ts(j,i) + 0.61 * pt(k+1,j,i) * qs(j,i) ) / &
    185                         ( vpt(k+1,j,i) * ( us(j,i)**2 + 1E-30 ) )
     184                        ( ts(j,i) + 0.61_wp * pt(k+1,j,i) * qs(j,i) ) / &
     185                        ( vpt(k+1,j,i) * ( us(j,i)**2 + 1E-30_wp ) )
    186186!
    187187!--          Limit the value range of the Richardson numbers.
     
    209209!--       Compute the absolute value of the horizontal velocity
    210210!--       (relative to the surface)
    211           uv_total = SQRT( ( 0.5 * ( u(k+1,j,i) + u(k+1,j,i+1)        &
    212                                    - u(k,j,i)   - u(k,j,i+1) ) )**2 + &
    213                            ( 0.5 * ( v(k+1,j,i) + v(k+1,j+1,i)        &
    214                                    - v(k,j,i)   - v(k,j+1,i) ) )**2 )   
    215 
    216 
    217           IF ( rif(j,i) >= 0.0 )  THEN
     211          uv_total = SQRT( ( 0.5_wp * ( u(k+1,j,i) + u(k+1,j,i+1)        &
     212                                      - u(k,j,i)   - u(k,j,i+1) ) )**2 + &
     213                           ( 0.5_wp * ( v(k+1,j,i) + v(k+1,j+1,i)        &
     214                                      - v(k,j,i)   - v(k,j+1,i) ) )**2 )   
     215
     216
     217          IF ( rif(j,i) >= 0.0_wp )  THEN
    218218!
    219219!--          Stable stratification
    220220             us(j,i) = kappa * uv_total / (                                &
    221221                                  LOG( z_p / z0(j,i) ) +                   &
    222                                   5.0 * rif(j,i) * ( z_p - z0(j,i) ) / z_p &
     222                                  5.0_wp * rif(j,i) * ( z_p - z0(j,i) ) / z_p &
    223223                                          )
    224224          ELSE
    225225!
    226226!--          Unstable stratification
    227              a = SQRT( SQRT( 1.0 - 16.0 * rif(j,i) ) )
    228              b = SQRT( SQRT( 1.0 - 16.0 * rif(j,i) / z_p * z0(j,i) ) )
     227             a = SQRT( SQRT( 1.0_wp - 16.0_wp * rif(j,i) ) )
     228             b = SQRT( SQRT( 1.0_wp - 16.0_wp * rif(j,i) / z_p * z0(j,i) ) )
    229229
    230230             us(j,i) = kappa * uv_total / (                                  &
    231231                       LOG( z_p / z0(j,i) ) -                                &
    232                        LOG( ( 1.0 + a )**2 * ( 1.0 + a**2 ) / (              &
    233                             ( 1.0 + b )**2 * ( 1.0 + b**2 )   ) ) +          &
    234                             2.0 * ( ATAN( a ) - ATAN( b ) )                  &
     232                       LOG( ( 1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / (              &
     233                            ( 1.0_wp + b )**2 * ( 1.0_wp + b**2 )   ) ) +          &
     234                               2.0_wp * ( ATAN( a ) - ATAN( b ) )                  &
    235235                                           )
    236236          ENDIF
     
    258258!
    259259!--       Compute Richardson-flux number for this point
    260           rifm = 0.5 * ( rif(j,i-1) + rif(j,i) )
    261           IF ( rifm >= 0.0 )  THEN
     260          rifm = 0.5_wp * ( rif(j,i-1) + rif(j,i) )
     261          IF ( rifm >= 0.0_wp )  THEN
    262262!
    263263!--          Stable stratification
    264264             usws(j,i) = kappa * ( u(k+1,j,i) - u(k,j,i) )/ (              &
    265265                                     LOG( z_p / z0(j,i) ) +               &
    266                                      5.0 * rifm * ( z_p - z0(j,i) ) / z_p &
    267                                               )
     266                                     5.0_wp * rifm * ( z_p - z0(j,i) ) / z_p &
     267                                                            )
    268268          ELSE
    269269!
    270270!--          Unstable stratification
    271              a = SQRT( SQRT( 1.0 - 16.0 * rifm ) )
    272              b = SQRT( SQRT( 1.0 - 16.0 * rifm / z_p * z0(j,i) ) )
     271             a = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm ) )
     272             b = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm / z_p * z0(j,i) ) )
    273273
    274274             usws(j,i) = kappa * ( u(k+1,j,i) - u(k,j,i) ) / (            &
    275275                         LOG( z_p / z0(j,i) ) -                           &
    276                          LOG( (1.0 + a )**2 * ( 1.0 + a**2 ) / (          &
    277                               (1.0 + b )**2 * ( 1.0 + b**2 )   ) ) +      &
    278                               2.0 * ( ATAN( a ) - ATAN( b ) )             &
     276                         LOG( (1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / (          &
     277                              (1.0_wp + b )**2 * ( 1.0_wp + b**2 )   ) ) +      &
     278                                 2.0_wp * ( ATAN( a ) - ATAN( b ) )             &
    279279                                                 )
    280280          ENDIF
    281           usws(j,i) = -usws(j,i) * 0.5 * ( us(j,i-1) + us(j,i) )
     281          usws(j,i) = -usws(j,i) * 0.5_wp * ( us(j,i-1) + us(j,i) )
    282282       ENDDO
    283283    ENDDO
     
    296296!
    297297!--       Compute Richardson-flux number for this point
    298           rifm = 0.5 * ( rif(j-1,i) + rif(j,i) )
    299           IF ( rifm >= 0.0 )  THEN
     298          rifm = 0.5_wp * ( rif(j-1,i) + rif(j,i) )
     299          IF ( rifm >= 0.0_wp )  THEN
    300300!
    301301!--          Stable stratification
    302302             vsws(j,i) = kappa * ( v(k+1,j,i) -  v(k,j,i) ) / (           &
    303303                                     LOG( z_p / z0(j,i) ) +               &
    304                                      5.0 * rifm * ( z_p - z0(j,i) ) / z_p &
    305                                               )
     304                                     5.0_wp * rifm * ( z_p - z0(j,i) ) / z_p &
     305                                                              )
    306306          ELSE
    307307!
    308308!--          Unstable stratification
    309              a = SQRT( SQRT( 1.0 - 16.0 * rifm ) )
    310              b = SQRT( SQRT( 1.0 - 16.0 * rifm / z_p * z0(j,i) ) )
     309             a = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm ) )
     310             b = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm / z_p * z0(j,i) ) )
    311311
    312312             vsws(j,i) = kappa * ( v(k+1,j,i) - v(k,j,i) ) / (            &
    313313                         LOG( z_p / z0(j,i) ) -                           &
    314                          LOG( (1.0 + a )**2 * ( 1.0 + a**2 ) / (          &
    315                               (1.0 + b )**2 * ( 1.0 + b**2 )   ) ) +      &
    316                               2.0 * ( ATAN( a ) - ATAN( b ) )             &
     314                         LOG( (1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / (          &
     315                              (1.0_wp + b )**2 * ( 1.0_wp + b**2 )   ) ) +      &
     316                                 2.0_wp * ( ATAN( a ) - ATAN( b ) )             &
    317317                                                 )
    318318          ENDIF
    319           vsws(j,i) = -vsws(j,i) * 0.5 * ( us(j-1,i) + us(j,i) )
     319          vsws(j,i) = -vsws(j,i) * 0.5_wp * ( us(j-1,i) + us(j,i) )
    320320       ENDDO
    321321    ENDDO
     
    331331          DO  i = nxlg, nxrg
    332332             DO  j = nysg, nyng
    333                 qs(j,i) = -qsws(j,i) / ( us(j,i) + 1E-30 )
     333                qs(j,i) = -qsws(j,i) / ( us(j,i) + 1E-30_wp )
    334334             ENDDO
    335335          ENDDO
     
    355355!--             in case of precursor runs)
    356356                IF ( coupled_run )  THEN
    357                    e_q = 6.1 * &
    358                         EXP( 0.07 * ( MIN(pt(0,j,i),pt(1,j,i)) - 273.15 ) )
    359                    q(k,j,i) = 0.622 * e_q / ( surface_pressure - e_q )
     357                   e_q = 6.1_wp * &
     358                           EXP( 0.07_wp * ( MIN(pt(0,j,i),pt(1,j,i)) - 273.15_wp ) )
     359                   q(k,j,i) = 0.622_wp * e_q / ( surface_pressure - e_q )
    360360                ENDIF
    361                 IF ( rif(j,i) >= 0.0 )  THEN
     361                IF ( rif(j,i) >= 0.0_wp )  THEN
    362362!
    363363!--                Stable stratification
    364364                   qs(j,i) = kappa * ( q(k+1,j,i) - q(k,j,i) ) / (          &
    365365                                  LOG( z_p / z0h(j,i) ) +                   &
    366                                   5.0 * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &
     366                                  5.0_wp * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &
    367367                                                                 )
    368368                ELSE
    369369!
    370370!--                Unstable stratification
    371                    a = SQRT( 1.0 - 16.0 * rif(j,i) )
    372                    b = SQRT( 1.0 - 16.0 * rif(j,i) * z0h(j,i) / z_p )
     371                   a = SQRT( 1.0_wp - 16.0_wp * rif(j,i) )
     372                   b = SQRT( 1.0_wp - 16.0_wp * rif(j,i) * z0h(j,i) / z_p )
    373373 
    374374                   qs(j,i) = kappa * ( q(k+1,j,i) - q(k,j,i) ) /   (        &
    375375                             LOG( z_p / z0h(j,i) ) -                        &
    376                               2.0 * LOG( (1.0 + a ) / ( 1.0 + b ) ) )
     376                              2.0_wp * LOG( (1.0_wp + a ) / ( 1.0_wp + b ) ) )
    377377                ENDIF
    378378
     
    429429          !$acc loop independent
    430430          DO  j = nysg, nyng
    431              e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.1 )**2
     431             e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.1_wp )**2
    432432!
    433433!--          As a test: cm = 0.4
    434 !            e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.4 )**2
     434!            e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.4_wp )**2
    435435             e(nzb_s_inner(j,i),j,i)   = e(nzb_s_inner(j,i)+1,j,i)
    436436          ENDDO
Note: See TracChangeset for help on using the changeset viewer.