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

Last change on this file since 1321 was 1321, checked in by raasch, 10 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 21.4 KB
RevLine 
[102]1 SUBROUTINE surface_coupler
2
[1036]3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later 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!
[1310]17! Copyright 1997-2014 Leibniz Universitaet Hannover
[1036]18!--------------------------------------------------------------------------------!
19!
[258]20! Current revisions:
[1092]21! ------------------
[1321]22!
23!
24! Former revisions:
25! -----------------
26! $Id: surface_coupler.f90 1321 2014-03-20 09:40:40Z raasch $
27!
28! 1320 2014-03-20 08:40:49Z raasch
[1320]29! ONLY-attribute added to USE-statements,
30! kind-parameters added to all INTEGER and REAL declaration statements,
31! kinds are defined in new module kinds,
32! old module precision_kind is removed,
33! revision history before 2012 removed,
34! comment fields (!:) to be used for variable explanations added to
35! all variable declaration statements
[102]36!
[1319]37! 1318 2014-03-17 13:35:16Z raasch
38! module interfaces removed
39!
[1093]40! 1092 2013-02-02 11:24:22Z raasch
41! unused variables removed
42!
[1037]43! 1036 2012-10-22 13:43:42Z raasch
44! code put under GPL (PALM 3.9)
45!
[881]46! 880 2012-04-13 06:28:59Z raasch
47! Bugfix: preprocessor statements for parallel execution added
48!
[110]49! 109 2007-08-28 15:26:47Z letzel
[102]50! Initial revision
51!
52! Description:
53! ------------
54! Data exchange at the interface between coupled models
55!------------------------------------------------------------------------------!
56
[1320]57    USE arrays_3d,                                                             &
58        ONLY:  pt, shf, qsws, qswst_remote, rho, sa, saswst, total_2d_a,       &
59               total_2d_o, tswst, u, usws, uswst, v, vsws, vswst
60
61    USE control_parameters,                                                    &
62        ONLY:  coupling_mode, coupling_mode_remote, coupling_topology,         &
63               humidity, humidity_remote, message_string, terminate_coupled,   &
64               terminate_coupled_remote, time_since_reference_point
65
66    USE cpulog,                                                                &
67        ONLY:  cpu_log, log_point
68
69    USE indices,                                                               &
70        ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, nx_a, nx_o, ny, nyn, nyng, nys, &
71               nysg, ny_a, ny_o, nzt
72
73    USE kinds
74
[102]75    USE pegrid
76
77    IMPLICIT NONE
78
[1320]79    REAL(wp)    ::  time_since_reference_point_rem        !:
80    REAL(wp)    ::  total_2d(-nbgp:ny+nbgp,-nbgp:nx+nbgp) !:
[102]81
[206]82#if defined( __parallel )
[102]83
[667]84    CALL cpu_log( log_point(39), 'surface_coupler', 'start' )
[102]85
[667]86
87
[102]88!
[108]89!-- In case of model termination initiated by the remote model
90!-- (terminate_coupled_remote > 0), initiate termination of the local model.
91!-- The rest of the coupler must then be skipped because it would cause an MPI
92!-- intercomminucation hang.
93!-- If necessary, the coupler will be called at the beginning of the next
94!-- restart run.
[667]95
96    IF ( coupling_topology == 0 ) THEN
[709]97       CALL MPI_SENDRECV( terminate_coupled,        1, MPI_INTEGER, target_id, &
98                          0,                                                   &
99                          terminate_coupled_remote, 1, MPI_INTEGER, target_id, &
[667]100                          0, comm_inter, status, ierr )
101    ELSE
102       IF ( myid == 0) THEN
103          CALL MPI_SENDRECV( terminate_coupled,        1, MPI_INTEGER, &
104                             target_id, 0,                             &
105                             terminate_coupled_remote, 1, MPI_INTEGER, & 
106                             target_id, 0,                             &
107                             comm_inter, status, ierr )
108       ENDIF
[709]109       CALL MPI_BCAST( terminate_coupled_remote, 1, MPI_INTEGER, 0, comm2d, &
110                       ierr )
[667]111
112       ALLOCATE( total_2d_a(-nbgp:ny_a+nbgp,-nbgp:nx_a+nbgp),       &
113                 total_2d_o(-nbgp:ny_o+nbgp,-nbgp:nx_o+nbgp) )
114
115    ENDIF
116
[108]117    IF ( terminate_coupled_remote > 0 )  THEN
[274]118       WRITE( message_string, * ) 'remote model "',                         &
119                                  TRIM( coupling_mode_remote ),             &
120                                  '" terminated',                           &
121                                  '&with terminate_coupled_remote = ',      &
122                                  terminate_coupled_remote,                 &
123                                  '&local model  "', TRIM( coupling_mode ), &
124                                  '" has',                                  &
125                                  '&terminate_coupled = ',                  &
[667]126                                   terminate_coupled
[258]127       CALL message( 'surface_coupler', 'PA0310', 1, 2, 0, 6, 0 )
[108]128       RETURN
129    ENDIF
[667]130 
[291]131
[108]132!
133!-- Exchange the current simulated time between the models,
[667]134!-- currently just for total_2ding
[709]135    IF ( coupling_topology == 0 ) THEN
136   
137       CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, 11, &
138                      comm_inter, ierr )
139       CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL, target_id, &
140                      11, comm_inter, status, ierr )
[667]141    ELSE
[709]142
[667]143       IF ( myid == 0 ) THEN
[709]144
145          CALL MPI_SEND( time_since_reference_point, 1, MPI_REAL, target_id, &
146                         11, comm_inter, ierr )
147          CALL MPI_RECV( time_since_reference_point_rem, 1, MPI_REAL,        &
[667]148                         target_id, 11, comm_inter, status, ierr )
[709]149
[667]150       ENDIF
[709]151
152       CALL MPI_BCAST( time_since_reference_point_rem, 1, MPI_REAL, 0, comm2d, &
153                       ierr )
154
[667]155    ENDIF
[102]156
157!
158!-- Exchange the interface data
159    IF ( coupling_mode == 'atmosphere_to_ocean' )  THEN
[667]160   
161!
[709]162!--    Horizontal grid size and number of processors is equal in ocean and
163!--    atmosphere
164       IF ( coupling_topology == 0 )  THEN
[102]165
166!
[709]167!--       Send heat flux at bottom surface to the ocean
168          CALL MPI_SEND( shf(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 12, &
169                         comm_inter, ierr )
[102]170!
[709]171!--       Send humidity flux at bottom surface to the ocean
[667]172          IF ( humidity )  THEN
[709]173             CALL MPI_SEND( qsws(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 13, &
174                            comm_inter, ierr )
[667]175          ENDIF
176!
[709]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 )
[108]180!
[709]181!--       Send the momentum flux (u) at bottom surface to the ocean
182          CALL MPI_SEND( usws(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 15, &
183                         comm_inter, ierr )
[102]184!
[709]185!--       Send the momentum flux (v) at bottom surface to the ocean
186          CALL MPI_SEND( vsws(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 16, &
187                         comm_inter, ierr )
[102]188!
[709]189!--       Receive u at the bottom surface from the ocean
190          CALL MPI_RECV( u(0,nysg,nxlg), 1, type_xy, target_id, 17, &
191                         comm_inter, status, ierr )
[667]192!
[709]193!--       Receive v at the bottom surface from the ocean
194          CALL MPI_RECV( v(0,nysg,nxlg), 1, type_xy, target_id, 18, &
195                         comm_inter, status, ierr )
[667]196!
197!--    Horizontal grid size or number of processors differs between
198!--    ocean and atmosphere
199       ELSE
200     
201!
[709]202!--       Send heat flux at bottom surface to the ocean
[667]203          total_2d_a = 0.0
[709]204          total_2d   = 0.0
[667]205          total_2d(nys:nyn,nxl:nxr) = shf(nys:nyn,nxl:nxr)
[709]206
207          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
208                           comm2d, ierr )
209          CALL interpolate_to_ocean( 12 )   
[667]210!
[709]211!--       Send humidity flux at bottom surface to the ocean
212          IF ( humidity )  THEN
[667]213             total_2d_a = 0.0
[709]214             total_2d   = 0.0
[667]215             total_2d(nys:nyn,nxl:nxr) = qsws(nys:nyn,nxl:nxr)
[709]216
217             CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, &
218                              0, comm2d, ierr )
219             CALL interpolate_to_ocean( 13 )
[667]220          ENDIF
221!
[709]222!--       Receive temperature at the bottom surface from the ocean
223          IF ( myid == 0 )  THEN
[667]224             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
225                            target_id, 14, comm_inter, status, ierr )   
226          ENDIF
227          CALL MPI_BARRIER( comm2d, ierr )
[709]228          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
229                          ierr )
[667]230          pt(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
231!
[709]232!--       Send momentum flux (u) at bottom surface to the ocean
[667]233          total_2d_a = 0.0 
[709]234          total_2d   = 0.0
[667]235          total_2d(nys:nyn,nxl:nxr) = usws(nys:nyn,nxl:nxr)
[709]236          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
237                           comm2d, ierr )
238          CALL interpolate_to_ocean( 15 )
[667]239!
[709]240!--       Send momentum flux (v) at bottom surface to the ocean
[667]241          total_2d_a = 0.0
[709]242          total_2d   = 0.0
[667]243          total_2d(nys:nyn,nxl:nxr) = vsws(nys:nyn,nxl:nxr)
[709]244          CALL MPI_REDUCE( total_2d, total_2d_a, ngp_a, MPI_REAL, MPI_SUM, 0, &
245                           comm2d, ierr )
246          CALL interpolate_to_ocean( 16 )
[667]247!
[709]248!--       Receive u at the bottom surface from the ocean
249          IF ( myid == 0 )  THEN
[667]250             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]251                            target_id, 17, comm_inter, status, ierr )
[667]252          ENDIF
253          CALL MPI_BARRIER( comm2d, ierr )
[709]254          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
255                          ierr )
[667]256          u(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
257!
[709]258!--       Receive v at the bottom surface from the ocean
259          IF ( myid == 0 )  THEN
[667]260             CALL MPI_RECV( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, &
[709]261                            target_id, 18, comm_inter, status, ierr )
[667]262          ENDIF
263          CALL MPI_BARRIER( comm2d, ierr )
[709]264          CALL MPI_BCAST( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, 0, comm2d, &
265                          ierr )
[667]266          v(0,nysg:nyng,nxlg:nxrg) = total_2d_a(nysg:nyng,nxlg:nxrg)
267
268       ENDIF
269
[102]270    ELSEIF ( coupling_mode == 'ocean_to_atmosphere' )  THEN
271
272!
[667]273!--    Horizontal grid size and number of processors is equal
274!--    in ocean and atmosphere
275       IF ( coupling_topology == 0 ) THEN
276!
[709]277!--       Receive heat flux at the sea surface (top) from the atmosphere
278          CALL MPI_RECV( tswst(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 12, &
279                         comm_inter, status, ierr )
[102]280!
[709]281!--       Receive humidity flux from the atmosphere (bottom)
[667]282!--       and add it to the heat flux at the sea surface (top)...
283          IF ( humidity_remote )  THEN
284             CALL MPI_RECV( qswst_remote(nysg,nxlg), ngp_xy, MPI_REAL, &
285                            target_id, 13, comm_inter, status, ierr )
286          ENDIF
287!
288!--       Send sea surface temperature to the atmosphere model
[709]289          CALL MPI_SEND( pt(nzt,nysg,nxlg), 1, type_xy, target_id, 14, &
290                         comm_inter, ierr )
[667]291!
292!--       Receive momentum flux (u) at the sea surface (top) from the atmosphere
[709]293          CALL MPI_RECV( uswst(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 15, &
294                         comm_inter, status, ierr )
[667]295!
296!--       Receive momentum flux (v) at the sea surface (top) from the atmosphere
[709]297          CALL MPI_RECV( vswst(nysg,nxlg), ngp_xy, MPI_REAL, target_id, 16, &
298                         comm_inter, status, ierr )
[667]299!
[709]300!--       Send u to the atmosphere
301          CALL MPI_SEND( u(nzt,nysg,nxlg), 1, type_xy, target_id, 17, &
302                         comm_inter, ierr )
[667]303!
[709]304!--       Send v to the atmosphere
305          CALL MPI_SEND( v(nzt,nysg,nxlg), 1, type_xy, target_id, 18, &
306                         comm_inter, ierr )
307!
[667]308!--    Horizontal gridsize or number of processors differs between
309!--    ocean and atmosphere
310       ELSE
311!
[709]312!--       Receive heat flux at the sea surface (top) from the atmosphere
313          IF ( myid == 0 )  THEN
[667]314             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]315                            target_id, 12, comm_inter, status, ierr )
[667]316          ENDIF
317          CALL MPI_BARRIER( comm2d, ierr )
[709]318          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
319                          ierr )
[667]320          tswst(nysg:nyng,nxlg:nxrg) = total_2d_o(nysg:nyng,nxlg:nxrg)
321!
[709]322!--       Receive humidity flux at the sea surface (top) from the atmosphere
323          IF ( humidity_remote )  THEN
324             IF ( myid == 0 )  THEN
[667]325                CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]326                               target_id, 13, comm_inter, status, ierr )
[667]327             ENDIF
328             CALL MPI_BARRIER( comm2d, ierr )
[709]329             CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, &
330                             comm2d, ierr)
[667]331             qswst_remote(nysg:nyng,nxlg:nxrg) = total_2d_o(nysg:nyng,nxlg:nxrg)
332          ENDIF
333!
334!--       Send surface temperature to atmosphere
335          total_2d_o = 0.0
[709]336          total_2d   = 0.0
[667]337          total_2d(nys:nyn,nxl:nxr) = pt(nzt,nys:nyn,nxl:nxr)
338
[709]339          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
340                           comm2d, ierr) 
341          CALL interpolate_to_atmos( 14 )
[667]342!
[709]343!--       Receive momentum flux (u) at the sea surface (top) from the atmosphere
344          IF ( myid == 0 )  THEN
[667]345             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]346                            target_id, 15, comm_inter, status, ierr )
[667]347          ENDIF
348          CALL MPI_BARRIER( comm2d, ierr )
349          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]350                          0, comm2d, ierr )
[667]351          uswst(nysg:nyng,nxlg:nxrg) = total_2d_o(nysg:nyng,nxlg:nxrg)
352!
[709]353!--       Receive momentum flux (v) at the sea surface (top) from the atmosphere
354          IF ( myid == 0 )  THEN
[667]355             CALL MPI_RECV( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
[709]356                            target_id, 16, comm_inter, status, ierr )
[667]357          ENDIF
358          CALL MPI_BARRIER( comm2d, ierr )
[709]359          CALL MPI_BCAST( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, 0, comm2d, &
360                          ierr )
[667]361          vswst(nysg:nyng,nxlg:nxrg) = total_2d_o(nysg:nyng,nxlg:nxrg)
362!
363!--       Send u to atmosphere
364          total_2d_o = 0.0 
[709]365          total_2d   = 0.0
[667]366          total_2d(nys:nyn,nxl:nxr) = u(nzt,nys:nyn,nxl:nxr)
[709]367          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
368                           comm2d, ierr )
369          CALL interpolate_to_atmos( 17 )
[667]370!
371!--       Send v to atmosphere
372          total_2d_o = 0.0
[709]373          total_2d   = 0.0
[667]374          total_2d(nys:nyn,nxl:nxr) = v(nzt,nys:nyn,nxl:nxr)
[709]375          CALL MPI_REDUCE( total_2d, total_2d_o, ngp_o, MPI_REAL, MPI_SUM, 0, &
376                           comm2d, ierr )
377          CALL interpolate_to_atmos( 18 )
[667]378       
379       ENDIF
380
381!
382!--    Conversions of fluxes received from atmosphere
383       IF ( humidity_remote )  THEN
[108]384!
[709]385!--       Here tswst is still the sum of atmospheric bottom heat fluxes,
386!--       * latent heat of vaporization in m2/s2, or 540 cal/g, or 40.65 kJ/mol
387!--       /(rho_atm(=1.0)*c_p)
388          tswst = tswst + qswst_remote * 2.2626108E6 / 1005.0
389!
[667]390!--        ...and convert it to a salinity flux at the sea surface (top)
[108]391!--       following Steinhorn (1991), JPO 21, pp. 1681-1683:
392!--       S'w' = -S * evaporation / ( rho_water * ( 1 - S ) )
393          saswst = -1.0 * sa(nzt,:,:) * qswst_remote /  &
[667]394                    ( rho(nzt,:,:) * ( 1.0 - sa(nzt,:,:) ) )
[108]395       ENDIF
396
397!
[102]398!--    Adjust the kinematic heat flux with respect to ocean density
399!--    (constants are the specific heat capacities for air and water)
[667]400!--    now tswst is the ocean top heat flux
[108]401       tswst = tswst / rho(nzt,:,:) * 1005.0 / 4218.0
[102]402
403!
[667]404!--    Adjust the momentum fluxes with respect to ocean density
405       uswst = uswst / rho(nzt,:,:)
406       vswst = vswst / rho(nzt,:,:)
[102]407
[667]408    ENDIF
409
[709]410    IF ( coupling_topology == 1 )  THEN
[667]411       DEALLOCATE( total_2d_o, total_2d_a )
412    ENDIF
413
414    CALL cpu_log( log_point(39), 'surface_coupler', 'stop' )
415
416#endif
417
418  END SUBROUTINE surface_coupler
419
420
421
[709]422  SUBROUTINE interpolate_to_atmos( tag )
[667]423
[880]424#if defined( __parallel )
425
[1320]426    USE arrays_3d,                                                             &
427        ONLY:  total_2d_a, total_2d_o
[667]428
[1320]429    USE indices,                                                               &
430        ONLY:  nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
431
432    USE kinds
433
434    USE pegrid,                                                                &
435        ONLY:  comm2d, comm_inter, ierr, MPI_DOUBLE_PRECISION, myid, ngp_a,    &
436               target_id
437
[667]438    IMPLICIT NONE
439
[1320]440    INTEGER(iwp) ::  dnx  !:
441    INTEGER(iwp) ::  dnx2 !:
442    INTEGER(iwp) ::  dny  !:
443    INTEGER(iwp) ::  dny2 !:
444    INTEGER(iwp) ::  i    !:
445    INTEGER(iwp) ::  ii   !:
446    INTEGER(iwp) ::  j    !:
447    INTEGER(iwp) ::  jj   !:
[667]448
[1320]449    INTEGER(iwp), intent(in) ::  tag !:
450
[667]451    CALL MPI_BARRIER( comm2d, ierr )
452
[709]453    IF ( myid == 0 )  THEN
454!
455!--    Cyclic boundary conditions for the total 2D-grid
[667]456       total_2d_o(-nbgp:-1,:) = total_2d_o(ny+1-nbgp:ny,:)
457       total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx+1-nbgp:nx)
458
459       total_2d_o(ny+1:ny+nbgp,:) = total_2d_o(0:nbgp-1,:)
460       total_2d_o(:,nx+1:nx+nbgp) = total_2d_o(:,0:nbgp-1)
461
[102]462!
[667]463!--    Number of gridpoints of the fine grid within one mesh of the coarse grid
464       dnx = (nx_o+1) / (nx_a+1) 
465       dny = (ny_o+1) / (ny_a+1) 
[102]466
467!
[709]468!--    Distance for interpolation around coarse grid points within the fine
469!--    grid (note: 2*dnx2 must not be equal with dnx)
[667]470       dnx2 = 2 * ( dnx / 2 )
471       dny2 = 2 * ( dny / 2 )
[102]472
[667]473       total_2d_a = 0.0
[102]474!
[667]475!--    Interpolation from ocean-grid-layer to atmosphere-grid-layer
476       DO  j = 0, ny_a
477          DO  i = 0, nx_a 
478             DO  jj = 0, dny2
479                DO  ii = 0, dnx2
480                   total_2d_a(j,i) = total_2d_a(j,i) &
481                                     + total_2d_o(j*dny+jj,i*dnx+ii)
482                ENDDO
483             ENDDO
484             total_2d_a(j,i) = total_2d_a(j,i) / ( ( dnx2 + 1 ) * ( dny2 + 1 ) )
485          ENDDO
486       ENDDO
487!
[709]488!--    Cyclic boundary conditions for atmosphere grid
[667]489       total_2d_a(-nbgp:-1,:) = total_2d_a(ny_a+1-nbgp:ny_a,:)
490       total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx_a+1-nbgp:nx_a)
491       
492       total_2d_a(ny_a+1:ny_a+nbgp,:) = total_2d_a(0:nbgp-1,:)
493       total_2d_a(:,nx_a+1:nx_a+nbgp) = total_2d_a(:,0:nbgp-1)
494!
495!--    Transfer of the atmosphere-grid-layer to the atmosphere
[709]496       CALL MPI_SEND( total_2d_a(-nbgp,-nbgp), ngp_a, MPI_REAL, target_id, &
497                      tag, comm_inter, ierr )
[102]498
499    ENDIF
500
[667]501    CALL MPI_BARRIER( comm2d, ierr )
[102]502
[880]503#endif
504
[667]505  END SUBROUTINE interpolate_to_atmos
[102]506
[667]507
[709]508  SUBROUTINE interpolate_to_ocean( tag )
[667]509
[880]510#if defined( __parallel )
511
[1320]512    USE arrays_3d,                                                             &
513        ONLY:  total_2d_a, total_2d_o
[667]514
[1320]515    USE indices,                                                               &
516        ONLY:  nbgp, nx, nx_a, nx_o, ny, ny_a, ny_o
517
518    USE kinds
519
520    USE pegrid,                                                                &
521        ONLY:  comm2d, comm_inter, ierr, MPI_DOUBLE_PRECISION, myid, ngp_o,    &
522               target_id
523
[667]524    IMPLICIT NONE
525
[1320]526    INTEGER(iwp)             ::  dnx !:
527    INTEGER(iwp)             ::  dny !:
528    INTEGER(iwp)             ::  i   !:
529    INTEGER(iwp)             ::  ii  !:
530    INTEGER(iwp)             ::  j   !:
531    INTEGER(iwp)             ::  jj  !:
532    INTEGER(iwp), intent(in) ::  tag !:
[667]533
[1320]534    REAL(wp)                 ::  fl  !:
535    REAL(wp)                 ::  fr  !:
536    REAL(wp)                 ::  myl !:
537    REAL(wp)                 ::  myr !:
[709]538
[667]539    CALL MPI_BARRIER( comm2d, ierr )
540
[709]541    IF ( myid == 0 )  THEN   
[667]542
543!
[709]544!--    Number of gridpoints of the fine grid within one mesh of the coarse grid
[667]545       dnx = ( nx_o + 1 ) / ( nx_a + 1 ) 
546       dny = ( ny_o + 1 ) / ( ny_a + 1 ) 
547
548!
[709]549!--    Cyclic boundary conditions for atmosphere grid
[667]550       total_2d_a(-nbgp:-1,:) = total_2d_a(ny+1-nbgp:ny,:)
551       total_2d_a(:,-nbgp:-1) = total_2d_a(:,nx+1-nbgp:nx)
552       
553       total_2d_a(ny+1:ny+nbgp,:) = total_2d_a(0:nbgp-1,:)
554       total_2d_a(:,nx+1:nx+nbgp) = total_2d_a(:,0:nbgp-1)
555!
[709]556!--    Bilinear Interpolation from atmosphere grid-layer to ocean grid-layer
[667]557       DO  j = 0, ny
558          DO  i = 0, nx
559             myl = ( total_2d_a(j+1,i)   - total_2d_a(j,i)   ) / dny
560             myr = ( total_2d_a(j+1,i+1) - total_2d_a(j,i+1) ) / dny
561             DO  jj = 0, dny-1
[709]562                fl = myl*jj + total_2d_a(j,i) 
563                fr = myr*jj + total_2d_a(j,i+1) 
[667]564                DO  ii = 0, dnx-1
565                   total_2d_o(j*dny+jj,i*dnx+ii) = ( fr - fl ) / dnx * ii + fl
566                ENDDO
567             ENDDO
568          ENDDO
569       ENDDO
570!
[709]571!--    Cyclic boundary conditions for ocean grid
[667]572       total_2d_o(-nbgp:-1,:) = total_2d_o(ny_o+1-nbgp:ny_o,:)
573       total_2d_o(:,-nbgp:-1) = total_2d_o(:,nx_o+1-nbgp:nx_o)
574
575       total_2d_o(ny_o+1:ny_o+nbgp,:) = total_2d_o(0:nbgp-1,:)
576       total_2d_o(:,nx_o+1:nx_o+nbgp) = total_2d_o(:,0:nbgp-1)
577
578       CALL MPI_SEND( total_2d_o(-nbgp,-nbgp), ngp_o, MPI_REAL, &
579                      target_id, tag, comm_inter, ierr )
580
581    ENDIF
582
583    CALL MPI_BARRIER( comm2d, ierr ) 
584
[880]585#endif
586
[667]587  END SUBROUTINE interpolate_to_ocean
Note: See TracBrowser for help on using the repository browser.