Ignore:
Timestamp:
Oct 22, 2018 5:32:49 PM (4 years ago)
Author:
eckhard
Message:

inifor: Added computation of geostrophic winds from COSMO input

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/UTIL/inifor/src/inifor.f90

    r3183 r3395  
    2626! -----------------
    2727! $Id$
     28! Added main loop support for computation of geostrophic winds and surface
     29!     pressure
     30! Cleaned up terminal output, show some meVssages only with --debug option
     31!
     32! 3183 2018-07-27 14:25:55Z suehring
    2833! Introduced new PALM grid stretching
    2934! Renamend initial-condition mode variable 'mode' to 'ic_mode'
     
    5661               fini_file_lists, preprocess, origin_lon, origin_lat,            &
    5762               output_file, io_group_list, output_var_table,                   &
    58                cosmo_grid, palm_grid, nx, ny, nz, ug, vg, p0, cfg,             &
    59                average_imin, average_imax, average_jmin, average_jmax
    60 
     63               cosmo_grid, palm_grid, nx, ny, nz, p0, cfg, f3,                 &
     64               averaging_width_ns, averaging_width_ew, phi_n, lambda_n,        &
     65               lam_centre, phi_centre
    6166    USE io
    6267    USE transform,                                                             &
    63         ONLY:  average_profile, average_2d, interpolate_2d, interpolate_3d
     68        ONLY:  average_profile, interpolate_2d, interpolate_3d,                &
     69               geostrophic_winds, extrapolate_density, extrapolate_pressure,   &
     70               get_surface_pressure
    6471    USE types
    6572   
     
    7077    INTEGER                                 ::  iter
    7178    REAL(dp), ALLOCATABLE, DIMENSION(:,:,:) ::  output_arr
     79    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  rho_centre
     80    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  ug_arr
     81    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  vg_arr
     82    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  rho_north
     83    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  rho_south
     84    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  rho_east
     85    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  rho_west
     86    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  p_north
     87    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  p_south
     88    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  p_east
     89    REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET ::  p_west
     90    REAL(dp), POINTER, DIMENSION(:)         ::  internal_arr
     91    REAL(dp), POINTER, DIMENSION(:)         ::  ug_vg_arr
    7292    TYPE(nc_var), POINTER                   ::  output_var
    7393    TYPE(io_group), POINTER                 ::  group
    7494    TYPE(container), ALLOCATABLE            ::  input_buffer(:)
     95    LOGICAL, SAVE                           ::  ug_vg_have_been_computed = .FALSE.
     96    LOGICAL, SAVE                           ::  debugging_output = .TRUE.
    7597   
    7698!> \mainpage About INIFOR
     
    199221 CALL run_control('time', 'alloc')
    200222                     
    201 
    202223                      CALL average_profile(                                    &
    203224                         input_buffer(output_var % input_id) % array(:,:,:),   &
    204                          output_arr(:,:,:), average_imin, average_imax,        &
    205                          average_jmin, average_jmax,                           &
    206                          output_var % intermediate_grid,                       &
    207                          output_var % grid)
    208  CALL run_control('time', 'comp')
     225                         output_arr(0,0,:),                                    &
     226                         output_var % averaging_grid)
     227
     228                      IF ( TRIM(output_var % name) ==                          &
     229                           'surface_forcing_surface_pressure' )  THEN
     230
     231                         IF ( cfg % p0_is_set )  THEN
     232                            output_arr(0,0,1) = p0
     233                         ELSE
     234                            CALL get_surface_pressure(                         &
     235                               output_arr(0,0,:), rho_centre,                  &
     236                               output_var % averaging_grid)
     237                         END IF
     238
     239                      END IF
     240 CALL run_control('time', 'comp')
     241
     242                   CASE ( 'internal profile' )
     243
     244                      message = "Averaging of internal profile for variable '" //&
     245                         TRIM(output_var % name) // "' is not supported."
     246
     247                      SELECT CASE (TRIM(output_var % name))
     248
     249                      CASE('internal_density_centre')
     250                         ALLOCATE( rho_centre( 1:output_var % grid % nz ) )
     251                         internal_arr => rho_centre
     252
     253                      CASE('internal_density_north')
     254                         ALLOCATE( rho_north( 1:output_var % grid % nz ) )
     255                         internal_arr => rho_north
     256
     257                      CASE('internal_density_south')
     258                         ALLOCATE( rho_south( 1:output_var % grid % nz ) )
     259                         internal_arr => rho_south
     260
     261                      CASE('internal_density_east')
     262                         ALLOCATE( rho_east( 1:output_var % grid % nz) )
     263                         internal_arr => rho_east
     264
     265                      CASE('internal_density_west')
     266                         ALLOCATE( rho_west( 1:output_var % grid % nz ) )
     267                         internal_arr => rho_west
     268
     269                      CASE('internal_pressure_north')
     270                         ALLOCATE( p_north( 1:output_var % grid % nz ) )
     271                         internal_arr => p_north
     272
     273                      CASE('internal_pressure_south')
     274                         ALLOCATE( p_south( 1:output_var % grid % nz ) )
     275                         internal_arr => p_south
     276
     277                      CASE('internal_pressure_east')
     278                         ALLOCATE( p_east( 1:output_var % grid % nz) )
     279                         internal_arr => p_east
     280
     281                      CASE('internal_pressure_west')
     282                         ALLOCATE( p_west( 1:output_var % grid % nz ) )
     283                         internal_arr => p_west
     284
     285                      CASE DEFAULT
     286                         CALL abort('main loop', message)
     287
     288                      END SELECT
     289 CALL run_control('time', 'alloc')
     290
     291
     292                      CALL average_profile(                                 &
     293                         input_buffer(output_var % input_id) % array(:,:,:),&
     294                         internal_arr(:),                                   &
     295                         output_var % averaging_grid)
     296
     297                      SELECT CASE (TRIM(output_var % name))
     298
     299                      CASE('internal_density_centre',                          &
     300                           'internal_density_north',                           &
     301                           'internal_density_south',                           &
     302                           'internal_density_east',                            &
     303                           'internal_density_west')
     304                         CALL extrapolate_density(internal_arr,                &
     305                                                  output_var % averaging_grid)
     306
     307                      CASE('internal_pressure_north')
     308                         CALL extrapolate_pressure(internal_arr, rho_north,    &
     309                                                   output_var % averaging_grid)
     310
     311                      CASE('internal_pressure_south')
     312                         CALL extrapolate_pressure(internal_arr, rho_south,    &
     313                                                   output_var % averaging_grid)
     314
     315                      CASE('internal_pressure_east')
     316                         CALL extrapolate_pressure(internal_arr, rho_east,     &
     317                                                   output_var % averaging_grid)
     318
     319                      CASE('internal_pressure_west')
     320                         CALL extrapolate_pressure(internal_arr, rho_west,     &
     321                                                   output_var % averaging_grid)
     322
     323                      CASE DEFAULT
     324                         CALL abort('main loop', message)
     325
     326                      END SELECT
     327
     328                      IF (.TRUE.)  THEN
     329                         ALLOCATE( output_arr(1,1,1:output_var % grid % nz) )
     330                         output_arr(1,1,:) = internal_arr(:)
     331                      END IF
     332 CALL run_control('time', 'comp')
     333
     334
     335                   ! This case gets called twice, the first time for ug, the
     336                   ! second time for vg. We compute ug and vg at the first call
     337                   ! and keep vg (and ug for that matter) around for the second
     338                   ! call.
     339                   CASE ( 'geostrophic winds' )
     340
     341                      IF (.NOT. ug_vg_have_been_computed )  THEN
     342                         ALLOCATE( ug_arr(output_var % grid % nz) )
     343                         ALLOCATE( vg_arr(output_var % grid % nz) )
     344
     345                         IF ( cfg % ug_is_set )  THEN
     346                            ug_arr = cfg % ug
     347                            vg_arr = cfg % vg
     348                         ELSE
     349                            CALL geostrophic_winds( p_north, p_south, p_east,     &
     350                                                    p_west, rho_centre, f3,      &
     351                                                    averaging_width_ew,           &
     352                                                    averaging_width_ns,           &
     353                                                    phi_n, lambda_n, phi_centre, lam_centre, &
     354                                                    ug_arr, vg_arr )
     355                         END IF
     356
     357                         ug_vg_have_been_computed = .TRUE.
     358
     359                      END IF
     360
     361                      ! Prepare output of geostrophic winds
     362                      SELECT CASE(TRIM(output_var % name))
     363                      CASE ('ls_forcing_ug')
     364                         ug_vg_arr => ug_arr
     365                      CASE ('ls_forcing_vg')
     366                         ug_vg_arr => vg_arr
     367                      END SELECT
     368
     369                      ALLOCATE( output_arr(1,1,output_var % grid % nz) )
     370                      output_arr(1,1,:) = ug_vg_arr(:)
    209371
    210372                   CASE ( 'average scalar' )
     
    222384                      SELECT CASE (TRIM(output_var % name))
    223385
    224                       CASE('ls_forcing_ug')
    225                           output_arr(1, 1, :) = ug
    226 
    227                       CASE('ls_forcing_vg')
    228                           output_arr(1, 1, :) = vg
     386                      !CASE('ls_forcing_ug')
     387                      !    output_arr(1, 1, :) = ug
     388
     389                      !CASE('ls_forcing_vg')
     390                      !    output_arr(1, 1, :) = vg
    229391
    230392                      CASE('nudging_tau')
     
    256418!- Section 2.3: Write current time step of current variable
    257419!------------------------------------------------------------------------------
    258                    message = "Writing variable '" // TRIM(output_var%name) // "'."
    259                    CALL report('main loop', message)
    260                    CALL update_output(output_var, output_arr, iter, output_file)
     420                   IF (.NOT. output_var % is_internal .OR. debugging_output)  THEN
     421                      message = "Writing variable '" // TRIM(output_var%name) // "'."
     422                      CALL report('main loop', message)
     423                      CALL update_output(output_var, output_arr, iter, output_file)
    261424 CALL run_control('time', 'write')
    262 
    263                    DEALLOCATE(output_arr)
     425                   END IF
     426
     427                   IF (ALLOCATED(output_arr))  DEALLOCATE(output_arr)
    264428 CALL run_control('time', 'alloc')
    265429
     
    267431
    268432             END DO ! ouput variables
     433
     434             IF ( group % kind == 'thermodynamics' )  THEN
     435                DEALLOCATE( rho_centre )
     436                IF ( .NOT. cfg % ug_is_set )  THEN
     437                   DEALLOCATE( rho_north )
     438                   DEALLOCATE( rho_south )
     439                   DEALLOCATE( rho_east )
     440                   DEALLOCATE( rho_west )
     441                   DEALLOCATE( p_north )
     442                   DEALLOCATE( p_south )
     443                   DEALLOCATE( p_east )
     444                   DEALLOCATE( p_west )
     445                END IF
     446             END IF
    269447
    270448             IF ( group % kind == 'running average' .OR. &
     
    273451                ! accumulated COSMO-DE quantities (precipitation).
    274452             ELSE
    275                 CALL report('main loop', 'Deallocating input buffer')
     453                CALL report('main loop', 'Deallocating input buffer', cfg % debug)
    276454                DEALLOCATE(input_buffer)
    277455             END IF
     
    281459
    282460          IF (ALLOCATED(input_buffer))  THEN
    283              CALL report('main loop', 'Deallocating input buffer')
     461             CALL report('main loop', 'Deallocating input buffer', cfg % debug)
    284462             DEALLOCATE(input_buffer)
    285463          END IF
     
    294472          END IF
    295473
    296           CALL report('main loop', message)
     474          CALL report('main loop', message, cfg % debug)
    297475
    298476       END IF ! IO group % to_be_processed
Note: See TracChangeset for help on using the changeset viewer.