Changeset 4196 for palm


Ignore:
Timestamp:
Aug 29, 2019 11:02:06 AM (5 years ago)
Author:
gronemeier
Message:

Consider rotation of model domain for claculation of Coriolis force:

  • check_parameters: Overwrite rotation_angle from namelist by value from static driver;
  • coriolis: Consider rotation of model domain;
  • header: Write information about rotation angle;
  • modules: Added rotation_angle;
  • netcdf_interface_mod: replaced rotation angle from input-netCDF file by namelist parameter 'rotation_angle';
  • ocean_mod: Consider rotation of model domain for calculating the Stokes drift;
  • parin: added rotation_angle to initialization_parameters;
Location:
palm/trunk/SOURCE
Files:
7 edited

Legend:

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

    r4183 r4196  
    2525! -----------------
    2626! 4172 2019-08-20 11:55:33Z oliver.maas
     27! Overwrite rotation_angle from namelist by value from static driver
     28!
     29! 11:55:33Z oliver.maas
    2730! removed conversion from recycle_absolute_quantities to raq, added check and
    2831! error message for correct input of recycling_method_for_thermodynamic_quantities
     
    10851088
    10861089!
    1087 !-- Overwrite latitude if necessary and compute Coriolis parameter.
     1090!-- Overwrite parameters from namelist if necessary and compute Coriolis parameter.
    10881091!-- @todo - move initialization of f and fs to coriolis_mod.
    10891092    IF ( input_pids_static )  THEN
    1090        latitude  = init_model%latitude
    1091        longitude = init_model%longitude
     1093       latitude       = init_model%latitude
     1094       longitude      = init_model%longitude
     1095       rotation_angle = init_model%rotation_angle
    10921096    ENDIF
    10931097
  • palm/trunk/SOURCE/coriolis.f90

    r4182 r4196  
    2525! -----------------
    2626! $Id$
     27! Consider rotation of model domain
     28!
     29! 4182 2019-08-22 15:20:23Z scharf
    2730! Corrected "Former revisions" section
    2831!
     
    6568       USE arrays_3d,                                                          &
    6669           ONLY:  tend, u, ug, v, vg, w
    67            
     70
     71       USE basic_constants_and_equations_mod,                                  &
     72           ONLY:  pi
     73
    6874       USE control_parameters,                                                 &
    69            ONLY:  f, fs, message_string
     75           ONLY:  f, fs, message_string, rotation_angle
    7076           
    7177       USE indices,                                                            &
     
    7682       IMPLICIT NONE
    7783
    78        INTEGER(iwp) ::  component  !<
    79        INTEGER(iwp) ::  i          !< running index x direction
    80        INTEGER(iwp) ::  j          !< running index y direction
    81        INTEGER(iwp) ::  k          !< running index z direction
    82 
     84       INTEGER(iwp) ::  component      !< component of momentum equation
     85       INTEGER(iwp) ::  i              !< running index x direction
     86       INTEGER(iwp) ::  j              !< running index y direction
     87       INTEGER(iwp) ::  k              !< running index z direction
     88
     89       REAL(wp)     ::  cos_rot_angle  !< cosine of model rotation angle
    8390       REAL(wp)     ::  flag           !< flag to mask topography
     91       REAL(wp)     ::  sin_rot_angle  !< sine of model rotation angle
     92
     93!
     94!--    Precalculate cosine and sine of rotation angle
     95       cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     96       sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    8497
    8598!
     
    99112!
    100113!--                   Predetermine flag to mask topography
    101                       flag = MERGE( 1.0_wp, 0.0_wp,                            &
    102                                     BTEST( wall_flags_0(k,j,i), 1 ) )
    103 
    104                       tend(k,j,i) = tend(k,j,i) + f  *    ( 0.25_wp *          &
    105                                    ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) +    &
    106                                      v(k,j+1,i) ) - vg(k) ) * flag             &
    107                                                 - fs *    ( 0.25_wp *          &
    108                                    ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) +  &
    109                                      w(k,j,i)   )                              &
    110                                                           ) * flag
     114                      flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 1 ) )
     115
     116                      tend(k,j,i) = tend(k,j,i) + flag *                                           &
     117                            ( f                                                                    &
     118                              * ( 0.25_wp * ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) + v(k,j+1,i) )  &
     119                                - vg(k) )                                                          &
     120                            - fs * cos_rot_angle                                                   &
     121                              * 0.25_wp * ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) + w(k,j,i) )    &
     122                            )
    111123                   ENDDO
    112124                ENDDO
     
    118130             !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k, flag) &
    119131             !$ACC PRESENT(wall_flags_0) &
    120              !$ACC PRESENT(u, ug) &
     132             !$ACC PRESENT(u, w, ug) &
    121133             !$ACC PRESENT(tend)
    122134             DO  i = nxl, nxr
     
    125137!
    126138!--                   Predetermine flag to mask topography
    127                       flag = MERGE( 1.0_wp, 0.0_wp,                            &
    128                                     BTEST( wall_flags_0(k,j,i), 2 ) )
    129 
    130                       tend(k,j,i) = tend(k,j,i) - f *     ( 0.25_wp *          &
    131                                    ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) +    &
    132                                      u(k,j,i+1) ) - ug(k) ) * flag
     139                      flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 2 ) )
     140
     141                      tend(k,j,i) = tend(k,j,i) - flag *                                           &
     142                            ( f                                                                    &
     143                              * ( 0.25_wp * ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) + u(k,j,i+1) )  &
     144                                - ug(k) )                                                          &
     145                            + fs * sin_rot_angle                                                   &
     146                              * 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) )    &
     147                            )
    133148                   ENDDO
    134149                ENDDO
     
    140155             !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k, flag) &
    141156             !$ACC PRESENT(wall_flags_0) &
    142              !$ACC PRESENT(u) &
     157             !$ACC PRESENT(u, v) &
    143158             !$ACC PRESENT(tend)
    144159             DO  i = nxl, nxr
     
    147162!
    148163!--                   Predetermine flag to mask topography
    149                       flag = MERGE( 1.0_wp, 0.0_wp,                            &
    150                                     BTEST( wall_flags_0(k,j,i), 3 ) )
    151 
    152                       tend(k,j,i) = tend(k,j,i) + fs * 0.25_wp *               &
    153                                    ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) +      &
    154                                      u(k+1,j,i+1) ) * flag
     164                      flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 3 ) )
     165
     166                      tend(k,j,i) = tend(k,j,i)                                                 &
     167                                  + fs * 0.25_wp * flag                                         &
     168                                    * ( cos_rot_angle                                           &
     169                                        * ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) + u(k+1,j,i+1) ) &
     170                                      + sin_rot_angle                                           &
     171                                        * ( v(k,j,i) + v(k+1,j,i) + v(k,j+1,i) + v(k+1,j+1,i) ) &
     172                                      )
    155173                   ENDDO
    156174                ENDDO
     
    176194       USE arrays_3d,                                                          &
    177195           ONLY:  tend, u, ug, v, vg, w
    178            
     196
     197       USE basic_constants_and_equations_mod,                                  &
     198           ONLY:  pi
     199
    179200       USE control_parameters,                                                 &
    180            ONLY:  f, fs, message_string
     201           ONLY:  f, fs, message_string, rotation_angle
    181202           
    182203       USE indices,                                                            &
     
    184205           
    185206       USE kinds
    186        
     207
    187208       IMPLICIT NONE
    188209
    189        INTEGER(iwp) ::  component  !<
     210       INTEGER(iwp) ::  component  !< component of momentum equation
    190211       INTEGER(iwp) ::  i          !< running index x direction
    191212       INTEGER(iwp) ::  j          !< running index y direction
    192213       INTEGER(iwp) ::  k          !< running index z direction
    193214
    194        REAL(wp)     ::  flag       !< flag to mask topography
     215       REAL(wp)     ::  cos_rot_angle  !< cosine of model rotation angle
     216       REAL(wp)     ::  flag           !< flag to mask topography
     217       REAL(wp)     ::  sin_rot_angle  !< sine of model rotation angle
     218
     219!
     220!--    Precalculate cosine and sine of rotation angle
     221       cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     222       sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    195223
    196224!
     
    206234                flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 1 ) )
    207235
    208                 tend(k,j,i) = tend(k,j,i) + f  *     ( 0.25_wp *               &
    209                                 ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) +       &
    210                                   v(k,j+1,i) ) - vg(k)                         &
    211                                                      ) * flag                  &
    212                                           - fs *     ( 0.25_wp *               &
    213                                 ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) +     &
    214                                   w(k,j,i)   )       ) * flag
     236                tend(k,j,i) = tend(k,j,i) + flag *                                                 &
     237                            ( f                                                                    &
     238                              * ( 0.25_wp * ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) + v(k,j+1,i) )  &
     239                                - vg(k) )                                                          &
     240                            - fs * cos_rot_angle                                                   &
     241                              * 0.25_wp * ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) + w(k,j,i) )    &
     242                            )
    215243             ENDDO
    216244
     
    223251                flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 2 ) )
    224252
    225                 tend(k,j,i) = tend(k,j,i) - f *        ( 0.25_wp *             &
    226                                 ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) +       &
    227                                   u(k,j,i+1) ) - ug(k) ) * flag
     253                tend(k,j,i) = tend(k,j,i) - flag *                                                 &
     254                            ( f                                                                    &
     255                              * ( 0.25_wp * ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) + u(k,j,i+1) )  &
     256                                - ug(k) )                                                          &
     257                            + fs * sin_rot_angle                                                   &
     258                              * 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) )    &
     259                            )
    228260             ENDDO
    229261
     
    236268                flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 3 ) )
    237269
    238                 tend(k,j,i) = tend(k,j,i) + fs * 0.25_wp *                     &
    239                                 ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) +         &
    240                                   u(k+1,j,i+1) ) * flag
     270                tend(k,j,i) = tend(k,j,i)                                                 &
     271                            + fs * 0.25_wp * flag                                         &
     272                              * ( cos_rot_angle                                           &
     273                                  * ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) + u(k+1,j,i+1) ) &
     274                                + sin_rot_angle                                           &
     275                                  * ( v(k,j,i) + v(k+1,j,i) + v(k,j+1,i) + v(k+1,j+1,i) ) &
     276                                )
    241277             ENDDO
    242278
  • palm/trunk/SOURCE/header.f90

    r4182 r4196  
    2525! -----------------
    2626! $Id$
     27! Write information about rotation angle
     28!
     29! 4182 2019-08-22 15:20:23Z scharf
    2730! Corrected "Former revisions" section
    2831!
     
    14661469!
    14671470!-- Geostrophic parameters
    1468     WRITE ( io, 410 )  latitude, longitude, omega, f, fs
     1471    WRITE ( io, 410 )  latitude, longitude, rotation_angle, omega, f, fs
    14691472
    14701473!
     
    17961799400 FORMAT (//' Physical quantities:'/ &
    17971800              ' -------------------'/)
    1798 410 FORMAT ('    Geograph. latitude  :   latitude  = ',F4.1,' degr'/   &
    1799             '    Geograph. longitude :   longitude = ',F4.1,' degr'/   &
     1801410 FORMAT ('    Geograph. latitude  :   latitude  = ',F5.1,' degr'/   &
     1802            '    Geograph. longitude :   longitude = ',F5.1,' degr'/   &
     1803            '    Rotation angle      :   rotation_angle = ',F5.1,' degr'/   &
    18001804            '    Angular velocity    :   omega  =',E10.3,' rad/s'/  &
    18011805            '    Coriolis parameter  :   f      = ',F9.6,' 1/s'/    &
  • palm/trunk/SOURCE/modules.f90

    r4184 r4196  
    2525! -----------------
    2626! $Id$
     27! Added rotation_angle
     28!
     29! 4184 2019-08-23 08:07:40Z oliver.maas
    2730! changed allocated length of recycling_method_for_thermodynamic_quantities
    2831! from 20 to 80 characters
     
    867870    REAL(wp) ::  rho_reference                                 !< reference state of density
    868871    REAL(wp) ::  rho_surface                                   !< surface value of density
     872    REAL(wp) ::  rotation_angle = 0.0_wp                       !< angle between real North and model North (clockwise)
    869873    REAL(wp) ::  roughness_length = 0.1_wp                     !< namelist parameter
    870874    REAL(wp) ::  simulated_time = 0.0_wp                       !< elapsed simulated time
  • palm/trunk/SOURCE/netcdf_interface_mod.f90

    r4182 r4196  
    2525! -----------------
    2626! $Id$
     27! replaced rotation angle from input-netCDF file
     28! by namelist parameter 'rotation_angle'
     29!
     30! 4182 2019-08-22 15:20:23Z scharf
    2731! Corrected "Former revisions" section
    2832!
     
    114118        ONLY:  biometeorology, fl_max,                                         &
    115119               max_masks, multi_agent_system_end,                              &
    116                multi_agent_system_start, var_fl_max, varnamelength
     120               multi_agent_system_start,                                       &
     121               rotation_angle,                                                 &
     122               var_fl_max, varnamelength
    117123    USE kinds
    118124#if defined( __netcdf )
     
    544550    LOGICAL, SAVE ::  init_netcdf = .FALSE.                  !<
    545551
    546     REAL(wp) ::  cos_ra                                      !< cosine of rotation_angle
     552    REAL(wp) ::  cos_rot_angle                               !< cosine of rotation_angle
    547553    REAL(wp) ::  eutm                                        !< easting (UTM)
    548554    REAL(wp) ::  nutm                                        !< northing (UTM)
    549555    REAL(wp) ::  shift_x                                     !< shift of x coordinate
    550556    REAL(wp) ::  shift_y                                     !< shift of y coordinate
    551     REAL(wp) ::  sin_ra                                      !< sine of rotation_angle
     557    REAL(wp) ::  sin_rot_angle                               !< sine of rotation_angle
    552558
    553559    REAL(wp), DIMENSION(1) ::  last_time_coordinate          !< last time value in file
     
    10261032!
    10271033!--       Write UTM coordinates
    1028           IF ( init_model%rotation_angle == 0.0_wp )  THEN
     1034          IF ( rotation_angle == 0.0_wp )  THEN
    10291035!
    10301036!--          1D in case of no rotation
    1031              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
     1037             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
    10321038!
    10331039!--          x coordinates
     
    10481054                ENDIF
    10491055
    1050                 netcdf_data = init_model%origin_x + cos_ra                     &
     1056                netcdf_data = init_model%origin_x + cos_rot_angle              &
    10511057                       * ( mask_i_global(mid,:mask_size(mid,1)) + shift_x ) * dx
    10521058
     
    10771083                ENDIF
    10781084
    1079                 netcdf_data = init_model%origin_y + cos_ra                     &
     1085                netcdf_data = init_model%origin_y + cos_rot_angle              &
    10801086                       * ( mask_j_global(mid,:mask_size(mid,2)) + shift_y ) * dy
    10811087
     
    10931099!--          2D in case of rotation
    10941100             ALLOCATE( netcdf_data_2d(mask_size(mid,1),mask_size(mid,2)) )
    1095              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    1096              sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     1101             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     1102             sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    10971103
    10981104             DO  k = 0, 2
     
    11131119               DO  j = 1, mask_size(mid,2)
    11141120                  DO  i = 1, mask_size(mid,1)
    1115                      netcdf_data_2d(i,j) = init_model%origin_x                 &
    1116                            + cos_ra * ( mask_i_global(mid,i) + shift_x ) * dx  &
    1117                            + sin_ra * ( mask_j_global(mid,j) + shift_y ) * dy
     1121                     netcdf_data_2d(i,j) = init_model%origin_x                        &
     1122                           + cos_rot_angle * ( mask_i_global(mid,i) + shift_x ) * dx  &
     1123                           + sin_rot_angle * ( mask_j_global(mid,j) + shift_y ) * dy
    11181124                  ENDDO
    11191125               ENDDO
     
    11281134               DO  j = 1, mask_size(mid,2)
    11291135                  DO  i = 1, mask_size(mid,1)
    1130                      netcdf_data_2d(i,j) = init_model%origin_y                 &
    1131                            - sin_ra * ( mask_i_global(mid,i) + shift_x ) * dx  &
    1132                            + cos_ra * ( mask_j_global(mid,j) + shift_y ) * dy
     1136                     netcdf_data_2d(i,j) = init_model%origin_y                        &
     1137                           - sin_rot_angle * ( mask_i_global(mid,i) + shift_x ) * dx  &
     1138                           + cos_rot_angle * ( mask_j_global(mid,j) + shift_y ) * dy
    11331139                  ENDDO
    11341140               ENDDO
     
    11481154          ALLOCATE( lat(mask_size(mid,1),mask_size(mid,2)) )
    11491155          ALLOCATE( lon(mask_size(mid,1),mask_size(mid,2)) )
    1150           cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    1151           sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     1156          cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     1157          sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    11521158
    11531159          DO  k = 0, 2
     
    11681174             DO  j = 1, mask_size(mid,2)
    11691175                DO  i = 1, mask_size(mid,1)
    1170                    eutm = init_model%origin_x                               &
    1171                         + cos_ra * ( mask_i_global(mid,i) + shift_x ) * dx  &
    1172                         + sin_ra * ( mask_j_global(mid,j) + shift_y ) * dy
    1173                    nutm = init_model%origin_y                               &
    1174                         - sin_ra * ( mask_i_global(mid,i) + shift_x ) * dx  &
    1175                         + cos_ra * ( mask_j_global(mid,j) + shift_y ) * dy
     1176                   eutm = init_model%origin_x                                      &
     1177                        + cos_rot_angle * ( mask_i_global(mid,i) + shift_x ) * dx  &
     1178                        + sin_rot_angle * ( mask_j_global(mid,j) + shift_y ) * dy
     1179                   nutm = init_model%origin_y                                      &
     1180                        - sin_rot_angle * ( mask_i_global(mid,i) + shift_x ) * dx  &
     1181                        + cos_rot_angle * ( mask_j_global(mid,j) + shift_y ) * dy
    11761182
    11771183                   CALL  convert_utm_to_geographic( crs_list,          &
     
    18451851!
    18461852!--          Write UTM coordinates
    1847              IF ( init_model%rotation_angle == 0.0_wp )  THEN
     1853             IF ( rotation_angle == 0.0_wp )  THEN
    18481854!
    18491855!--             1D in case of no rotation
    1850                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
     1856                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
    18511857!
    18521858!--             x coordinates
     
    18681874
    18691875                   DO  i = 0, nx
    1870                      netcdf_data(i) = init_model%origin_x            &
    1871                                     + cos_ra * ( i + shift_x ) * dx
     1876                     netcdf_data(i) = init_model%origin_x                      &
     1877                                    + cos_rot_angle * ( i + shift_x ) * dx
    18721878                   ENDDO
    18731879
    18741880                   nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_eutm_3d(k,av),&
    1875                                            netcdf_data, start = (/ 1 /),   &
     1881                                           netcdf_data, start = (/ 1 /),       &
    18761882                                           count = (/ nx+1 /) )
    18771883                   CALL netcdf_handle_error( 'netcdf_define_header', 555 )
     
    18981904
    18991905                   DO  j = 0, ny
    1900                       netcdf_data(j) = init_model%origin_y            &
    1901                                      + cos_ra * ( j + shift_y ) * dy
     1906                      netcdf_data(j) = init_model%origin_y                     &
     1907                                     + cos_rot_angle * ( j + shift_y ) * dy
    19021908                   ENDDO
    19031909
    19041910                   nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_nutm_3d(k,av),&
    1905                                            netcdf_data, start = (/ 1 /),   &
     1911                                           netcdf_data, start = (/ 1 /),       &
    19061912                                           count = (/ ny+1 /) )
    19071913                   CALL netcdf_handle_error( 'netcdf_define_header', 556 )
     
    19141920!--             2D in case of rotation
    19151921                ALLOCATE( netcdf_data_2d(0:nx,0:ny) )
    1916                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    1917                 sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     1922                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     1923                sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    19181924
    19191925                DO  k = 0, 2
     
    19341940                  DO  j = 0, ny
    19351941                     DO  i = 0, nx
    1936                         netcdf_data_2d(i,j) = init_model%origin_x            &
    1937                                             + cos_ra * ( i + shift_x ) * dx  &
    1938                                             + sin_ra * ( j + shift_y ) * dy
     1942                        netcdf_data_2d(i,j) = init_model%origin_x                   &
     1943                                            + cos_rot_angle * ( i + shift_x ) * dx  &
     1944                                            + sin_rot_angle * ( j + shift_y ) * dy
    19391945                     ENDDO
    19401946                  ENDDO
     
    19471953                  DO  j = 0, ny
    19481954                     DO  i = 0, nx
    1949                         netcdf_data_2d(i,j) = init_model%origin_y            &
    1950                                             - sin_ra * ( i + shift_x ) * dx  &
    1951                                             + cos_ra * ( j + shift_y ) * dy
     1955                        netcdf_data_2d(i,j) = init_model%origin_y                   &
     1956                                            - sin_rot_angle * ( i + shift_x ) * dx  &
     1957                                            + cos_rot_angle * ( j + shift_y ) * dy
    19521958                     ENDDO
    19531959                  ENDDO
     
    19992005             ALLOCATE( lat(nxl:nxr,nys:nyn) )
    20002006             ALLOCATE( lon(nxl:nxr,nys:nyn) )
    2001              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    2002              sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     2007             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     2008             sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    20032009
    20042010             DO  k = 0, 2
     
    20192025                DO  j = nys, nyn
    20202026                   DO  i = nxl, nxr
    2021                       eutm = init_model%origin_x            &
    2022                            + cos_ra * ( i + shift_x ) * dx  &
    2023                            + sin_ra * ( j + shift_y ) * dy
    2024                       nutm = init_model%origin_y            &
    2025                            - sin_ra * ( i + shift_x ) * dx  &
    2026                            + cos_ra * ( j + shift_y ) * dy
     2027                      eutm = init_model%origin_x                   &
     2028                           + cos_rot_angle * ( i + shift_x ) * dx  &
     2029                           + sin_rot_angle * ( j + shift_y ) * dy
     2030                      nutm = init_model%origin_y                   &
     2031                           - sin_rot_angle * ( i + shift_x ) * dx  &
     2032                           + cos_rot_angle * ( j + shift_y ) * dy
    20272033
    20282034                      CALL  convert_utm_to_geographic( crs_list,          &
     
    29132919!
    29142920!--          Write UTM coordinates
    2915              IF ( init_model%rotation_angle == 0.0_wp )  THEN
     2921             IF ( rotation_angle == 0.0_wp )  THEN
    29162922!
    29172923!--             1D in case of no rotation
    2918                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
     2924                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
    29192925!
    29202926!--             x coordinates
     
    29362942
    29372943                   DO  i = 0, nx
    2938                      netcdf_data(i) = init_model%origin_x            &
    2939                                     + cos_ra * ( i + shift_x ) * dx
     2944                     netcdf_data(i) = init_model%origin_x                      &
     2945                                    + cos_rot_angle * ( i + shift_x ) * dx
    29402946                   ENDDO
    29412947
    29422948                   nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_eutm_xy(k,av),&
    2943                                            netcdf_data, start = (/ 1 /),   &
     2949                                           netcdf_data, start = (/ 1 /),       &
    29442950                                           count = (/ nx+1 /) )
    29452951                   CALL netcdf_handle_error( 'netcdf_define_header', 555 )
     
    29662972
    29672973                   DO  j = 0, ny
    2968                       netcdf_data(j) = init_model%origin_y            &
    2969                                      + cos_ra * ( j + shift_y ) * dy
     2974                      netcdf_data(j) = init_model%origin_y                     &
     2975                                     + cos_rot_angle * ( j + shift_y ) * dy
    29702976                   ENDDO
    29712977
    29722978                   nc_stat = NF90_PUT_VAR( id_set_xy(av), id_var_nutm_xy(k,av),&
    2973                                            netcdf_data, start = (/ 1 /),   &
     2979                                           netcdf_data, start = (/ 1 /),       &
    29742980                                           count = (/ ny+1 /) )
    29752981                   CALL netcdf_handle_error( 'netcdf_define_header', 556 )
     
    29822988!--             2D in case of rotation
    29832989                ALLOCATE( netcdf_data_2d(0:nx,0:ny) )
    2984                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    2985                 sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     2990                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     2991                sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    29862992
    29872993                DO  k = 0, 2
     
    30023008                  DO  j = 0, ny
    30033009                     DO  i = 0, nx
    3004                         netcdf_data_2d(i,j) = init_model%origin_x            &
    3005                                             + cos_ra * ( i + shift_x ) * dx  &
    3006                                             + sin_ra * ( j + shift_y ) * dy
     3010                        netcdf_data_2d(i,j) = init_model%origin_x                   &
     3011                                            + cos_rot_angle * ( i + shift_x ) * dx  &
     3012                                            + sin_rot_angle * ( j + shift_y ) * dy
    30073013                     ENDDO
    30083014                  ENDDO
     
    30153021                  DO  j = 0, ny
    30163022                     DO  i = 0, nx
    3017                         netcdf_data_2d(i,j) = init_model%origin_y            &
    3018                                             - sin_ra * ( i + shift_x ) * dx  &
    3019                                             + cos_ra * ( j + shift_y ) * dy
     3023                        netcdf_data_2d(i,j) = init_model%origin_y                   &
     3024                                            - sin_rot_angle * ( i + shift_x ) * dx  &
     3025                                            + cos_rot_angle * ( j + shift_y ) * dy
    30203026                     ENDDO
    30213027                  ENDDO
     
    30373043             ALLOCATE( lat(nxl:nxr,nys:nyn) )
    30383044             ALLOCATE( lon(nxl:nxr,nys:nyn) )
    3039              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    3040              sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     3045             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     3046             sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    30413047
    30423048             DO  k = 0, 2
     
    30573063                DO  j = nys, nyn
    30583064                   DO  i = nxl, nxr
    3059                       eutm = init_model%origin_x            &
    3060                            + cos_ra * ( i + shift_x ) * dx  &
    3061                            + sin_ra * ( j + shift_y ) * dy
    3062                       nutm = init_model%origin_y            &
    3063                            - sin_ra * ( i + shift_x ) * dx  &
    3064                            + cos_ra * ( j + shift_y ) * dy
     3065                      eutm = init_model%origin_x                   &
     3066                           + cos_rot_angle * ( i + shift_x ) * dx  &
     3067                           + sin_rot_angle * ( j + shift_y ) * dy
     3068                      nutm = init_model%origin_y                   &
     3069                           - sin_rot_angle * ( i + shift_x ) * dx  &
     3070                           + cos_rot_angle * ( j + shift_y ) * dy
    30653071
    30663072                      CALL  convert_utm_to_geographic( crs_list,          &
     
    37983804!
    37993805!--          Write UTM coordinates
    3800              IF ( init_model%rotation_angle == 0.0_wp )  THEN
     3806             IF ( rotation_angle == 0.0_wp )  THEN
    38013807!
    38023808!--             1D in case of no rotation
    3803                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
     3809                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
    38043810!
    38053811!--             x coordinates
     
    38213827
    38223828                   DO  i = 0, nx
    3823                      netcdf_data(i) = init_model%origin_x            &
    3824                                     + cos_ra * ( i + shift_x ) * dx
     3829                     netcdf_data(i) = init_model%origin_x                      &
     3830                                    + cos_rot_angle * ( i + shift_x ) * dx
    38253831                   ENDDO
    38263832
    38273833                   nc_stat = NF90_PUT_VAR( id_set_xz(av), id_var_eutm_xz(k,av),&
    3828                                            netcdf_data, start = (/ 1 /),   &
     3834                                           netcdf_data, start = (/ 1 /),       &
    38293835                                           count = (/ nx+1 /) )
    38303836                   CALL netcdf_handle_error( 'netcdf_define_header', 555 )
     
    38553861                      ELSE
    38563862                         netcdf_data(i) = init_model%origin_y &
    3857                                      + cos_ra * ( section(i,2) + shift_y ) * dy
     3863                                     + cos_rot_angle * ( section(i,2) + shift_y ) * dy
    38583864                      ENDIF
    38593865                   ENDDO
     
    38713877!--             2D in case of rotation
    38723878                ALLOCATE( netcdf_data_2d(0:nx,1:ns) )
    3873                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    3874                 sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     3879                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     3880                sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    38753881
    38763882                DO  k = 0, 2
     
    38943900                      ELSE
    38953901                         DO  i = 0, nx
    3896                             netcdf_data_2d(i,j) = init_model%origin_x          &
    3897                                     + cos_ra * ( i + shift_x ) * dx            &
    3898                                     + sin_ra * ( section(j,2) + shift_y ) * dy
     3902                            netcdf_data_2d(i,j) = init_model%origin_x                 &
     3903                                    + cos_rot_angle * ( i + shift_x ) * dx            &
     3904                                    + sin_rot_angle * ( section(j,2) + shift_y ) * dy
    38993905                         ENDDO
    39003906                      ENDIF
     
    39113917                      ELSE
    39123918                         DO  i = 0, nx
    3913                             netcdf_data_2d(i,j) = init_model%origin_y          &
    3914                                     - sin_ra * ( i + shift_x ) * dx            &
    3915                                     + cos_ra * ( section(j,2) + shift_y ) * dy
     3919                            netcdf_data_2d(i,j) = init_model%origin_y                 &
     3920                                    - sin_rot_angle * ( i + shift_x ) * dx            &
     3921                                    + cos_rot_angle * ( section(j,2) + shift_y ) * dy
    39163922                         ENDDO
    39173923                      ENDIF
     
    39303936             ALLOCATE( lat(0:nx,1:ns) )
    39313937             ALLOCATE( lon(0:nx,1:ns) )
    3932              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    3933              sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     3938             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     3939             sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    39343940
    39353941             DO  k = 0, 2
     
    39543960                   ELSE
    39553961                      DO  i = 0, nx
    3956                          eutm = init_model%origin_x            &
    3957                               + cos_ra * ( i + shift_x ) * dx  &
    3958                               + sin_ra * ( section(j,2) + shift_y ) * dy
    3959                          nutm = init_model%origin_y            &
    3960                               - sin_ra * ( i + shift_x ) * dx  &
    3961                               + cos_ra * ( section(j,2) + shift_y ) * dy
     3962                         eutm = init_model%origin_x                   &
     3963                              + cos_rot_angle * ( i + shift_x ) * dx  &
     3964                              + sin_rot_angle * ( section(j,2) + shift_y ) * dy
     3965                         nutm = init_model%origin_y                   &
     3966                              - sin_rot_angle * ( i + shift_x ) * dx  &
     3967                              + cos_rot_angle * ( section(j,2) + shift_y ) * dy
    39623968
    39633969                         CALL  convert_utm_to_geographic( crs_list,          &
     
    46484654!
    46494655!--          Write UTM coordinates
    4650              IF ( init_model%rotation_angle == 0.0_wp )  THEN
     4656             IF ( rotation_angle == 0.0_wp )  THEN
    46514657!
    46524658!--             1D in case of no rotation
    4653                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
     4659                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
    46544660!
    46554661!--             x coordinates
     
    46754681                      ELSE
    46764682                         netcdf_data(i) = init_model%origin_x &
    4677                                      + cos_ra * ( section(i,3) + shift_x ) * dx
     4683                                     + cos_rot_angle * ( section(i,3) + shift_x ) * dx
    46784684                      ENDIF
    46794685                   ENDDO
     
    47054711
    47064712                   DO  i = 0, ny
    4707                      netcdf_data(i) = init_model%origin_y            &
    4708                                     + cos_ra * ( i + shift_y ) * dy
     4713                     netcdf_data(i) = init_model%origin_y                      &
     4714                                    + cos_rot_angle * ( i + shift_y ) * dy
    47094715                   ENDDO
    47104716
     
    47214727!--             2D in case of rotation
    47224728                ALLOCATE( netcdf_data_2d(1:ns,0:ny) )
    4723                 cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    4724                 sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     4729                cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     4730                sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    47254731
    47264732                DO  k = 0, 2
     
    47444750                            netcdf_data_2d(i,:) = -1.0_wp !section averaged along x
    47454751                         ELSE
    4746                             netcdf_data_2d(i,j) = init_model%origin_x          &
    4747                                     + cos_ra * ( section(i,3) + shift_x ) * dx &
    4748                                     + sin_ra * ( j + shift_y ) * dy
     4752                            netcdf_data_2d(i,j) = init_model%origin_x                 &
     4753                                    + cos_rot_angle * ( section(i,3) + shift_x ) * dx &
     4754                                    + sin_rot_angle * ( j + shift_y ) * dy
    47494755                         ENDIF
    47504756                      ENDDO
     
    47614767                            netcdf_data_2d(i,:) = -1.0_wp !section averaged along x
    47624768                         ELSE
    4763                             netcdf_data_2d(i,j) = init_model%origin_y          &
    4764                                     - sin_ra * ( section(i,3) + shift_x ) * dx &
    4765                                     + cos_ra * ( j + shift_y ) * dy
     4769                            netcdf_data_2d(i,j) = init_model%origin_y                 &
     4770                                    - sin_rot_angle * ( section(i,3) + shift_x ) * dx &
     4771                                    + cos_rot_angle * ( j + shift_y ) * dy
    47664772                         ENDIF
    47674773                      ENDDO
     
    47804786             ALLOCATE( lat(1:ns,0:ny) )
    47814787             ALLOCATE( lon(1:ns,0:ny) )
    4782              cos_ra = COS( init_model%rotation_angle * pi / 180.0_wp )
    4783              sin_ra = SIN( init_model%rotation_angle * pi / 180.0_wp )
     4788             cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     4789             sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
    47844790
    47854791             DO  k = 0, 2
     
    48044810                         lon(i,:) = -180.0_wp  ! section averaged along x
    48054811                      ELSE
    4806                          eutm = init_model%origin_x            &
    4807                               + cos_ra * ( section(i,3) + shift_x ) * dx  &
    4808                               + sin_ra * ( j + shift_y ) * dy
    4809                          nutm = init_model%origin_y            &
    4810                               - sin_ra * ( section(i,3) + shift_x ) * dx  &
    4811                               + cos_ra * ( j + shift_y ) * dy
     4812                         eutm = init_model%origin_x                              &
     4813                              + cos_rot_angle * ( section(i,3) + shift_x ) * dx  &
     4814                              + sin_rot_angle * ( j + shift_y ) * dy
     4815                         nutm = init_model%origin_y                              &
     4816                              - sin_rot_angle * ( section(i,3) + shift_x ) * dx  &
     4817                              + cos_rot_angle * ( j + shift_y ) * dy
    48124818
    48134819                         CALL  convert_utm_to_geographic( crs_list,          &
     
    69066912    nc_stat = NF90_PUT_ATT( ncid, NF90_GLOBAL, 'origin_z', init_model%origin_z )
    69076913    CALL netcdf_handle_error( 'netcdf_create_global_atts 11', error_no )
    6908     nc_stat = NF90_PUT_ATT( ncid, NF90_GLOBAL, 'rotation_angle', init_model%rotation_angle )
     6914    nc_stat = NF90_PUT_ATT( ncid, NF90_GLOBAL, 'rotation_angle', rotation_angle )
    69096915    CALL netcdf_handle_error( 'netcdf_create_global_atts 12', error_no )
    69106916
     
    70477053!
    70487054!-- Define UTM coordinates
    7049     IF ( init_model%rotation_angle == 0.0_wp )  THEN
     7055    IF ( rotation_angle == 0.0_wp )  THEN
    70507056       CALL netcdf_create_var( id_set, (/ id_dim_x(0) /), 'E_UTM', NF90_DOUBLE, id_var_eutm(0), 'm', 'easting', 000, 000, 000 )
    70517057       CALL netcdf_create_var( id_set, (/ id_dim_y(0) /), 'N_UTM', NF90_DOUBLE, id_var_nutm(0), 'm', 'northing', 000, 000, 000 )
  • palm/trunk/SOURCE/ocean_mod.f90

    r4182 r4196  
    2525! -----------------
    2626! $Id$
     27! Consider rotation of model domain for calculating the Stokes drift
     28!
     29! 4182 2019-08-22 15:20:23Z scharf
    2730! Corrected "Former revisions" section
    2831!
     
    21262129               v_stokes_zw, w, tend
    21272130
     2131    USE basic_constants_and_equations_mod,                                     &
     2132        ONLY:  pi
     2133
    21282134    USE control_parameters,                                                    &
    2129         ONLY:  f, fs, message_string
     2135        ONLY:  f, fs, message_string, rotation_angle
    21302136
    21312137    USE grid_variables,                                                        &
     
    21372143    IMPLICIT NONE
    21382144
    2139     INTEGER(iwp) ::  component  !< component of momentum equation
    2140     INTEGER(iwp) ::  i          !< loop index along x
    2141     INTEGER(iwp) ::  j          !< loop index along y
    2142     INTEGER(iwp) ::  k          !< loop index along z
    2143 
     2145    INTEGER(iwp) ::  component      !< component of momentum equation
     2146    INTEGER(iwp) ::  i              !< loop index along x
     2147    INTEGER(iwp) ::  j              !< loop index along y
     2148    INTEGER(iwp) ::  k              !< loop index along z
     2149
     2150    REAL(wp)     ::  cos_rot_angle  !< cosine of model rotation angle
     2151    REAL(wp)     ::  sin_rot_angle  !< sine of model rotation angle
    21442152
    21452153!
     
    21822190!--    w-component
    21832191       CASE ( 3 )
     2192
     2193!
     2194!--       Precalculate cosine and sine of rotation angle
     2195          cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     2196          sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
     2197
    21842198          DO  i = nxl, nxr
    21852199             DO  j = nys, nyn
     
    21972211                                                   + v(k+1,j+1,i) - v(k,j+1,i) &
    21982212                                                   ) * ddzu(k)  )              &
    2199                                            + fs * u_stokes_zw(k)
     2213                                           + fs * (                            &
     2214                                               sin_rot_angle * v_stokes_zw(k)  &
     2215                                             + cos_rot_angle * u_stokes_zw(k)  &
     2216                                                  )
    22002217                ENDDO
    22012218             ENDDO
     
    22262243               v_stokes_zw, w, tend
    22272244
     2245    USE basic_constants_and_equations_mod,                                     &
     2246        ONLY:  pi
     2247
    22282248    USE control_parameters,                                                    &
    2229         ONLY:  f, fs, message_string
     2249        ONLY:  f, fs, message_string, rotation_angle
    22302250
    22312251    USE grid_variables,                                                        &
     
    22372257    IMPLICIT NONE
    22382258
    2239     INTEGER(iwp) ::  component  !< component of momentum equation
    2240     INTEGER(iwp) ::  i          !< loop index along x
    2241     INTEGER(iwp) ::  j          !< loop index along y
    2242     INTEGER(iwp) ::  k          !< loop incex along z
    2243 
     2259    INTEGER(iwp) ::  component      !< component of momentum equation
     2260    INTEGER(iwp) ::  i              !< loop index along x
     2261    INTEGER(iwp) ::  j              !< loop index along y
     2262    INTEGER(iwp) ::  k              !< loop incex along z
     2263
     2264    REAL(wp)     ::  cos_rot_angle  !< cosine of model rotation angle
     2265    REAL(wp)     ::  sin_rot_angle  !< sine of model rotation angle
    22442266
    22452267!
     
    22732295!--    w-component
    22742296       CASE ( 3 )
     2297
     2298!
     2299!--       Precalculate cosine and sine of rotation angle
     2300          cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
     2301          sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
     2302
    22752303          DO  k = nzb+1, nzt
    22762304             tend(k,j,i) = tend(k,j,i) + u_stokes_zw(k) * (              &
     
    22862314                                                   + v(k+1,j+1,i) - v(k,j+1,i) &
    22872315                                                   ) * ddzu(k)  )              &
    2288                                        + fs * u_stokes_zw(k)
     2316                                       + fs * ( sin_rot_angle * v_stokes_zw(k) &
     2317                                              + cos_rot_angle * u_stokes_zw(k) &
     2318                                              )
    22892319          ENDDO
    22902320
  • palm/trunk/SOURCE/parin.f90

    r4191 r4196  
    2525! -----------------
    2626! $Id$
     27! added rotation_angle to initialization_parameters
     28!
     29! 4191 2019-08-27 15:45:07Z gronemeier
    2730! bugfix: add recycling_method_for_thermodynamic_quantities to inipar namelist
    2831!
     
    191194             recycling_width, recycling_yshift,                                &
    192195             reference_state, residual_limit,                                  &
     196             rotation_angle,                                                   &
    193197             roughness_length,                                                 &
    194198             scalar_advec,   &
     
    264268             recycling_width, recycling_yshift,                                &
    265269             reference_state, residual_limit,                                  &
     270             rotation_angle,                                                   &
    266271             roughness_length, scalar_advec,                                   &
    267272             scalar_rayleigh_damping,                                          &
Note: See TracChangeset for help on using the changeset viewer.