source: palm/trunk/SOURCE/surface_coupler.f90 @ 4759

Last change on this file since 4759 was 4671, checked in by pavelkrc, 4 years ago

Radiative transfer model RTM version 4.1

  • Property svn:keywords set to Id
File size: 30.5 KB
RevLine 
[1682]1!> @file surface_coupler.f90
[2000]2!------------------------------------------------------------------------------!
[2696]3! This file is part of the PALM model system.
[1036]4!
[2000]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.
[1036]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/>.
16!
[4360]17! Copyright 1997-2020 Leibniz Universitaet Hannover
[2000]18!------------------------------------------------------------------------------!
[1036]19!
[258]20! Current revisions:
[1092]21! ------------------
[4671]22!
23!
[1321]24! Former revisions:
25! -----------------
26! $Id: surface_coupler.f90 4671 2020-09-09 20:27:58Z suehring $
[4671]27! Implementation of downward facing USM and LSM surfaces
28!
29! 4429 2020-02-27 15:24:30Z raasch
[4429]30! bugfix: preprocessor directives rearranged for serial mode
31!
32! 4360 2020-01-07 11:25:50Z suehring
[4182]33! Corrected "Former revisions" section
34!
35! 3655 2019-01-07 16:51:22Z knoop
[3274]36! Modularization of all bulk cloud physics code components
[1321]37!
[4182]38! 109 2007-08-28 15:26:47Z letzel
39! Initial revision
40!
[102]41! Description:
42! ------------
[1682]43!> Data exchange at the interface between coupled models
[102]44!------------------------------------------------------------------------------!
[1682]45 SUBROUTINE surface_coupler
[4429]46#if defined( __parallel )
[1682]47 
[102]48
[1320]49    USE arrays_3d,                                                             &
[2232]50        ONLY:  pt, rho_ocean, sa, total_2d_a, total_2d_o, u, v
[1320]51
[3274]52    USE basic_constants_and_equations_mod,                                     &
53        ONLY:  c_p, l_v
[1427]54
[1320]55    USE control_parameters,                                                    &
56        ONLY:  coupling_mode, coupling_mode_remote, coupling_topology,         &
[2232]57               humidity, humidity_remote, land_surface, message_string,        &
58               terminate_coupled, terminate_coupled_remote,                    &
59               time_since_reference_point, urban_surface
[1320]60
61    USE cpulog,                                                                &
62        ONLY:  cpu_log, log_point
63
64    USE indices,                                                               &
65        ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, nx_a, nx_o, ny, nyn, nyng, nys, &
66               nysg, ny_a, ny_o, nzt
67
68    USE kinds
69
[102]70    USE pegrid
71
[2232]72    USE surface_mod,                                                           &
73        ONLY :  surf_def_h, surf_lsm_h, surf_type, surf_usm_h
74
[102]75    IMPLICIT NONE
76
[2232]77    INTEGER(iwp) ::  i                                    !< index variable x-direction
78    INTEGER(iwp) ::  j                                    !< index variable y-direction
79    INTEGER(iwp) ::  m                                    !< running index for surface elements
80
81    REAL(wp)    ::  cpw = 4218.0_wp                       !< heat capacity of water at constant pressure
[1682]82    REAL(wp)    ::  time_since_reference_point_rem        !<
83    REAL(wp)    ::  total_2d(-nbgp:ny+nbgp,-nbgp:nx+nbgp) !<
[102]84
[2232]85    REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) ::  surface_flux !< dummy array for surface fluxes on 2D grid
[1427]86
[2232]87
[667]88    CALL cpu_log( log_point(39), 'surface_coupler', 'start' )
[102]89
[667]90
91
[102]92!
[108]93!-- In case of model termination initiated by the remote model
94!-- (terminate_coupled_remote > 0), initiate termination of the local model.
95!-- The rest of the coupler must then be skipped because it would cause an MPI
96!-- intercomminucation hang.
97!-- If necessary, the coupler will be called at the beginning of the next
98!-- restart run.
[667]99
100    IF ( coupling_topology == 0 ) THEN
[709]101       CALL MPI_SENDRECV( terminate_coupled,        1, MPI_INTEGER, target_id, &
102                          0,                                                   &
103                          terminate_coupled_remote, 1, MPI_INTEGER, target_id, &
[667]104                          0, comm_inter, status, ierr )
105    ELSE
106       IF ( myid == 0) THEN
107          CALL MPI_SENDRECV( terminate_coupled,        1, MPI_INTEGER, &
108                             target_id, 0,                             &
109                             terminate_coupled_remote, 1, MPI_INTEGER, & 
110                             target_id, 0,                             &
111                             comm_inter, status, ierr )
112       ENDIF
[709]113       CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
114                       ierr )
[667]115
116       ALLOCATE( total_2d_a(-nbgp:ny_a+nbgp,-nbgp:nx_a+nbgp),       &
117                 total_2d_o(-nbgp:ny_o+nbgp,-nbgp:nx_o+nbgp) )
118
119    ENDIF
120
[108]121    IF ( terminate_coupled_remote > 0 )  THEN
[3045]122       WRITE( message_string, * ) 'remote model "',                            &
123                                  TRIM( coupling_mode_remote ),                &
124                                  '" terminated',                              &
[3046]125                                  '&with terminate_coupled_remote = ',         &
[3045]126                                  terminate_coupled_remote,                    &
[3046]127                                  '&local model  "', TRIM( coupling_mode ),    &
[3045]128                                  '" has',                                     &
[3046]129                                  '&terminate_coupled = ',                     &
[667]130                                   terminate_coupled
[258]131       CALL message( 'surface_coupler', 'PA0310', 1, 2, 0, 6, 0 )
[108]132       RETURN
133    ENDIF
[667]134 
[291]135
[108]136!
137!-- Exchange the current simulated time between the models,
[2232]138!-- currently just for total_2d
[709]139    IF ( coupling_topology == 0 ) THEN
140   
141       CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, 11, &
142                      comm_inter, ierr )
143       CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, target_id, &
144                      11, comm_inter, status, ierr )
[667]145    ELSE
[709]146
[667]147       IF ( myid == 0 ) THEN
[709]148
149          CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
150                         11, comm_inter, ierr )
151          CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL,        &
[667]152                         target_id, 11, comm_inter, status, ierr )
[709]153
[667]154       ENDIF
[709]155
156       CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
157                       ierr )
158
[667]159    ENDIF
[102]160
161!
162!-- Exchange the interface data
163    IF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
[667]164   
165!
[709]166!--    Horizontal grid size and number of processors is equal in ocean and
167!--    atmosphere
168       IF ( coupling_topology == 0 )  THEN
[102]169
170!
[2232]171!--       Send heat flux at bottom surface to the ocean. First, transfer from
172!--       1D surface type to 2D grid.
[4671]173          CALL transfer_1D_to_2D_equal( surf_def_h(0)%shf, surf_lsm_h(0)%shf,     &
174                                        surf_usm_h(0)%shf )
[2232]175          CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
176                         12, comm_inter, ierr )
[102]177!
[2232]178!--       Send humidity flux at bottom surface to the ocean. First, transfer
179!--       from 1D surface type to 2D grid.
[4671]180          CALL transfer_1D_to_2D_equal( surf_def_h(0)%qsws, surf_lsm_h(0)%qsws,   &
181                                        surf_usm_h(0)%qsws )
[667]182          IF ( humidity )  THEN
[2232]183             CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL,         &
184                            target_id, 13, comm_inter, ierr )
[667]185          ENDIF
186!
[709]187!--       Receive temperature at the bottom surface from the ocean
[2232]188          CALL MPI_RECV( pt(0,nysg,nxlg), 1, type_xy, target_id, 14,           &
[709]189                         comm_inter, status, ierr )
[108]190!
[2232]191!--       Send the momentum flux (u) at bottom surface to the ocean. First,
192!--       transfer from 1D surface type to 2D grid.
[4671]193          CALL transfer_1D_to_2D_equal( surf_def_h(0)%usws, surf_lsm_h(0)%usws,   &
194                                        surf_usm_h(0)%usws )
[2232]195          CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
196                         15, comm_inter, ierr )
[102]197!
[2232]198!--       Send the momentum flux (v) at bottom surface to the ocean. First,
199!--       transfer from 1D surface type to 2D grid.
[4671]200          CALL transfer_1D_to_2D_equal( surf_def_h(0)%vsws, surf_lsm_h(0)%vsws,   &
201                                        surf_usm_h(0)%vsws )
[2232]202          CALL MPI_SEND( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, &
203                         16, comm_inter, ierr )
[102]204!
[709]205!--       Receive u at the bottom surface from the ocean
[2232]206          CALL MPI_RECV( u(0,nysg,nxlg), 1, type_xy, target_id, 17,            &
[709]207                         comm_inter, status, ierr )
[667]208!
[709]209!--       Receive v at the bottom surface from the ocean
[2232]210          CALL MPI_RECV( v(0,nysg,nxlg), 1, type_xy, target_id, 18,            &
[709]211                         comm_inter, status, ierr )
[667]212!
213!--    Horizontal grid size or number of processors differs between
214!--    ocean and atmosphere
215       ELSE
216     
217!
[709]218!--       Send heat flux at bottom surface to the ocean
[1353]219          total_2d_a = 0.0_wp
220          total_2d   = 0.0_wp
[2232]221!
222!--       Transfer from 1D surface type to 2D grid.
[4671]223          CALL transfer_1D_to_2D_unequal( surf_def_h(0)%shf, surf_lsm_h(0)%shf,   &
224                                          surf_usm_h(0)%shf )
[709]225
[2232]226          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0,  &
[709]227                           comm2d, ierr )
228          CALL interpolate_to_ocean( 12 )   
[667]229!
[709]230!--       Send humidity flux at bottom surface to the ocean
231          IF ( humidity )  THEN
[1353]232             total_2d_a = 0.0_wp
233             total_2d   = 0.0_wp
[2232]234!
235!--          Transfer from 1D surface type to 2D grid.
236             CALL transfer_1D_to_2D_unequal( surf_def_h(0)%qsws,              &
[4671]237                                             surf_lsm_h(0)%qsws,                 &
238                                             surf_usm_h(0)%qsws )
[709]239
240             CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, &
241                              0, comm2d, ierr )
242             CALL interpolate_to_ocean( 13 )
[667]243          ENDIF
244!
[709]245!--       Receive temperature at the bottom surface from the ocean
246          IF ( myid == 0 )  THEN
[2232]247             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL,          &
[667]248                            target_id, 14, comm_inter, status, ierr )   
249          ENDIF
250          CALL MPI_BARRIER( comm2d, ierr )
[709]251          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
252                          ierr )
[667]253          pt(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
254!
[709]255!--       Send momentum flux (u) at bottom surface to the ocean
[1353]256          total_2d_a = 0.0_wp 
257          total_2d   = 0.0_wp
[2232]258!
259!--       Transfer from 1D surface type to 2D grid.
[4671]260          CALL transfer_1D_to_2D_unequal( surf_def_h(0)%usws, surf_lsm_h(0)%usws, &
261                                          surf_usm_h(0)%usws )
[709]262          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
263                           comm2d, ierr )
264          CALL interpolate_to_ocean( 15 )
[667]265!
[709]266!--       Send momentum flux (v) at bottom surface to the ocean
[1353]267          total_2d_a = 0.0_wp
268          total_2d   = 0.0_wp
[2232]269!
270!--       Transfer from 1D surface type to 2D grid.
[4671]271          CALL transfer_1D_to_2D_unequal( surf_def_h(0)%usws, surf_lsm_h(0)%usws, &
272                                          surf_usm_h(0)%usws )
[709]273          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
274                           comm2d, ierr )
275          CALL interpolate_to_ocean( 16 )
[667]276!
[709]277!--       Receive u at the bottom surface from the ocean
278          IF ( myid == 0 )  THEN
[667]279             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]280                            target_id, 17, comm_inter, status, ierr )
[667]281          ENDIF
282          CALL MPI_BARRIER( comm2d, ierr )
[709]283          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
284                          ierr )
[667]285          u(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
286!
[709]287!--       Receive v at the bottom surface from the ocean
288          IF ( myid == 0 )  THEN
[667]289             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]290                            target_id, 18, comm_inter, status, ierr )
[667]291          ENDIF
292          CALL MPI_BARRIER( comm2d, ierr )
[709]293          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
294                          ierr )
[667]295          v(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
296
297       ENDIF
298
[102]299    ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
300
301!
[667]302!--    Horizontal grid size and number of processors is equal
303!--    in ocean and atmosphere
304       IF ( coupling_topology == 0 ) THEN
305!
[709]306!--       Receive heat flux at the sea surface (top) from the atmosphere
[2232]307          CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 12, &
[709]308                         comm_inter, status, ierr )
[2232]309          CALL transfer_2D_to_1D_equal( surf_def_h(2)%shf )
[102]310!
[709]311!--       Receive humidity flux from the atmosphere (bottom)
[667]312!--       and add it to the heat flux at the sea surface (top)...
313          IF ( humidity_remote )  THEN
[2232]314             CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, &
[667]315                            target_id, 13, comm_inter, status, ierr )
[2232]316             CALL transfer_2D_to_1D_equal( surf_def_h(2)%qsws )
[667]317          ENDIF
318!
319!--       Send sea surface temperature to the atmosphere model
[709]320          CALL MPI_SEND( pt(nzt,nysg,nxlg), 1, type_xy, target_id, 14, &
321                         comm_inter, ierr )
[667]322!
323!--       Receive momentum flux (u) at the sea surface (top) from the atmosphere
[2232]324          CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 15, &
[709]325                         comm_inter, status, ierr )
[2232]326          CALL transfer_2D_to_1D_equal( surf_def_h(2)%usws )
[667]327!
328!--       Receive momentum flux (v) at the sea surface (top) from the atmosphere
[2232]329          CALL MPI_RECV( surface_flux(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 16, &
[709]330                         comm_inter, status, ierr )
[2232]331          CALL transfer_2D_to_1D_equal( surf_def_h(2)%vsws )
[667]332!
[709]333!--       Send u to the atmosphere
334          CALL MPI_SEND( u(nzt,nysg,nxlg), 1, type_xy, target_id, 17, &
335                         comm_inter, ierr )
[667]336!
[709]337!--       Send v to the atmosphere
338          CALL MPI_SEND( v(nzt,nysg,nxlg), 1, type_xy, target_id, 18, &
339                         comm_inter, ierr )
340!
[667]341!--    Horizontal gridsize or number of processors differs between
342!--    ocean and atmosphere
343       ELSE
344!
[709]345!--       Receive heat flux at the sea surface (top) from the atmosphere
346          IF ( myid == 0 )  THEN
[667]347             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]348                            target_id, 12, comm_inter, status, ierr )
[667]349          ENDIF
350          CALL MPI_BARRIER( comm2d, ierr )
[709]351          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
352                          ierr )
[2232]353          CALL transfer_2D_to_1D_unequal( surf_def_h(2)%shf )
[667]354!
[709]355!--       Receive humidity flux at the sea surface (top) from the atmosphere
356          IF ( humidity_remote )  THEN
357             IF ( myid == 0 )  THEN
[667]358                CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]359                               target_id, 13, comm_inter, status, ierr )
[667]360             ENDIF
361             CALL MPI_BARRIER( comm2d, ierr )
[709]362             CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, &
363                             comm2d, ierr)
[2232]364             CALL transfer_2D_to_1D_unequal( surf_def_h(2)%qsws )
[667]365          ENDIF
366!
367!--       Send surface temperature to atmosphere
[1353]368          total_2d_o = 0.0_wp
369          total_2d   = 0.0_wp
[667]370          total_2d(nys:nyn,nxl:nxr) = pt(nzt,nys:nyn,nxl:nxr)
371
[709]372          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
373                           comm2d, ierr) 
374          CALL interpolate_to_atmos( 14 )
[667]375!
[709]376!--       Receive momentum flux (u) at the sea surface (top) from the atmosphere
377          IF ( myid == 0 )  THEN
[667]378             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]379                            target_id, 15, comm_inter, status, ierr )
[667]380          ENDIF
381          CALL MPI_BARRIER( comm2d, ierr )
382          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]383                          0, comm2d, ierr )
[2232]384          CALL transfer_2D_to_1D_unequal( surf_def_h(2)%usws )
[667]385!
[709]386!--       Receive momentum flux (v) at the sea surface (top) from the atmosphere
387          IF ( myid == 0 )  THEN
[667]388             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]389                            target_id, 16, comm_inter, status, ierr )
[667]390          ENDIF
391          CALL MPI_BARRIER( comm2d, ierr )
[709]392          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
393                          ierr )
[2232]394          CALL transfer_2D_to_1D_unequal( surf_def_h(2)%vsws )
[667]395!
396!--       Send u to atmosphere
[1353]397          total_2d_o = 0.0_wp 
398          total_2d   = 0.0_wp
[667]399          total_2d(nys:nyn,nxl:nxr) = u(nzt,nys:nyn,nxl:nxr)
[709]400          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
401                           comm2d, ierr )
402          CALL interpolate_to_atmos( 17 )
[667]403!
404!--       Send v to atmosphere
[1353]405          total_2d_o = 0.0_wp
406          total_2d   = 0.0_wp
[667]407          total_2d(nys:nyn,nxl:nxr) = v(nzt,nys:nyn,nxl:nxr)
[709]408          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
409                           comm2d, ierr )
410          CALL interpolate_to_atmos( 18 )
[667]411       
412       ENDIF
413
414!
415!--    Conversions of fluxes received from atmosphere
416       IF ( humidity_remote )  THEN
[108]417!
[2232]418!--       Here top heat flux is still the sum of atmospheric bottom heat fluxes,
[709]419!--       * latent heat of vaporization in m2/s2, or 540 cal/g, or 40.65 kJ/mol
420!--       /(rho_atm(=1.0)*c_p)
[2232]421          DO  m = 1, surf_def_h(2)%ns
422             i = surf_def_h(2)%i(m)
423             j = surf_def_h(2)%j(m)
424             
425             surf_def_h(2)%shf(m) = surf_def_h(2)%shf(m) +                     &
[3274]426                                    surf_def_h(2)%qsws(m) * l_v / c_p
[709]427!
[2232]428!--          ...and convert it to a salinity flux at the sea surface (top)
429!--          following Steinhorn (1991), JPO 21, pp. 1681-1683:
430!--          S'w' = -S * evaporation / ( rho_water * ( 1 - S ) )
431             surf_def_h(2)%sasws(m) = -1.0_wp * sa(nzt,j,i) * 0.001_wp *       &
432                                      surf_def_h(2)%qsws(m) /                  &
433                                    ( rho_ocean(nzt,j,i) *                     &
434                                      ( 1.0_wp - sa(nzt,j,i) * 0.001_wp )      &
435                                    )
436          ENDDO
[108]437       ENDIF
438
439!
[102]440!--    Adjust the kinematic heat flux with respect to ocean density
[2232]441!--    (constants are the specific heat capacities for air and water), as well
442!--    as momentum fluxes
443       DO  m = 1, surf_def_h(2)%ns
444          i = surf_def_h(2)%i(m)
445          j = surf_def_h(2)%j(m)
446          surf_def_h(2)%shf(m) = surf_def_h(2)%shf(m) / rho_ocean(nzt,j,i) *   &
[3274]447                                 c_p / cpw
[102]448
[2232]449          surf_def_h(2)%usws(m) = surf_def_h(2)%usws(m) / rho_ocean(nzt,j,i)
450          surf_def_h(2)%vsws(m) = surf_def_h(2)%vsws(m) / rho_ocean(nzt,j,i)
451       ENDDO
[102]452
[667]453    ENDIF
454
[709]455    IF ( coupling_topology == 1 )  THEN
[667]456       DEALLOCATE( total_2d_o, total_2d_a )
457    ENDIF
458
459    CALL cpu_log( log_point(39), 'surface_coupler', 'stop' )
460
461
[2232]462     CONTAINS 
463
464!       Description:
465!------------------------------------------------------------------------------!
466!>      Data transfer from 1D surface-data type to 2D dummy array for equal
467!>      grids in atmosphere and ocean.
468!------------------------------------------------------------------------------!
469        SUBROUTINE transfer_1D_to_2D_equal( def_1d, lsm_1d, usm_1d )
470
471           IMPLICIT NONE
472
473            INTEGER(iwp) ::  i   !< running index x
474            INTEGER(iwp) ::  j   !< running index y
475            INTEGER(iwp) ::  m   !< running index surface type
476
477            REAL(wp), DIMENSION(1:surf_def_h(0)%ns) ::  def_1d !< 1D surface flux, default surfaces
[4671]478            REAL(wp), DIMENSION(1:surf_lsm_h(0)%ns) ::  lsm_1d !< 1D surface flux, natural surfaces
479            REAL(wp), DIMENSION(1:surf_usm_h(0)%ns) ::  usm_1d !< 1D surface flux, urban surfaces
[2232]480!
481!--         Transfer surface flux at default surfaces to 2D grid
482            DO  m = 1, surf_def_h(0)%ns
483               i = surf_def_h(0)%i(m)
484               j = surf_def_h(0)%j(m)
485               surface_flux(j,i) = def_1d(m)
486            ENDDO
487!
488!--         Transfer surface flux at natural surfaces to 2D grid
489            IF ( land_surface )  THEN
490               DO  m = 1, SIZE(lsm_1d)
[4671]491                  i = surf_lsm_h(0)%i(m)
492                  j = surf_lsm_h(0)%j(m)
[2232]493                  surface_flux(j,i) = lsm_1d(m)
494               ENDDO
495            ENDIF
496!
497!--         Transfer surface flux at natural surfaces to 2D grid
498            IF ( urban_surface )  THEN
499               DO  m = 1, SIZE(usm_1d)
[4671]500                  i = surf_usm_h(0)%i(m)
501                  j = surf_usm_h(0)%j(m)
[2232]502                  surface_flux(j,i) = usm_1d(m)
503               ENDDO
504            ENDIF
505
506        END SUBROUTINE transfer_1D_to_2D_equal
507
508!       Description:
509!------------------------------------------------------------------------------!
510!>      Data transfer from 2D array for equal grids onto 1D surface-data type
511!>      array.
512!------------------------------------------------------------------------------!
513        SUBROUTINE transfer_2D_to_1D_equal( def_1d )
514
515           IMPLICIT NONE
516
517            INTEGER(iwp) ::  i   !< running index x
518            INTEGER(iwp) ::  j   !< running index y
519            INTEGER(iwp) ::  m   !< running index surface type
520
521            REAL(wp), DIMENSION(1:surf_def_h(2)%ns) ::  def_1d !< 1D surface flux, default surfaces
522!
523!--         Transfer surface flux to 1D surface type, only for default surfaces
524            DO  m = 1, surf_def_h(2)%ns
525               i = surf_def_h(2)%i(m)
526               j = surf_def_h(2)%j(m)
527               def_1d(m) = surface_flux(j,i)
528            ENDDO
529
530        END SUBROUTINE transfer_2D_to_1D_equal
531
532!       Description:
533!------------------------------------------------------------------------------!
534!>      Data transfer from 1D surface-data type to 2D dummy array from unequal
535!>      grids in atmosphere and ocean.
536!------------------------------------------------------------------------------!
537        SUBROUTINE transfer_1D_to_2D_unequal( def_1d, lsm_1d, usm_1d )
538
539           IMPLICIT NONE
540
541            INTEGER(iwp) ::  i   !< running index x
542            INTEGER(iwp) ::  j   !< running index y
543            INTEGER(iwp) ::  m   !< running index surface type
544
545            REAL(wp), DIMENSION(1:surf_def_h(0)%ns) ::  def_1d !< 1D surface flux, default surfaces
[4671]546            REAL(wp), DIMENSION(1:surf_lsm_h(0)%ns)    ::  lsm_1d !< 1D surface flux, natural surfaces
547            REAL(wp), DIMENSION(1:surf_usm_h(0)%ns)    ::  usm_1d !< 1D surface flux, urban surfaces
[2232]548!
549!--         Transfer surface flux at default surfaces to 2D grid. Transfer no
550!--         ghost-grid points since total_2d is a global array.
551            DO  m = 1, SIZE(def_1d)
552               i = surf_def_h(0)%i(m)
553               j = surf_def_h(0)%j(m)
554
555               IF ( i >= nxl  .AND.  i <= nxr  .AND.                           &
556                    j >= nys  .AND.  j <= nyn )  THEN
557                  total_2d(j,i) = def_1d(m)
558               ENDIF
559            ENDDO
560!
561!--         Transfer surface flux at natural surfaces to 2D grid
562            IF ( land_surface )  THEN
563               DO  m = 1, SIZE(lsm_1d)
[4671]564                  i = surf_lsm_h(0)%i(m)
565                  j = surf_lsm_h(0)%j(m)
[2232]566
567                  IF ( i >= nxl  .AND.  i <= nxr  .AND.                        &
568                       j >= nys  .AND.  j <= nyn )  THEN
569                     total_2d(j,i) = lsm_1d(m)
570                  ENDIF
571               ENDDO
572            ENDIF
573!
574!--         Transfer surface flux at natural surfaces to 2D grid
575            IF ( urban_surface )  THEN
576               DO  m = 1, SIZE(usm_1d)
[4671]577                  i = surf_usm_h(0)%i(m)
578                  j = surf_usm_h(0)%j(m)
[2232]579
580                  IF ( i >= nxl  .AND.  i <= nxr  .AND.                        &
581                       j >= nys  .AND.  j <= nyn )  THEN
582                     total_2d(j,i) = usm_1d(m)
583                  ENDIF
584               ENDDO
585            ENDIF
586
587        END SUBROUTINE transfer_1D_to_2D_unequal
588
589!       Description:
590!------------------------------------------------------------------------------!
591!>      Data transfer from 2D dummy array from unequal grids to 1D surface-data
592!>      type.
593!------------------------------------------------------------------------------!
594        SUBROUTINE transfer_2D_to_1D_unequal( def_1d )
595
596           IMPLICIT NONE
597
598            INTEGER(iwp) ::  i   !< running index x
599            INTEGER(iwp) ::  j   !< running index y
600            INTEGER(iwp) ::  m   !< running index surface type
601
602            REAL(wp), DIMENSION(1:surf_def_h(2)%ns) ::  def_1d !< 1D surface flux, default surfaces
603!
604!--         Transfer 2D surface flux to default surfaces data type. Transfer no
605!--         ghost-grid points since total_2d is a global array.
606            DO  m = 1, SIZE(def_1d)
607               i = surf_def_h(2)%i(m)
608               j = surf_def_h(2)%j(m)
609
610               IF ( i >= nxl  .AND.  i <= nxr  .AND.                           &
611                    j >= nys  .AND.  j <= nyn )  THEN
612                  def_1d(m) = total_2d_o(j,i)
613               ENDIF
614            ENDDO
615
616
617        END SUBROUTINE transfer_2D_to_1D_unequal
618
[4429]619#endif
[667]620  END SUBROUTINE surface_coupler
621
622
623
[1682]624!------------------------------------------------------------------------------!
625! Description:
626! ------------
627!> @todo Missing subroutine description.
628!------------------------------------------------------------------------------!
[4429]629#if defined( __parallel )
630
[709]631  SUBROUTINE interpolate_to_atmos( tag )
[667]632
[1320]633    USE arrays_3d,                                                             &
634        ONLY:  total_2d_a, total_2d_o
[667]635
[1320]636    USE indices,                                                               &
637        ONLY:  nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
638
639    USE kinds
640
[1324]641    USE pegrid
[1320]642
[667]643    IMPLICIT NONE
644
[1682]645    INTEGER(iwp) ::  dnx  !<
646    INTEGER(iwp) ::  dnx2 !<
647    INTEGER(iwp) ::  dny  !<
648    INTEGER(iwp) ::  dny2 !<
649    INTEGER(iwp) ::  i    !<
650    INTEGER(iwp) ::  ii   !<
651    INTEGER(iwp) ::  j    !<
652    INTEGER(iwp) ::  jj   !<
[667]653
[1682]654    INTEGER(iwp), intent(in) ::  tag !<
[1320]655
[667]656    CALL MPI_BARRIER( comm2d, ierr )
657
[709]658    IF ( myid == 0 )  THEN
659!
660!--    Cyclic boundary conditions for the total 2D-grid
[667]661       total_2d_o(-nbgp:-1,:) = total_2d_o(ny+1-nbgp:ny,:)
662       total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx+1-nbgp:nx)
663
664       total_2d_o(ny+1:ny+nbgp,:) = total_2d_o(0:nbgp-1,:)
665       total_2d_o(:,nx+1:nx+nbgp) = total_2d_o(:,0:nbgp-1)
666
[102]667!
[667]668!--    Number of gridpoints of the fine grid within one mesh of the coarse grid
669       dnx = (nx_o+1) / (nx_a+1) 
670       dny = (ny_o+1) / (ny_a+1) 
[102]671
672!
[709]673!--    Distance for interpolation around coarse grid points within the fine
674!--    grid (note: 2*dnx2 must not be equal with dnx)
[667]675       dnx2 = 2 * ( dnx / 2 )
676       dny2 = 2 * ( dny / 2 )
[102]677
[1353]678       total_2d_a = 0.0_wp
[102]679!
[667]680!--    Interpolation from ocean-grid-layer to atmosphere-grid-layer
681       DO  j = 0, ny_a
682          DO  i = 0, nx_a 
683             DO  jj = 0, dny2
684                DO  ii = 0, dnx2
685                   total_2d_a(j,i) = total_2d_a(j,i) &
686                                     + total_2d_o(j*dny+jj,i*dnx+ii)
687                ENDDO
688             ENDDO
689             total_2d_a(j,i) = total_2d_a(j,i) / ( ( dnx2 + 1 ) * ( dny2 + 1 ) )
690          ENDDO
691       ENDDO
692!
[709]693!--    Cyclic boundary conditions for atmosphere grid
[667]694       total_2d_a(-nbgp:-1,:) = total_2d_a(ny_a+1-nbgp:ny_a,:)
695       total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx_a+1-nbgp:nx_a)
696       
697       total_2d_a(ny_a+1:ny_a+nbgp,:) = total_2d_a(0:nbgp-1,:)
698       total_2d_a(:,nx_a+1:nx_a+nbgp) = total_2d_a(:,0:nbgp-1)
699!
700!--    Transfer of the atmosphere-grid-layer to the atmosphere
[709]701       CALL MPI_SEND( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, target_id, &
702                      tag, comm_inter, ierr )
[102]703
704    ENDIF
705
[667]706    CALL MPI_BARRIER( comm2d, ierr )
[102]707
[4429]708  END SUBROUTINE interpolate_to_atmos
709
[880]710#endif
711
[102]712
[1682]713!------------------------------------------------------------------------------!
714! Description:
715! ------------
716!> @todo Missing subroutine description.
717!------------------------------------------------------------------------------!
[4429]718#if defined( __parallel )
719
[709]720  SUBROUTINE interpolate_to_ocean( tag )
[667]721
[1320]722    USE arrays_3d,                                                             &
723        ONLY:  total_2d_a, total_2d_o
[667]724
[1320]725    USE indices,                                                               &
726        ONLY:  nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
727
728    USE kinds
729
[1324]730    USE pegrid
[1320]731
[667]732    IMPLICIT NONE
733
[1682]734    INTEGER(iwp)             ::  dnx !<
735    INTEGER(iwp)             ::  dny !<
736    INTEGER(iwp)             ::  i   !<
737    INTEGER(iwp)             ::  ii  !<
738    INTEGER(iwp)             ::  j   !<
739    INTEGER(iwp)             ::  jj  !<
740    INTEGER(iwp), intent(in) ::  tag !<
[667]741
[1682]742    REAL(wp)                 ::  fl  !<
743    REAL(wp)                 ::  fr  !<
744    REAL(wp)                 ::  myl !<
745    REAL(wp)                 ::  myr !<
[709]746
[667]747    CALL MPI_BARRIER( comm2d, ierr )
748
[709]749    IF ( myid == 0 )  THEN   
[667]750
751!
[709]752!--    Number of gridpoints of the fine grid within one mesh of the coarse grid
[667]753       dnx = ( nx_o + 1 ) / ( nx_a + 1 ) 
754       dny = ( ny_o + 1 ) / ( ny_a + 1 ) 
755
756!
[709]757!--    Cyclic boundary conditions for atmosphere grid
[667]758       total_2d_a(-nbgp:-1,:) = total_2d_a(ny+1-nbgp:ny,:)
759       total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx+1-nbgp:nx)
760       
761       total_2d_a(ny+1:ny+nbgp,:) = total_2d_a(0:nbgp-1,:)
762       total_2d_a(:,nx+1:nx+nbgp) = total_2d_a(:,0:nbgp-1)
763!
[709]764!--    Bilinear Interpolation from atmosphere grid-layer to ocean grid-layer
[667]765       DO  j = 0, ny
766          DO  i = 0, nx
767             myl = ( total_2d_a(j+1,i)   - total_2d_a(j,i)   ) / dny
768             myr = ( total_2d_a(j+1,i+1) - total_2d_a(j,i+1) ) / dny
769             DO  jj = 0, dny-1
[709]770                fl = myl*jj + total_2d_a(j,i) 
771                fr = myr*jj + total_2d_a(j,i+1) 
[667]772                DO  ii = 0, dnx-1
773                   total_2d_o(j*dny+jj,i*dnx+ii) = ( fr - fl ) / dnx * ii + fl
774                ENDDO
775             ENDDO
776          ENDDO
777       ENDDO
778!
[709]779!--    Cyclic boundary conditions for ocean grid
[667]780       total_2d_o(-nbgp:-1,:) = total_2d_o(ny_o+1-nbgp:ny_o,:)
781       total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx_o+1-nbgp:nx_o)
782
783       total_2d_o(ny_o+1:ny_o+nbgp,:) = total_2d_o(0:nbgp-1,:)
784       total_2d_o(:,nx_o+1:nx_o+nbgp) = total_2d_o(:,0:nbgp-1)
785
786       CALL MPI_SEND( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
787                      target_id, tag, comm_inter, ierr )
788
789    ENDIF
790
791    CALL MPI_BARRIER( comm2d, ierr ) 
792
[4429]793  END SUBROUTINE interpolate_to_ocean
794
[880]795#endif
Note: See TracBrowser for help on using the repository browser.