Changeset 4797 for palm/trunk/SOURCE
- Timestamp:
- Nov 26, 2020 4:02:39 PM (4 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/ocean_mod.f90
r4768 r4797 1 1 !> @file ocean_mod.f90 2 !------------------------------------------------------------------------------ !2 !--------------------------------------------------------------------------------------------------! 3 3 ! This file is part of the PALM model system. 4 4 ! 5 ! PALM is free software: you can redistribute it and/or modify it under the 6 ! terms of the GNU General Public License as published by the Free Software 7 ! Foundation, either version 3 of the License, or (at your option) any later 8 ! version. 9 ! 10 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY 11 ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR 12 ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. 13 ! 14 ! You should have received a copy of the GNU General Public License along with 15 ! PALM. If not, see <http://www.gnu.org/licenses/>. 5 ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General 6 ! Public License as published by the Free Software Foundation, either version 3 of the License, or 7 ! (at your option) any later version. 8 ! 9 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the 10 ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 ! Public License for more details. 12 ! 13 ! You should have received a copy of the GNU General Public License along with PALM. If not, see 14 ! <http://www.gnu.org/licenses/>. 16 15 ! 17 16 ! Copyright 2017-2020 Leibniz Universitaet Hannover 18 !-------------------------------------------------------------------------------- !17 !--------------------------------------------------------------------------------------------------! 19 18 ! 20 19 ! Current revisions: … … 25 24 ! ----------------- 26 25 ! $Id$ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4768 2020-11-02 19:11:23Z suehring 27 29 ! Enable 3D data output also with 64-bit precision 28 30 ! … … 32 34 ! 4671 2020-09-09 20:27:58Z pavelkrc 33 35 ! Implementation of downward facing USM and LSM surfaces 34 ! 36 ! 35 37 ! 4535 2020-05-15 12:07:23Z raasch 36 38 ! bugfix for restart data format query 37 ! 39 ! 38 40 ! 4517 2020-05-03 14:29:30Z raasch 39 41 ! added restart with MPI-IO for reading local arrays 40 ! 42 ! 41 43 ! 4495 2020-04-13 20:11:20Z raasch 42 44 ! restart data handling with MPI-IO added 43 ! 45 ! 44 46 ! 4481 2020-03-31 18:55:54Z maronga 45 47 ! vector directives added to force vectorization on Intel19 compiler 46 ! 48 ! 47 49 ! 4346 2019-12-18 11:55:56Z motisi 48 ! Introduction of wall_flags_total_0, which currently sets bits based on static 49 ! topographyinformation used in wall_flags_static_050 ! 50 ! Introduction of wall_flags_total_0, which currently sets bits based on static topography 51 ! information used in wall_flags_static_0 52 ! 51 53 ! 4329 2019-12-10 15:46:36Z motisi 52 54 ! Renamed wall_flags_0 to wall_flags_static_0 53 ! 55 ! 54 56 ! 4272 2019-10-23 15:18:57Z schwenkel 55 ! Further modularization of boundary conditions: moved boundary conditions to 56 ! respective modules 57 ! Further modularization of boundary conditions: moved boundary conditions to respective modules 57 58 ! 58 59 ! 4196 2019-08-29 11:02:06Z gronemeier 59 60 ! Consider rotation of model domain for calculating the Stokes drift 60 ! 61 ! 61 62 ! 4182 2019-08-22 15:20:23Z scharf 62 63 ! Corrected "Former revisions" section 63 ! 64 ! 64 65 ! 4110 2019-07-22 17:05:21Z suehring 65 ! Pass integer flag array as well as boundary flags to WS scalar advection 66 ! routine 67 ! 66 ! Pass integer flag array as well as boundary flags to WS scalar advection routine 67 ! 68 68 ! 4109 2019-07-22 17:00:34Z suehring 69 69 ! implemented ocean_actions 70 ! 70 ! 71 71 ! 3767 2019-02-27 08:18:02Z raasch 72 ! unused variable for file index and tmp_2d removed from rrd-subroutine parameter 73 ! list 74 ! 72 ! unused variable for file index and tmp_2d removed from rrd-subroutine parameter list 73 ! 75 74 ! 3719 2019-02-06 13:10:18Z kanani 76 ! Changed log_point to log_point_s, otherwise this overlaps with 77 ! 'all progn.equations'cpu measurement.78 ! 75 ! Changed log_point to log_point_s, otherwise this overlaps with 'all progn.equations' 76 ! cpu measurement. 77 ! 79 78 ! 3684 2019-01-20 20:20:58Z knoop 80 79 ! nopointer option removed 80 ! 81 81 ! 3294 2018-10-01 02:37:10Z raasch 82 82 ! initial revision … … 89 89 ! Description: 90 90 ! ------------ 91 !> This module contains relevant code for PALM's ocean mode, e.g. the 92 !> prognostic equation for salinity, the equation of state for seawater,93 !> the Craik Leibovich force (Stokes force), and wavebreaking effects94 !------------------------------------------------------------------------------ !91 !> This module contains relevant code for PALM's ocean mode, e.g. the prognostic equation for 92 !> salinity, the equation of state for seawater, the Craik Leibovich force (Stokes force), and wave 93 !> breaking effects 94 !--------------------------------------------------------------------------------------------------! 95 95 MODULE ocean_mod 96 97 98 USE arrays_3d, &99 ONLY: prho, prho_1, rho_ocean, rho_1, sa, sa_init, sa_1, sa_2, sa_3,&100 sa_ p, tsa_m, flux_l_sa, flux_s_sa, diss_l_sa, diss_s_sa101 102 USE control_parameters, &103 ONLY: atmos_ocean_sign, bottom_salinityflux, 104 constant_top_salinityflux, restart_data_format_output, ocean_mode, top_salinityflux,&105 w all_salinityflux, loop_optimization, ws_scheme_sca96 97 98 USE arrays_3d, & 99 ONLY: diss_l_sa, diss_s_sa, flux_l_sa, flux_s_sa, prho, prho_1, rho_1, rho_ocean, sa, & 100 sa_1, sa_2, sa_3, sa_init, sa_p, tsa_m 101 102 USE control_parameters, & 103 ONLY: atmos_ocean_sign, bottom_salinityflux, constant_top_salinityflux, loop_optimization,& 104 ocean_mode, restart_data_format_output, top_salinityflux, wall_salinityflux, & 105 ws_scheme_sca 106 106 107 107 USE kinds 108 108 109 USE pegrid, &109 USE pegrid, & 110 110 ONLY: threads_per_task 111 111 112 USE statistics, &112 USE statistics, & 113 113 ONLY: sums_wssas_ws_l 114 114 115 USE indices, &116 ONLY: advc_flags_s, n xl, nxr, nyn, nys, nzb, nzt, wall_flags_total_0, nbgp115 USE indices, & 116 ONLY: advc_flags_s, nbgp, nxl, nxr, nyn, nys, nzb, nzt, wall_flags_total_0 117 117 118 118 USE restart_data_mpi_io_mod, & … … 120 120 wrd_mpi_io_global_array 121 121 122 USE surface_mod, & 123 ONLY: bc_h, surf_def_v, surf_def_h, surf_lsm_h, surf_lsm_v, & 124 surf_usm_h, surf_usm_v 122 USE surface_mod, & 123 ONLY: bc_h, surf_def_v, surf_def_h, surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v 125 124 126 125 USE exchange_horiz_mod, & … … 132 131 CHARACTER (LEN=20) :: bc_sa_t = 'neumann' !< namelist parameter 133 132 134 INTEGER(iwp) :: ibc_sa_t !< integer flag for bc_sa_t 135 INTEGER(iwp) :: iran_ocean = -1234567 !< random number used for wave breaking effects 136 133 INTEGER(iwp) :: ibc_sa_t !< integer flag for bc_sa_t 134 INTEGER(iwp) :: iran_ocean = -1234567 !< random number used for wave breaking effects 137 135 INTEGER(iwp) :: sa_vertical_gradient_level_ind(10) = -9999 !< grid index values of sa_vertical_gradient_level(s) 138 136 139 137 LOGICAL :: salinity = .TRUE. !< switch for using salinity 140 138 LOGICAL :: stokes_force = .FALSE. !< switch to switch on the Stokes force 139 LOGICAL :: surface_cooling_switched_off = .FALSE. !< variable to check if surface heat flux has been switched off 141 140 LOGICAL :: wave_breaking = .FALSE. !< switch to switch on wave breaking effects 142 LOGICAL :: surface_cooling_switched_off = .FALSE. !< variable to check if surface heat flux has been switched off143 141 144 142 REAL(wp) :: alpha_wave_breaking = 3.0_wp !< coefficient for wave breaking generated turbulence from Noh et al. (2004), JPO … … 153 151 REAL(wp) :: u_star_wave_breaking !< to store the absolute value of friction velocity at the ocean surface 154 152 155 REAL(wp), DIMENSION(12), PARAMETER :: nom = &156 (/ 9.9984085444849347D2, 7.3471625860981584D0, &157 -5.3211231792841769D-2, 3.6492439109814549D-4, &158 2.5880571023991390D0, -6.7168282786692354D-3, &159 1.9203202055760151D-3, 1.1798263740430364D-2, &160 9.8920219266399117D-8, 4.6996642771754730D-6, &161 -2.5862187075154352D-8, -3.2921414007960662D-12 /)162 !< constants used in equation of state for seawater163 164 REAL(wp), DIMENSION(13), PARAMETER :: den = &165 (/ 1.0D0, 7.2815210113327091D-3, &166 -4.4787265461983921D-5, 3.3851002965802430D-7, &167 1.3651202389758572D-10, 1.7632126669040377D-3, &168 -8.8066583251206474D-6, -1.8832689434804897D-10, &169 5.7463776745432097D-6, 1.4716275472242334D-9, &170 6.7103246285651894D-6, -2.4461698007024582D-17, &171 -9.1534417604289062D-18 /)172 !< constants used in equation of state for seawater153 REAL(wp), DIMENSION(12), PARAMETER :: nom = & 154 (/ 9.9984085444849347D2, 7.3471625860981584D0, & 155 -5.3211231792841769D-2, 3.6492439109814549D-4, & 156 2.5880571023991390D0, -6.7168282786692354D-3, & 157 1.9203202055760151D-3, 1.1798263740430364D-2, & 158 9.8920219266399117D-8, 4.6996642771754730D-6, & 159 -2.5862187075154352D-8, -3.2921414007960662D-12 /) 160 !< constants used in equation of state for seawater 161 162 REAL(wp), DIMENSION(13), PARAMETER :: den = & 163 (/ 1.0D0, 7.2815210113327091D-3, & 164 -4.4787265461983921D-5, 3.3851002965802430D-7, & 165 1.3651202389758572D-10, 1.7632126669040377D-3, & 166 -8.8066583251206474D-6, -1.8832689434804897D-10, & 167 5.7463776745432097D-6, 1.4716275472242334D-9, & 168 6.7103246285651894D-6, -2.4461698007024582D-17, & 169 -9.1534417604289062D-18 /) 170 !< constants used in equation of state for seawater 173 171 174 172 SAVE 175 173 176 PUBLIC :: bc_sa_t, ibc_sa_t, prho_reference, sa_surface, 177 sa_vertical_gradient , sa_vertical_gradient_level,&178 sa_vertical_gradient_level_ind, stokes_force,wave_breaking174 PUBLIC :: bc_sa_t, ibc_sa_t, prho_reference, sa_surface, sa_vertical_gradient, & 175 sa_vertical_gradient_level, sa_vertical_gradient_level_ind, stokes_force, & 176 wave_breaking 179 177 180 178 … … 285 283 ! 286 284 !-- Add INTERFACES that must be available to other modules (alphabetical order) 287 PUBLIC eqn_state_seawater, ocean_actions, ocean_check_data_output, & 288 ocean_check_data_output_pr, ocean_check_parameters, & 289 ocean_data_output_2d, ocean_data_output_3d, ocean_exchange_horiz, & 290 ocean_define_netcdf_grid, ocean_header, ocean_init, & 291 ocean_init_arrays, ocean_parin, ocean_prognostic_equations, & 292 ocean_swap_timelevel, ocean_rrd_global, ocean_rrd_local, & 293 ocean_wrd_global, ocean_wrd_local, ocean_3d_data_averaging, & 294 ocean_boundary_conditions, stokes_drift_terms, wave_breaking_term 285 PUBLIC eqn_state_seawater, ocean_actions, ocean_check_data_output, ocean_check_data_output_pr, & 286 ocean_check_parameters, ocean_data_output_2d, ocean_data_output_3d, & 287 ocean_exchange_horiz, ocean_define_netcdf_grid, ocean_header, ocean_init, & 288 ocean_init_arrays, ocean_parin, ocean_prognostic_equations, ocean_swap_timelevel, & 289 ocean_rrd_global, ocean_rrd_local, ocean_wrd_global, ocean_wrd_local, & 290 ocean_3d_data_averaging, ocean_boundary_conditions, stokes_drift_terms, & 291 wave_breaking_term 295 292 296 293 297 294 CONTAINS 298 295 299 !------------------------------------------------------------------------------! 300 ! Description: 301 ! ------------ 302 !> Equation of state for seawater as a function of potential temperature, 303 !> salinity, and pressure. 296 !--------------------------------------------------------------------------------------------------! 297 ! Description: 298 ! ------------ 299 !> Equation of state for seawater as a function of potential temperature, salinity, and pressure. 304 300 !> For coefficients see Jackett et al., 2006: J. Atm. Ocean Tech. 305 301 !> eqn_state_seawater calculates the potential density referred at hyp(0). 306 302 !> eqn_state_seawater_func calculates density. 307 303 !> TODO: so far, routine is not adjusted to use topography 308 !------------------------------------------------------------------------------ !304 !--------------------------------------------------------------------------------------------------! 309 305 SUBROUTINE eqn_state_seawater 310 306 311 USE arrays_3d, &307 USE arrays_3d, & 312 308 ONLY: hyp, prho, pt_p, rho_ocean, sa_p 313 USE indices, &309 USE indices, & 314 310 ONLY: nxl, nxr, nyn, nys, nzb, nzt 315 311 … … 355 351 sa2 = sa1 * sa1 356 352 357 pnom = nom(1) + nom(2)*pt1 + nom(3)*pt2 + &358 nom(4)*pt3 + nom(5)*sa1 + nom(6)*sa1*pt1 + &353 pnom = nom(1) + nom(2)*pt1 + nom(3)*pt2 + & 354 nom(4)*pt3 + nom(5)*sa1 + nom(6)*sa1*pt1 + & 359 355 nom(7)*sa2 360 356 361 pden = den(1) + den(2)*pt1 + den(3)*pt2 + &362 den(4)*pt3 + den(5)*pt4 + den(6)*sa1 + &363 den(7)*sa1*pt1 + den(8)*sa1*pt3 + den(9)*sa15 + &357 pden = den(1) + den(2)*pt1 + den(3)*pt2 + & 358 den(4)*pt3 + den(5)*pt4 + den(6)*sa1 + & 359 den(7)*sa1*pt1 + den(8)*sa1*pt3 + den(9)*sa15 + & 364 360 den(10)*sa15*pt2 365 361 ! … … 367 363 prho(k,j,i) = pnom / pden 368 364 369 pnom = pnom + nom(8)*p1 + nom(9)*p1*pt2 + &365 pnom = pnom + nom(8)*p1 + nom(9)*p1*pt2 + & 370 366 nom(10)*p1*sa1 + nom(11)*p2 + nom(12)*p2*pt2 371 367 372 pden = pden + den(11)*p1 + den(12)*p2*pt3 + &368 pden = pden + den(11)*p1 + den(12)*p2*pt3 + & 373 369 den(13)*p3*pt1 374 370 … … 380 376 381 377 ! 382 !-- Neumann conditions are assumed at top boundary and currently also at 383 !-- bottom boundary(see comment lines below)378 !-- Neumann conditions are assumed at top boundary and currently also at bottom boundary 379 !-- (see comment lines below) 384 380 prho(nzt+1,j,i) = prho(nzt,j,i) 385 381 rho_ocean(nzt+1,j,i) = rho_ocean(nzt,j,i) … … 411 407 412 408 413 !------------------------------------------------------------------------------ !409 !--------------------------------------------------------------------------------------------------! 414 410 ! Description: 415 411 ! ------------ 416 412 !> Same as above, but for grid point i,j 417 !------------------------------------------------------------------------------ !413 !--------------------------------------------------------------------------------------------------! 418 414 SUBROUTINE eqn_state_seawater_ij( i, j ) 419 415 420 USE arrays_3d, &416 USE arrays_3d, & 421 417 ONLY: hyp, prho, pt_p, rho_ocean, sa_p 422 418 423 USE indices, &419 USE indices, & 424 420 ONLY: nzb, nzt 425 421 … … 446 442 REAL(wp) :: sa2 !< temporary scalar user for calculating density 447 443 444 448 445 DO k = nzb+1, nzt 449 446 ! … … 463 460 sa2 = sa1 * sa1 464 461 465 pnom = nom(1) + nom(2)*pt1 + nom(3)*pt2 + &462 pnom = nom(1) + nom(2)*pt1 + nom(3)*pt2 + & 466 463 nom(4)*pt3 + nom(5)*sa1 + nom(6)*sa1*pt1 + nom(7)*sa2 467 464 468 pden = den(1) + den(2)*pt1 + den(3)*pt2 + &469 den(4)*pt3 + den(5)*pt4 + den(6)*sa1 + &470 den(7)*sa1*pt1 + den(8)*sa1*pt3 + den(9)*sa15 + &465 pden = den(1) + den(2)*pt1 + den(3)*pt2 + & 466 den(4)*pt3 + den(5)*pt4 + den(6)*sa1 + & 467 den(7)*sa1*pt1 + den(8)*sa1*pt3 + den(9)*sa15 + & 471 468 den(10)*sa15*pt2 472 469 ! … … 474 471 prho(k,j,i) = pnom / pden 475 472 476 pnom = pnom + nom(8)*p1 + nom(9)*p1*pt2 +&473 pnom = pnom + nom(8)*p1 + nom(9)*p1*pt2 + & 477 474 nom(10)*p1*sa1 + nom(11)*p2 + nom(12)*p2*pt2 478 pden = pden + den(11)*p1 + den(12)*p2*pt3 + & 475 476 pden = pden + den(11)*p1 + den(12)*p2*pt3 + & 479 477 den(13)*p3*pt1 480 478 … … 510 508 511 509 512 !------------------------------------------------------------------------------ !510 !--------------------------------------------------------------------------------------------------! 513 511 ! Description: 514 512 ! ------------ 515 513 !> Equation of state for seawater as a function 516 !------------------------------------------------------------------------------ !514 !--------------------------------------------------------------------------------------------------! 517 515 REAL(wp) FUNCTION eqn_state_seawater_func( p, pt, sa ) 518 516 … … 532 530 REAL(wp) :: sa2 !< temporary scalar user for calculating density 533 531 532 534 533 ! 535 534 !-- Pressure is needed in dbar … … 549 548 550 549 551 eqn_state_seawater_func = &552 ( nom(1) + nom(2)*pt1 + nom(3)*pt2 + nom(4)*pt3 + &553 nom(5)*sa + nom(6)*sa*pt1 + nom(7)*sa2 + nom(8)*p1 + &554 nom(9)*p1*pt2 + nom(10)*p1*sa + nom(11)*p2 + nom(12)*p2*pt2 &555 ) / &556 ( den(1) + den(2)*pt1 + den(3)*pt2 + den(4)*pt3 + &557 den(5)*pt4 + den(6)*sa + den(7)*sa*pt1 + den(8)*sa*pt3 + &558 den(9)*sa15 + den(10)*sa15*pt2 + den(11)*p1 + den(12)*p2*pt3 + &559 den(13)*p3*pt1 &560 )550 eqn_state_seawater_func = & 551 ( nom(1) + nom(2)*pt1 + nom(3)*pt2 + nom(4)*pt3 + & 552 nom(5)*sa + nom(6)*sa*pt1 + nom(7)*sa2 + nom(8)*p1 + & 553 nom(9)*p1*pt2 + nom(10)*p1*sa + nom(11)*p2 + nom(12)*p2*pt2 & 554 ) / & 555 ( den(1) + den(2)*pt1 + den(3)*pt2 + den(4)*pt3 + & 556 den(5)*pt4 + den(6)*sa + den(7)*sa*pt1 + den(8)*sa*pt3 + & 557 den(9)*sa15 + den(10)*sa15*pt2 + den(11)*p1 + den(12)*p2*pt3 + & 558 den(13)*p3*pt1 & 559 ) 561 560 562 561 … … 564 563 565 564 566 !------------------------------------------------------------------------------ !565 !--------------------------------------------------------------------------------------------------! 567 566 ! Description: 568 567 ! ------------ 569 568 !> Reads the ocean parameters namelist 570 !------------------------------------------------------------------------------ !569 !--------------------------------------------------------------------------------------------------! 571 570 SUBROUTINE ocean_parin 572 571 … … 576 575 577 576 578 NAMELIST /ocean_parameters/ bc_sa_t, bottom_salinityflux, salinity, &579 sa_surface, sa_vertical_gradient, sa_vertical_gradient_level,&580 stokes_waveheight, stokes_wavelength, surface_cooling_spinup_time,&581 top_salinityflux, wall_salinityflux, wave_breaking577 NAMELIST /ocean_parameters/ bc_sa_t, bottom_salinityflux, salinity, sa_surface, & 578 sa_vertical_gradient, sa_vertical_gradient_level, & 579 stokes_waveheight, stokes_wavelength, surface_cooling_spinup_time,& 580 top_salinityflux, wall_salinityflux, wave_breaking 582 581 583 582 ! … … 607 606 END SUBROUTINE ocean_parin 608 607 609 !------------------------------------------------------------------------------ !608 !--------------------------------------------------------------------------------------------------! 610 609 ! Description: 611 610 ! ------------ 612 611 !> Check parameters routine for the ocean mode 613 !------------------------------------------------------------------------------ !612 !--------------------------------------------------------------------------------------------------! 614 613 SUBROUTINE ocean_check_parameters 615 614 616 USE control_parameters, & 617 ONLY: coupling_char, coupling_mode, initializing_actions, & 618 message_string, use_top_fluxes 619 620 USE pmc_interface, & 615 USE control_parameters, & 616 ONLY: coupling_char, coupling_mode, initializing_actions, message_string, use_top_fluxes 617 618 USE pmc_interface, & 621 619 ONLY: nested_run 622 620 … … 638 636 ! 639 637 !-- Check ocean setting 640 IF ( TRIM( coupling_mode ) == 'uncoupled' .AND. & 641 TRIM( coupling_char ) == '_O' .AND. & 638 IF ( TRIM( coupling_mode ) == 'uncoupled' .AND. TRIM( coupling_char ) == '_O' .AND. & 642 639 .NOT. ocean_mode ) THEN 643 640 644 641 ! 645 !-- Check whether an (uncoupled) atmospheric run has been declared as an 646 !-- ocean run (this settingis done via palmrun-option -y)647 message_string = 'ocean mode does not allow coupling_char = "' // 648 TRIM( coupling_char ) //'" set by palmrun-option "-y"'642 !-- Check whether an (uncoupled) atmospheric run has been declared as an ocean run (this setting 643 !-- is done via palmrun-option -y) 644 message_string = 'ocean mode does not allow coupling_char = "' // TRIM( coupling_char ) // & 645 '" set by palmrun-option "-y"' 649 646 CALL message( 'ocean_check_parameters', 'PA0317', 1, 2, 0, 6, 0 ) 650 647 … … 665 662 ibc_sa_t = 1 666 663 ELSE 667 message_string = 'unknown boundary condition: bc_sa_t = "' // & 668 TRIM( bc_sa_t ) // '"' 664 message_string = 'unknown boundary condition: bc_sa_t = "' // TRIM( bc_sa_t ) // '"' 669 665 CALL message( 'ocean_check_parameters', 'PA0068', 1, 2, 0, 6, 0 ) 670 666 ENDIF … … 673 669 674 670 IF ( .NOT. salinity ) THEN 675 IF ( ( bottom_salinityflux /= 0.0_wp .AND. & 676 bottom_salinityflux /= 9999999.9_wp ) .OR. & 677 ( top_salinityflux /= 0.0_wp .AND. & 678 top_salinityflux /= 9999999.9_wp ) ) & 671 IF ( ( bottom_salinityflux /= 0.0_wp .AND. bottom_salinityflux /= 9999999.9_wp ) .OR. & 672 ( top_salinityflux /= 0.0_wp .AND. top_salinityflux /= 9999999.9_wp ) ) & 679 673 THEN 680 message_string = 'salinityflux must not be set for ocean run ' // & 681 'without salinity' 674 message_string = 'salinityflux must not be set for ocean run ' // 'without salinity' 682 675 CALL message( 'ocean_check_parameters', 'PA0509', 1, 2, 0, 6, 0 ) 683 676 ENDIF … … 685 678 686 679 IF ( ibc_sa_t == 1 .AND. top_salinityflux == 9999999.9_wp ) THEN 687 message_string = 'boundary condition: bc_sa_t = "' // &688 TRIM( bc_sa_t ) //'" requires to set top_salinityflux'680 message_string = 'boundary condition: bc_sa_t = "' // TRIM( bc_sa_t ) // & 681 '" requires to set top_salinityflux' 689 682 CALL message( 'ocean_check_parameters', 'PA0069', 1, 2, 0, 6, 0 ) 690 683 ENDIF 691 684 692 685 ! 693 !-- A fixed salinity at the top implies Dirichlet boundary condition for 694 !-- salinity. In this case specification of a constant salinity flux is 695 !-- forbidden. 696 IF ( ibc_sa_t == 0 .AND. constant_top_salinityflux .AND. & 697 top_salinityflux /= 0.0_wp ) THEN 698 message_string = 'boundary condition: bc_sa_t = "' // & 699 TRIM( bc_sa_t ) // '" is not allowed with ' // & 700 'top_salinityflux /= 0.0' 686 !-- A fixed salinity at the top implies Dirichlet boundary condition for salinity. In this case 687 !-- specification of a constant salinity flux is forbidden. 688 IF ( ibc_sa_t == 0 .AND. constant_top_salinityflux .AND. top_salinityflux /= 0.0_wp ) THEN 689 message_string = 'boundary condition: bc_sa_t = "' // TRIM( bc_sa_t ) // & 690 '" is not allowed with top_salinityflux /= 0.0' 701 691 CALL message( 'ocean_check_parameters', 'PA0070', 1, 2, 0, 6, 0 ) 702 692 ENDIF … … 707 697 stokes_force = .TRUE. 708 698 ELSE 709 IF ( ( stokes_waveheight <= 0.0_wp .AND. stokes_wavelength > 0.0_wp ) & 710 .OR. & 711 ( stokes_waveheight > 0.0_wp .AND. stokes_wavelength <= 0.0_wp ) & 712 .OR. & 713 ( stokes_waveheight < 0.0_wp .AND. stokes_wavelength < 0.0_wp ) ) & 699 IF ( ( stokes_waveheight <= 0.0_wp .AND. stokes_wavelength > 0.0_wp ) .OR. & 700 ( stokes_waveheight > 0.0_wp .AND. stokes_wavelength <= 0.0_wp ) .OR. & 701 ( stokes_waveheight < 0.0_wp .AND. stokes_wavelength < 0.0_wp ) ) & 714 702 THEN 715 message_string = 'wrong settings for stokes_wavelength and/or ' // & 716 'stokes_waveheight' 703 message_string = 'wrong settings for stokes_wavelength and/or stokes_waveheight' 717 704 CALL message( 'ocean_check_parameters', 'PA0460', 1, 2, 0, 6, 0 ) 718 705 ENDIF … … 722 709 723 710 724 !------------------------------------------------------------------------------ !711 !--------------------------------------------------------------------------------------------------! 725 712 ! Description: 726 713 ! ------------ 727 714 !> Check data output. 728 !------------------------------------------------------------------------------ !715 !--------------------------------------------------------------------------------------------------! 729 716 SUBROUTINE ocean_check_data_output( var, unit ) 730 717 731 718 IMPLICIT NONE 732 719 … … 751 738 752 739 753 !------------------------------------------------------------------------------ !740 !--------------------------------------------------------------------------------------------------! 754 741 ! Description: 755 742 ! ------------ 756 743 !> Check data output of profiles 757 !------------------------------------------------------------------------------ !744 !--------------------------------------------------------------------------------------------------! 758 745 SUBROUTINE ocean_check_data_output_pr( variable, var_count, unit, dopr_unit ) 759 746 760 USE arrays_3d, &747 USE arrays_3d, & 761 748 ONLY: zu, zw 762 749 763 USE control_parameters, &750 USE control_parameters, & 764 751 ONLY: data_output_pr 765 752 … … 772 759 IMPLICIT NONE 773 760 761 CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit 774 762 CHARACTER (LEN=*) :: unit !< 775 763 CHARACTER (LEN=*) :: variable !< 776 CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit777 764 778 765 INTEGER(iwp) :: var_count !< 766 779 767 780 768 SELECT CASE ( TRIM( variable ) ) … … 837 825 838 826 839 !------------------------------------------------------------------------------ !827 !--------------------------------------------------------------------------------------------------! 840 828 ! Description: 841 829 ! ------------ 842 830 !> Define appropriate grid for netcdf variables. 843 831 !> It is called out from subroutine netcdf. 844 !------------------------------------------------------------------------------ !832 !--------------------------------------------------------------------------------------------------! 845 833 SUBROUTINE ocean_define_netcdf_grid( var, found, grid_x, grid_y, grid_z ) 846 834 847 835 IMPLICIT NONE 848 836 … … 854 842 LOGICAL, INTENT(OUT) :: found !< flag if output variable is found 855 843 844 856 845 found = .TRUE. 857 846 … … 860 849 SELECT CASE ( TRIM( var ) ) 861 850 862 CASE ( 'rho_sea_water', 'rho_sea_water_xy', 'rho_sea_water_xz', &863 ' rho_sea_water_yz', 'sa', 'sa_xy', 'sa_xz', 'sa_yz' )851 CASE ( 'rho_sea_water', 'rho_sea_water_xy', 'rho_sea_water_xz', 'rho_sea_water_yz', 'sa', & 852 'sa_xy', 'sa_xz', 'sa_yz' ) 864 853 grid_x = 'x' 865 854 grid_y = 'y' … … 877 866 878 867 879 !------------------------------------------------------------------------------ !868 !--------------------------------------------------------------------------------------------------! 880 869 ! Description: 881 870 ! ------------ 882 871 !> Average 3D data. 883 !------------------------------------------------------------------------------ !872 !--------------------------------------------------------------------------------------------------! 884 873 SUBROUTINE ocean_3d_data_averaging( mode, variable ) 885 886 887 USE arrays_3d, &874 875 876 USE arrays_3d, & 888 877 ONLY: rho_ocean, sa 889 878 890 USE averaging, &879 USE averaging, & 891 880 ONLY: rho_ocean_av, sa_av 892 881 893 USE control_parameters, &882 USE control_parameters, & 894 883 ONLY: average_count_3d 895 884 896 USE indices, &885 USE indices, & 897 886 ONLY: nxlg, nxrg, nyng, nysg, nzb, nzt 898 887 … … 905 894 INTEGER(iwp) :: j !< loop index 906 895 INTEGER(iwp) :: k !< loop index 896 907 897 908 898 IF ( mode == 'allocate' ) THEN … … 936 926 DO j = nysg, nyng 937 927 DO k = nzb, nzt+1 938 rho_ocean_av(k,j,i) = rho_ocean_av(k,j,i) + & 939 rho_ocean(k,j,i) 928 rho_ocean_av(k,j,i) = rho_ocean_av(k,j,i) + rho_ocean(k,j,i) 940 929 ENDDO 941 930 ENDDO … … 964 953 965 954 CASE ( 'rho_sea_water' ) 966 IF ( ALLOCATED( rho_ocean_av ) ) THEN955 IF ( ALLOCATED( rho_ocean_av ) ) THEN 967 956 DO i = nxlg, nxrg 968 957 DO j = nysg, nyng 969 958 DO k = nzb, nzt+1 970 rho_ocean_av(k,j,i) = rho_ocean_av(k,j,i) / &959 rho_ocean_av(k,j,i) = rho_ocean_av(k,j,i) / & 971 960 REAL( average_count_3d, KIND=wp ) 972 961 ENDDO … … 976 965 977 966 CASE ( 'sa' ) 978 IF ( ALLOCATED( sa_av ) ) THEN967 IF ( ALLOCATED( sa_av ) ) THEN 979 968 DO i = nxlg, nxrg 980 969 DO j = nysg, nyng 981 970 DO k = nzb, nzt+1 982 sa_av(k,j,i) = sa_av(k,j,i) / & 983 REAL( average_count_3d, KIND=wp ) 971 sa_av(k,j,i) = sa_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 984 972 ENDDO 985 973 ENDDO … … 994 982 995 983 996 !------------------------------------------------------------------------------ !984 !--------------------------------------------------------------------------------------------------! 997 985 ! Description: 998 986 ! ------------ 999 987 !> Define 2D output variables. 1000 !------------------------------------------------------------------------------! 1001 SUBROUTINE ocean_data_output_2d( av, variable, found, grid, mode, local_pf, & 1002 nzb_do, nzt_do ) 1003 1004 USE arrays_3d, & 988 !--------------------------------------------------------------------------------------------------! 989 SUBROUTINE ocean_data_output_2d( av, variable, found, grid, mode, local_pf, nzb_do, nzt_do ) 990 991 USE arrays_3d, & 1005 992 ONLY: rho_ocean, sa 1006 993 1007 USE averaging, &994 USE averaging, & 1008 995 ONLY: rho_ocean_av, sa_av 1009 996 1010 USE indices, & 1011 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, & 1012 wall_flags_total_0 997 USE indices, & 998 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, wall_flags_total_0 1013 999 1014 1000 IMPLICIT NONE … … 1031 1017 REAL(wp) :: fill_value = -999.0_wp !< value for the _FillValue attribute 1032 1018 1033 REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< local 1034 !< array to which output data is resorted to 1019 REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< local array to which output data is resorted to 1035 1020 1036 1021 REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to selected output variable 1037 1022 1023 1038 1024 found = .TRUE. 1039 1025 resorted = .FALSE. … … 1048 1034 to_be_resorted => rho_ocean 1049 1035 ELSE 1050 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) THEN1036 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) THEN 1051 1037 ALLOCATE( rho_ocean_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1052 1038 rho_ocean_av = REAL( fill_value, KIND = wp ) … … 1059 1045 to_be_resorted => sa 1060 1046 ELSE 1061 IF ( .NOT. ALLOCATED( sa_av ) ) THEN1047 IF ( .NOT. ALLOCATED( sa_av ) ) THEN 1062 1048 ALLOCATE( sa_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1063 1049 sa_av = REAL( fill_value, KIND = wp ) … … 1077 1063 DO j = nys, nyn 1078 1064 DO k = nzb_do, nzt_do 1079 local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), & 1080 REAL( fill_value, KIND = wp ), & 1081 BTEST( wall_flags_total_0(k,j,i), flag_nr ) ) 1065 local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), REAL( fill_value, KIND = wp ), & 1066 BTEST( wall_flags_total_0(k,j,i), flag_nr ) ) 1082 1067 ENDDO 1083 1068 ENDDO … … 1085 1070 resorted = .TRUE. 1086 1071 ENDIF 1087 1072 1088 1073 END SUBROUTINE ocean_data_output_2d 1089 1074 1090 1091 !------------------------------------------------------------------------------ !1075 1076 !--------------------------------------------------------------------------------------------------! 1092 1077 ! Description: 1093 1078 ! ------------ 1094 1079 !> Define 3D output variables. 1095 !------------------------------------------------------------------------------ !1080 !--------------------------------------------------------------------------------------------------! 1096 1081 SUBROUTINE ocean_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do ) 1097 1098 1099 USE arrays_3d, &1082 1083 1084 USE arrays_3d, & 1100 1085 ONLY: rho_ocean, sa 1101 1086 1102 USE averaging, &1087 USE averaging, & 1103 1088 ONLY: rho_ocean_av, sa_av 1104 1089 1105 USE indices, & 1106 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, & 1107 wall_flags_total_0 1090 USE indices, & 1091 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, wall_flags_total_0 1108 1092 1109 1093 IMPLICIT NONE … … 1129 1113 REAL(wp), DIMENSION(:,:,:), POINTER :: to_be_resorted !< points to selected output variable 1130 1114 1115 1131 1116 found = .TRUE. 1132 1117 resorted = .FALSE. … … 1141 1126 to_be_resorted => rho_ocean 1142 1127 ELSE 1143 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) THEN1128 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) THEN 1144 1129 ALLOCATE( rho_ocean_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1145 1130 rho_ocean_av = REAL( fill_value, KIND = wp ) … … 1152 1137 to_be_resorted => sa 1153 1138 ELSE 1154 IF ( .NOT. ALLOCATED( sa_av ) ) THEN1139 IF ( .NOT. ALLOCATED( sa_av ) ) THEN 1155 1140 ALLOCATE( sa_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1156 1141 sa_av = REAL( fill_value, KIND = wp ) … … 1169 1154 DO j = nys, nyn 1170 1155 DO k = nzb_do, nzt_do 1171 local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), & 1172 REAL( fill_value, KIND = wp ), & 1173 BTEST( wall_flags_total_0(k,j,i), flag_nr ) ) 1156 local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i), REAL( fill_value, KIND = wp ), & 1157 BTEST( wall_flags_total_0(k,j,i), flag_nr ) ) 1174 1158 ENDDO 1175 1159 ENDDO … … 1206 1190 1207 1191 1208 !------------------------------------------------------------------------------ !1192 !--------------------------------------------------------------------------------------------------! 1209 1193 ! Description: 1210 1194 ! ------------ 1211 1195 !> Header output for ocean parameters 1212 !------------------------------------------------------------------------------ !1196 !--------------------------------------------------------------------------------------------------! 1213 1197 SUBROUTINE ocean_header( io ) 1214 1198 1215 1216 1199 IMPLICIT NONE 1217 1200 1218 1201 INTEGER(iwp), INTENT(IN) :: io !< Unit of the output file 1202 1219 1203 1220 1204 ! … … 1230 1214 ENDIF 1231 1215 1232 1 FORMAT (//' Ocean settings:'/ & 1233 ' ------------------------------------------'/) 1234 2 FORMAT (' --> Craik-Leibovich vortex force and Stokes drift switched', & 1235 ' on'/ & 1216 1 FORMAT (//' Ocean settings:'/' ------------------------------------------'/) 1217 2 FORMAT (' --> Craik-Leibovich vortex force and Stokes drift switched',' on'/ & 1236 1218 ' waveheight: ',F4.1,' m wavelength: ',F6.1,' m') 1237 3 FORMAT (' --> wave breaking generated turbulence switched on'/ & 1238 ' alpha: ',F4.1/ & 1239 ' timescale:',F5.1,' s') 1219 3 FORMAT (' --> wave breaking generated turbulence switched on'/ & 1220 ' alpha: ',F4.1/' timescale:',F5.1,' s') 1240 1221 4 FORMAT (' --> prognostic salinity equation is switched off' ) 1241 1222 5 FORMAT (' --> surface heat flux is switched off after ',F8.1,' s') … … 1244 1225 1245 1226 1246 !------------------------------------------------------------------------------ !1227 !--------------------------------------------------------------------------------------------------! 1247 1228 ! Description: 1248 1229 ! ------------ 1249 1230 !> Allocate arrays and assign pointers. 1250 !------------------------------------------------------------------------------ !1231 !--------------------------------------------------------------------------------------------------! 1251 1232 SUBROUTINE ocean_init_arrays 1252 1233 1253 USE indices, &1234 USE indices, & 1254 1235 ONLY: nxlg, nxrg, nyn, nyng, nys, nysg, nzb, nzt 1255 1236 1256 1237 IMPLICIT NONE 1257 1238 1258 ALLOCATE( prho_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &1259 rho_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &1260 sa_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &1239 ALLOCATE( prho_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), & 1240 rho_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), & 1241 sa_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), & 1261 1242 sa_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 1243 1262 1244 1263 1245 IF ( salinity ) THEN … … 1282 1264 1283 1265 prho => prho_1 1284 rho_ocean => rho_1 ! routines calc_mean_profile and diffusion_e require 1285 ! density to be a pointer 1266 rho_ocean => rho_1 ! routines calc_mean_profile and diffusion_e require density to be a pointer 1286 1267 1287 1268 ! … … 1296 1277 1297 1278 1298 !------------------------------------------------------------------------------ !1279 !--------------------------------------------------------------------------------------------------! 1299 1280 ! Description: 1300 1281 ! ------------ 1301 1282 !> Initialization of quantities needed for the ocean mode 1302 !------------------------------------------------------------------------------ !1283 !--------------------------------------------------------------------------------------------------! 1303 1284 SUBROUTINE ocean_init 1304 1285 1305 1286 1306 USE arrays_3d, &1307 ONLY: dzu, dzw, hyp, pt_init, ref_state, u_stokes_zu, u_stokes_zw, &1308 v_stokes_z u, v_stokes_zw, zu, zw1309 1310 USE basic_constants_and_equations_mod, &1287 USE arrays_3d, & 1288 ONLY: dzu, dzw, hyp, pt_init, ref_state, u_stokes_zu, u_stokes_zw, v_stokes_zu, & 1289 v_stokes_zw, zu, zw 1290 1291 USE basic_constants_and_equations_mod, & 1311 1292 ONLY: g 1312 1293 1313 USE basic_constants_and_equations_mod, &1294 USE basic_constants_and_equations_mod, & 1314 1295 ONLY: pi 1315 1296 1316 USE control_parameters, & 1317 ONLY: initializing_actions, molecular_viscosity, rho_surface, & 1318 rho_reference, surface_pressure, top_momentumflux_u, & 1319 top_momentumflux_v, use_single_reference_value 1320 1321 USE indices, & 1297 USE control_parameters, & 1298 ONLY: initializing_actions, molecular_viscosity, rho_surface, rho_reference, & 1299 surface_pressure, top_momentumflux_u, top_momentumflux_v, use_single_reference_value 1300 1301 USE indices, & 1322 1302 ONLY: nxl, nxlg, nxrg, nyng, nys, nysg, nzb, nzt 1323 1303 1324 1304 USE kinds 1325 1305 1326 USE pegrid, &1306 USE pegrid, & 1327 1307 ONLY: myid 1328 1308 1329 USE statistics, &1309 USE statistics, & 1330 1310 ONLY: hom, statistic_regions 1331 1311 … … 1337 1317 INTEGER(iwp) :: n !< loop index 1338 1318 1339 REAL(wp) :: alpha !< angle of surface stress1340 REAL(wp) :: dum !< dummy argument1341 REAL(wp) :: pt_l !< local scalar for pt used in equation of state function1342 REAL(wp) :: sa_l !< local scalar for sa used in equation of state function1319 REAL(wp) :: alpha !< angle of surface stress 1320 REAL(wp) :: dum !< dummy argument 1321 REAL(wp) :: pt_l !< local scalar for pt used in equation of state function 1322 REAL(wp) :: sa_l !< local scalar for sa used in equation of state function 1343 1323 REAL(wp) :: velocity_amplitude !< local scalar for amplitude of Stokes drift velocity 1344 REAL(wp) :: x !< temporary variable to store surface stress along x1345 REAL(wp) :: y !< temporary variable to store surface stress along y1324 REAL(wp) :: x !< temporary variable to store surface stress along x 1325 REAL(wp) :: y !< temporary variable to store surface stress along y 1346 1326 1347 1327 REAL(wp), DIMENSION(nzb:nzt+1) :: rho_ocean_init !< local array for initial density … … 1351 1331 1352 1332 ! 1353 !-- In case of no restart run, calculate the inital salinity profile vcusing the1354 !-- g iven salinity gradients1333 !-- In case of no restart run, calculate the inital salinity profile using the given salinity 1334 !-- gradients. 1355 1335 IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 1356 1336 1357 1337 sa_init = sa_surface 1358 1338 ! 1359 !-- Last arguments gives back the gradient at top level to be used as 1360 !-- possible Neumann boundary condition. This is not realized for the ocean 1361 !-- mode, therefore a dummy argument is used. 1339 !-- Last arguments gives back the gradient at top level to be used as possible Neumann boundary 1340 !-- condition. This is not realized for the ocean mode, therefore a dummy argument is used. 1362 1341 IF ( salinity ) THEN 1363 CALL init_vertical_profiles( sa_vertical_gradient_level_ind, &1364 sa_vertical_gradient_level, &1365 sa_vertical_gradient, sa_init, &1342 CALL init_vertical_profiles( sa_vertical_gradient_level_ind, & 1343 sa_vertical_gradient_level, & 1344 sa_vertical_gradient, sa_init, & 1366 1345 sa_surface, dum ) 1367 1346 ENDIF … … 1370 1349 ! 1371 1350 !-- Initialize required 3d-arrays 1372 IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. &1351 IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. & 1373 1352 TRIM( initializing_actions ) /= 'cyclic_fill' ) THEN 1374 1353 ! … … 1395 1374 1396 1375 ! 1397 !-- Initialize new time levels (only done in order to set boundary values 1398 !-- including ghost points) 1376 !-- Initialize new time levels (only done in order to set boundary values including ghost points) 1399 1377 sa_p = sa 1400 1378 ! 1401 !-- Allthough tendency arrays are set in prognostic_equations, they have 1402 !-- have to be predefined here because they are used (but multiplied with 0) 1403 !-- there before they are set. 1379 !-- Allthough tendency arrays are set in prognostic_equations, they have to be predefined here 1380 !-- because they are used (but multiplied with 0) there before they are set. 1404 1381 tsa_m = 0.0_wp 1405 1382 … … 1416 1393 1417 1394 ! 1418 !-- Change sign of buoyancy/stability terms because density gradient is used 1419 !-- instead of the potential temperature gradient to calculate the buoyancy1395 !-- Change sign of buoyancy/stability terms because density gradient is used instead of the 1396 !-- potential temperature gradient to calculate the buoyancy. 1420 1397 atmos_ocean_sign = -1.0_wp 1421 1398 1422 1399 ! 1423 !-- Calculate initial vertical profile of hydrostatic pressure (in Pa) 1424 !-- and the reference density (used later in buoyancy term)1425 !-- First step: Calculate pressure using reference density 1400 !-- Calculate initial vertical profile of hydrostatic pressure (in Pa) and the reference density 1401 !-- (used later in buoyancy term). 1402 !-- First step: Calculate pressure using reference density. 1426 1403 hyp(nzt+1) = surface_pressure * 100.0_wp 1427 1404 hyp(nzt) = hyp(nzt+1) + rho_surface * g * 0.5_wp * dzu(nzt+1) … … 1435 1412 1436 1413 ! 1437 !-- Second step: Iteratively calculate in situ density (based on presssure) 1438 !-- and pressure(based on in situ density)1414 !-- Second step: Iteratively calculate in situ density (based on presssure) and pressure 1415 !-- (based on in situ density) 1439 1416 DO n = 1, 5 1440 1417 … … 1456 1433 hyp(nzt) = hyp(nzt+1) + rho_surface * g * 0.5_wp * dzu(nzt+1) 1457 1434 DO k = nzt-1, 0, -1 1458 hyp(k) = hyp(k+1) + g * 0.5_wp * ( rho_ocean_init(k) & 1459 + rho_ocean_init(k+1) ) * dzu(k+1) 1435 hyp(k) = hyp(k+1) + g * 0.5_wp * ( rho_ocean_init(k) + rho_ocean_init(k+1) ) * dzu(k+1) 1460 1436 ENDDO 1461 1437 … … 1470 1446 pt_l = 0.5_wp * ( pt_init(k) + pt_init(k+1) ) 1471 1447 1472 prho_reference = prho_reference + dzu(k+1) * & 1473 eqn_state_seawater_func( 0.0_wp, pt_l, sa_l ) 1448 prho_reference = prho_reference + dzu(k+1) * eqn_state_seawater_func( 0.0_wp, pt_l, sa_l ) 1474 1449 1475 1450 ENDDO … … 1478 1453 1479 1454 ! 1480 !-- Calculate the 3d array of initial in situ and potential density, 1481 !-- based on the initial temperature and salinity profile1455 !-- Calculate the 3d array of initial in situ and potential density, based on the initial 1456 !-- temperature and salinity profile. 1482 1457 CALL eqn_state_seawater 1483 1458 … … 1520 1495 ENDIF 1521 1496 1522 velocity_amplitude = ( pi * stokes_waveheight / stokes_wavelength )**2 * &1497 velocity_amplitude = ( pi * stokes_waveheight / stokes_wavelength )**2 * & 1523 1498 SQRT( g * stokes_wavelength / ( 2.0_wp * pi ) ) 1524 1499 1525 1500 DO k = nzb, nzt 1526 u_stokes_zu(k) = velocity_amplitude * COS( alpha ) * &1501 u_stokes_zu(k) = velocity_amplitude * COS( alpha ) * & 1527 1502 EXP( 4.0_wp * pi * zu(k) / stokes_wavelength ) 1528 u_stokes_zw(k) = velocity_amplitude * COS( alpha ) * &1503 u_stokes_zw(k) = velocity_amplitude * COS( alpha ) * & 1529 1504 EXP( 4.0_wp * pi * zw(k) / stokes_wavelength ) 1530 v_stokes_zu(k) = velocity_amplitude * SIN( alpha ) * &1505 v_stokes_zu(k) = velocity_amplitude * SIN( alpha ) * & 1531 1506 EXP( 4.0_wp * pi * zu(k) / stokes_wavelength ) 1532 1507 v_stokes_zw(k) = velocity_amplitude * SIN( alpha ) * & … … 1545 1520 ! 1546 1521 !-- Calculate friction velocity at ocean surface 1547 u_star_wave_breaking = SQRT( SQRT( top_momentumflux_u**2 + & 1548 top_momentumflux_v**2 ) ) 1549 ! 1550 !-- Set the time scale of random forcing. The vertical grid spacing at the 1551 !-- ocean surface is assumed as the length scale of turbulence. 1522 u_star_wave_breaking = SQRT( SQRT( top_momentumflux_u**2 + top_momentumflux_v**2 ) ) 1523 ! 1524 !-- Set the time scale of random forcing. The vertical grid spacing at the ocean surface is 1525 !-- assumed as the length scale of turbulence. 1552 1526 !-- Formula follows Noh et al. (2004), JPO 1553 timescale_wave_breaking = 0.1_wp * dzw(nzt) / alpha_wave_breaking / & 1554 u_star_wave_breaking 1555 ! 1556 !-- Set random number seeds differently on the processor cores in order to 1557 !-- create different random number sequences 1527 timescale_wave_breaking = 0.1_wp * dzw(nzt) / alpha_wave_breaking / u_star_wave_breaking 1528 ! 1529 !-- Set random number seeds differently on the processor cores in order to create different 1530 !-- random number sequences 1558 1531 iran_ocean = iran_ocean + myid 1559 1532 ENDIF … … 1562 1535 1563 1536 1564 !------------------------------------------------------------------------------ !1537 !--------------------------------------------------------------------------------------------------! 1565 1538 ! Description: 1566 1539 ! ------------ 1567 1540 !> Call for all grid points 1568 !------------------------------------------------------------------------------ !1541 !--------------------------------------------------------------------------------------------------! 1569 1542 SUBROUTINE ocean_actions( location ) 1570 1543 … … 1588 1561 1589 1562 1590 !------------------------------------------------------------------------------ !1563 !--------------------------------------------------------------------------------------------------! 1591 1564 ! Description: 1592 1565 ! ------------ 1593 1566 !> Call for grid points i,j 1594 !------------------------------------------------------------------------------ !1567 !--------------------------------------------------------------------------------------------------! 1595 1568 SUBROUTINE ocean_actions_ij( i, j, location ) 1596 1569 1597 1598 INTEGER(iwp), INTENT(IN) :: i !< grid index in x-direction1599 INTEGER(iwp), INTENT(IN) :: j !< grid index in y-direction1600 1570 CHARACTER (LEN=*), INTENT(IN) :: location !< call location string 1601 INTEGER(iwp) :: dummy !< call location string 1571 1572 INTEGER(iwp) :: dummy !< call location string 1573 INTEGER(iwp), INTENT(IN) :: i !< grid index in x-direction 1574 INTEGER(iwp), INTENT(IN) :: j !< grid index in y-direction 1575 1602 1576 1603 1577 IF ( ocean_mode ) dummy = i + j … … 1619 1593 1620 1594 1621 !------------------------------------------------------------------------------ !1595 !--------------------------------------------------------------------------------------------------! 1622 1596 ! Description: 1623 1597 ! ------------ 1624 1598 !> Prognostic equation for salinity. 1625 1599 !> Vector-optimized version 1626 !------------------------------------------------------------------------------ !1600 !--------------------------------------------------------------------------------------------------! 1627 1601 SUBROUTINE ocean_prognostic_equations 1628 1602 1629 USE advec_s_bc_mod, &1603 USE advec_s_bc_mod, & 1630 1604 ONLY: advec_s_bc 1631 1605 1632 USE advec_s_pw_mod, &1606 USE advec_s_pw_mod, & 1633 1607 ONLY: advec_s_pw 1634 1608 1635 USE advec_s_up_mod, &1609 USE advec_s_up_mod, & 1636 1610 ONLY: advec_s_up 1637 1611 1638 USE advec_ws, &1612 USE advec_ws, & 1639 1613 ONLY: advec_s_ws 1640 1614 1641 USE arrays_3d, &1615 USE arrays_3d, & 1642 1616 ONLY: rdf_sc, tend, tsa_m 1643 1617 1644 USE control_parameters, & 1645 ONLY: bc_dirichlet_l, & 1646 bc_dirichlet_n, & 1647 bc_dirichlet_r, & 1648 bc_dirichlet_s, & 1649 bc_radiation_l, & 1650 bc_radiation_n, & 1651 bc_radiation_r, & 1652 bc_radiation_s, & 1653 dt_3d, intermediate_timestep_count, & 1654 intermediate_timestep_count_max, scalar_advec, simulated_time, & 1655 timestep_scheme, tsc, ws_scheme_sca 1656 1657 USE cpulog, & 1618 USE control_parameters, & 1619 ONLY: bc_dirichlet_l, & 1620 bc_dirichlet_n, & 1621 bc_dirichlet_r, & 1622 bc_dirichlet_s, & 1623 bc_radiation_l, & 1624 bc_radiation_n, & 1625 bc_radiation_r, & 1626 bc_radiation_s, & 1627 dt_3d, & 1628 intermediate_timestep_count, & 1629 intermediate_timestep_count_max, & 1630 scalar_advec, & 1631 simulated_time, & 1632 timestep_scheme, & 1633 tsc, & 1634 ws_scheme_sca 1635 1636 USE cpulog, & 1658 1637 ONLY: cpu_log, log_point_s 1659 1638 1660 USE diffusion_s_mod, &1639 USE diffusion_s_mod, & 1661 1640 ONLY: diffusion_s 1662 1641 … … 1669 1648 REAL(wp) :: sbt !< weighting factor for sub-time step 1670 1649 1650 1671 1651 ! 1672 1652 !-- Switch of the surface heat flux, if requested 1673 1653 IF ( surface_cooling_spinup_time /= 999999.9_wp ) THEN 1674 IF ( .NOT. surface_cooling_switched_off .AND. &1654 IF ( .NOT. surface_cooling_switched_off .AND. & 1675 1655 simulated_time >= surface_cooling_spinup_time ) THEN 1676 1656 … … 1682 1662 1683 1663 ! 1684 !-- Compute prognostic equations for the ocean mode 1685 !-- First, start with salinity 1664 !-- Compute prognostic equations for the ocean mode. 1665 !-- First, start with salinity. 1686 1666 IF ( salinity ) THEN 1687 1667 … … 1709 1689 IF ( timestep_scheme(1:5) == 'runge' ) THEN 1710 1690 IF ( ws_scheme_sca ) THEN 1711 CALL advec_s_ws( advc_flags_s, sa, 'sa', &1712 bc_dirichlet_l .OR. bc_radiation_l, &1713 bc_dirichlet_n .OR. bc_radiation_n, &1714 bc_dirichlet_r .OR. bc_radiation_r, &1691 CALL advec_s_ws( advc_flags_s, sa, 'sa', & 1692 bc_dirichlet_l .OR. bc_radiation_l, & 1693 bc_dirichlet_n .OR. bc_radiation_n, & 1694 bc_dirichlet_r .OR. bc_radiation_r, & 1715 1695 bc_dirichlet_s .OR. bc_radiation_s ) 1716 1696 ELSE … … 1722 1702 ENDIF 1723 1703 1724 CALL diffusion_s( sa, &1725 surf_def_h(0)%sasws, surf_def_h(1)%sasws, &1726 surf_def_h(2)%sasws, &1727 surf_lsm_h(0)%sasws, surf_lsm_h(1)%sasws, &1728 surf_usm_h(0)%sasws, surf_usm_h(1)%sasws, &1729 surf_def_v(0)%sasws, surf_def_v(1)%sasws, &1730 surf_def_v(2)%sasws, surf_def_v(3)%sasws, &1731 surf_lsm_v(0)%sasws, surf_lsm_v(1)%sasws, &1732 surf_lsm_v(2)%sasws, surf_lsm_v(3)%sasws, &1733 surf_usm_v(0)%sasws, surf_usm_v(1)%sasws, &1704 CALL diffusion_s( sa, & 1705 surf_def_h(0)%sasws, surf_def_h(1)%sasws, & 1706 surf_def_h(2)%sasws, & 1707 surf_lsm_h(0)%sasws, surf_lsm_h(1)%sasws, & 1708 surf_usm_h(0)%sasws, surf_usm_h(1)%sasws, & 1709 surf_def_v(0)%sasws, surf_def_v(1)%sasws, & 1710 surf_def_v(2)%sasws, surf_def_v(3)%sasws, & 1711 surf_lsm_v(0)%sasws, surf_lsm_v(1)%sasws, & 1712 surf_lsm_v(2)%sasws, surf_lsm_v(3)%sasws, & 1713 surf_usm_v(0)%sasws, surf_usm_v(1)%sasws, & 1734 1714 surf_usm_v(2)%sasws, surf_usm_v(3)%sasws ) 1735 1715 … … 1743 1723 !DIR$ IVDEP 1744 1724 DO k = nzb+1, nzt 1745 sa_p(k,j,i) = sa(k,j,i) + ( dt_3d * ( sbt * tend(k,j,i) + & 1746 tsc(3) * tsa_m(k,j,i) ) & 1747 - tsc(5) * rdf_sc(k) * & 1748 ( sa(k,j,i) - sa_init(k) ) & 1749 ) & 1750 * MERGE( 1.0_wp, 0.0_wp, & 1751 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1725 sa_p(k,j,i) = sa(k,j,i) + ( dt_3d * ( sbt * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) ) & 1726 - tsc(5) * rdf_sc(k) * ( sa(k,j,i) - sa_init(k) ) & 1727 ) & 1728 * MERGE( 1.0_wp, 0.0_wp, & 1729 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1752 1730 ) 1753 1731 IF ( sa_p(k,j,i) < 0.0_wp ) sa_p(k,j,i) = 0.1_wp * sa(k,j,i) … … 1767 1745 ENDDO 1768 1746 ENDDO 1769 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) & 1770 THEN 1747 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1771 1748 DO i = nxl, nxr 1772 1749 DO j = nys, nyn 1773 1750 DO k = nzb+1, nzt 1774 tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 1775 5.3125_wp * tsa_m(k,j,i) 1751 tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tsa_m(k,j,i) 1776 1752 ENDDO 1777 1753 ENDDO … … 1793 1769 1794 1770 1795 !------------------------------------------------------------------------------ !1771 !--------------------------------------------------------------------------------------------------! 1796 1772 ! Description: 1797 1773 ! ------------ 1798 1774 !> Prognostic equations for ocean mode (so far, salinity only) 1799 1775 !> Cache-optimized version 1800 !------------------------------------------------------------------------------ !1776 !--------------------------------------------------------------------------------------------------! 1801 1777 SUBROUTINE ocean_prognostic_equations_ij( i, j, i_omp_start, tn ) 1802 1778 1803 USE advec_s_pw_mod, &1779 USE advec_s_pw_mod, & 1804 1780 ONLY: advec_s_pw 1805 1781 1806 USE advec_s_up_mod, &1782 USE advec_s_up_mod, & 1807 1783 ONLY: advec_s_up 1808 1784 1809 USE advec_ws, &1785 USE advec_ws, & 1810 1786 ONLY: advec_s_ws 1811 1787 1812 USE arrays_3d, &1788 USE arrays_3d, & 1813 1789 ONLY: diss_l_sa, diss_s_sa, flux_l_sa, flux_s_sa, rdf_sc, tend, tsa_m 1814 1790 1815 USE control_parameters, & 1816 ONLY: bc_dirichlet_l, & 1817 bc_dirichlet_n, & 1818 bc_dirichlet_r, & 1819 bc_dirichlet_s, & 1820 bc_radiation_l, & 1821 bc_radiation_n, & 1822 bc_radiation_r, & 1823 bc_radiation_s, & 1824 dt_3d, intermediate_timestep_count, & 1825 intermediate_timestep_count_max, simulated_time, & 1826 timestep_scheme, tsc, ws_scheme_sca 1827 1828 USE diffusion_s_mod, & 1791 USE control_parameters, & 1792 ONLY: bc_dirichlet_l, & 1793 bc_dirichlet_n, & 1794 bc_dirichlet_r, & 1795 bc_dirichlet_s, & 1796 bc_radiation_l, & 1797 bc_radiation_n, & 1798 bc_radiation_r, & 1799 bc_radiation_s, & 1800 dt_3d, & 1801 intermediate_timestep_count, & 1802 intermediate_timestep_count_max, & 1803 simulated_time, & 1804 timestep_scheme, & 1805 tsc, & 1806 ws_scheme_sca 1807 1808 USE diffusion_s_mod, & 1829 1809 ONLY: diffusion_s 1830 1810 … … 1832 1812 1833 1813 INTEGER(iwp) :: i !< loop index x direction 1834 INTEGER(iwp) :: i_omp_start !< first loop index of i-loop in calling & 1835 !< routine prognostic_equations 1814 INTEGER(iwp) :: i_omp_start !< first loop index of i-loop in calling routine prognostic_equations 1836 1815 INTEGER(iwp) :: j !< loop index y direction 1837 1816 INTEGER(iwp) :: k !< loop index z direction … … 1842 1821 !-- Switch of the surface heat flux, if requested 1843 1822 IF ( surface_cooling_spinup_time /= 999999.9_wp ) THEN 1844 IF ( .NOT. surface_cooling_switched_off .AND. &1823 IF ( .NOT. surface_cooling_switched_off .AND. & 1845 1824 simulated_time >= surface_cooling_spinup_time ) THEN 1846 1825 … … 1852 1831 1853 1832 ! 1854 !-- Compute prognostic equations for the ocean mode 1855 !-- First, start with tendency-terms for salinity 1833 !-- Compute prognostic equations for the ocean mode. 1834 !-- First, start with tendency-terms for salinity. 1856 1835 IF ( salinity ) THEN 1857 1836 1858 1837 tend(:,j,i) = 0.0_wp 1859 IF ( timestep_scheme(1:5) == 'runge' ) & 1860 THEN 1838 IF ( timestep_scheme(1:5) == 'runge' ) THEN 1861 1839 IF ( ws_scheme_sca ) THEN 1862 CALL advec_s_ws( advc_flags_s, &1863 i, j, sa, 'sa', flux_s_sa, diss_s_sa, flux_l_sa, &1864 diss_l_sa, i_omp_start, tn,&1865 bc_dirichlet_l .OR. bc_radiation_l, &1866 bc_dirichlet_n .OR. bc_radiation_n, &1867 bc_dirichlet_r .OR. bc_radiation_r, &1840 CALL advec_s_ws( advc_flags_s, & 1841 i, j, sa, 'sa', flux_s_sa, diss_s_sa, flux_l_sa, diss_l_sa, & 1842 i_omp_start, tn, & 1843 bc_dirichlet_l .OR. bc_radiation_l, & 1844 bc_dirichlet_n .OR. bc_radiation_n, & 1845 bc_dirichlet_r .OR. bc_radiation_r, & 1868 1846 bc_dirichlet_s .OR. bc_radiation_s ) 1869 1847 ELSE … … 1873 1851 CALL advec_s_up( i, j, sa ) 1874 1852 ENDIF 1875 CALL diffusion_s( i, j, sa, & 1876 surf_def_h(0)%sasws, surf_def_h(1)%sasws, & 1877 surf_def_h(2)%sasws, & 1878 surf_lsm_h(0)%sasws, surf_lsm_h(1)%sasws, & 1879 surf_usm_h(0)%sasws, surf_usm_h(1)%sasws, & 1880 surf_def_v(0)%sasws, surf_def_v(1)%sasws, & 1881 surf_def_v(2)%sasws, surf_def_v(3)%sasws, & 1882 surf_lsm_v(0)%sasws, surf_lsm_v(1)%sasws, & 1883 surf_lsm_v(2)%sasws, surf_lsm_v(3)%sasws, & 1884 surf_usm_v(0)%sasws, surf_usm_v(1)%sasws, & 1885 surf_usm_v(2)%sasws, surf_usm_v(3)%sasws ) 1853 CALL diffusion_s( i, j, sa, & 1854 surf_def_h(0)%sasws, surf_def_h(1)%sasws, surf_def_h(2)%sasws, & 1855 surf_lsm_h(0)%sasws, surf_lsm_h(1)%sasws, & 1856 surf_usm_h(0)%sasws, surf_usm_h(1)%sasws, & 1857 surf_def_v(0)%sasws, surf_def_v(1)%sasws, surf_def_v(2)%sasws, & 1858 surf_def_v(3)%sasws, & 1859 surf_lsm_v(0)%sasws, surf_lsm_v(1)%sasws, surf_lsm_v(2)%sasws, & 1860 surf_lsm_v(3)%sasws, & 1861 surf_usm_v(0)%sasws, surf_usm_v(1)%sasws, surf_usm_v(2)%sasws, & 1862 surf_usm_v(3)%sasws ) 1886 1863 1887 1864 ! CALL user_actions( i, j, 'sa-tendency' ) ToDo: find general solution for dependency between modules … … 1891 1868 DO k = nzb+1, nzt 1892 1869 1893 sa_p(k,j,i) = sa(k,j,i) + ( dt_3d * & 1894 ( tsc(2) * tend(k,j,i) + & 1895 tsc(3) * tsa_m(k,j,i) ) & 1896 - tsc(5) * rdf_sc(k) & 1897 * ( sa(k,j,i) - sa_init(k) ) & 1898 ) * MERGE( 1.0_wp, 0.0_wp, & 1899 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1870 sa_p(k,j,i) = sa(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) ) & 1871 - tsc(5) * rdf_sc(k) * ( sa(k,j,i) - sa_init(k) ) & 1872 ) * MERGE( 1.0_wp, 0.0_wp, & 1873 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1874 ) 1900 1875 1901 1876 IF ( sa_p(k,j,i) < 0.0_wp ) sa_p(k,j,i) = 0.1_wp * sa(k,j,i) … … 1910 1885 tsa_m(k,j,i) = tend(k,j,i) 1911 1886 ENDDO 1912 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) & 1913 THEN 1887 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1914 1888 DO k = nzb+1, nzt 1915 tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 1916 5.3125_wp * tsa_m(k,j,i) 1889 tsa_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tsa_m(k,j,i) 1917 1890 ENDDO 1918 1891 ENDIF … … 1927 1900 END SUBROUTINE ocean_prognostic_equations_ij 1928 1901 1929 !------------------------------------------------------------------------------ !1902 !--------------------------------------------------------------------------------------------------! 1930 1903 ! Description: 1931 1904 ! ------------ 1932 1905 !> Boundary conditions for ocean model 1933 !------------------------------------------------------------------------------ !1906 !--------------------------------------------------------------------------------------------------! 1934 1907 SUBROUTINE ocean_boundary_conditions 1935 1908 … … 1943 1916 1944 1917 ! 1945 !-- Boundary conditions for salinity 1946 !-- Bottom boundary: Neumann condition because salinity flux is always 1947 !-- given. 1918 !-- Boundary conditions for salinity. 1919 !-- Bottom boundary: Neumann condition because salinity flux is always given. 1948 1920 DO l = 0, 1 1949 1921 !$OMP PARALLEL DO PRIVATE( i, j, k ) … … 1965 1937 END SUBROUTINE ocean_boundary_conditions 1966 1938 1967 !------------------------------------------------------------------------------ !1939 !--------------------------------------------------------------------------------------------------! 1968 1940 ! Description: 1969 1941 ! ------------ 1970 1942 !> Swapping of timelevels. 1971 !------------------------------------------------------------------------------ !1943 !--------------------------------------------------------------------------------------------------! 1972 1944 SUBROUTINE ocean_swap_timelevel( mod_count ) 1973 1945 … … 1994 1966 1995 1967 1996 !------------------------------------------------------------------------------ !1968 !--------------------------------------------------------------------------------------------------! 1997 1969 ! Description: 1998 1970 ! ------------ 1999 1971 !> Read module-specific global restart data (Fortran binary format). 2000 !------------------------------------------------------------------------------ !1972 !--------------------------------------------------------------------------------------------------! 2001 1973 SUBROUTINE ocean_rrd_global_ftn( found ) 2002 1974 2003 1975 2004 USE control_parameters, &1976 USE control_parameters, & 2005 1977 ONLY: length, restart_string 2006 1978 … … 2063 2035 2064 2036 2065 !------------------------------------------------------------------------------ !2037 !--------------------------------------------------------------------------------------------------! 2066 2038 ! Description: 2067 2039 ! ------------ 2068 2040 !> Read module-specific global restart data (MPI-IO). 2069 !------------------------------------------------------------------------------ !2041 !--------------------------------------------------------------------------------------------------! 2070 2042 SUBROUTINE ocean_rrd_global_mpi 2071 2043 … … 2088 2060 2089 2061 2090 !------------------------------------------------------------------------------ !2062 !--------------------------------------------------------------------------------------------------! 2091 2063 ! Description: 2092 2064 ! ------------ 2093 2065 !> Read module-specific local restart data arrays (Fortran binary format). 2094 !------------------------------------------------------------------------------! 2095 SUBROUTINE ocean_rrd_local_ftn( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, & 2096 nxr_on_file, nynf, nync, nyn_on_file, nysf, & 2097 nysc, nys_on_file, tmp_3d, found ) 2098 2099 USE averaging, & 2066 !--------------------------------------------------------------------------------------------------! 2067 SUBROUTINE ocean_rrd_local_ftn( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, nxr_on_file, nynf, nync, & 2068 nyn_on_file, nysf, nysc, nys_on_file, tmp_3d, found ) 2069 2070 USE averaging, & 2100 2071 ONLY: rho_ocean_av, sa_av 2101 2072 2102 USE control_parameters, &2073 USE control_parameters, & 2103 2074 ONLY: length, restart_string 2104 2075 2105 USE indices, &2076 USE indices, & 2106 2077 ONLY: nbgp, nxlg, nxrg, nyng, nysg, nzb, nzt 2107 2078 … … 2127 2098 LOGICAL, INTENT(OUT) :: found 2128 2099 2129 REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d !< 2100 REAL(wp), & 2101 DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) & 2102 :: tmp_3d !< 2130 2103 2131 2104 … … 2139 2112 ENDIF 2140 2113 IF ( k == 1 ) READ ( 13 ) tmp_3d 2141 rho_ocean_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &2142 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)2114 rho_ocean_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 2115 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 2143 2116 2144 2117 CASE ( 'sa' ) 2145 2118 IF ( k == 1 ) READ ( 13 ) tmp_3d 2146 sa(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &2147 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)2119 sa(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 2120 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 2148 2121 2149 2122 CASE ( 'sa_av' ) … … 2152 2125 ENDIF 2153 2126 IF ( k == 1 ) READ ( 13 ) tmp_3d 2154 sa_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &2155 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)2127 sa_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 2128 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 2156 2129 2157 2130 CASE DEFAULT … … 2163 2136 2164 2137 2165 !------------------------------------------------------------------------------ !2138 !--------------------------------------------------------------------------------------------------! 2166 2139 ! Description: 2167 2140 ! ------------ 2168 2141 !> Read module-specific local restart data arrays (MPI-IO). 2169 !------------------------------------------------------------------------------ !2142 !--------------------------------------------------------------------------------------------------! 2170 2143 SUBROUTINE ocean_rrd_local_mpi 2171 2144 2172 USE averaging, &2145 USE averaging, & 2173 2146 ONLY: rho_ocean_av, sa_av 2174 2147 2175 USE indices, &2148 USE indices, & 2176 2149 ONLY: nxlg, nxrg, nyng, nysg, nzb, nzt 2177 2150 … … 2183 2156 CALL rd_mpi_io_check_array( 'rho_ocean_av' , found = array_found ) 2184 2157 IF ( array_found ) THEN 2185 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) ALLOCATE( rho_ocean_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 2158 IF ( .NOT. ALLOCATED( rho_ocean_av ) ) & 2159 ALLOCATE( rho_ocean_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 2186 2160 CALL rrd_mpi_io( 'rho_ocean_av', rho_ocean_av ) 2187 2161 ENDIF … … 2198 2172 2199 2173 2200 !------------------------------------------------------------------------------ !2174 !--------------------------------------------------------------------------------------------------! 2201 2175 ! Description: 2202 2176 ! ------------ 2203 2177 !> This routine writes the respective restart data for the ocean module. 2204 !------------------------------------------------------------------------------ !2178 !--------------------------------------------------------------------------------------------------! 2205 2179 SUBROUTINE ocean_wrd_global 2206 2180 … … 2258 2232 CALL wrd_mpi_io_global_array( 'sa_vertical_gradient', sa_vertical_gradient ) 2259 2233 CALL wrd_mpi_io_global_array( 'sa_vertical_gradient_level', sa_vertical_gradient_level ) 2260 CALL wrd_mpi_io_global_array( 'sa_vertical_gradient_level_ind', sa_vertical_gradient_level_ind ) 2234 CALL wrd_mpi_io_global_array( 'sa_vertical_gradient_level_ind', & 2235 sa_vertical_gradient_level_ind ) 2261 2236 CALL wrd_mpi_io( 'stokes_waveheight', stokes_waveheight ) 2262 2237 CALL wrd_mpi_io( 'stokes_wavelength', stokes_wavelength ) … … 2271 2246 2272 2247 2273 !------------------------------------------------------------------------------ !2248 !--------------------------------------------------------------------------------------------------! 2274 2249 ! Description: 2275 2250 ! ------------ 2276 2251 !> This routine writes the respective restart data for the ocean module. 2277 !------------------------------------------------------------------------------ !2252 !--------------------------------------------------------------------------------------------------! 2278 2253 SUBROUTINE ocean_wrd_local 2279 2254 2280 USE averaging, &2255 USE averaging, & 2281 2256 ONLY: rho_ocean_av, sa_av 2282 2257 … … 2310 2285 2311 2286 2312 !------------------------------------------------------------------------------ !2313 ! Description: 2314 ! ------------ 2315 !> This routine calculates the Craik Leibovich vortex force and the additional 2316 !> effect of the Stokes drift on the Coriolis force2287 !--------------------------------------------------------------------------------------------------! 2288 ! Description: 2289 ! ------------ 2290 !> This routine calculates the Craik Leibovich vortex force and the additional effect of the Stokes 2291 !> drift on the Coriolis force. 2317 2292 !> Call for all gridpoints. 2318 !------------------------------------------------------------------------------ !2293 !--------------------------------------------------------------------------------------------------! 2319 2294 SUBROUTINE stokes_drift_terms( component ) 2320 2295 2321 USE arrays_3d, & 2322 ONLY: ddzu, u, u_stokes_zu, u_stokes_zw, v, v_stokes_zu, & 2323 v_stokes_zw, w, tend 2324 2325 USE basic_constants_and_equations_mod, & 2296 USE arrays_3d, & 2297 ONLY: ddzu, u, u_stokes_zu, u_stokes_zw, v, v_stokes_zu, v_stokes_zw, w, tend 2298 2299 USE basic_constants_and_equations_mod, & 2326 2300 ONLY: pi 2327 2301 2328 USE control_parameters, &2302 USE control_parameters, & 2329 2303 ONLY: f, fs, message_string, rotation_angle 2330 2304 2331 USE grid_variables, &2305 USE grid_variables, & 2332 2306 ONLY: ddx, ddy 2333 2307 2334 USE indices, &2308 USE indices, & 2335 2309 ONLY: nxl, nxr, nys, nysv, nyn, nzb, nzt 2336 2310 … … 2344 2318 REAL(wp) :: cos_rot_angle !< cosine of model rotation angle 2345 2319 REAL(wp) :: sin_rot_angle !< sine of model rotation angle 2320 2346 2321 2347 2322 ! … … 2355 2330 DO j = nysv, nyn 2356 2331 DO k = nzb+1, nzt 2357 tend(k,j,i) = tend(k,j,i) + v_stokes_zu(k) * ( &2358 0.5 * ( v(k,j+1,i) - v(k,j+1,i-1)&2359 + v(k,j,i) - v(k,j,i-1) ) * ddx&2360 - 0.5 * ( u(k,j+1,i) - u(k,j-1,i) ) * ddy&2361 ) &2332 tend(k,j,i) = tend(k,j,i) + v_stokes_zu(k) * ( & 2333 0.5 * ( v(k,j+1,i) - v(k,j+1,i-1) & 2334 + v(k,j,i) - v(k,j,i-1) ) * ddx & 2335 - 0.5 * ( u(k,j+1,i) - u(k,j-1,i) ) * ddy & 2336 ) & 2362 2337 + f * v_stokes_zu(k) 2363 2338 ENDDO … … 2371 2346 DO j = nysv, nyn 2372 2347 DO k = nzb+1, nzt 2373 tend(k,j,i) = tend(k,j,i) - u_stokes_zu(k) * ( &2374 0.5 * ( v(k,j,i+1) - v(k,j,i-1) ) * ddx&2375 - 0.5 * ( u(k,j,i) - u(k,j-1,i)&2376 + u(k,j,i+1) - u(k,j-1,i+1) ) * ddy&2348 tend(k,j,i) = tend(k,j,i) - u_stokes_zu(k) * ( & 2349 0.5 * ( v(k,j,i+1) - v(k,j,i-1) ) * ddx & 2350 - 0.5 * ( u(k,j,i) - u(k,j-1,i) & 2351 + u(k,j,i+1) - u(k,j-1,i+1) ) * ddy & 2377 2352 ) & 2378 2353 - f * u_stokes_zu(k) … … 2393 2368 DO j = nys, nyn 2394 2369 DO k = nzb+1, nzt 2395 tend(k,j,i) = tend(k,j,i) + u_stokes_zw(k) * ( & 2396 0.5 * ( u(k+1,j,i) - u(k,j,i) & 2397 + u(k+1,j,i+1) - u(k,j,i+1) & 2398 ) * ddzu(k+1) & 2399 - 0.5 * ( w(k,j,i+1) - w(k,j,i-1) & 2400 ) * ddx ) & 2401 - v_stokes_zw(k) * ( & 2402 0.5 * ( w(k,j+1,i) - w(k,j-1,i) & 2403 ) * ddy & 2404 - 0.5 * ( v(k+1,j,i) - v(k,j,i) & 2405 + v(k+1,j+1,i) - v(k,j+1,i) & 2406 ) * ddzu(k) ) & 2407 + fs * ( & 2408 sin_rot_angle * v_stokes_zw(k) & 2409 + cos_rot_angle * u_stokes_zw(k) & 2410 ) 2370 tend(k,j,i) = tend(k,j,i) + u_stokes_zw(k) * ( & 2371 0.5 * ( u(k+1,j,i) - u(k,j,i) & 2372 + u(k+1,j,i+1) - u(k,j,i+1) & 2373 ) * ddzu(k+1) & 2374 - 0.5 * ( w(k,j,i+1) - w(k,j,i-1) & 2375 ) * ddx ) & 2376 - v_stokes_zw(k) * ( & 2377 0.5 * ( w(k,j+1,i) - w(k,j-1,i) & 2378 ) * ddy & 2379 - 0.5 * ( v(k+1,j,i) - v(k,j,i) & 2380 + v(k+1,j+1,i) - v(k,j+1,i) & 2381 ) * ddzu(k) ) & 2382 + fs * ( sin_rot_angle * v_stokes_zw(k) & 2383 + cos_rot_angle * u_stokes_zw(k) & 2384 ) 2411 2385 ENDDO 2412 2386 ENDDO … … 2414 2388 2415 2389 CASE DEFAULT 2416 WRITE( message_string, * ) 'wrong component of Stokes force: ', & 2417 component 2390 WRITE( message_string, * ) 'wrong component of Stokes force: ', component 2418 2391 CALL message( 'stokes_drift_terms', 'PA0091', 1, 2, 0, 6, 0 ) 2419 2392 … … 2423 2396 2424 2397 2425 !------------------------------------------------------------------------------ !2426 ! Description: 2427 ! ------------ 2428 !> This routine calculates the Craik Leibovich vortex force and the additional 2429 !> effect of the Stokes drift on the Coriolis force2398 !--------------------------------------------------------------------------------------------------! 2399 ! Description: 2400 ! ------------ 2401 !> This routine calculates the Craik Leibovich vortex force and the additional effect of the Stokes 2402 !> drift on the Coriolis force. 2430 2403 !> Call for gridpoints i,j. 2431 !------------------------------------------------------------------------------ !2404 !--------------------------------------------------------------------------------------------------! 2432 2405 2433 2406 SUBROUTINE stokes_drift_terms_ij( i, j, component ) 2434 2407 2435 USE arrays_3d, & 2436 ONLY: ddzu, u, u_stokes_zu, u_stokes_zw, v, v_stokes_zu, & 2437 v_stokes_zw, w, tend 2438 2439 USE basic_constants_and_equations_mod, & 2408 USE arrays_3d, & 2409 ONLY: ddzu, u, u_stokes_zu, u_stokes_zw, v, v_stokes_zu, v_stokes_zw, w, tend 2410 2411 USE basic_constants_and_equations_mod, & 2440 2412 ONLY: pi 2441 2413 2442 USE control_parameters, &2414 USE control_parameters, & 2443 2415 ONLY: f, fs, message_string, rotation_angle 2444 2416 2445 USE grid_variables, &2417 USE grid_variables, & 2446 2418 ONLY: ddx, ddy 2447 2419 2448 USE indices, &2420 USE indices, & 2449 2421 ONLY: nzb, nzt 2450 2422 … … 2459 2431 REAL(wp) :: sin_rot_angle !< sine of model rotation angle 2460 2432 2433 2461 2434 ! 2462 2435 !-- Compute Stokes terms for the respective velocity components … … 2467 2440 CASE ( 1 ) 2468 2441 DO k = nzb+1, nzt 2469 tend(k,j,i) = tend(k,j,i) + v_stokes_zu(k) * ( &2470 0.5 * ( v(k,j+1,i) - v(k,j+1,i-1)&2471 + v(k,j,i) - v(k,j,i-1) ) * ddx&2472 - 0.5 * ( u(k,j+1,i) - u(k,j-1,i) ) * ddy&2473 ) &2442 tend(k,j,i) = tend(k,j,i) + v_stokes_zu(k) * ( & 2443 0.5 * ( v(k,j+1,i) - v(k,j+1,i-1) & 2444 + v(k,j,i) - v(k,j,i-1) ) * ddx & 2445 - 0.5 * ( u(k,j+1,i) - u(k,j-1,i) ) * ddy & 2446 ) & 2474 2447 + f * v_stokes_zu(k) 2475 2448 ENDDO … … 2478 2451 CASE ( 2 ) 2479 2452 DO k = nzb+1, nzt 2480 tend(k,j,i) = tend(k,j,i) - u_stokes_zu(k) * ( &2481 0.5 * ( v(k,j,i+1) - v(k,j,i-1) ) * ddx&2482 - 0.5 * ( u(k,j,i) - u(k,j-1,i)&2483 + u(k,j,i+1) - u(k,j-1,i+1) ) * ddy&2484 ) &2453 tend(k,j,i) = tend(k,j,i) - u_stokes_zu(k) * ( & 2454 0.5 * ( v(k,j,i+1) - v(k,j,i-1) ) * ddx & 2455 - 0.5 * ( u(k,j,i) - u(k,j-1,i) & 2456 + u(k,j,i+1) - u(k,j-1,i+1) ) * ddy & 2457 ) & 2485 2458 - f * u_stokes_zu(k) 2486 2459 ENDDO … … 2496 2469 2497 2470 DO k = nzb+1, nzt 2498 tend(k,j,i) = tend(k,j,i) + u_stokes_zw(k) * ( & 2499 0.5 * ( u(k+1,j,i) - u(k,j,i) & 2500 + u(k+1,j,i+1) - u(k,j,i+1) & 2501 ) * ddzu(k+1) & 2502 - 0.5 * ( w(k,j,i+1) - w(k,j,i-1) & 2503 ) * ddx ) & 2504 - v_stokes_zw(k) * ( & 2505 0.5 * ( w(k,j+1,i) - w(k,j-1,i) & 2506 ) * ddy & 2507 - 0.5 * ( v(k+1,j,i) - v(k,j,i) & 2508 + v(k+1,j+1,i) - v(k,j+1,i) & 2509 ) * ddzu(k) ) & 2510 + fs * ( sin_rot_angle * v_stokes_zw(k) & 2511 + cos_rot_angle * u_stokes_zw(k) & 2471 tend(k,j,i) = tend(k,j,i) + u_stokes_zw(k) * ( 0.5 * ( u(k+1,j,i) - u(k,j,i) & 2472 + u(k+1,j,i+1) - u(k,j,i+1) & 2473 ) * ddzu(k+1) & 2474 - 0.5 * ( w(k,j,i+1) - w(k,j,i-1) & 2475 ) * ddx ) & 2476 - v_stokes_zw(k) * ( 0.5 * ( w(k,j+1,i) - w(k,j-1,i) & 2477 ) * ddy & 2478 - 0.5 * ( v(k+1,j,i) - v(k,j,i) & 2479 + v(k+1,j+1,i) - v(k,j+1,i) & 2480 ) * ddzu(k) ) & 2481 + fs * ( sin_rot_angle * v_stokes_zw(k) & 2482 + cos_rot_angle * u_stokes_zw(k) & 2512 2483 ) 2513 2484 ENDDO … … 2522 2493 2523 2494 2524 !------------------------------------------------------------------------------ !2525 ! Description: 2526 ! ------------ 2527 !> This routine calculates turbulence generated by wave breaking near the ocean 2528 !> surface, followinga parameterization given in Noh et al. (2004), JPO2495 !--------------------------------------------------------------------------------------------------! 2496 ! Description: 2497 ! ------------ 2498 !> This routine calculates turbulence generated by wave breaking near the ocean surface, following 2499 !> a parameterization given in Noh et al. (2004), JPO 2529 2500 !> Call for all gridpoints. 2530 !> TODO: so far, this routine only works if the model time step has about the 2531 !> s ame value as the time scale of wave breaking!2532 !------------------------------------------------------------------------------ !2501 !> TODO: so far, this routine only works if the model time step has about the same value as the time 2502 !> scale of wave breaking! 2503 !--------------------------------------------------------------------------------------------------! 2533 2504 SUBROUTINE wave_breaking_term( component ) 2534 2505 2535 USE arrays_3d, &2506 USE arrays_3d, & 2536 2507 ONLY: u_p, v_p 2537 2508 2538 USE control_parameters, &2509 USE control_parameters, & 2539 2510 ONLY: dt_3d, message_string 2540 2511 2541 USE indices, &2512 USE indices, & 2542 2513 ONLY: nxl, nxlu, nxr, nys, nysv, nyn, nzt 2543 2514 … … 2548 2519 INTEGER(iwp) :: j !< loop index along y 2549 2520 2550 REAL(wp) :: random_gauss !< function that creates a random number with a 2551 !< Gaussian distribution 2521 REAL(wp) :: random_gauss !< function that creates a random number with a Gaussian distribution 2552 2522 2553 2523 … … 2562 2532 DO i = nxlu, nxr 2563 2533 DO j = nys, nyn 2564 u_p(nzt,j,i) = u_p(nzt,j,i) + &2565 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) &2566 * alpha_wave_breaking * u_star_wave_breaking &2534 u_p(nzt,j,i) = u_p(nzt,j,i) + & 2535 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) & 2536 * alpha_wave_breaking * u_star_wave_breaking & 2567 2537 / timescale_wave_breaking * dt_3d 2568 2538 ENDDO … … 2573 2543 DO i = nxl, nxr 2574 2544 DO j = nysv, nyn 2575 v_p(nzt,j,i) = v_p(nzt,j,i) + &2576 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) &2577 * alpha_wave_breaking * u_star_wave_breaking &2545 v_p(nzt,j,i) = v_p(nzt,j,i) + & 2546 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) & 2547 * alpha_wave_breaking * u_star_wave_breaking & 2578 2548 / timescale_wave_breaking * dt_3d 2579 2549 ENDDO … … 2581 2551 2582 2552 CASE DEFAULT 2583 WRITE( message_string, * ) 'wrong component of wave breaking: ', & 2584 component 2553 WRITE( message_string, * ) 'wrong component of wave breaking: ', component 2585 2554 CALL message( 'stokes_drift_terms', 'PA0466', 1, 2, 0, 6, 0 ) 2586 2555 … … 2590 2559 2591 2560 2592 !------------------------------------------------------------------------------ !2593 ! Description: 2594 ! ------------ 2595 !> This routine calculates turbulence generated by wave breaking near the ocean 2596 !> surface, following aparameterization given in Noh et al. (2004), JPO2561 !--------------------------------------------------------------------------------------------------! 2562 ! Description: 2563 ! ------------ 2564 !> This routine calculates turbulence generated by wave breaking near the ocean surface, following a 2565 !> parameterization given in Noh et al. (2004), JPO 2597 2566 !> Call for gridpoint i,j. 2598 !> TODO: so far, this routine only works if the model time step has about the 2599 !> same value as the timescale of wave breaking!2600 !------------------------------------------------------------------------------ !2567 !> TODO: so far, this routine only works if the model time step has about the same value as the time 2568 !> scale of wave breaking! 2569 !--------------------------------------------------------------------------------------------------! 2601 2570 SUBROUTINE wave_breaking_term_ij( i, j, component ) 2602 2571 2603 USE arrays_3d, &2572 USE arrays_3d, & 2604 2573 ONLY: u_p, v_p 2605 2574 2606 USE control_parameters, &2575 USE control_parameters, & 2607 2576 ONLY: dt_3d, message_string 2608 2577 2609 USE indices, &2578 USE indices, & 2610 2579 ONLY: nzt 2611 2580 … … 2616 2585 INTEGER(iwp) :: j !< loop index along y 2617 2586 2618 REAL(wp) :: random_gauss !< function that creates a random number with a 2619 !< Gaussian distribution 2587 REAL(wp) :: random_gauss !< function that creates a random number with a Gaussian distribution 2620 2588 2621 2589 ! … … 2626 2594 !-- u-/v-component 2627 2595 CASE ( 1 ) 2628 u_p(nzt,j,i) = u_p(nzt,j,i) + &2629 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) &2630 * alpha_wave_breaking * u_star_wave_breaking &2596 u_p(nzt,j,i) = u_p(nzt,j,i) + & 2597 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) & 2598 * alpha_wave_breaking * u_star_wave_breaking & 2631 2599 / timescale_wave_breaking * dt_3d 2632 2600 2633 2601 CASE ( 2 ) 2634 v_p(nzt,j,i) = v_p(nzt,j,i) + &2635 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) &2636 * alpha_wave_breaking * u_star_wave_breaking &2602 v_p(nzt,j,i) = v_p(nzt,j,i) + & 2603 ( random_gauss( iran_ocean, 1.0_wp ) - 1.0_wp ) & 2604 * alpha_wave_breaking * u_star_wave_breaking & 2637 2605 / timescale_wave_breaking * dt_3d 2638 2606 2639 2607 CASE DEFAULT 2640 WRITE( message_string, * ) 'wrong component of wave breaking: ', & 2641 component 2608 WRITE( message_string, * ) 'wrong component of wave breaking: ', component 2642 2609 CALL message( 'stokes_drift_terms', 'PA0466', 1, 2, 0, 6, 0 ) 2643 2610 -
palm/trunk/SOURCE/outflow_turbulence.f90
r4360 r4797 1 1 !> @file outflow_turbulence.f90 2 !------------------------------------------------------------------------------ !2 !--------------------------------------------------------------------------------------------------! 3 3 ! This file is part of the PALM model system. 4 4 ! 5 ! PALM is free software: you can redistribute it and/or modify it under the 6 ! terms of the GNU General Public License as published by the Free Software 7 ! Foundation, either version 3 of the License, or (at your option) any later 8 ! version. 5 ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General 6 ! Public License as published by the Free Software Foundation, either version 3 of the License, or 7 ! (at your option) any later version. 9 8 ! 10 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY 11 ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR12 ! A PARTICULAR PURPOSE. See the GNU GeneralPublic License for more details.9 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the 10 ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 ! Public License for more details. 13 12 ! 14 ! You should have received a copy of the GNU General Public License along with 15 ! PALM. If not, see<http://www.gnu.org/licenses/>.13 ! You should have received a copy of the GNU General Public License along with PALM. If not, see 14 ! <http://www.gnu.org/licenses/>. 16 15 ! 17 16 ! Copyright 1997-2020 Leibniz Universitaet Hannover 18 !------------------------------------------------------------------------------ !17 !--------------------------------------------------------------------------------------------------! 19 18 ! 20 19 ! Current revisions: 21 20 ! ----------------- 22 ! 23 ! 21 ! 22 ! 24 23 ! Former revisions: 25 24 ! ----------------- 26 25 ! $Id: outflow_turbulence.f90 3241 2018-09-12 15:02:00Z raasch $ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 3241 2018-09-12 15:02:00Z raasch 27 29 ! Corrected "Former revisions" section 28 ! 30 ! 29 31 ! 3241 2018-09-12 15:02:00Z raasch 30 32 ! unused variables removed … … 32 34 ! 2050 2016-11-08 15:00:55Z gronemeier 33 35 ! Initial version 34 ! 36 ! 35 37 ! 36 38 ! Description: 37 39 ! ------------ 38 !> Routine based on inflow_turbulence.f90. Copies values of 3d data from a 2d 39 !> vertical source plane (defined by outflow_source_plane) to the outflow 40 !> boundary. 41 !------------------------------------------------------------------------------! 40 !> Routine based on inflow_turbulence.f90. Copies values of 3d data from a 2d vertical source plane 41 !> (defined by outflow_source_plane) to the outflow boundary. 42 !--------------------------------------------------------------------------------------------------! 42 43 SUBROUTINE outflow_turbulence 43 44 44 USE arrays_3d, &45 USE arrays_3d, & 45 46 ONLY: e, pt, q, s, u, v, w 46 47 47 USE control_parameters, &48 USE control_parameters, & 48 49 ONLY: humidity, passive_scalar, outflow_source_plane 49 50 50 USE cpulog, &51 USE cpulog, & 51 52 ONLY: cpu_log, log_point 52 53 53 USE grid_variables, &54 USE grid_variables, & 54 55 ONLY: ddx 55 56 56 USE indices, &57 USE indices, & 57 58 ONLY: nbgp, nx, nxr, nyn, nys, nyng, nysg, nzb, nzt 58 59 59 60 USE kinds 60 61 61 USE pegrid!, &62 USE pegrid!, & 62 63 !ONLY: comm1dx, id_outflow, id_outflow_source, ierr, myidx, status 63 64 … … 70 71 INTEGER(iwp) :: ngp_ofv !< number of grid points stored in outflow_val 71 72 72 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,7,nbgp) :: & 73 outflow_val !< values to be copied to the outflow boundary 73 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,7,nbgp) :: outflow_val !< values to be copied to the outflow boundary 74 74 75 75 … … 94 94 outflow_val(k,j,4,l) = pt(k,j,i) 95 95 outflow_val(k,j,5,l) = e(k,j,i) 96 IF ( humidity ) &96 IF ( humidity ) & 97 97 outflow_val(k,j,6,l) = q(k,j,i) 98 IF ( passive_scalar ) &98 IF ( passive_scalar ) & 99 99 outflow_val(k,j,7,l) = s(k,j,i) 100 100 … … 115 115 outflow_val(k,j,4,l) = pt(k,j,i) 116 116 outflow_val(k,j,5,l) = e(k,j,i) 117 IF ( humidity ) &117 IF ( humidity ) & 118 118 outflow_val(k,j,6,l) = q(k,j,i) 119 IF ( passive_scalar ) &119 IF ( passive_scalar ) & 120 120 outflow_val(k,j,7,l) = s(k,j,i) 121 121 … … 131 131 IF ( myidx == id_outflow_source .AND. myidx /= id_outflow ) THEN 132 132 133 CALL MPI_SEND( outflow_val(nzb,nysg,1,1), ngp_ofv, MPI_REAL, & 134 id_outflow, 1, comm1dx, ierr ) 133 CALL MPI_SEND( outflow_val(nzb,nysg,1,1), ngp_ofv, MPI_REAL, id_outflow, 1, comm1dx, ierr ) 135 134 136 135 ELSEIF ( myidx /= id_outflow_source .AND. myidx == id_outflow ) THEN 137 136 138 137 outflow_val = 0.0_wp 139 CALL MPI_RECV( outflow_val(nzb,nysg,1,1), ngp_ofv, MPI_REAL, 140 id_outflow_source, 1, comm1dx,status, ierr )138 CALL MPI_RECV( outflow_val(nzb,nysg,1,1), ngp_ofv, MPI_REAL, id_outflow_source, 1, comm1dx, & 139 status, ierr ) 141 140 142 141 ENDIF … … 157 156 e(k,j,nx+1:nx+nbgp) = MAX( e(k,j,nx+1:nx+nbgp), 0.0_wp ) 158 157 159 IF ( humidity ) &158 IF ( humidity ) & 160 159 q(k,j,nx+1:nx+nbgp) = outflow_val(k,j,6,1:nbgp) 161 IF ( passive_scalar ) &160 IF ( passive_scalar ) & 162 161 s(k,j,nx+1:nx+nbgp) = outflow_val(k,j,7,1:nbgp) 163 162 -
palm/trunk/SOURCE/palm.f90
r4539 r4797 1 1 ! !> @file palm.f90 2 !------------------------------------------------------------------------------ !2 !--------------------------------------------------------------------------------------------------! 3 3 ! This file is part of the PALM model system. 4 4 ! 5 ! PALM is free software: you can redistribute it and/or modify it under the 6 ! terms of the GNU General Public License as published by the Free Software 7 ! Foundation, either version 3 of the License, or (at your option) any later 8 ! version. 9 ! 10 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY 11 ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR 12 ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. 13 ! 14 ! You should have received a copy of the GNU General Public License along with 15 ! PALM. If not, see <http://www.gnu.org/licenses/>. 5 ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General 6 ! Public License as published by the Free Software Foundation, either version 3 of the License, or 7 ! (at your option) any later version. 8 ! 9 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the 10 ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 ! Public License for more details. 12 ! 13 ! You should have received a copy of the GNU General Public License along with PALM. If not, see 14 ! <http://www.gnu.org/licenses/>. 16 15 ! 17 16 ! Copyright 1997-2020 Leibniz Universitaet Hannover 18 !------------------------------------------------------------------------------ !17 !--------------------------------------------------------------------------------------------------! 19 18 ! 20 19 ! Current revisions: 21 20 ! ----------------- 22 ! 21 ! 23 22 ! 24 23 ! Former revisions: 25 24 ! ----------------- 26 25 ! $Id$ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4539 2020-05-18 14:05:17Z raasch 27 29 ! log point name changed 28 ! 30 ! 29 31 ! 4535 2020-05-15 12:07:23Z raasch 30 32 ! bugfix for restart data format query 31 ! 33 ! 32 34 ! 4496 2020-04-15 08:37:26Z raasch 33 35 ! bugfix: coupling character added to restart output filename 34 ! 36 ! 35 37 ! 4495 2020-04-13 20:11:20Z raasch 36 38 ! restart data handling with MPI-IO added 37 ! 39 ! 38 40 ! 4457 2020-03-11 14:20:43Z raasch 39 41 ! use statement for exchange horiz added 40 ! 42 ! 41 43 ! 4444 2020-03-05 15:59:50Z raasch 42 44 ! bugfix: cpp-directives for serial mode added 43 ! 45 ! 44 46 ! 4414 2020-02-19 20:16:04Z suehring 45 47 ! Call to module_interface_init_numerics 46 ! 48 ! 47 49 ! 4400 2020-02-10 20:32:41Z suehring 48 50 ! Add interface to initialize data output with dom 49 ! 51 ! 50 52 ! 4360 2020-01-07 11:25:50Z suehring 51 53 ! implement new palm_date_time_mod 52 ! 54 ! 53 55 ! 4094 2019-07-12 09:24:21Z gronemeier 54 56 ! Corrected "Former revisions" section 55 ! 57 ! 56 58 ! 4039 2019-06-18 10:32:41Z suehring 57 59 ! Rename subroutines in module for diagnostic quantities 58 ! 60 ! 59 61 ! 4017 2019-06-06 12:16:46Z schwenkel 60 62 ! new module for calculation and output of diagnostic quantities added 61 ! 63 ! 62 64 ! 3885 2019-04-11 11:29:34Z kanani 63 ! Changes related to global restructuring of location messages and introduction 65 ! Changes related to global restructuring of location messages and introduction 64 66 ! of additional debug messages 65 ! 67 ! 66 68 ! 3761 2019-02-25 15:31:42Z raasch 67 69 ! unused variable removed 68 ! 70 ! 69 71 ! 3719 2019-02-06 13:10:18Z kanani 70 72 ! Included cpu measurement for wall/soil spinup 71 ! 73 ! 72 74 ! 3703 2019-01-29 16:43:53Z knoop 73 75 ! Some interface calls moved to module_interface + cleanup 74 ! 76 ! 75 77 ! 3648 2019-01-02 16:35:46Z suehring 76 78 ! Rename subroutines for surface-data output … … 82 84 ! Description: 83 85 ! ------------ 84 !> Large-Eddy Simulation (LES) model for atmospheric and oceanic boundary-layer 85 !> flows 86 !> Large-Eddy Simulation (LES) model for atmospheric and oceanic boundary-layer flows, 86 87 !> see the PALM homepage https://palm-model.org for further information 87 !------------------------------------------------------------------------------ !88 !--------------------------------------------------------------------------------------------------! 88 89 PROGRAM palm 89 90 90 91 91 92 USE arrays_3d 92 93 93 94 #if defined( __parallel ) 94 USE bulk_cloud_model_mod, &95 USE bulk_cloud_model_mod, & 95 96 ONLY: bulk_cloud_model, microphysics_morrison, microphysics_seifert 96 97 #endif 97 98 98 USE control_parameters, &99 ONLY: coupling_char, do2d_at_begin, do3d_at_begin, io_blocks, 100 io_group, message_string, restart_data_format_output, runnr, simulated_time_chr, spinup,&101 time_since_reference_point, user_interface_current_revision, &99 USE control_parameters, & 100 ONLY: coupling_char, do2d_at_begin, do3d_at_begin, io_blocks, io_group, message_string, & 101 restart_data_format_output, runnr, simulated_time_chr, spinup, & 102 time_since_reference_point, user_interface_current_revision, & 102 103 user_interface_required_revision, version, write_binary 103 104 104 105 #if defined( __parallel ) 105 USE control_parameters, &106 ONLY: child_domain, constant_diffusion, humidity, 107 initializing_actions, neutral,passive_scalar108 #endif 109 110 USE cpulog, &111 ONLY: cpu_log, log_point, cpu_statistics112 113 #if defined( __parallel ) 114 USE cpulog, &106 USE control_parameters, & 107 ONLY: child_domain, constant_diffusion, humidity, initializing_actions, neutral, & 108 passive_scalar 109 #endif 110 111 USE cpulog, & 112 ONLY: cpu_log, cpu_statistics, log_point 113 114 #if defined( __parallel ) 115 USE cpulog, & 115 116 ONLY: log_point_s 116 117 #endif 117 118 118 USE diagnostic_output_quantities_mod, &119 USE diagnostic_output_quantities_mod, & 119 120 ONLY: doq_calculate 120 121 121 122 #if defined( __parallel ) 122 USE exchange_horiz_mod, &123 USE exchange_horiz_mod, & 123 124 ONLY: exchange_horiz 124 125 125 USE indices, &126 USE indices, & 126 127 ONLY: nbgp 127 128 #endif … … 129 130 USE kinds 130 131 131 USE module_interface, &132 ONLY: module_interface_init_numerics, &133 module_interface_init_output, &132 USE module_interface, & 133 ONLY: module_interface_init_numerics, & 134 module_interface_init_output, & 134 135 module_interface_last_actions 135 136 136 137 137 USE multi_agent_system_mod, &138 USE multi_agent_system_mod, & 138 139 ONLY: agents_active, mas_last_actions 139 140 140 USE netcdf_data_input_mod, &141 ONLY: netcdf_data_input_inquire_file, netcdf_data_input_init, &141 USE netcdf_data_input_mod, & 142 ONLY: netcdf_data_input_inquire_file, netcdf_data_input_init, & 142 143 netcdf_data_input_surface_data, netcdf_data_input_topo 143 144 … … 145 146 146 147 #if defined( __parallel ) 147 USE pmc_particle_interface, &148 USE pmc_particle_interface, & 148 149 ONLY: pmcp_g_alloc_win 149 150 150 USE pmc_interface, &151 ONLY: nested_run, pmci_child_initialize, pmci_init, 152 pmci_ modelconfiguration, pmci_parent_initialize153 #endif 154 155 USE restart_data_mpi_io_mod, &151 USE pmc_interface, & 152 ONLY: nested_run, pmci_child_initialize, pmci_init, pmci_modelconfiguration, & 153 pmci_parent_initialize 154 #endif 155 156 USE restart_data_mpi_io_mod, & 156 157 ONLY: rd_mpi_io_close, rd_mpi_io_open 157 158 158 USE surface_data_output_mod, &159 USE surface_data_output_mod, & 159 160 ONLY: surface_data_output_last_action 160 161 161 USE write_restart_data_mod, &162 USE write_restart_data_mod, & 162 163 ONLY: wrd_global, wrd_local 163 164 … … 172 173 !-- Local variables 173 174 CHARACTER(LEN=9) :: time_to_string !< 174 INTEGER(iwp) :: i !< loop counter for blocked I/O 175 176 INTEGER(iwp) :: i !< loop counter for blocked I/O 175 177 #if defined( __parallel) && defined( _OPENACC ) 176 INTEGER( iwp) :: local_comm !< local communicator (shared memory)177 INTEGER(iwp) :: local_num_procs !< local number of processes178 INTEGER(iwp) :: local_id !< local id179 INTEGER( acc_device_kind) :: device_type !< device type for OpenACC180 INTEGER(iwp) :: num_devices !< number of devices visible to OpenACC181 INTEGER(iwp) :: my_device !< device used by this process178 INTEGER(acc_device_kind) :: device_type !< device type for OpenACC 179 INTEGER(iwp) :: local_comm !< local communicator (shared memory) 180 INTEGER(iwp) :: local_num_procs !< local number of processes 181 INTEGER(iwp) :: local_id !< local id 182 INTEGER(iwp) :: num_devices !< number of devices visible to OpenACC 183 INTEGER(iwp) :: my_device !< device used by this process 182 184 #endif 183 185 … … 187 189 #if defined( __parallel ) 188 190 ! 189 !-- MPI initialisation. comm2d is preliminary set, because 190 !-- it will be defined in init_pegrid but isused before in cpu_log.191 !-- MPI initialisation. comm2d is preliminary set, because it will be defined in init_pegrid but is 192 !-- used before in cpu_log. 191 193 CALL MPI_INIT( ierr ) 192 194 193 195 ! 194 !-- Initialize the coupling for nested-domain runs 195 !-- comm_palm is the communicator which includes all PEs (MPI processes) 196 !-- available for this (nested) model. If it is not a nested run, comm_palm 197 !-- is returned as MPI_COMM_WORLD 196 !-- Initialize the coupling for nested-domain runs comm_palm is the communicator which includes all 197 !-- PEs (MPI processes) available for this (nested) model. If it is not a nested run, comm_palm is 198 !-- returned as MPI_COMM_WORLD. 198 199 CALL cpu_log( log_point_s(70), 'pmci_init', 'start' ) 199 200 CALL pmci_init( comm_palm ) … … 201 202 comm2d = comm_palm 202 203 ! 203 !-- Get the (preliminary) number of MPI processes and the local PE-id (in case 204 !-- of a further communicator splitting in init_coupling, these numbers will 205 !-- be changed in init_pegrid). 204 !-- Get the (preliminary) number of MPI processes and the local PE-id (in case of a further 205 !-- communicator splitting in init_coupling, these numbers will be changed in init_pegrid). 206 206 IF ( nested_run ) THEN 207 207 … … 214 214 CALL MPI_COMM_RANK( MPI_COMM_WORLD, myid, ierr ) 215 215 ! 216 !-- Initialize PE topology in case of coupled atmosphere-ocean runs (comm_palm 217 !-- will be splittedin init_coupling)216 !-- Initialize PE topology in case of coupled atmosphere-ocean runs (comm_palm will be splitted 217 !-- in init_coupling) 218 218 CALL init_coupling 219 219 ENDIF … … 221 221 #ifdef _OPENACC 222 222 ! 223 !-- Select OpenACC device to use in this process. For this find out how many 224 !-- neighbors there arerunning on the same node and which id this process is.223 !-- Select OpenACC device to use in this process. For this find out how many neighbors there are 224 !-- running on the same node and which id this process is. 225 225 IF ( nested_run ) THEN 226 CALL MPI_COMM_SPLIT_TYPE( comm_palm, MPI_COMM_TYPE_SHARED, 0, 227 MPI_INFO_NULL, local_comm,ierr )226 CALL MPI_COMM_SPLIT_TYPE( comm_palm, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, local_comm, & 227 ierr ) 228 228 ELSE 229 CALL MPI_COMM_SPLIT_TYPE( MPI_COMM_WORLD, MPI_COMM_TYPE_SHARED, 0, &230 MPI_INFO_NULL,local_comm, ierr )229 CALL MPI_COMM_SPLIT_TYPE( MPI_COMM_WORLD, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, & 230 local_comm, ierr ) 231 231 ENDIF 232 232 CALL MPI_COMM_SIZE( local_comm, local_num_procs, ierr ) … … 234 234 235 235 ! 236 !-- This loop including the barrier is a workaround for PGI compiler versions 237 !-- up to and including 18.4. Later releases are able to select their GPUs in238 !-- parallel, without running into spuriouserrors.236 !-- This loop including the barrier is a workaround for PGI compiler versions up to and including 237 !-- 18.4. Later releases are able to select their GPUs in parallel, without running into spurious 238 !-- errors. 239 239 DO i = 0, local_num_procs-1 240 240 CALL MPI_BARRIER( local_comm, ierr ) … … 267 267 268 268 ! 269 !-- Initialize dvrp logging. Also, one PE maybe split from the global 270 !-- communicator for doing the dvrp output. In that case, the number of 271 !-- PEs available for PALM is reduced by one and communicator comm_palm 272 !-- is changed respectively. 269 !-- Initialize dvrp logging. Also, one PE maybe split from the global communicator for doing the 270 !-- dvrp output. In that case, the number of PEs available for PALM is reduced by one and 271 !-- communicator comm_palm is changed respectively. 273 272 #if defined( __parallel ) 274 273 CALL MPI_COMM_RANK( comm_palm, myid, ierr ) … … 281 280 ! 282 281 !-- Check for the user's interface version 283 IF ( user_interface_current_revision /= user_interface_required_revision ) & 284 THEN 285 message_string = 'current user-interface revision "' // & 286 TRIM( user_interface_current_revision ) // '" does ' // & 287 'not match the required revision ' // & 282 IF ( user_interface_current_revision /= user_interface_required_revision ) THEN 283 message_string = 'current user-interface revision "' // & 284 TRIM( user_interface_current_revision ) // '" does ' // & 285 'not match the required revision ' // & 288 286 TRIM( user_interface_required_revision ) 289 287 CALL message( 'palm', 'PA0169', 1, 2, 0, 6, 0 ) … … 297 295 CALL netcdf_data_input_inquire_file 298 296 ! 299 !-- Read topography input data if required. This is required before the 300 !-- numerical grid is finally created in init_grid301 CALL netcdf_data_input_topo 302 ! 303 !-- Generate grid parameters, initialize generic topography and further process 304 !-- topography information if required297 !-- Read topography input data if required. This is required before the numerical grid is finally 298 !-- created in init_grid. 299 CALL netcdf_data_input_topo 300 ! 301 !-- Generate grid parameters, initialize generic topography and further process topography 302 !-- information if required. 305 303 CALL init_grid 306 304 ! 307 !-- Initialize boundary conditions and numerics such as the multigrid solver or 308 !-- the advectionroutine305 !-- Initialize boundary conditions and numerics such as the multigrid solver or the advection 306 !-- routine 309 307 CALL module_interface_init_numerics 310 308 ! 311 !-- Read global attributes if available. 312 CALL netcdf_data_input_init 313 ! 314 !-- Read surface classification data, e.g. vegetation and soil types, water 315 !-- surfaces, etc., if available. Some of these data is required before 316 !-- check parameters is invoked. 309 !-- Read global attributes if available. 310 CALL netcdf_data_input_init 311 ! 312 !-- Read surface classification data, e.g. vegetation and soil types, water surfaces, etc., if 313 !-- available. Some of these data is required before check parameters is invoked. 317 314 CALL netcdf_data_input_surface_data 318 315 ! … … 331 328 ! 332 329 !-- Receive and interpolate initial data on children. 333 !-- Child initialization must be made first if the model is both child and 334 !-- parent if necessary 330 !-- Child initialization must be made first if the model is both child and parent if necessary. 335 331 IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 336 332 CALL pmci_child_initialize … … 355 351 ENDIF 356 352 IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN 357 CALL exchange_horiz( qr, nbgp ) 353 CALL exchange_horiz( qr, nbgp ) 358 354 CALL exchange_horiz( nr, nbgp ) 359 355 ENDIF … … 374 370 375 371 ! 376 !-- Integration of the non-atmospheric equations (land surface model, urban 377 !-- surface model) 372 !-- Integration of the non-atmospheric equations (land surface model, urban surface model) 378 373 IF ( spinup ) THEN 379 374 CALL cpu_log( log_point(41), 'wall/soil spinup', 'start' ) … … 457 452 458 453 CALL cpu_log( log_point(22), 'write-restart-data', 'stop' ) 459 454 460 455 ENDIF 461 456 ! … … 471 466 CALL cpu_log( log_point(4), 'last actions', 'start' ) 472 467 473 IF ( myid == 0 .AND. agents_active ) CALL mas_last_actions ! ToDo: move to module_interface468 IF ( myid == 0 .AND. agents_active ) CALL mas_last_actions ! ToDo: move to module_interface 474 469 475 470 CALL module_interface_last_actions … … 482 477 483 478 ! 484 !-- Write run number to file (used by palmrun to create unified cycle numbers 485 !-- for output files 479 !-- Write run number to file (used by palmrun to create unified cycle numbers for output files). 486 480 IF ( myid == 0 .AND. runnr > 0 ) THEN 487 481 OPEN( 90, FILE='RUN_NUMBER', FORM='FORMATTED' ) -
palm/trunk/SOURCE/palm_date_time_mod.f90
r4680 r4797 24 24 ! ----------------- 25 25 ! $Id$ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4680 2020-09-16 10:20:34Z gronemeier 26 29 ! Add option to fix date and time; renamed set_reference_date_time to init_date_time 27 30 ! … … 31 34 ! 4227 2019-09-10 18:04:34Z gronemeier 32 35 ! Complete rework of module date_and_time_mod: 33 ! - renamed module to prevent confusion with 34 ! FORTRAN Standard routine date_and_time 36 ! - renamed module to prevent confusion with FORTRAN Standard routine date_and_time 35 37 ! - introduce date_time_type 36 38 ! - add set_reference_date_time … … 46 48 ! Description: 47 49 ! ------------ 48 !> This routine calculates all needed information on date and time used by 49 !> other modules 50 !> This routine calculates all needed information on date and time used by other modules 50 51 !> 51 52 !> @todo Consider leap seconds … … 74 75 INTEGER(iwp), PARAMETER :: southward_equinox = 264_iwp !< Sep 21 (leap year: Sep 20) 75 76 77 CHARACTER(LEN=3), DIMENSION(days_per_week), PARAMETER :: & 78 weekdays = (/"Mon", "Tue", "Wed", "Thu", "Fri", "Sat", "Sun"/) !< names of weekdays 79 80 INTEGER, DIMENSION(months_per_year), PARAMETER :: & 81 days_per_month_noleapyear = (/31,28,31,30,31,30,31,31,30,31,30,31/) !< days for each month (no leap year) 82 83 INTEGER, DIMENSION(months_per_year), PARAMETER :: & 84 days_per_month_leapyear = (/31,29,31,30,31,30,31,31,30,31,30,31/) !< days for each month (leap year) 85 86 INTEGER, DIMENSION(121), PARAMETER :: leap_years = & !< list of leap years 87 (/ 1804_iwp, 1808_iwp, 1812_iwp, 1816_iwp, 1820_iwp, 1824_iwp, 1828_iwp, 1832_iwp, & 88 1836_iwp, 1840_iwp, 1844_iwp, 1848_iwp, 1852_iwp, 1856_iwp, 1860_iwp, 1864_iwp, & 89 1868_iwp, 1872_iwp, 1876_iwp, 1880_iwp, 1884_iwp, 1888_iwp, 1892_iwp, 1896_iwp, & 90 1904_iwp, 1908_iwp, 1912_iwp, 1916_iwp, 1920_iwp, 1924_iwp, 1928_iwp, 1932_iwp, & 91 1936_iwp, 1940_iwp, 1944_iwp, 1948_iwp, 1952_iwp, 1956_iwp, 1960_iwp, 1964_iwp, & 92 1968_iwp, 1972_iwp, 1976_iwp, 1980_iwp, 1984_iwp, 1988_iwp, 1992_iwp, 1996_iwp, & 93 2000_iwp, 2004_iwp, 2008_iwp, 2012_iwp, 2016_iwp, 2020_iwp, 2024_iwp, 2028_iwp, & 94 2032_iwp, 2036_iwp, 2040_iwp, 2044_iwp, 2048_iwp, 2052_iwp, 2056_iwp, 2060_iwp, & 95 2064_iwp, 2068_iwp, 2072_iwp, 2076_iwp, 2080_iwp, 2084_iwp, 2088_iwp, 2092_iwp, & 96 2096_iwp, 2104_iwp, 2108_iwp, 2112_iwp, 2116_iwp, 2120_iwp, 2124_iwp, 2128_iwp, & 97 2132_iwp, 2136_iwp, 2140_iwp, 2144_iwp, 2148_iwp, 2152_iwp, 2156_iwp, 2160_iwp, & 98 2164_iwp, 2168_iwp, 2172_iwp, 2176_iwp, 2180_iwp, 2184_iwp, 2188_iwp, 2192_iwp, & 99 2196_iwp, 2204_iwp, 2208_iwp, 2212_iwp, 2216_iwp, 2220_iwp, 2224_iwp, 2228_iwp, & 100 2232_iwp, 2236_iwp, 2240_iwp, 2244_iwp, 2248_iwp, 2252_iwp, 2256_iwp, 2260_iwp, & 101 2264_iwp, 2268_iwp, 2272_iwp, 2276_iwp, 2280_iwp, 2284_iwp, 2288_iwp, 2292_iwp, & 102 2296_iwp /) 103 76 104 REAL(wp), PARAMETER :: seconds_per_minute = 60.0_wp !< seconds in a minute 77 105 REAL(wp), PARAMETER :: seconds_per_hour = seconds_per_minute * minutes_per_hour !< seconds in an hour 78 106 REAL(wp), PARAMETER :: seconds_per_day = seconds_per_hour * hours_per_day !< seconds in a day 79 80 CHARACTER(LEN=3), DIMENSION(days_per_week), PARAMETER :: &81 weekdays = (/"Mon", "Tue", "Wed", "Thu", "Fri", "Sat", "Sun"/) !< names of weekdays82 83 INTEGER, DIMENSION(months_per_year), PARAMETER :: &84 days_per_month_noleapyear = (/31,28,31,30,31,30,31,31,30,31,30,31/) !< days for each month (no leap year)85 86 INTEGER, DIMENSION(months_per_year), PARAMETER :: &87 days_per_month_leapyear = (/31,29,31,30,31,30,31,31,30,31,30,31/) !< days for each month (leap year)88 89 INTEGER, DIMENSION(121), PARAMETER :: leap_years = & !< list of leap years90 (/1804_iwp, 1808_iwp, 1812_iwp, 1816_iwp, 1820_iwp, 1824_iwp, 1828_iwp, 1832_iwp, &91 1836_iwp, 1840_iwp, 1844_iwp, 1848_iwp, 1852_iwp, 1856_iwp, 1860_iwp, 1864_iwp, &92 1868_iwp, 1872_iwp, 1876_iwp, 1880_iwp, 1884_iwp, 1888_iwp, 1892_iwp, 1896_iwp, &93 1904_iwp, 1908_iwp, 1912_iwp, 1916_iwp, 1920_iwp, 1924_iwp, 1928_iwp, 1932_iwp, &94 1936_iwp, 1940_iwp, 1944_iwp, 1948_iwp, 1952_iwp, 1956_iwp, 1960_iwp, 1964_iwp, &95 1968_iwp, 1972_iwp, 1976_iwp, 1980_iwp, 1984_iwp, 1988_iwp, 1992_iwp, 1996_iwp, &96 2000_iwp, 2004_iwp, 2008_iwp, 2012_iwp, 2016_iwp, 2020_iwp, 2024_iwp, 2028_iwp, &97 2032_iwp, 2036_iwp, 2040_iwp, 2044_iwp, 2048_iwp, 2052_iwp, 2056_iwp, 2060_iwp, &98 2064_iwp, 2068_iwp, 2072_iwp, 2076_iwp, 2080_iwp, 2084_iwp, 2088_iwp, 2092_iwp, &99 2096_iwp, 2104_iwp, 2108_iwp, 2112_iwp, 2116_iwp, 2120_iwp, 2124_iwp, 2128_iwp, &100 2132_iwp, 2136_iwp, 2140_iwp, 2144_iwp, 2148_iwp, 2152_iwp, 2156_iwp, 2160_iwp, &101 2164_iwp, 2168_iwp, 2172_iwp, 2176_iwp, 2180_iwp, 2184_iwp, 2188_iwp, 2192_iwp, &102 2196_iwp, 2204_iwp, 2208_iwp, 2212_iwp, 2216_iwp, 2220_iwp, 2224_iwp, 2228_iwp, &103 2232_iwp, 2236_iwp, 2240_iwp, 2244_iwp, 2248_iwp, 2252_iwp, 2256_iwp, 2260_iwp, &104 2264_iwp, 2268_iwp, 2272_iwp, 2276_iwp, 2280_iwp, 2284_iwp, 2288_iwp, 2292_iwp, &105 2296_iwp /)106 107 107 108 ! … … 146 147 ! 147 148 !-- Public Interfaces 148 PUBLIC &149 get_date_time, &149 PUBLIC & 150 get_date_time, & 150 151 init_date_time 151 152 ! 152 153 !-- Public variables 153 PUBLIC &154 date_time_str_len, &155 days_per_week, &156 hours_per_day, &157 minutes_per_hour, &158 months_per_year, &159 northward_equinox, &160 seconds_per_minute, &161 seconds_per_hour, &162 seconds_per_day, &163 southward_equinox, &154 PUBLIC & 155 date_time_str_len, & 156 days_per_week, & 157 hours_per_day, & 158 minutes_per_hour, & 159 months_per_year, & 160 northward_equinox, & 161 seconds_per_minute, & 162 seconds_per_hour, & 163 seconds_per_day, & 164 southward_equinox, & 164 165 weekdays 165 166 … … 180 181 LOGICAL, INTENT(IN), OPTIONAL :: use_fixed_date !< flag to fix date 181 182 LOGICAL, INTENT(IN), OPTIONAL :: use_fixed_time !< flag to fix time 183 182 184 ! 183 185 !-- Check if date and time are already set … … 185 187 !> @note This error should never be observed by a user. 186 188 !> It can only appear if the code was modified. 187 WRITE( message_string, * ) 'Multiple calls to init_date_time detected.&' // &189 WRITE( message_string, * ) 'Multiple calls to init_date_time detected.&' // & 188 190 'This routine must not be called more than once.' 189 191 CALL message( 'init_date_time', 'PA0680', 2, 2, 0, 6, 0 ) … … 212 214 !> via 'reference_date_time_str' or previously set by calling routine 'init_date_time'. 213 215 !--------------------------------------------------------------------------------------------------! 214 SUBROUTINE get_date_time( time_since_reference, reference_date_time_str, & 215 year, month, day, hour, minute, second, zone, & 216 second_of_day, second_of_year, & 217 day_of_year, day_of_week, weekday, date_time_str, & 218 days_per_month, days_per_year ) 216 SUBROUTINE get_date_time( time_since_reference, reference_date_time_str, year, month, day, hour, & 217 minute, second, zone, second_of_day, second_of_year, day_of_year, & 218 day_of_week, weekday, date_time_str, days_per_month, days_per_year ) 219 219 220 220 CHARACTER(LEN=date_time_str_len), INTENT(OUT), OPTIONAL :: date_time_str !< date-time as string … … 223 223 CHARACTER(LEN=3), INTENT(OUT), OPTIONAL :: weekday !< weekday 224 224 225 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day !< day of month 226 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day_of_week !< number of weekday 227 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day_of_year !< day of the year 228 INTEGER(iwp), INTENT(OUT), OPTIONAL :: hour !< hour of day 229 INTEGER(iwp), INTENT(OUT), OPTIONAL :: minute !< minute of hour 230 INTEGER(iwp), INTENT(OUT), OPTIONAL :: month !< month of year 231 INTEGER(iwp), INTENT(OUT), OPTIONAL :: year !< year 232 INTEGER(iwp), INTENT(OUT), OPTIONAL :: zone !< time zone 233 INTEGER(iwp), INTENT(OUT), OPTIONAL :: days_per_year !< days per year 225 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day !< day of month 226 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day_of_week !< number of weekday 227 INTEGER(iwp), INTENT(OUT), OPTIONAL :: day_of_year !< day of the year 228 INTEGER(iwp), INTENT(OUT), OPTIONAL :: days_per_year !< days per year 229 INTEGER(iwp), INTENT(OUT), OPTIONAL :: hour !< hour of day 230 INTEGER(iwp), INTENT(OUT), OPTIONAL :: minute !< minute of hour 231 INTEGER(iwp), INTENT(OUT), OPTIONAL :: month !< month of year 232 INTEGER(iwp), INTENT(OUT), OPTIONAL :: year !< year 233 INTEGER(iwp), INTENT(OUT), OPTIONAL :: zone !< time zone 234 234 235 INTEGER(iwp), DIMENSION(months_per_year), INTENT(OUT), OPTIONAL :: days_per_month !< days per year 235 236 … … 247 248 !> @note This error should never be observed by a user. 248 249 !> It can only appear if the code was modified. 249 WRITE( message_string, * ) 'No reference date-time defined. '// &250 'Returning date-time information is not possible. ' // &251 'Either specify reference_date_time_str ' // &250 WRITE( message_string, * ) 'No reference date-time defined. '// & 251 'Returning date-time information is not possible. ' // & 252 'Either specify reference_date_time_str ' // & 252 253 'or set a reference via init_date_time.' 253 254 CALL message( 'get_date_time', 'PA0677', 2, 2, 0, 6, 0 ) … … 273 274 !-- If date shall be fixed, revert it to the reference date if changed 274 275 IF ( date_is_fixed ) THEN 275 IF ( date_time%year /= internal_reference_date_time%year .OR. &276 get_day_of_year( date_time ) /= get_day_of_year( internal_reference_date_time ) ) &276 IF ( date_time%year /= internal_reference_date_time%year .OR. & 277 get_day_of_year( date_time ) /= get_day_of_year( internal_reference_date_time ) ) & 277 278 THEN 278 279 279 date_time%year = internal_reference_date_time%year 280 date_time%month = internal_reference_date_time%month 281 date_time%day = internal_reference_date_time%day 282 283 date_time = update_leapyear_setting( date_time ) 284 280 date_time%year = internal_reference_date_time%year 281 date_time%month = internal_reference_date_time%month 282 date_time%day = internal_reference_date_time%day 283 date_time = update_leapyear_setting( date_time ) 285 284 date_time%second_of_year = get_second_of_year( date_time ) 286 285 … … 311 310 plus_minus = '+' 312 311 ENDIF 313 WRITE( UNIT = date_time_str, &314 FMT = '(I4,"-",I2.2,"-",I2.2,1X,I2.2,":",I2.2,":",I2.2,1X,A1,I2.2)' ) &315 date_time%year, date_time%month, date_time%day, &316 date_time%hour, date_time%minute, INT( date_time%second ), &312 WRITE( UNIT = date_time_str, & 313 FMT = '(I4,"-",I2.2,"-",I2.2,1X,I2.2,":",I2.2,":",I2.2,1X,A1,I2.2)' ) & 314 date_time%year, date_time%month, date_time%day, & 315 date_time%hour, date_time%minute, INT( date_time%second ), & 317 316 plus_minus, ABS( date_time%zone ) 318 317 ENDIF … … 326 325 !> Convert a date-time string into a date_time object. 327 326 !--------------------------------------------------------------------------------------------------! 328 FUNCTION convert_string_to_date_time( date_time_str ) RESULT( date_time )327 FUNCTION convert_string_to_date_time( date_time_str ) RESULT( date_time ) 329 328 330 329 CHARACTER(LEN=date_time_str_len), INTENT(IN) :: date_time_str !< date-time as string … … 337 336 !-- Decompose string to date-time information 338 337 READ( UNIT = date_time_str( 1: 4), IOSTAT = read_status, FMT = '(I4)' ) date_time%year 339 IF ( read_status == 0 ) &338 IF ( read_status == 0 ) THEN 340 339 READ( UNIT = date_time_str( 6: 7), IOSTAT = read_status, FMT = '(I2)' ) date_time%month 341 IF ( read_status == 0 ) & 340 ENDIF 341 IF ( read_status == 0 ) THEN 342 342 READ( UNIT = date_time_str( 9:10), IOSTAT = read_status, FMT = '(I2)' ) date_time%day 343 IF ( read_status == 0 ) & 343 ENDIF 344 IF ( read_status == 0 ) THEN 344 345 READ( UNIT = date_time_str(12:13), IOSTAT = read_status, FMT = '(I2)' ) date_time%hour 345 IF ( read_status == 0 ) & 346 ENDIF 347 IF ( read_status == 0 ) THEN 346 348 READ( UNIT = date_time_str(15:16), IOSTAT = read_status, FMT = '(I2)' ) date_time%minute 347 IF ( read_status == 0 ) & 349 ENDIF 350 IF ( read_status == 0 ) THEN 348 351 READ( UNIT = date_time_str(18:19), IOSTAT = read_status, FMT = '(F2.0)' ) date_time%second 349 IF ( read_status == 0 ) & 352 ENDIF 353 IF ( read_status == 0 ) THEN 350 354 READ( UNIT = date_time_str(21:23), IOSTAT = read_status, FMT = '(I3)' ) date_time%zone 355 ENDIF 351 356 352 357 IF ( read_status /= 0 ) THEN 353 WRITE( message_string, * ) 'Error while converting date-time string. ' // &354 'Please check format of date-time string: "' // &355 TRIM( date_time_str ) // '". ' // &358 WRITE( message_string, * ) 'Error while converting date-time string. ' // & 359 'Please check format of date-time string: "' // & 360 TRIM( date_time_str ) // '". ' // & 356 361 'Format must be "YYYY-MM-DD hh:mm:ss ZZZ".' 357 362 CALL message( 'convert_string_to_date_time', 'PA0678', 2, 2, 0, 6, 0 ) … … 367 372 ! 368 373 !-- Shift time to UTC and set zone to UTC 369 date_time = add_date_time( REAL( -1 * date_time%zone, KIND = wp ) & 370 * seconds_per_hour, & 374 date_time = add_date_time( REAL( -1 * date_time%zone, KIND = wp ) * seconds_per_hour, & 371 375 date_time ) 372 376 date_time%zone = 0_iwp … … 401 405 date_time%second_of_year = date_time%second_of_year + inc_seconds 402 406 ! 403 !-- Check if year changes 407 !-- Check if year changes. 404 408 !-- First, if year is reduced 405 409 DO WHILE ( date_time%second_of_year < 0.0_wp ) 406 date_time%year = date_time%year - 1_iwp407 date_time = update_leapyear_setting( date_time )408 seconds_per_year = REAL( date_time%days_per_year * seconds_per_day, KIND = wp )410 date_time%year = date_time%year - 1_iwp 411 date_time = update_leapyear_setting( date_time ) 412 seconds_per_year = REAL( date_time%days_per_year * seconds_per_day, KIND = wp ) 409 413 date_time%second_of_year = date_time%second_of_year + seconds_per_year 410 414 ENDDO … … 412 416 !-- Now, if year is increased 413 417 DO WHILE ( date_time%second_of_year > seconds_per_year ) 414 date_time%year = date_time%year + 1_iwp415 date_time = update_leapyear_setting( date_time )418 date_time%year = date_time%year + 1_iwp 419 date_time = update_leapyear_setting( date_time ) 416 420 date_time%second_of_year = date_time%second_of_year - seconds_per_year 417 seconds_per_year = REAL( date_time%days_per_year * seconds_per_day, KIND = wp )421 seconds_per_year = REAL( date_time%days_per_year * seconds_per_day, KIND = wp ) 418 422 ENDDO 419 423 ! 420 424 !-- Based on updated year and second_of_year, update month, day, hour, minute, second 421 425 DO i = 1, months_per_year 422 IF ( date_time%second_of_year < SUM( date_time%days_per_month(1:i) ) * seconds_per_day ) &426 IF ( date_time%second_of_year < SUM( date_time%days_per_month(1:i) ) * seconds_per_day ) & 423 427 THEN 424 428 date_time%month = i 425 seconds_left = date_time%second_of_year &426 - REAL( SUM( date_time%days_per_month(1:i-1) ), KIND=wp )&427 * seconds_per_day429 seconds_left = date_time%second_of_year & 430 - REAL( SUM( date_time%days_per_month(1:i-1) ), KIND=wp ) & 431 * seconds_per_day 428 432 date_time%day = INT( seconds_left / seconds_per_day, KIND = iwp ) + 1_iwp 429 seconds_left = seconds_left &430 - REAL( date_time%day - 1_iwp, KIND = wp ) * seconds_per_day433 seconds_left = seconds_left & 434 - REAL( date_time%day - 1_iwp, KIND = wp ) * seconds_per_day 431 435 date_time%hour = INT( seconds_left / seconds_per_hour, KIND = iwp ) 432 436 seconds_left = seconds_left - REAL( date_time%hour, KIND = wp ) * seconds_per_hour … … 445 449 !> Return day of year of given date. 446 450 !--------------------------------------------------------------------------------------------------! 447 FUNCTION get_day_of_year( date_time ) RESULT( day_of_year )451 FUNCTION get_day_of_year( date_time ) RESULT( day_of_year ) 448 452 449 453 INTEGER(iwp) :: day_of_year !< day of the year … … 466 470 !> Return second of day of given time. 467 471 !--------------------------------------------------------------------------------------------------! 468 FUNCTION get_second_of_day( date_time ) RESULT( second_of_day )472 FUNCTION get_second_of_day( date_time ) RESULT( second_of_day ) 469 473 470 474 REAL(wp) :: second_of_day !< second of the day … … 473 477 474 478 475 second_of_day = date_time%second &476 + REAL( ( date_time%hour * minutes_per_hour ) + date_time%minute, KIND = wp ) &477 * seconds_per_minute479 second_of_day = date_time%second & 480 + REAL( ( date_time%hour * minutes_per_hour ) + date_time%minute, KIND = wp ) & 481 * seconds_per_minute 478 482 479 483 END FUNCTION get_second_of_day … … 485 489 !> Return second of year of given date-time. 486 490 !--------------------------------------------------------------------------------------------------! 487 FUNCTION get_second_of_year( date_time ) RESULT( second_of_year )491 FUNCTION get_second_of_year( date_time ) RESULT( second_of_year ) 488 492 489 493 REAL(wp) :: second_of_year !< second of the year … … 492 496 493 497 494 second_of_year = get_second_of_day( date_time ) &495 + REAL( get_day_of_year( date_time ) - 1_iwp, KIND = wp ) * seconds_per_day498 second_of_year = get_second_of_day( date_time ) & 499 + REAL( get_day_of_year( date_time ) - 1_iwp, KIND = wp ) * seconds_per_day 496 500 497 501 END FUNCTION get_second_of_year … … 503 507 !> Return index of the day of the week of the given date-time. 504 508 !--------------------------------------------------------------------------------------------------! 505 FUNCTION get_day_of_week( date_time_in ) RESULT( day_of_week ) 506 507 INTEGER(iwp) :: date_time_internal_reference_day_of_week !< day of week of reference date 509 FUNCTION get_day_of_week( date_time_in ) RESULT( day_of_week ) 510 508 511 INTEGER(iwp) :: day_difference !< day between given date and reference 509 512 INTEGER(iwp) :: day_of_week !< day of the week 513 INTEGER(iwp) :: date_time_internal_reference_day_of_week !< day of week of reference date 510 514 511 515 TYPE(date_time_type), INTENT(IN) :: date_time_in !< date of which to get the weekday … … 548 552 ENDIF 549 553 ! 550 !-- Remove full weeks of day_difference and shift day_of_week of reference by 551 !-- remaining days 554 !-- Remove full weeks of day_difference and shift day_of_week of reference by remaining days. 552 555 day_of_week = date_time_internal_reference_day_of_week + MOD( day_difference, days_per_week ) 553 556 … … 570 573 !> Check if given year is a leap year and update days per month accordingly. 571 574 !--------------------------------------------------------------------------------------------------! 572 FUNCTION update_leapyear_setting( date_time_in ) RESULT( date_time_out )575 FUNCTION update_leapyear_setting( date_time_in ) RESULT( date_time_out ) 573 576 574 577 TYPE(date_time_type), INTENT(IN) :: date_time_in !< input date-time … … 593 596 ! ------------ 594 597 !> Check if given date and time are valid. Returns 0 if all checks are passed. 595 !> @todo Revise error message. ATM, gives only last errorneous value even if 596 !> multiple values violatethe bounds.597 !--------------------------------------------------------------------------------------------------! 598 FUNCTION check_date( date_time, date_time_str ) RESULT( error_code )598 !> @todo Revise error message. ATM, gives only last errorneous value even if multiple values violate 599 !> the bounds. 600 !--------------------------------------------------------------------------------------------------! 601 FUNCTION check_date( date_time, date_time_str ) RESULT( error_code ) 599 602 600 603 … … 612 615 ! 613 616 !-- Check date 614 IF ( date_time%month < 1_iwp .OR. & 615 date_time%month > months_per_year ) THEN 617 IF ( date_time%month < 1_iwp .OR. date_time%month > months_per_year ) THEN 616 618 error_code = 2 617 619 ELSE 618 IF ( date_time%day < 1_iwp .OR. &620 IF ( date_time%day < 1_iwp .OR. & 619 621 date_time%day > date_time%days_per_month(date_time%month) ) THEN 620 622 error_code = 3 … … 623 625 ! 624 626 !-- Check time 625 IF ( date_time%hour < 0_iwp .OR. & 626 date_time%hour > hours_per_day ) THEN 627 IF ( date_time%hour < 0_iwp .OR. date_time%hour > hours_per_day ) THEN 627 628 error_code = 4 628 629 ELSE 629 IF ( date_time%minute < 0_iwp .OR. & 630 date_time%minute > minutes_per_hour ) THEN 630 IF ( date_time%minute < 0_iwp .OR. date_time%minute > minutes_per_hour ) THEN 631 631 error_code = 5 632 632 ELSE 633 IF ( date_time%second < 0.0_wp .OR. & 634 date_time%second >= seconds_per_minute ) THEN 633 IF ( date_time%second < 0.0_wp .OR. date_time%second >= seconds_per_minute ) THEN 635 634 error_code = 6 636 635 ENDIF … … 638 637 ENDIF 639 638 ! 640 !-- Check time zone 639 !-- Check time zone. 641 640 !-- Bounds defined by maximum and minimum time zone present on earth 642 IF ( date_time%zone < -12_iwp .OR. & 643 date_time%zone > 14_iwp ) THEN 641 IF ( date_time%zone < -12_iwp .OR. date_time%zone > 14_iwp ) THEN 644 642 error_code = 7 645 643 ENDIF … … 647 645 !-- Raise error if any check is marked invalid 648 646 IF ( error_code /= 0 ) THEN 649 WRITE( message_string, * ) 'Date-time values out of bounds: "' // &650 TRIM( error_str_list(error_code) ) // &651 '" is out of bounds. Please check date-time string: "' // &652 TRIM( date_time_str ) // '". ' // &647 WRITE( message_string, * ) 'Date-time values out of bounds: "' // & 648 TRIM( error_str_list(error_code) ) // & 649 '" is out of bounds. Please check date-time string: "' // & 650 TRIM( date_time_str ) // '". ' // & 653 651 'Format must be "YYYY-MM-DD hh:mm:ss ZZZ".' 654 652 CALL message( 'check_date', 'PA0679', 2, 2, 0, 6, 0 ) -
palm/trunk/SOURCE/parin.f90
r4783 r4797 1 1 !> @file parin.f90 2 !------------------------------------------------------------------------------ !2 !--------------------------------------------------------------------------------------------------! 3 3 ! This file is part of the PALM model system. 4 4 ! 5 ! PALM is free software: you can redistribute it and/or modify it under the 6 ! terms of the GNU General Public License as published by the Free Software 7 ! Foundation, either version 3 of the License, or (at your option) any later 8 ! version. 9 ! 10 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY 11 ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR 12 ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. 13 ! 14 ! You should have received a copy of the GNU General Public License along with 15 ! PALM. If not, see <http://www.gnu.org/licenses/>. 5 ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General 6 ! Public License as published by the Free Software Foundation, either version 3 of the License, or 7 ! (at your option) any later version. 8 ! 9 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the 10 ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 ! Public License for more details. 12 ! 13 ! You should have received a copy of the GNU General Public License along with PALM. If not, see 14 ! <http://www.gnu.org/licenses/>. 16 15 ! 17 16 ! Copyright 1997-2020 Leibniz Universitaet Hannover 18 !------------------------------------------------------------------------------ !17 !--------------------------------------------------------------------------------------------------! 19 18 ! 20 19 ! Current revisions: … … 25 24 ! ----------------- 26 25 ! $Id$ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4783 2020-11-13 13:58:45Z raasch 27 29 ! bugfix for reading restart data with MPI-I/O (does not work with blockwise I/O) 28 30 ! … … 76 78 ! 77 79 ! 4022 2019-06-12 11:52:39Z suehring 78 ! Change default top boundary condition for pressure to Neumann in offline 79 ! nesting case 80 ! Change default top boundary condition for pressure to Neumann in offline nesting case 80 81 ! 81 82 ! 4017 2019-06-06 12:16:46Z schwenkel … … 83 84 ! 84 85 ! 3885 2019-04-11 11:29:34Z kanani 85 ! Changes related to global restructuring of location messages and introduction 86 ! of additional debugmessages86 ! Changes related to global restructuring of location messages and introduction of additional debug 87 ! messages 87 88 ! 88 89 ! 3806 2019-03-21 12:45:50Z raasch … … 107 108 !> 108 109 !> @todo: Revise max_pr_cs (profiles for chemistry) 109 !------------------------------------------------------------------------------ !110 !--------------------------------------------------------------------------------------------------! 110 111 SUBROUTINE parin 111 112 112 113 113 USE arrays_3d, & 114 ONLY: pt_init, q_init, ref_state, s_init, sa_init, & 115 ug, u_init, v_init, vg 114 USE arrays_3d, & 115 ONLY: pt_init, q_init, ref_state, s_init, sa_init, ug, u_init, v_init, vg 116 116 117 117 USE chem_modules … … 119 119 USE control_parameters 120 120 121 USE cpulog, &121 USE cpulog, & 122 122 ONLY: cpu_log_barrierwait 123 123 124 USE grid_variables, &124 USE grid_variables, & 125 125 ONLY: dx, dy 126 126 127 USE indices, &127 USE indices, & 128 128 ONLY: nx, ny, nz 129 129 130 130 USE kinds 131 131 132 USE model_1d_mod, &132 USE model_1d_mod, & 133 133 ONLY: damp_level_1d, dt_pr_1d, dt_run_control_1d, end_time_1d 134 134 135 USE module_interface, &135 USE module_interface, & 136 136 ONLY: module_interface_parin 137 137 138 USE netcdf_interface, &138 USE netcdf_interface, & 139 139 ONLY: netcdf_data_format, netcdf_deflate, netcdf_precision 140 140 141 141 USE pegrid 142 142 143 USE pmc_interface, &143 USE pmc_interface, & 144 144 ONLY: nested_run, nesting_mode 145 145 146 USE profil_parameter, &146 USE profil_parameter, & 147 147 ONLY: cross_profiles, profile_columns, profile_rows 148 148 149 USE progress_bar, &149 USE progress_bar, & 150 150 ONLY : progress_bar_disabled 151 151 152 USE read_restart_data_mod, &152 USE read_restart_data_mod, & 153 153 ONLY: rrd_global 154 154 155 USE statistics, &155 USE statistics, & 156 156 ONLY: hom, hom_sum, pr_palm, statistic_regions 157 157 158 USE turbulence_closure_mod, &158 USE turbulence_closure_mod, & 159 159 ONLY: rans_const_c, rans_const_sigma 160 160 … … 169 169 INTEGER(iwp) :: ioerr !< error flag for open/read/write 170 170 171 NAMELIST /inipar/ alpha_surface, approximation, bc_e_b, & 172 bc_lr, bc_ns, bc_p_b, bc_p_t, bc_pt_b, bc_pt_t, bc_q_b, & 173 bc_q_t,bc_s_b, bc_s_t, bc_uv_b, bc_uv_t, & 174 building_height, building_length_x, & 175 building_length_y, building_wall_left, building_wall_south, & 176 calc_soil_moisture_during_spinup, & 177 call_psolver_at_all_substeps, & 178 canyon_height, & 179 canyon_width_x, canyon_width_y, canyon_wall_left, & 180 canyon_wall_south, cfl_factor, check_realistic_q, cloud_droplets, & 181 collective_wait, complex_terrain, & 182 conserve_volume_flow, & 183 conserve_volume_flow_mode, constant_flux_layer, & 184 coupling_start_time, & 185 cycle_mg, damp_level_1d, & 186 data_output_during_spinup, & 187 origin_date_time, & 188 dissipation_1d, & 189 dp_external, dp_level_b, dp_smooth, dpdxy, & 190 dt, dt_pr_1d, dt_run_control_1d, dt_spinup, dx, dy, dz, dz_max, & 191 dz_stretch_factor, dz_stretch_level, dz_stretch_level_start, & 192 dz_stretch_level_end, end_time_1d, ensemble_member_nr, e_init, & 193 e_min, fft_method, flux_input_mode, flux_output_mode, & 194 galilei_transformation, humidity, & 195 inflow_damping_height, inflow_damping_width, & 196 inflow_disturbance_begin, inflow_disturbance_end, & 197 initializing_actions, km_constant, & 198 large_scale_forcing, large_scale_subsidence, latitude, & 199 longitude, & 200 loop_optimization, lsf_exception, masking_method, mg_cycles, & 201 mg_switch_to_pe0_level, mixing_length_1d, momentum_advec, & 202 monotonic_limiter_z, & 203 netcdf_precision, neutral, ngsrb, & 204 nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega, & 205 omega_sor, outflow_source_plane, passive_scalar, & 206 prandtl_number, psolver, pt_damping_factor, pt_damping_width, & 207 pt_reference, pt_surface, pt_surface_heating_rate, & 208 pt_surface_initial_change, pt_vertical_gradient, & 209 pt_vertical_gradient_level, q_surface, q_surface_initial_change, & 210 q_vertical_gradient, q_vertical_gradient_level, & 211 random_generator, random_heatflux, rans_const_c, rans_const_sigma,& 212 rayleigh_damping_factor, rayleigh_damping_height, & 213 recycling_method_for_thermodynamic_quantities, recycling_width, & 214 reference_state, residual_limit, & 215 rotation_angle, & 216 roughness_length, & 217 scalar_advec, & 218 scalar_rayleigh_damping, & 219 spinup_time, spinup_pt_amplitude, spinup_pt_mean, & 220 statistic_regions, subs_vertical_gradient, & 221 subs_vertical_gradient_level, surface_heatflux, surface_pressure, & 222 surface_scalarflux, surface_waterflux, & 223 s_surface, s_surface_initial_change, s_vertical_gradient, & 224 s_vertical_gradient_level, timestep_scheme, & 225 topography, topography_grid_convention, top_heatflux, & 226 top_momentumflux_u, top_momentumflux_v, & 227 top_scalarflux, transpose_compute_overlap, & 228 tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y, & 229 tunnel_wall_depth, turbulence_closure, & 230 turbulent_inflow, turbulent_outflow, & 231 use_subsidence_tendencies, ug_surface, ug_vertical_gradient, & 232 use_fixed_date, use_fixed_time, & 233 use_free_convection_scaling, & 234 ug_vertical_gradient_level, use_surface_fluxes, use_cmax, & 235 use_top_fluxes, use_ug_for_galilei_tr, use_upstream_for_tke, & 236 uv_heights, u_bulk, u_profile, vdi_checks, vg_surface, & 237 vg_vertical_gradient, & 238 vg_vertical_gradient_level, v_bulk, v_profile,& 239 wall_adjustment, wall_heatflux, wall_humidityflux, & 240 wall_scalarflux, y_shift, zeta_max, zeta_min, & 241 z0h_factor 242 243 NAMELIST /initialization_parameters/ alpha_surface, & 244 approximation, bc_e_b, & 245 bc_lr, bc_ns, bc_p_b, bc_p_t, bc_pt_b, bc_pt_t, bc_q_b, & 246 bc_q_t,bc_s_b, bc_s_t, bc_uv_b, bc_uv_t, & 247 building_height, building_length_x, & 248 building_length_y, building_wall_left, building_wall_south, & 249 calc_soil_moisture_during_spinup, & 250 call_psolver_at_all_substeps, & 251 canyon_height, & 252 canyon_width_x, canyon_width_y, canyon_wall_left, & 253 canyon_wall_south, cfl_factor, check_realistic_q, cloud_droplets, & 254 collective_wait, complex_terrain, & 255 conserve_volume_flow, & 256 conserve_volume_flow_mode, constant_flux_layer, & 257 coupling_start_time, & 258 cycle_mg, damp_level_1d, & 259 data_output_during_spinup, & 260 origin_date_time, & 261 dissipation_1d, & 262 dp_external, dp_level_b, dp_smooth, dpdxy, & 263 dt, dt_pr_1d, dt_run_control_1d, dt_spinup, dx, dy, dz, dz_max, & 264 dz_stretch_factor, dz_stretch_level, dz_stretch_level_start, & 265 dz_stretch_level_end, end_time_1d, ensemble_member_nr, e_init, & 266 e_min, fft_method, flux_input_mode, flux_output_mode, & 267 galilei_transformation, humidity, & 268 inflow_damping_height, inflow_damping_width, & 269 inflow_disturbance_begin, inflow_disturbance_end, & 270 initializing_actions, km_constant, & 271 large_scale_forcing, large_scale_subsidence, latitude, & 272 longitude, & 273 loop_optimization, lsf_exception, masking_method, mg_cycles, & 274 mg_switch_to_pe0_level, mixing_length_1d, momentum_advec, & 275 monotonic_limiter_z, & 276 netcdf_precision, neutral, ngsrb, & 277 nsor, nsor_ini, nudging, nx, ny, nz, ocean_mode, omega, & 278 omega_sor, outflow_source_plane, passive_scalar, & 279 prandtl_number, psolver, pt_damping_factor, pt_damping_width, & 280 pt_surface_heating_rate, pt_reference, pt_surface, & 281 pt_surface_initial_change, pt_vertical_gradient, & 282 pt_vertical_gradient_level, q_surface, q_surface_initial_change, & 283 q_vertical_gradient, q_vertical_gradient_level, & 284 random_generator, random_heatflux, rans_const_c, rans_const_sigma,& 285 rayleigh_damping_factor, rayleigh_damping_height, & 286 recycling_method_for_thermodynamic_quantities, recycling_width, & 287 reference_state, residual_limit, & 288 restart_data_format, restart_data_format_input, restart_data_format_output, & 289 rotation_angle, & 290 roughness_length, scalar_advec, & 291 scalar_rayleigh_damping, & 292 spinup_time, spinup_pt_amplitude, spinup_pt_mean, & 293 statistic_regions, subs_vertical_gradient, & 294 subs_vertical_gradient_level, surface_heatflux, surface_pressure, & 295 surface_scalarflux, surface_waterflux, & 296 s_surface, s_surface_initial_change, s_vertical_gradient, & 297 s_vertical_gradient_level, timestep_scheme, & 298 topography, topography_grid_convention, top_heatflux, & 299 top_momentumflux_u, top_momentumflux_v, & 300 top_scalarflux, transpose_compute_overlap, & 301 tunnel_height, tunnel_length, tunnel_width_x, tunnel_width_y, & 302 tunnel_wall_depth, turbulence_closure, & 303 turbulent_inflow, turbulent_outflow, & 304 use_subsidence_tendencies, ug_surface, ug_vertical_gradient, & 305 ug_vertical_gradient_level, use_surface_fluxes, use_cmax, & 306 use_top_fluxes, use_ug_for_galilei_tr, use_upstream_for_tke, & 307 use_fixed_date, use_fixed_time, & 308 use_free_convection_scaling, & 309 uv_heights, u_bulk, u_profile, vdi_checks, & 310 vg_surface, vg_vertical_gradient, & 311 vg_vertical_gradient_level, v_bulk, v_profile, & 312 wall_adjustment, wall_heatflux, wall_humidityflux, & 313 wall_scalarflux, y_shift, zeta_max, zeta_min, z0h_factor 314 315 NAMELIST /d3par/ averaging_interval, averaging_interval_pr, & 316 cpu_log_barrierwait, create_disturbances, & 317 cross_profiles, data_output, data_output_masks, & 318 data_output_pr, data_output_2d_on_each_pe, & 319 debug_output, & 320 debug_output_timestep, & 321 disturbance_amplitude, & 322 disturbance_energy_limit, disturbance_level_b, & 323 disturbance_level_t, do2d_at_begin, do3d_at_begin, & 324 dt, dt_averaging_input, dt_averaging_input_pr, & 325 dt_coupling, dt_data_output, dt_data_output_av, dt_disturb, & 326 dt_domask, dt_dopr, dt_dopr_listing, dt_dots, dt_do2d_xy, & 327 dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart, & 328 dt_run_control,end_time, force_print_header, mask_k_over_surface, & 329 mask_scale_x, & 330 mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop, & 331 mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate, & 332 normalizing_region, npex, npey, nz_do3d, & 333 profile_columns, profile_rows, & 334 restart_time, section_xy, section_xz, section_yz, & 335 skip_time_data_output, skip_time_data_output_av, skip_time_dopr, & 336 skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz, & 337 skip_time_do3d, skip_time_domask, synchronous_exchange, & 338 termination_time_needed 339 340 NAMELIST /runtime_parameters/ averaging_interval, averaging_interval_pr, & 341 cpu_log_barrierwait, create_disturbances, & 342 cross_profiles, data_output, data_output_masks, & 343 data_output_pr, data_output_2d_on_each_pe, & 344 debug_output, & 345 debug_output_timestep, & 346 disturbance_amplitude, & 347 disturbance_energy_limit, disturbance_level_b, & 348 disturbance_level_t, do2d_at_begin, do3d_at_begin, & 349 dt, dt_averaging_input, dt_averaging_input_pr, & 350 dt_coupling, dt_data_output, dt_data_output_av, dt_disturb, & 351 dt_domask, dt_dopr, dt_dopr_listing, dt_dots, dt_do2d_xy, & 352 dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_max, dt_restart, & 353 dt_run_control,end_time, force_print_header, mask_k_over_surface, & 354 mask_scale_x, & 355 mask_scale_y, mask_scale_z, mask_x, mask_y, mask_z, mask_x_loop, & 356 mask_y_loop, mask_z_loop, netcdf_data_format, netcdf_deflate, & 357 normalizing_region, npex, npey, nz_do3d, & 358 profile_columns, profile_rows, & 359 restart_time, section_xy, section_xz, section_yz, & 360 restart_data_format, restart_data_format_input, restart_data_format_output, & 361 skip_time_data_output, skip_time_data_output_av, skip_time_dopr, & 362 skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz, & 363 skip_time_do3d, skip_time_domask, synchronous_exchange, & 364 termination_time_needed 365 366 NAMELIST /envpar/ progress_bar_disabled, host, & 367 maximum_cpu_time_allowed, maximum_parallel_io_streams, & 368 read_svf, revision, run_identifier, tasks_per_node, & 369 write_binary, write_svf 370 371 ! 372 !-- First read values of environment variables (this NAMELIST file is 373 !-- generated by palmrun) 171 NAMELIST /inipar/ alpha_surface, & 172 approximation, & 173 bc_e_b, & 174 bc_lr, & 175 bc_ns, & 176 bc_p_b, & 177 bc_p_t, & 178 bc_pt_b, & 179 bc_pt_t, & 180 bc_q_b, & 181 bc_q_t, & 182 bc_s_b, & 183 bc_s_t, & 184 bc_uv_b, & 185 bc_uv_t, & 186 building_height, & 187 building_length_x, & 188 building_length_y, & 189 building_wall_left, & 190 building_wall_south, & 191 calc_soil_moisture_during_spinup, & 192 call_psolver_at_all_substeps, & 193 canyon_height, & 194 canyon_wall_left, & 195 canyon_wall_south, & 196 canyon_width_x, & 197 canyon_width_y, & 198 cfl_factor, & 199 check_realistic_q, & 200 cloud_droplets, & 201 collective_wait, & 202 complex_terrain, & 203 conserve_volume_flow, & 204 conserve_volume_flow_mode, & 205 constant_flux_layer, & 206 coupling_start_time, & 207 cycle_mg, & 208 damp_level_1d, & 209 data_output_during_spinup, & 210 dissipation_1d, & 211 dp_external, & 212 dp_level_b, & 213 dp_smooth, & 214 dpdxy, & 215 dt, & 216 dt_pr_1d, & 217 dt_run_control_1d, & 218 dt_spinup, & 219 dx, & 220 dy, & 221 dz, & 222 dz_max, & 223 dz_stretch_factor, & 224 dz_stretch_level, & 225 dz_stretch_level_end, & 226 dz_stretch_level_start, & 227 e_init, & 228 e_min, & 229 end_time_1d, & 230 ensemble_member_nr, & 231 fft_method, & 232 flux_input_mode, & 233 flux_output_mode, & 234 galilei_transformation, & 235 humidity, & 236 inflow_damping_height, & 237 inflow_damping_width, & 238 inflow_disturbance_begin, & 239 inflow_disturbance_end, & 240 initializing_actions, & 241 km_constant, & 242 large_scale_forcing, & 243 large_scale_subsidence, & 244 latitude, & 245 longitude, & 246 loop_optimization, & 247 lsf_exception, & 248 masking_method, & 249 mg_cycles, & 250 mg_switch_to_pe0_level, & 251 mixing_length_1d, & 252 momentum_advec, & 253 monotonic_limiter_z, & 254 netcdf_precision, & 255 neutral, & 256 ngsrb, & 257 nsor, & 258 nsor_ini, & 259 nudging, & 260 nx, & 261 ny, & 262 nz, & 263 ocean_mode, & 264 omega, & 265 omega_sor, & 266 origin_date_time, & 267 outflow_source_plane, & 268 passive_scalar, & 269 prandtl_number, & 270 psolver, & 271 pt_damping_factor, & 272 pt_damping_width, & 273 pt_reference, & 274 pt_surface, & 275 pt_surface_heating_rate, & 276 pt_surface_initial_change, & 277 pt_vertical_gradient, & 278 pt_vertical_gradient_level, & 279 q_surface, & 280 q_surface_initial_change, & 281 q_vertical_gradient, & 282 q_vertical_gradient_level, & 283 random_generator, & 284 random_heatflux, & 285 rans_const_c, & 286 rans_const_sigma, & 287 rayleigh_damping_factor, & 288 rayleigh_damping_height, & 289 recycling_method_for_thermodynamic_quantities, & 290 recycling_width, & 291 reference_state, & 292 residual_limit, & 293 rotation_angle, & 294 roughness_length, & 295 scalar_advec, & 296 scalar_rayleigh_damping, & 297 spinup_time, & 298 spinup_pt_amplitude, & 299 spinup_pt_mean, & 300 statistic_regions, & 301 subs_vertical_gradient, & 302 subs_vertical_gradient_level, & 303 surface_heatflux, & 304 surface_pressure, & 305 surface_scalarflux, & 306 surface_waterflux, & 307 s_surface, & 308 s_surface_initial_change, & 309 s_vertical_gradient, & 310 s_vertical_gradient_level, & 311 timestep_scheme, & 312 topography, & 313 topography_grid_convention, & 314 top_heatflux, & 315 top_momentumflux_u, & 316 top_momentumflux_v, & 317 top_scalarflux, & 318 transpose_compute_overlap, & 319 tunnel_height, & 320 tunnel_length, & 321 tunnel_wall_depth, & 322 tunnel_width_x, & 323 tunnel_width_y, & 324 turbulence_closure, & 325 turbulent_inflow, & 326 turbulent_outflow, & 327 u_bulk, & 328 u_profile, & 329 ug_surface, & 330 ug_vertical_gradient, & 331 ug_vertical_gradient_level, & 332 use_cmax, & 333 use_fixed_date, & 334 use_fixed_time, & 335 use_free_convection_scaling, & 336 use_ug_for_galilei_tr, & 337 use_subsidence_tendencies, & 338 use_surface_fluxes, & 339 use_top_fluxes, & 340 use_upstream_for_tke, & 341 uv_heights, & 342 v_bulk, & 343 v_profile, & 344 vdi_checks, & 345 vg_surface, & 346 vg_vertical_gradient, & 347 vg_vertical_gradient_level, & 348 wall_adjustment, & 349 wall_heatflux, & 350 wall_humidityflux, & 351 wall_scalarflux, & 352 y_shift, & 353 zeta_max, & 354 zeta_min, & 355 z0h_factor 356 357 NAMELIST /initialization_parameters/ alpha_surface, & 358 approximation, & 359 bc_e_b, & 360 bc_lr, & 361 bc_ns, & 362 bc_p_b, & 363 bc_p_t, & 364 bc_pt_b, & 365 bc_pt_t, & 366 bc_q_b, & 367 bc_q_t, & 368 bc_s_b, & 369 bc_s_t, & 370 bc_uv_b, & 371 bc_uv_t, & 372 building_height, & 373 building_length_x, & 374 building_length_y, & 375 building_wall_left, & 376 building_wall_south, & 377 calc_soil_moisture_during_spinup, & 378 call_psolver_at_all_substeps, & 379 canyon_height, & 380 canyon_wall_left, & 381 canyon_wall_south, & 382 canyon_width_x, & 383 canyon_width_y, & 384 cfl_factor, & 385 check_realistic_q, & 386 cloud_droplets, & 387 collective_wait, & 388 complex_terrain, & 389 conserve_volume_flow, & 390 conserve_volume_flow_mode, & 391 constant_flux_layer, & 392 coupling_start_time, & 393 cycle_mg, & 394 damp_level_1d, & 395 data_output_during_spinup, & 396 dissipation_1d, & 397 dp_external, & 398 dp_level_b, & 399 dp_smooth, dpdxy, & 400 dt, & 401 dt_pr_1d, & 402 dt_run_control_1d, & 403 dt_spinup, & 404 dx, & 405 dy, & 406 dz, & 407 dz_max, & 408 dz_stretch_factor, & 409 dz_stretch_level, & 410 dz_stretch_level_start, & 411 dz_stretch_level_end, & 412 e_init, & 413 e_min, & 414 end_time_1d, & 415 ensemble_member_nr, & 416 fft_method, & 417 flux_input_mode, & 418 flux_output_mode, & 419 galilei_transformation, & 420 humidity, & 421 inflow_damping_height, & 422 inflow_damping_width, & 423 inflow_disturbance_begin, & 424 inflow_disturbance_end, & 425 initializing_actions, & 426 km_constant, & 427 large_scale_forcing, & 428 large_scale_subsidence, & 429 latitude, & 430 longitude, & 431 loop_optimization, & 432 lsf_exception, & 433 masking_method, & 434 mg_cycles, & 435 mg_switch_to_pe0_level, & 436 mixing_length_1d, & 437 momentum_advec, & 438 monotonic_limiter_z, & 439 netcdf_precision, & 440 neutral, & 441 ngsrb, & 442 nsor, & 443 nsor_ini, & 444 nudging, & 445 nx, & 446 ny, & 447 nz, & 448 ocean_mode, & 449 omega, & 450 omega_sor, & 451 origin_date_time, & 452 outflow_source_plane, & 453 passive_scalar, & 454 prandtl_number, & 455 psolver, & 456 pt_damping_factor, & 457 pt_damping_width, & 458 pt_reference, pt_surface, & 459 pt_surface_heating_rate, & 460 pt_surface_initial_change, & 461 pt_vertical_gradient, & 462 pt_vertical_gradient_level, & 463 q_surface, & 464 q_surface_initial_change, & 465 q_vertical_gradient, & 466 q_vertical_gradient_level, & 467 random_generator, & 468 random_heatflux, & 469 rans_const_c, & 470 rans_const_sigma, & 471 rayleigh_damping_factor, & 472 rayleigh_damping_height, & 473 recycling_method_for_thermodynamic_quantities, & 474 recycling_width, & 475 reference_state, & 476 residual_limit, & 477 restart_data_format, & 478 restart_data_format_input, & 479 restart_data_format_output, & 480 rotation_angle, & 481 roughness_length, & 482 s_surface, & 483 s_surface_initial_change, & 484 s_vertical_gradient, & 485 s_vertical_gradient_level, & 486 scalar_advec, & 487 scalar_rayleigh_damping, & 488 spinup_time, & 489 spinup_pt_amplitude, & 490 spinup_pt_mean, & 491 statistic_regions, & 492 subs_vertical_gradient, & 493 subs_vertical_gradient_level, & 494 surface_heatflux, & 495 surface_pressure, & 496 surface_scalarflux, & 497 surface_waterflux, & 498 timestep_scheme, & 499 topography, & 500 topography_grid_convention, & 501 top_heatflux, & 502 top_momentumflux_u, & 503 top_momentumflux_v, & 504 top_scalarflux, & 505 transpose_compute_overlap, & 506 tunnel_height, & 507 tunnel_length, & 508 tunnel_wall_depth, & 509 tunnel_width_x, & 510 tunnel_width_y, & 511 turbulence_closure, & 512 turbulent_inflow, & 513 turbulent_outflow, & 514 u_bulk, & 515 u_profile, & 516 ug_surface, & 517 ug_vertical_gradient, & 518 ug_vertical_gradient_level, & 519 use_fixed_date, & 520 use_fixed_time, & 521 use_free_convection_scaling, & 522 use_ug_for_galilei_tr, & 523 use_subsidence_tendencies, & 524 use_surface_fluxes, use_cmax, & 525 use_top_fluxes, & 526 use_upstream_for_tke, & 527 uv_heights, & 528 v_bulk, & 529 v_profile, & 530 vdi_checks, & 531 vg_surface, & 532 vg_vertical_gradient, & 533 vg_vertical_gradient_level, & 534 wall_adjustment, & 535 wall_heatflux, & 536 wall_humidityflux, & 537 wall_scalarflux, & 538 y_shift, & 539 zeta_max, & 540 zeta_min, & 541 z0h_factor 542 543 NAMELIST /d3par/ averaging_interval, & 544 averaging_interval_pr, & 545 cpu_log_barrierwait, & 546 create_disturbances, & 547 cross_profiles, & 548 data_output, & 549 data_output_2d_on_each_pe, & 550 data_output_masks, & 551 data_output_pr, & 552 debug_output, & 553 debug_output_timestep, & 554 disturbance_amplitude, & 555 disturbance_energy_limit, & 556 disturbance_level_b, & 557 disturbance_level_t, & 558 do2d_at_begin, & 559 do3d_at_begin, & 560 dt, & 561 dt_averaging_input, & 562 dt_averaging_input_pr, & 563 dt_coupling, & 564 dt_data_output, & 565 dt_data_output_av, & 566 dt_disturb, & 567 dt_domask, & 568 dt_dopr, & 569 dt_dopr_listing, & 570 dt_dots, & 571 dt_do2d_xy, & 572 dt_do2d_xz, & 573 dt_do2d_yz, & 574 dt_do3d, & 575 dt_max, & 576 dt_restart, & 577 dt_run_control, & 578 end_time, & 579 force_print_header, & 580 mask_k_over_surface, & 581 mask_scale_x, & 582 mask_scale_y, & 583 mask_scale_z, & 584 mask_x, & 585 mask_y, & 586 mask_z, & 587 mask_x_loop, & 588 mask_y_loop, & 589 mask_z_loop, & 590 netcdf_data_format, & 591 netcdf_deflate, & 592 normalizing_region, & 593 npex, & 594 npey, & 595 nz_do3d, & 596 profile_columns, & 597 profile_rows, & 598 restart_time, & 599 section_xy, & 600 section_xz, & 601 section_yz, & 602 skip_time_data_output, & 603 skip_time_data_output_av, & 604 skip_time_dopr, & 605 skip_time_do2d_xy, & 606 skip_time_do2d_xz, & 607 skip_time_do2d_yz, & 608 skip_time_do3d, & 609 skip_time_domask, & 610 synchronous_exchange, & 611 termination_time_needed 612 613 NAMELIST /runtime_parameters/ averaging_interval, & 614 averaging_interval_pr, & 615 cpu_log_barrierwait, & 616 create_disturbances, & 617 cross_profiles, & 618 data_output, & 619 data_output_2d_on_each_pe, & 620 data_output_masks, & 621 data_output_pr, & 622 debug_output, & 623 debug_output_timestep, & 624 disturbance_amplitude, & 625 disturbance_energy_limit, & 626 disturbance_level_b, & 627 disturbance_level_t, & 628 do2d_at_begin, & 629 do3d_at_begin, & 630 dt, & 631 dt_averaging_input, & 632 dt_averaging_input_pr, & 633 dt_coupling, & 634 dt_data_output, & 635 dt_data_output_av, & 636 dt_disturb, & 637 dt_domask, & 638 dt_dopr, & 639 dt_dopr_listing, & 640 dt_dots, & 641 dt_do2d_xy, & 642 dt_do2d_xz, & 643 dt_do2d_yz, & 644 dt_do3d, & 645 dt_max, & 646 dt_restart, & 647 dt_run_control, & 648 end_time, & 649 force_print_header, & 650 mask_k_over_surface, & 651 mask_scale_x, & 652 mask_scale_y, & 653 mask_scale_z, & 654 mask_x, & 655 mask_y, & 656 mask_z, & 657 mask_x_loop, & 658 mask_y_loop, & 659 mask_z_loop, & 660 netcdf_data_format, & 661 netcdf_deflate, & 662 normalizing_region, & 663 npex, & 664 npey, & 665 nz_do3d, & 666 profile_columns, & 667 profile_rows, & 668 restart_time, & 669 section_xy, & 670 section_xz, & 671 section_yz, & 672 restart_data_format, & 673 restart_data_format_input, & 674 restart_data_format_output, & 675 skip_time_data_output, & 676 skip_time_data_output_av, & 677 skip_time_dopr, & 678 skip_time_do2d_xy, & 679 skip_time_do2d_xz, & 680 skip_time_do2d_yz, & 681 skip_time_do3d, & 682 skip_time_domask, & 683 synchronous_exchange, & 684 termination_time_needed 685 686 NAMELIST /envpar/ host, & 687 maximum_cpu_time_allowed, & 688 maximum_parallel_io_streams, & 689 progress_bar_disabled, & 690 read_svf, & 691 revision, & 692 run_identifier, & 693 tasks_per_node, & 694 write_binary, & 695 write_svf 696 697 ! 698 !-- First read values of environment variables (this NAMELIST file is generated by palmrun) 374 699 CALL location_message( 'reading environment parameters from ENVPAR', 'start' ) 375 700 376 OPEN 701 OPEN( 90, FILE='ENVPAR', STATUS='OLD', FORM='FORMATTED', IOSTAT=ioerr ) 377 702 378 703 IF ( ioerr /= 0 ) THEN 379 message_string = 'local file ENVPAR not found' // &704 message_string = 'local file ENVPAR not found' // & 380 705 '&some variables for steering may not be properly set' 381 706 CALL message( 'parin', 'PA0276', 0, 1, 0, 6, 0 ) 382 707 ELSE 383 READ 708 READ( 90, envpar, IOSTAT=ioerr ) 384 709 IF ( ioerr < 0 ) THEN 385 message_string = 'no envpar-NAMELIST found in local file ' // &710 message_string = 'no envpar-NAMELIST found in local file ' // & 386 711 'ENVPAR& or some variables for steering may ' // & 387 712 'not be properly set' 388 713 CALL message( 'parin', 'PA0278', 0, 1, 0, 6, 0 ) 389 714 ELSEIF ( ioerr > 0 ) THEN 390 message_string = 'errors in local file ENVPAR' // &715 message_string = 'errors in local file ENVPAR' // & 391 716 '&some variables for steering may not be properly set' 392 717 CALL message( 'parin', 'PA0277', 0, 1, 0, 6, 0 ) 393 718 ENDIF 394 CLOSE 719 CLOSE( 90 ) 395 720 ENDIF 396 721 … … 414 739 415 740 11 REWIND ( 11 ) 416 READ 417 418 message_string = 'namelist inipar is deprecated and will be ' // &419 'removed in near future. & Please use namelist ' // &741 READ( 11, inipar, ERR=12, END=13 ) 742 743 message_string = 'namelist inipar is deprecated and will be ' // & 744 'removed in near future. & Please use namelist ' // & 420 745 'initialization_parameters instead' 421 746 CALL message( 'parin', 'PA0017', 0, 1, 0, 6, 0 ) … … 431 756 432 757 ! 433 !-- Try to read runtime parameters given by the user for this run 434 !-- (namelist "runtime_parameters"). The namelist "runtime_parmeters" 435 !-- can be omitted. In that case default values are used for the 758 !-- Try to read runtime parameters given by the user for this run (namelist "runtime_parameters"). 759 !-- The namelist "runtime_parmeters" can be omitted. In that case default values are used for the 436 760 !-- parameters. 437 761 14 line = ' ' 438 762 439 REWIND 763 REWIND( 11 ) 440 764 line = ' ' 441 765 DO WHILE ( INDEX( line, '&runtime_parameters' ) == 0 ) 442 READ 766 READ( 11, '(A)', END=16 ) line 443 767 ENDDO 444 BACKSPACE 768 BACKSPACE( 11 ) 445 769 446 770 ! 447 771 !-- Read namelist 448 READ 772 READ( 11, runtime_parameters, ERR = 15 ) 449 773 GOTO 18 450 774 … … 453 777 CALL parin_fail_message( 'runtime_parameters', line ) 454 778 455 16 REWIND 779 16 REWIND( 11 ) 456 780 line = ' ' 457 781 DO WHILE ( INDEX( line, '&d3par' ) == 0 ) 458 READ 782 READ( 11, '(A)', END=18 ) line 459 783 ENDDO 460 784 BACKSPACE ( 11 ) … … 462 786 ! 463 787 !-- Read namelist 464 READ 465 466 message_string = 'namelist d3par is deprecated and will be ' // &467 'removed in near future. &Please use namelist ' // &788 READ( 11, d3par, ERR = 17, END = 18 ) 789 790 message_string = 'namelist d3par is deprecated and will be ' // & 791 'removed in near future. &Please use namelist ' // & 468 792 'runtime_parameters instead' 469 793 CALL message( 'parin', 'PA0487', 0, 1, 0, 6, 0 ) … … 483 807 ! 484 808 !-- Calculate the number of groups into which parallel I/O is split. 485 !-- The default for files which are opened by all PEs (or where each 486 !-- PE opens his own independent file) is, that all PEs are doing input/output 487 !-- in parallel at the same time. This might cause performance or even more 488 !-- severe problems depending on the configuration of the underlying file 809 !-- The default for files which are opened by all PEs (or where each PE opens its own independent 810 !-- file) is, that all PEs are doing input/output in parallel at the same time. This might cause 811 !-- performance or even more severe problems depending on the configuration of the underlying file 489 812 !-- system. 490 !-- Calculation of the number of blocks and the I/O group must be based on all 491 !-- PEs involved in this run. Since myid and numprocs are related to the 492 !-- comm2d communicator, which gives only a subset of all PEs in case of 493 !-- nested runs, that information must be inquired again from the global 813 !-- Calculation of the number of blocks and the I/O group must be based on all PEs involved in this 814 !-- run. Since myid and numprocs are related to the comm2d communicator, which gives only a subset 815 !-- of all PEs in case of nested runs, that information must be inquired again from the global 494 816 !-- communicator. 495 817 !-- First, set the default: … … 501 823 global_procs = 1 502 824 #endif 503 IF ( maximum_parallel_io_streams == -1 .OR. & 504 maximum_parallel_io_streams > global_procs ) THEN 825 IF ( maximum_parallel_io_streams == -1 .OR. maximum_parallel_io_streams > global_procs ) THEN 505 826 maximum_parallel_io_streams = global_procs 506 827 ENDIF 507 828 ! 508 !-- Now calculate the number of io_blocks and the io_group to which the 509 !-- respective PE belongs. I/O of the groups is done in serial, but in parallel 510 !-- for all PEs belonging to the same group. 829 !-- Now calculate the number of io_blocks and the io_group to which the respective PE belongs. I/O 830 !-- of the groups is done in serial, but in parallel for all PEs belonging to the same group. 511 831 io_blocks = global_procs / maximum_parallel_io_streams 512 832 io_group = MOD( global_id+1, io_blocks ) 513 833 514 834 ! 515 !-- If required, read control parameters from restart file (produced by 516 !-- a prior run). All PEs are reading from file created by PE0 (see 517 !-- check_open) 835 !-- If required, read control parameters from restart file (produced by a prior run). All PEs are 836 !-- reading from file created by PE0 (see check_open) 518 837 IF ( TRIM( initializing_actions ) == 'read_restart_data' ) THEN 519 838 520 839 ! 521 !-- If not set by the user in the runtime parameters, the data format for restart 522 !-- input needs tobe set now! This is normally done later in check parameters.523 IF ( TRIM( restart_data_format ) /= 'fortran_binary' .AND. &524 TRIM( restart_data_format ) /= 'mpi' .AND. &840 !-- If not set by the user in the runtime parameters, the data format for restart input needs to 841 !-- be set now! This is normally done later in check parameters. 842 IF ( TRIM( restart_data_format ) /= 'fortran_binary' .AND. & 843 TRIM( restart_data_format ) /= 'mpi' .AND. & 525 844 TRIM( restart_data_format ) /= 'mpi_shared_memory' ) THEN 526 845 message_string = 'illegal restart data format "' // TRIM( restart_data_format ) // '"' … … 553 872 runnr = runnr + 1 554 873 ! 555 !-- In case of a restart run, the number of user-defined profiles on 556 !-- the restart file (already stored in max_pr_user) has to match the 557 !-- one given for the current run. max_pr_user_tmp is calculated in 558 !-- user_parin and max_pr_user is read in via rrd_global. 874 !-- In case of a restart run, the number of user-defined profiles on the restart file (already 875 !-- stored in max_pr_user) has to match the one given for the current run. max_pr_user_tmp is 876 !-- calculated in user_parin and max_pr_user is read in via rrd_global. 559 877 IF ( max_pr_user /= max_pr_user_tmp ) THEN 560 WRITE( message_string, * ) 'the number of user-defined ', &561 'profiles given in data_output_pr (', max_pr_user_tmp,&562 ') does not match the one ',&563 'found in the restart file (', max_pr_user, ')'878 WRITE( message_string, * ) 'the number of user-defined ', & 879 'profiles given in data_output_pr (', max_pr_user_tmp, & 880 ') does not match the one ', & 881 'found in the restart file (', max_pr_user, ')' 564 882 CALL message( 'user_parin', 'UI0009', 1, 2, 0, 6, 0 ) 565 883 ENDIF … … 569 887 ! 570 888 !-- Activate spinup 571 IF ( land_surface .OR.urban_surface ) THEN889 IF ( land_surface .OR. urban_surface ) THEN 572 890 IF ( spinup_time > 0.0_wp ) THEN 573 891 coupling_start_time = spinup_time … … 579 897 IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 580 898 spinup = .TRUE. 581 ELSEIF ( TRIM( initializing_actions ) == 'read_restart_data'&582 .AND.time_since_reference_point > 0.0_wp ) THEN899 ELSEIF ( TRIM( initializing_actions ) == 'read_restart_data' .AND. & 900 time_since_reference_point > 0.0_wp ) THEN 583 901 data_output_during_spinup = .FALSE. !< required for correct ntdim calculation 584 902 !< in check_parameters for restart run … … 593 911 IF ( nested_run ) THEN 594 912 IF ( nesting_mode == 'vertical' ) THEN 595 IF ( bc_lr /= 'cyclic' .OR.bc_ns /= 'cyclic' ) THEN596 WRITE ( message_string, *) 'bc_lr and bc_ns were set to ,', &597 'cyclic for vertical nesting'913 IF ( bc_lr /= 'cyclic' .OR. bc_ns /= 'cyclic' ) THEN 914 WRITE ( message_string, *) 'bc_lr and bc_ns were set to ,', & 915 'cyclic for vertical nesting' 598 916 CALL message( 'parin', 'PA0428', 0, 0, 0, 6, 0 ) 599 917 bc_lr = 'cyclic' … … 609 927 ENDIF 610 928 ! 611 !-- For other nesting modes only set boundary conditions for 612 !-- nested domains. 929 !-- For other nesting modes only set boundary conditions for nested domains. 613 930 ELSE 614 931 IF ( child_domain ) THEN … … 625 942 ENDIF 626 943 ! 627 !-- Set boundary conditions also in case the model is offline-nested in 628 !-- larger-scale models. 944 !-- Set boundary conditions also in case the model is offline-nested in larger-scale models. 629 945 IF ( nesting_offline ) THEN 630 946 bc_lr = 'nesting_offline' … … 639 955 640 956 ! 641 !-- In case of nested runs, make sure that initializing_actions = 642 !-- 'set_constant_profiles' even though the constant-profiles 643 !-- initializations for the prognostic variables will be overwritten 644 !-- by pmci_child_initialize and pmci_parent_initialize. This is, 645 !-- however, important e.g. to make sure that diagnostic variables 646 !-- are set properly. An exception is made in case of restart runs and 647 !-- if user decides to do everything by its own. 648 IF ( child_domain .AND. .NOT. ( & 649 TRIM( initializing_actions ) == 'read_restart_data' .OR. & 650 TRIM( initializing_actions ) == 'set_constant_profiles' .OR. & 651 TRIM( initializing_actions ) == 'by_user' ) ) THEN 652 message_string = 'initializing_actions = ' // & 653 TRIM( initializing_actions ) // ' has been ' // & 654 'changed to set_constant_profiles in child ' // & 655 'domain.' 957 !-- In case of nested runs, make sure that initializing_actions = 'set_constant_profiles' even 958 !-- though the constant-profiles initializations for the prognostic variables will be overwritten by 959 !-- pmci_child_initialize and pmci_parent_initialize. This is, however, important e.g. to make sure 960 !-- that diagnostic variables are set properly. An exception is made in case of restart runs and if 961 !-- user decides to do everything by its own. 962 IF ( child_domain .AND. .NOT. ( TRIM( initializing_actions ) == 'read_restart_data' .OR. & 963 TRIM( initializing_actions ) == 'set_constant_profiles' .OR. & 964 TRIM( initializing_actions ) == 'by_user' ) & 965 ) THEN 966 message_string = 'initializing_actions = ' // TRIM( initializing_actions ) // & 967 ' has been ' // 'changed to set_constant_profiles in child ' // 'domain.' 656 968 CALL message( 'parin', 'PA0492', 0, 0, 0, 6, 0 ) 657 969 … … 659 971 ENDIF 660 972 ! 661 !-- Check validity of lateral boundary conditions. This has to be done 662 !-- here because they are already used in init_pegrid and init_grid and 663 !-- therefore cannot be check in check_parameters 664 IF ( bc_lr /= 'cyclic' .AND. bc_lr /= 'dirichlet/radiation' .AND. & 665 bc_lr /= 'radiation/dirichlet' .AND. bc_lr /= 'nested' .AND. & 973 !-- Check validity of lateral boundary conditions. This has to be done here because they are already 974 !-- used in init_pegrid and init_grid and therefore cannot be check in check_parameters 975 IF ( bc_lr /= 'cyclic' .AND. bc_lr /= 'dirichlet/radiation' .AND. & 976 bc_lr /= 'radiation/dirichlet' .AND. bc_lr /= 'nested' .AND. & 666 977 bc_lr /= 'nesting_offline' ) THEN 667 message_string = 'unknown boundary condition: bc_lr = "' // & 668 TRIM( bc_lr ) // '"' 978 message_string = 'unknown boundary condition: bc_lr = "' // TRIM( bc_lr ) // '"' 669 979 CALL message( 'parin', 'PA0049', 1, 2, 0, 6, 0 ) 670 980 ENDIF 671 IF ( bc_ns /= 'cyclic' .AND. bc_ns /= 'dirichlet/radiation' .AND.&672 bc_ns /= 'radiation/dirichlet' .AND. bc_ns /= 'nested' .AND.&981 IF ( bc_ns /= 'cyclic' .AND. bc_ns /= 'dirichlet/radiation' .AND. & 982 bc_ns /= 'radiation/dirichlet' .AND. bc_ns /= 'nested' .AND. & 673 983 bc_ns /= 'nesting_offline' ) THEN 674 message_string = 'unknown boundary condition: bc_ns = "' // & 675 TRIM( bc_ns ) // '"' 984 message_string = 'unknown boundary condition: bc_ns = "' // TRIM( bc_ns ) // '"' 676 985 CALL message( 'parin', 'PA0050', 1, 2, 0, 6, 0 ) 677 986 ENDIF … … 686 995 ! 687 996 !-- Radiation-Dirichlet conditions are allowed along one of the horizontal directions only. 688 !-- In general, such conditions along x and y may work, but require a) some code changes 689 !-- (e.g. concerning allocation of c_u, c_v ... arrays), and b) a careful model setup by the690 !-- user, inorder to guarantee that there is a clearly defined outflow at two sides.997 !-- In general, such conditions along x and y may work, but require a) some code changes (e.g. 998 !-- concerning allocation of c_u, c_v ... arrays), and b) a careful model setup by the user, in 999 !-- order to guarantee that there is a clearly defined outflow at two sides. 691 1000 !-- Otherwise, the radiation condition may produce wrong results. 692 1001 IF ( ( bc_lr_dirrad .OR. bc_lr_raddir ) .AND. ( bc_ns_dirrad .OR. bc_ns_raddir ) ) THEN 693 message_string = 'bc_lr = "' // TRIM( bc_lr ) // '" and bc_ns = "' // 694 TRIM( bc_ns ) //'" are not allowed to be set at the same time'1002 message_string = 'bc_lr = "' // TRIM( bc_lr ) // '" and bc_ns = "' // TRIM( bc_ns ) // & 1003 '" are not allowed to be set at the same time' 695 1004 CALL message( 'parin', 'PA0589', 1, 2, 0, 6, 0 ) 696 1005 ENDIF 697 1006 ! 698 !-- Check in case of initial run, if the grid point numbers are well 699 !-- defined and allocate some arrays which are already needed in 700 !-- init_pegrid or check_parameters. During restart jobs, these arrays 701 !-- will be allocated in rrd_global. All other arrays are allocated 702 !-- in init_3d_model. 1007 !-- Check in case of initial run, if the grid point numbers are well defined and allocate some 1008 !-- arrays which are already needed in init_pegrid or check_parameters. During restart jobs, these 1009 !-- arrays will be allocated in rrd_global. All other arrays are allocated in init_3d_model. 703 1010 IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 704 1011 705 1012 IF ( nx <= 0 ) THEN 706 WRITE( message_string, * ) 'no value or wrong value given', & 707 ' for nx: nx=', nx 1013 WRITE( message_string, * ) 'no value or wrong value given', ' for nx: nx=', nx 708 1014 CALL message( 'parin', 'PA0273', 1, 2, 0, 6, 0 ) 709 1015 ENDIF 710 1016 IF ( ny <= 0 ) THEN 711 WRITE( message_string, * ) 'no value or wrong value given', & 712 ' for ny: ny=', ny 1017 WRITE( message_string, * ) 'no value or wrong value given', ' for ny: ny=', ny 713 1018 CALL message( 'parin', 'PA0274', 1, 2, 0, 6, 0 ) 714 1019 ENDIF 715 1020 IF ( nz <= 0 ) THEN 716 WRITE( message_string, * ) 'no value or wrong value given', & 717 ' for nz: nz=', nz 1021 WRITE( message_string, * ) 'no value or wrong value given', ' for nz: nz=', nz 718 1022 CALL message( 'parin', 'PA0275', 1, 2, 0, 6, 0 ) 719 1023 ENDIF 720 1024 721 1025 ! 722 !-- As a side condition, routine flow_statistics require at least 14 723 !-- vertical grid levels (theyare used to store time-series data)1026 !-- As a side condition, routine flow_statistics require at least 14 vertical grid levels (they 1027 !-- are used to store time-series data) 724 1028 !> @todo Remove this restriction 725 1029 IF ( nz < 14 ) THEN … … 729 1033 730 1034 ! 731 !-- ATTENTION: in case of changes to the following statement please 732 !-- also check the allocatestatement in routine rrd_global733 ALLOCATE( pt_init(0:nz+1), q_init(0:nz+1), s_init(0:nz+1), &734 ref_state(0:nz+1), sa_init(0:nz+1), ug(0:nz+1), &735 u_init(0:nz+1), v_init(0:nz+1), vg(0:nz+1), &736 hom(0:nz+1,2,pr_palm+max_pr_user+max_pr_cs+max_pr_salsa,0:statistic_regions), &1035 !-- ATTENTION: in case of changes to the following statement please also check the allocate 1036 !-- statement in routine rrd_global 1037 ALLOCATE( pt_init(0:nz+1), q_init(0:nz+1), s_init(0:nz+1), & 1038 ref_state(0:nz+1), sa_init(0:nz+1), ug(0:nz+1), & 1039 u_init(0:nz+1), v_init(0:nz+1), vg(0:nz+1), & 1040 hom(0:nz+1,2,pr_palm+max_pr_user+max_pr_cs+max_pr_salsa,0:statistic_regions), & 737 1041 hom_sum(0:nz+1,pr_palm+max_pr_user+max_pr_cs+max_pr_salsa,0:statistic_regions) ) 738 1042
Note: See TracChangeset
for help on using the changeset viewer.