Changeset 4351 for palm/trunk/SOURCE


Ignore:
Timestamp:
Dec 19, 2019 10:43:19 AM (5 years ago)
Author:
oliver.maas
Message:

added output of wu, wv, wtheta and wq to enable covariance calculation according to temporal EC method

File:
1 edited

Legend:

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

    r4346 r4351  
    2525! -----------------
    2626! $Id$
     27! added output of wu, wv, wtheta and wq to enable covariance calculation
     28! according to temporal EC method
     29!
     30! 4346 2019-12-18 11:55:56Z motisi
    2731! Introduction of wall_flags_total_0, which currently sets bits based on static
    2832! topography information used in wall_flags_static_0
     
    8589        ONLY:  ddzu,                                                           &
    8690               pt,                                                             &
     91               q,                                                              &
    8792               u,                                                              &
    8893               v,                                                              &
     
    145150    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti     !< rotation(u,v,w) aka turbulence intensity
    146151    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti_av  !< avg. rotation(u,v,w) aka turbulence intensity   
    147     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu     !< uu
    148     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu_av  !< mean of uu   
    149     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv     !< vv
    150     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv_av  !< mean of vv   
    151     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww     !< ww
    152     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww_av  !< mean of ww
     152    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu         !< uu
     153    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu_av      !< mean of uu
     154    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv         !< vv
     155    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv_av      !< mean of vv
     156    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww         !< ww
     157    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww_av      !< mean of ww
     158    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wu         !< wu
     159    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wu_av      !< mean of wu
     160    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wv         !< wv
     161    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wv_av      !< mean of wv
     162    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wtheta     !< wtheta
     163    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wtheta_av  !< mean of wtheta
     164    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wq         !< wq
     165    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wq_av      !< mean of wq
    153166
    154167
     
    275288             ENDIF
    276289             ww_av = 0.0_wp
     290           
     291           CASE ( 'wu' )
     292             IF ( .NOT. ALLOCATED( wu_av ) )  THEN
     293                ALLOCATE( wu_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     294             ENDIF
     295             wu_av = 0.0_wp
     296           
     297           CASE ( 'wv' )
     298             IF ( .NOT. ALLOCATED( wv_av ) )  THEN
     299                ALLOCATE( wv_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     300             ENDIF
     301             wv_av = 0.0_wp
     302               
     303           CASE ( 'wtheta' )
     304             IF ( .NOT. ALLOCATED( wtheta_av ) )  THEN
     305                ALLOCATE( wtheta_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     306             ENDIF
     307             wtheta_av = 0.0_wp
     308               
     309           CASE ( 'wq' )
     310             IF ( .NOT. ALLOCATED( wq_av ) )  THEN
     311                ALLOCATE( wq_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     312             ENDIF
     313             wq_av = 0.0_wp
    277314             
    278315          CASE ( 'theta_2m*' )
     
    341378             ENDIF
    342379             
     380          CASE ( 'wu' )
     381             IF ( ALLOCATED( wu_av ) ) THEN
     382                DO  i = nxl, nxr
     383                   DO  j = nys, nyn
     384                      DO  k = nzb, nzt+1
     385                         wu_av(k,j,i) = wu_av(k,j,i) + wu(k,j,i)
     386                      ENDDO
     387                   ENDDO
     388                ENDDO
     389             ENDIF
     390             
     391          CASE ( 'wv' )
     392             IF ( ALLOCATED( wv_av ) ) THEN
     393                DO  i = nxl, nxr
     394                   DO  j = nys, nyn
     395                      DO  k = nzb, nzt+1
     396                         wv_av(k,j,i) = wv_av(k,j,i) + wv(k,j,i)
     397                      ENDDO
     398                   ENDDO
     399                ENDDO
     400             ENDIF
     401             
     402          CASE ( 'wtheta' )
     403             IF ( ALLOCATED( wtheta_av ) ) THEN
     404                DO  i = nxl, nxr
     405                   DO  j = nys, nyn
     406                      DO  k = nzb, nzt+1
     407                         wtheta_av(k,j,i) = wtheta_av(k,j,i) + wtheta(k,j,i)
     408                      ENDDO
     409                   ENDDO
     410                ENDDO
     411             ENDIF
     412                         
     413          CASE ( 'wq' )
     414             IF ( ALLOCATED( wq_av ) ) THEN
     415                DO  i = nxl, nxr
     416                   DO  j = nys, nyn
     417                      DO  k = nzb, nzt+1
     418                         wq_av(k,j,i) = wq_av(k,j,i) + wq(k,j,i)
     419                      ENDDO
     420                   ENDDO
     421                ENDDO
     422             ENDIF
     423             
    343424          CASE ( 'theta_2m*' )
    344425             IF ( ALLOCATED( pt_2m_av ) ) THEN
     
    407488                      DO  k = nzb, nzt+1
    408489                         ww_av(k,j,i) = ww_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     490                      ENDDO
     491                   ENDDO
     492                ENDDO
     493             ENDIF
     494
     495          CASE ( 'wu' )
     496             IF ( ALLOCATED( wu_av ) ) THEN
     497                DO  i = nxl, nxr
     498                   DO  j = nys, nyn
     499                      DO  k = nzb, nzt+1
     500                         wu_av(k,j,i) = wu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     501                      ENDDO
     502                   ENDDO
     503                ENDDO
     504             ENDIF
     505             
     506          CASE ( 'wv' )
     507             IF ( ALLOCATED( wv_av ) ) THEN
     508                DO  i = nxl, nxr
     509                   DO  j = nys, nyn
     510                      DO  k = nzb, nzt+1
     511                         wv_av(k,j,i) = wv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     512                      ENDDO
     513                   ENDDO
     514                ENDDO
     515             ENDIF
     516             
     517          CASE ( 'wtheta' )
     518             IF ( ALLOCATED( wtheta_av ) ) THEN
     519                DO  i = nxl, nxr
     520                   DO  j = nys, nyn
     521                      DO  k = nzb, nzt+1
     522                         wtheta_av(k,j,i) = wtheta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
     523                      ENDDO
     524                   ENDDO
     525                ENDDO
     526             ENDIF
     527             
     528          CASE ( 'wq' )
     529             IF ( ALLOCATED( wq_av ) ) THEN
     530                DO  i = nxl, nxr
     531                   DO  j = nys, nyn
     532                      DO  k = nzb, nzt+1
     533                         wq_av(k,j,i) = wq_av(k,j,i) / REAL( average_count_3d, KIND=wp )
    409534                      ENDDO
    410535                   ENDDO
     
    468593       CASE ( 'ww' )
    469594          unit = 'm2/s2'
     595         
     596       CASE ( 'wu' )
     597          unit = 'm2/s2'
     598             
     599       CASE ( 'wv' )
     600          unit = 'm2/s2'
     601             
     602       CASE ( 'wtheta' )
     603          unit = 'Km/s'
     604             
     605       CASE ( 'wq' )
     606          unit = 'm/s'
    470607!
    471608!--    Treat horizotal cross-section output quanatities
     
    512649!
    513650!--    s grid
    514        CASE ( 'ti', 'ti_xy', 'ti_xz', 'ti_yz'  )
     651       CASE ( 'ti', 'ti_xy', 'ti_xz', 'ti_yz',                                 &
     652              'wu', 'wu_xy', 'wu_xz', 'wu_yz',                                 &
     653              'wv', 'wv_xy', 'wv_xz', 'wv_yz',                                 &
     654              'wtheta', 'wtheta_xy', 'wtheta_xz', 'wtheta_yz',                 &
     655              'wq', 'wq_xy', 'wq_xz', 'wq_yz')
    515656
    516657          grid_x = 'x'
     
    650791          ENDIF
    651792          flag_nr = 3
     793         
     794          IF ( mode == 'xy' )  grid = 'zw'
     795         
     796       CASE ( 'wu_xy', 'wu_xz', 'wu_yz' )
     797          IF ( av == 0 )  THEN
     798             to_be_resorted => wu
     799          ELSE
     800             IF ( .NOT. ALLOCATED( wu_av ) ) THEN
     801                ALLOCATE( wu_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     802                wu_av = REAL( fill_value, KIND = wp )
     803             ENDIF
     804             to_be_resorted => wu_av
     805          ENDIF
     806          flag_nr = 0
     807         
     808          IF ( mode == 'xy' )  grid = 'zw'
     809         
     810       CASE ( 'wv_xy', 'wv_xz', 'wv_yz' )
     811          IF ( av == 0 )  THEN
     812             to_be_resorted => wv
     813          ELSE
     814             IF ( .NOT. ALLOCATED( wv_av ) ) THEN
     815                ALLOCATE( wv_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     816                wv_av = REAL( fill_value, KIND = wp )
     817             ENDIF
     818             to_be_resorted => wv_av
     819          ENDIF
     820          flag_nr = 0
     821         
     822          IF ( mode == 'xy' )  grid = 'zw'
     823         
     824       CASE ( 'wtheta_xy', 'wtheta_xz', 'wtheta_yz' )
     825          IF ( av == 0 )  THEN
     826             to_be_resorted => wtheta
     827          ELSE
     828             IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
     829                ALLOCATE( wtheta_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     830                wtheta_av = REAL( fill_value, KIND = wp )
     831             ENDIF
     832             to_be_resorted => wtheta_av
     833          ENDIF
     834          flag_nr = 0
     835         
     836          IF ( mode == 'xy' )  grid = 'zw'
     837         
     838       CASE ( 'wq_xy', 'wq_xz', 'wq_yz' )
     839          IF ( av == 0 )  THEN
     840             to_be_resorted => wq
     841          ELSE
     842             IF ( .NOT. ALLOCATED( wq_av ) ) THEN
     843                ALLOCATE( wq_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     844                wq_av = REAL( fill_value, KIND = wp )
     845             ENDIF
     846             to_be_resorted => wq_av
     847          ENDIF
     848          flag_nr = 0
    652849         
    653850          IF ( mode == 'xy' )  grid = 'zw'
     
    801998          flag_nr = 3
    802999
     1000       CASE ( 'wu' )
     1001          IF ( av == 0 )  THEN
     1002             to_be_resorted => wu
     1003          ELSE
     1004             IF ( .NOT. ALLOCATED( wu_av ) ) THEN
     1005                ALLOCATE( wu_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1006                wu_av = REAL( fill_value, KIND = wp )
     1007             ENDIF
     1008             to_be_resorted => wu_av
     1009          ENDIF
     1010          flag_nr = 0
     1011
     1012       CASE ( 'wv' )
     1013          IF ( av == 0 )  THEN
     1014             to_be_resorted => wv
     1015          ELSE
     1016             IF ( .NOT. ALLOCATED( wv_av ) ) THEN
     1017                ALLOCATE( wv_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1018                wv_av = REAL( fill_value, KIND = wp )
     1019             ENDIF
     1020             to_be_resorted => wv_av
     1021          ENDIF
     1022          flag_nr = 0
     1023         
     1024       CASE ( 'wtheta' )
     1025          IF ( av == 0 )  THEN
     1026             to_be_resorted => wtheta
     1027          ELSE
     1028             IF ( .NOT. ALLOCATED( wtheta_av ) ) THEN
     1029                ALLOCATE( wtheta_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1030                wtheta_av = REAL( fill_value, KIND = wp )
     1031             ENDIF
     1032             to_be_resorted => wtheta_av
     1033          ENDIF
     1034          flag_nr = 0
     1035         
     1036       CASE ( 'wq' )
     1037          IF ( av == 0 )  THEN
     1038             to_be_resorted => wq
     1039          ELSE
     1040             IF ( .NOT. ALLOCATED( wq_av ) ) THEN
     1041                ALLOCATE( wq_av(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1042                wq_av = REAL( fill_value, KIND = wp )
     1043             ENDIF
     1044             to_be_resorted => wq_av
     1045          ENDIF
     1046          flag_nr = 0
     1047
    8031048       CASE DEFAULT
    8041049          found = .FALSE.
     
    8991144          grid = 'w'
    9001145          flag_nr = 3
     1146   
     1147       CASE ( 'wu' )
     1148          IF ( av == 0 )  THEN
     1149             to_be_resorted => wu
     1150          ELSE
     1151             to_be_resorted => wu_av
     1152          ENDIF
     1153          grid = 's'
     1154          flag_nr = 0
     1155   
     1156       CASE ( 'wv' )
     1157          IF ( av == 0 )  THEN
     1158             to_be_resorted => wv
     1159          ELSE
     1160             to_be_resorted => wv_av
     1161          ENDIF
     1162          grid = 's'
     1163          flag_nr = 0
     1164             
     1165       CASE ( 'wtheta' )
     1166          IF ( av == 0 )  THEN
     1167             to_be_resorted => wtheta
     1168          ELSE
     1169             to_be_resorted => wtheta_av
     1170          ENDIF
     1171          grid = 's'
     1172          flag_nr = 0
     1173               
     1174       CASE ( 'wq' )
     1175          IF ( av == 0 )  THEN
     1176             to_be_resorted => wq
     1177          ELSE
     1178             to_be_resorted => wq_av
     1179          ENDIF
     1180          grid = 's'
     1181          flag_nr = 0
    9011182
    9021183       CASE DEFAULT
     
    10001281                ALLOCATE( ww(nzb:nzt+1,nys:nyn,nxl:nxr) )
    10011282                ww = 0.0_wp
     1283             ENDIF
     1284!
     1285!--       Allocate array for wu
     1286          CASE ( 'wu' )
     1287             IF ( .NOT. ALLOCATED( wu ) )  THEN
     1288                ALLOCATE( wu(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1289                wu = 0.0_wp
     1290             ENDIF
     1291!
     1292!--       Allocate array for wv
     1293          CASE ( 'wv' )
     1294             IF ( .NOT. ALLOCATED( wv ) )  THEN
     1295                ALLOCATE( wv(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1296                wv = 0.0_wp
     1297             ENDIF
     1298!
     1299!--       Allocate array for wtheta
     1300          CASE ( 'wtheta' )
     1301             IF ( .NOT. ALLOCATED( wtheta ) )  THEN
     1302                ALLOCATE( wtheta(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1303                wtheta = 0.0_wp
     1304             ENDIF
     1305!
     1306!--       Allocate array for wq
     1307          CASE ( 'wq' )
     1308             IF ( .NOT. ALLOCATED( wq ) )  THEN
     1309                ALLOCATE( wq(nzb:nzt+1,nys:nyn,nxl:nxr) )
     1310                wq = 0.0_wp
    10021311             ENDIF
    10031312!
     
    10861395                   DO  k = nzb+1, nzt
    10871396                      uu(k,j,i) = u(k,j,i) * u(k,j,i)                          &
    1088                        * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1) )
     1397                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1398                       BTEST( wall_flags_total_0(k,j,i), 1) )
    10891399                   ENDDO
    10901400                ENDDO
     
    10971407                   DO  k = nzb+1, nzt
    10981408                      vv(k,j,i) = v(k,j,i) * v(k,j,i)                          &
    1099                        * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2) )
     1409                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1410                       BTEST( wall_flags_total_0(k,j,i), 2) )
    11001411                   ENDDO
    11011412                ENDDO
     
    11081419                   DO  k = nzb+1, nzt-1
    11091420                      ww(k,j,i) = w(k,j,i) * w(k,j,i)                          &
    1110                        * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3) )
     1421                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1422                       BTEST( wall_flags_total_0(k,j,i), 3) )
     1423                   ENDDO
     1424                ENDDO
     1425             ENDDO
     1426!
     1427!--       wu
     1428          CASE ( 'wu' )
     1429             DO  i = nxl, nxr
     1430                DO  j = nys, nyn
     1431                   DO  k = nzb+1, nzt-1
     1432                      wu(k,j,i) = w(k,j,i) * u(k,j,i)                          &
     1433                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1434                       BTEST( wall_flags_total_0(k,j,i), 0) )
     1435                   ENDDO
     1436                ENDDO
     1437             ENDDO
     1438!
     1439!--       wv
     1440          CASE ( 'wv' )
     1441             DO  i = nxl, nxr
     1442                DO  j = nys, nyn
     1443                   DO  k = nzb+1, nzt-1
     1444                      wv(k,j,i) = w(k,j,i) * v(k,j,i)                          &
     1445                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1446                       BTEST( wall_flags_total_0(k,j,i), 0) )
     1447                   ENDDO
     1448                ENDDO
     1449             ENDDO
     1450!
     1451!--       wtheta
     1452          CASE ( 'wtheta' )
     1453             DO  i = nxl, nxr
     1454                DO  j = nys, nyn
     1455                   DO  k = nzb+1, nzt-1
     1456                      wtheta(k,j,i) = w(k,j,i) * pt(k,j,i)                     &
     1457                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1458                       BTEST( wall_flags_total_0(k,j,i), 0) )
     1459                   ENDDO
     1460                ENDDO
     1461             ENDDO
     1462!
     1463!--       wq
     1464          CASE ( 'wq' )
     1465             DO  i = nxl, nxr
     1466                DO  j = nys, nyn
     1467                   DO  k = nzb+1, nzt-1
     1468                      wq(k,j,i) = w(k,j,i) * q(k,j,i)                          &
     1469                       * MERGE( 1.0_wp, 0.0_wp,                                &
     1470                       BTEST( wall_flags_total_0(k,j,i), 0) )
    11111471                   ENDDO
    11121472                ENDDO
     
    14921852    ENDIF
    14931853
     1854    IF ( ALLOCATED( wu_av ) )  THEN
     1855       CALL wrd_write_string( 'wu_av' )
     1856       WRITE ( 14 )  wu_av
     1857    ENDIF
     1858
     1859    IF ( ALLOCATED( wv_av ) )  THEN
     1860       CALL wrd_write_string( 'wv_av' )
     1861       WRITE ( 14 )  wv_av
     1862    ENDIF
     1863   
     1864    IF ( ALLOCATED( wtheta_av ) )  THEN
     1865       CALL wrd_write_string( 'wtheta_av' )
     1866       WRITE ( 14 )  wtheta_av
     1867    ENDIF
     1868   
     1869    IF ( ALLOCATED( wq_av ) )  THEN
     1870       CALL wrd_write_string( 'wq_av' )
     1871       WRITE ( 14 )  wq_av
     1872    ENDIF
    14941873
    14951874 END SUBROUTINE doq_wrd_local
Note: See TracChangeset for help on using the changeset viewer.