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

Last change on this file since 4180 was 4180, checked in by scharf, 5 years ago

removed comments in 'Former revisions' section that are older than 01.01.2019

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