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

Last change on this file since 1353 was 1353, checked in by heinze, 7 years ago

REAL constants provided with KIND-attribute

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