New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
icblbc.F90 in NEMO/trunk/src/OCE/ICB – NEMO

source: NEMO/trunk/src/OCE/ICB/icblbc.F90 @ 14433

Last change on this file since 14433 was 14433, checked in by smasson, 3 years ago

trunk: merge dev_r14312_MPI_Interface into the trunk, #2598

  • Property svn:keywords set to Id
File size: 35.8 KB
Line 
1MODULE icblbc
2   !!======================================================================
3   !!                       ***  MODULE  icblbc  ***
4   !! Ocean physics:  routines to handle boundary exchanges for icebergs
5   !!======================================================================
6   !! History :  3.3  !  2010-01  (Martin&Adcroft) Original code
7   !!             -   !  2011-03  (Madec)          Part conversion to NEMO form
8   !!             -   !                            Removal of mapping from another grid
9   !!             -   !  2011-04  (Alderson)       Split into separate modules
10   !!             -   !  2011-05  (Alderson)       MPP exchanges written based on lib_mpp
11   !!             -   !  2011-05  (Alderson)       MPP and single processor boundary conditions added
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   icb_lbc       : -  Pass icebergs across cyclic boundaries
16   !!   icb_lbc_mpp   : -  In MPP pass icebergs from linked list between processors
17   !!                      as they advect around
18   !!                   -  Lagrangian processes cannot be handled by existing NEMO MPP
19   !!                      routines because they do not lie on regular jpi,jpj grids
20   !!                   -  Processor exchanges are handled as in lib_mpp whenever icebergs step
21   !!                      across boundary of interior domain (nicbdi-nicbei, nicbdj-nicbej)
22   !!                      so that iceberg does not exist in more than one processor
23   !!                   -  North fold exchanges controlled by three arrays:
24   !!                         nicbflddest - unique processor numbers that current one exchanges with
25   !!                         nicbfldproc - processor number that current grid point exchanges with
26   !!                         nicbfldpts  - packed i,j point in exchanging processor
27   !!----------------------------------------------------------------------
28   USE par_oce                             ! ocean parameters
29   USE dom_oce                             ! ocean domain
30   USE in_out_manager                      ! IO parameters
31   USE lib_mpp                             ! MPI code and lk_mpp in particular
32   USE icb_oce                             ! define iceberg arrays
33   USE icbutl                              ! iceberg utility routines
34
35   IMPLICIT NONE
36   PRIVATE
37
38#if ! defined key_mpi_off
39
40!$AGRIF_DO_NOT_TREAT
41   INCLUDE 'mpif.h'
42!$AGRIF_END_DO_NOT_TREAT
43
44   TYPE, PUBLIC :: buffer
45      INTEGER :: size = 0
46      REAL(wp), DIMENSION(:,:), POINTER ::   data
47   END TYPE buffer
48
49   TYPE(buffer), POINTER       ::   obuffer_n=>NULL() , ibuffer_n=>NULL()
50   TYPE(buffer), POINTER       ::   obuffer_s=>NULL() , ibuffer_s=>NULL()
51   TYPE(buffer), POINTER       ::   obuffer_e=>NULL() , ibuffer_e=>NULL()
52   TYPE(buffer), POINTER       ::   obuffer_w=>NULL() , ibuffer_w=>NULL()
53
54   ! north fold exchange buffers
55   TYPE(buffer), POINTER       ::   obuffer_f=>NULL() , ibuffer_f=>NULL()
56
57   INTEGER, PARAMETER, PRIVATE ::   jp_delta_buf = 25             ! Size by which to increment buffers
58   INTEGER, PARAMETER, PRIVATE ::   jp_buffer_width = 15+nkounts  ! items to store for each berg
59
60#endif
61
62   PUBLIC   icb_lbc
63   PUBLIC   icb_lbc_mpp
64
65   !! * Substitutions
66#  include "do_loop_substitute.h90"
67   !!----------------------------------------------------------------------
68   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
69   !! $Id$
70   !! Software governed by the CeCILL license (see ./LICENSE)
71   !!----------------------------------------------------------------------
72CONTAINS
73
74   SUBROUTINE icb_lbc()
75      !!----------------------------------------------------------------------
76      !!                 ***  SUBROUTINE icb_lbc  ***
77      !!
78      !! ** Purpose :   in non-mpp case need to deal with cyclic conditions
79      !!                including north-fold
80      !!----------------------------------------------------------------------
81      TYPE(iceberg), POINTER ::   this
82      TYPE(point)  , POINTER ::   pt
83      !!----------------------------------------------------------------------
84
85      !! periodic east/west boundaries
86      !! =============================
87
88      IF( l_Iperio ) THEN
89
90         this => first_berg
91         DO WHILE( ASSOCIATED(this) )
92            pt => this%current_point
93            IF( pt%xi > REAL(mig(nicbei),wp) + 0.5_wp ) THEN
94               pt%xi = ricb_right + MOD(pt%xi, 1._wp ) - 1._wp
95            ELSE IF( pt%xi < REAL(mig(nicbdi),wp) - 0.5_wp ) THEN
96               pt%xi = ricb_left + MOD(pt%xi, 1._wp )
97            ENDIF
98            this => this%next
99         END DO
100         !
101      ENDIF
102
103      !! north/south boundaries
104      !! ======================
105      IF( l_Jperio)      CALL ctl_stop(' north-south periodicity not implemented for icebergs')
106      ! north fold
107      IF( l_IdoNFold )   CALL icb_lbc_nfld()
108      !
109   END SUBROUTINE icb_lbc
110
111
112   SUBROUTINE icb_lbc_nfld()
113      !!----------------------------------------------------------------------
114      !!                 ***  SUBROUTINE icb_lbc_nfld  ***
115      !!
116      !! ** Purpose :   single processor north fold exchange
117      !!----------------------------------------------------------------------
118      TYPE(iceberg), POINTER ::   this
119      TYPE(point)  , POINTER ::   pt
120      INTEGER                ::   iine, ijne, ipts
121      INTEGER                ::   iiglo, ijglo
122      !!----------------------------------------------------------------------
123      !
124      this => first_berg
125      DO WHILE( ASSOCIATED(this) )
126         pt => this%current_point
127         ijne = INT( pt%yj + 0.5 )
128         IF( pt%yj > REAL(mjg(nicbej),wp) + 0.5_wp ) THEN
129            !
130            iine = INT( pt%xi + 0.5 )
131            ipts  = nicbfldpts (mi1(iine))
132            !
133            ! moving across the cut line means both position and
134            ! velocity must change
135            ijglo = INT( ipts/nicbpack )
136            iiglo = ipts - nicbpack*ijglo
137            pt%xi = iiglo - ( pt%xi - REAL(iine,wp) )
138            pt%yj = ijglo - ( pt%yj - REAL(ijne,wp) )
139            pt%uvel = -1._wp * pt%uvel
140            pt%vvel = -1._wp * pt%vvel
141         ENDIF
142         this => this%next
143      END DO
144      !
145   END SUBROUTINE icb_lbc_nfld
146
147#if ! defined key_mpi_off
148   !!----------------------------------------------------------------------
149   !!            MPI massively parallel processing library
150   !!----------------------------------------------------------------------
151
152   SUBROUTINE icb_lbc_mpp()
153      !!----------------------------------------------------------------------
154      !!                 ***  SUBROUTINE icb_lbc_mpp  ***
155      !!
156      !! ** Purpose :   multi processor exchange
157      !!
158      !! ** Method  :   identify direction for exchange, pack into a buffer
159      !!                which is basically a real array and delete from linked list
160      !!                length of buffer is exchanged first with receiving processor
161      !!                then buffer is sent if necessary
162      !!----------------------------------------------------------------------
163      TYPE(iceberg)         , POINTER     ::   tmpberg, this
164      TYPE(point)           , POINTER     ::   pt
165      INTEGER                             ::   ibergs_to_send_e, ibergs_to_send_w
166      INTEGER                             ::   ibergs_to_send_n, ibergs_to_send_s
167      INTEGER                             ::   ibergs_rcvd_from_e, ibergs_rcvd_from_w
168      INTEGER                             ::   ibergs_rcvd_from_n, ibergs_rcvd_from_s
169      INTEGER                             ::   i, ibergs_start, ibergs_end
170      INTEGER                             ::   ipe_N, ipe_S, ipe_W, ipe_E
171      REAL(wp), DIMENSION(2)              ::   zewbergs, zwebergs, znsbergs, zsnbergs
172      INTEGER                             ::   iml_req1, iml_req2, iml_req3, iml_req4
173      INTEGER                             ::   iml_req5, iml_req6, iml_req7, iml_req8, iml_err
174      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   iml_stat
175
176      ! set up indices of neighbouring processors
177      ipe_N = -1
178      ipe_S = -1
179      ipe_W = -1
180      ipe_E = -1
181      IF( mpinei(jpwe) >= 0 ) ipe_W = mpinei(jpwe)
182      IF( mpinei(jpea) >= 0 ) ipe_E = mpinei(jpea)
183      IF( mpinei(jpso) >= 0 ) ipe_S = mpinei(jpso)
184      IF( mpinei(jpno) >= 0 ) ipe_N = mpinei(jpno)
185      !
186      ! at northern line of processors with north fold handle bergs differently
187      IF( l_IdoNFold )   ipe_N = -1
188
189      ! if there's only one processor in x direction then don't let mpp try to handle periodicity
190      IF( jpni == 1 ) THEN
191         ipe_E = -1
192         ipe_W = -1
193      ENDIF
194
195      IF( nn_verbose_level >= 2 ) THEN
196         WRITE(numicb,*) 'processor west  : ', ipe_W
197         WRITE(numicb,*) 'processor east  : ', ipe_E
198         WRITE(numicb,*) 'processor north : ', ipe_N
199         WRITE(numicb,*) 'processor south : ', ipe_S
200         WRITE(numicb,*) 'processor nimpp : ', nimpp
201         WRITE(numicb,*) 'processor njmpp : ', njmpp
202         CALL flush( numicb )
203      ENDIF
204
205      ! periodicity is handled here when using mpp when there is more than one processor in
206      ! the i direction, but it also has to happen when jpni=1 case so this is dealt with
207      ! in icb_lbc and called here
208
209      IF( jpni == 1 ) CALL icb_lbc()
210
211      ! Note that xi is adjusted when swapping because of periodic condition
212
213      IF( nn_verbose_level > 0 ) THEN
214         ! store the number of icebergs on this processor at start
215         ibergs_start = icb_utl_count()
216      ENDIF
217
218      ibergs_to_send_e   = 0
219      ibergs_to_send_w   = 0
220      ibergs_to_send_n   = 0
221      ibergs_to_send_s   = 0
222      ibergs_rcvd_from_e = 0
223      ibergs_rcvd_from_w = 0
224      ibergs_rcvd_from_n = 0
225      ibergs_rcvd_from_s = 0
226
227      IF( ASSOCIATED(first_berg) ) THEN      ! Find number of bergs that headed east/west
228         this => first_berg
229         DO WHILE (ASSOCIATED(this))
230            pt => this%current_point
231            IF( ipe_E >= 0 .AND. pt%xi > REAL(mig(nicbei),wp) + 0.5_wp ) THEN
232               tmpberg => this
233               this => this%next
234               ibergs_to_send_e = ibergs_to_send_e + 1
235               IF( nn_verbose_level >= 4 ) THEN
236                  WRITE(numicb,*) 'bergstep ',nktberg,' packing berg ',tmpberg%number(:),' for transfer to east'
237                  CALL flush( numicb )
238               ENDIF
239               ! deal with periodic case
240               tmpberg%current_point%xi = ricb_right + MOD(tmpberg%current_point%xi, 1._wp ) - 1._wp
241               ! now pack it into buffer and delete from list
242               CALL icb_pack_into_buffer( tmpberg, obuffer_e, ibergs_to_send_e)
243               CALL icb_utl_delete(first_berg, tmpberg)
244            ELSE IF( ipe_W >= 0 .AND. pt%xi < REAL(mig(nicbdi),wp) - 0.5_wp ) THEN
245               tmpberg => this
246               this => this%next
247               ibergs_to_send_w = ibergs_to_send_w + 1
248               IF( nn_verbose_level >= 4 ) THEN
249                  WRITE(numicb,*) 'bergstep ',nktberg,' packing berg ',tmpberg%number(:),' for transfer to west'
250                  CALL flush( numicb )
251               ENDIF
252               ! deal with periodic case
253               tmpberg%current_point%xi = ricb_left + MOD(tmpberg%current_point%xi, 1._wp )
254               ! now pack it into buffer and delete from list
255               CALL icb_pack_into_buffer( tmpberg, obuffer_w, ibergs_to_send_w)
256               CALL icb_utl_delete(first_berg, tmpberg)
257            ELSE
258               this => this%next
259            ENDIF
260         END DO
261      ENDIF
262      IF( nn_verbose_level >= 3) THEN
263         WRITE(numicb,*) 'bergstep ',nktberg,' send ew: ', ibergs_to_send_e, ibergs_to_send_w
264         CALL flush(numicb)
265      ENDIF
266
267      ! send bergs east and receive bergs from west (ie ones that were sent east) and vice versa
268
269      ! pattern here is copied from lib_mpp code
270
271      IF( mpinei(jpwe) >= 0  )   zewbergs(1) = ibergs_to_send_w
272      IF( mpinei(jpea) >= 0  )   zwebergs(1) = ibergs_to_send_e
273      IF( mpinei(jpwe) >= 0  )   CALL mppsend( 11, zewbergs(1), 1, ipe_W, iml_req2)
274      IF( mpinei(jpea) >= 0  )   CALL mppsend( 12, zwebergs(1), 1, ipe_E, iml_req3)
275      IF( mpinei(jpea) >= 0  )   CALL mpprecv( 11, zewbergs(2), 1, ipe_E )
276      IF( mpinei(jpwe) >= 0  )   CALL mpprecv( 12, zwebergs(2), 1, ipe_W )
277      IF( mpinei(jpwe) >= 0  )   CALL mpi_wait( iml_req2, iml_stat, iml_err )
278      IF( mpinei(jpea) >= 0  )   CALL mpi_wait( iml_req3, iml_stat, iml_err )
279      IF( mpinei(jpea) >= 0  )   ibergs_rcvd_from_e = INT( zewbergs(2) )
280      IF( mpinei(jpwe) >= 0  )   ibergs_rcvd_from_w = INT( zwebergs(2) )
281     
282      IF( nn_verbose_level >= 3) THEN
283         WRITE(numicb,*) 'bergstep ',nktberg,' recv ew: ', ibergs_rcvd_from_w, ibergs_rcvd_from_e
284         CALL flush(numicb)
285      ENDIF
286     
287      IF( ibergs_to_send_w > 0 ) CALL mppsend( 13, obuffer_w%data, ibergs_to_send_w*jp_buffer_width, ipe_W, iml_req2 )
288      IF( ibergs_to_send_e > 0 ) CALL mppsend( 14, obuffer_e%data, ibergs_to_send_e*jp_buffer_width, ipe_E, iml_req3 )
289      IF( ibergs_rcvd_from_e > 0 ) THEN
290         CALL icb_increase_ibuffer(ibuffer_e, ibergs_rcvd_from_e)
291         CALL mpprecv( 13, ibuffer_e%data, ibergs_rcvd_from_e*jp_buffer_width )
292      ENDIF
293      IF( ibergs_rcvd_from_w > 0 ) THEN
294         CALL icb_increase_ibuffer(ibuffer_w, ibergs_rcvd_from_w)
295         CALL mpprecv( 14, ibuffer_w%data, ibergs_rcvd_from_w*jp_buffer_width )
296      ENDIF
297      IF( ibergs_to_send_w > 0 ) CALL mpi_wait( iml_req2, iml_stat, iml_err )
298      IF( ibergs_to_send_e > 0 ) CALL mpi_wait( iml_req3, iml_stat, iml_err )
299      DO i = 1, ibergs_rcvd_from_e
300         IF( nn_verbose_level >= 4 ) THEN
301            WRITE(numicb,*) 'bergstep ',nktberg,' unpacking berg ',INT(ibuffer_e%data(16,i)),' from east'
302            CALL FLUSH( numicb )
303         ENDIF
304         CALL icb_unpack_from_buffer(first_berg, ibuffer_e, i)
305      END DO
306      DO i = 1, ibergs_rcvd_from_w
307         IF( nn_verbose_level >= 4 ) THEN
308            WRITE(numicb,*) 'bergstep ',nktberg,' unpacking berg ',INT(ibuffer_w%data(16,i)),' from west'
309            CALL FLUSH( numicb )
310         ENDIF
311         CALL icb_unpack_from_buffer(first_berg, ibuffer_w, i)
312      END DO
313
314      ! Find number of bergs that headed north/south
315      ! (note: this block should technically go ahead of the E/W recv block above
316      !  to handle arbitrary orientation of PEs. But for simplicity, it is
317      !  here to accomodate diagonal transfer of bergs between PEs -AJA)
318
319      IF( ASSOCIATED(first_berg) ) THEN
320         this => first_berg
321         DO WHILE (ASSOCIATED(this))
322            pt => this%current_point
323            IF( ipe_N >= 0 .AND. pt%yj > REAL(mjg(nicbej),wp) + 0.5_wp ) THEN
324               tmpberg => this
325               this => this%next
326               ibergs_to_send_n = ibergs_to_send_n + 1
327               IF( nn_verbose_level >= 4 ) THEN
328                  WRITE(numicb,*) 'bergstep ',nktberg,' packing berg ',tmpberg%number(:),' for transfer to north'
329                  CALL flush( numicb )
330               ENDIF
331               CALL icb_pack_into_buffer( tmpberg, obuffer_n, ibergs_to_send_n)
332               CALL icb_utl_delete(first_berg, tmpberg)
333            ELSE IF( ipe_S >= 0 .AND. pt%yj < REAL(mjg(nicbdj),wp) - 0.5_wp ) THEN
334               tmpberg => this
335               this => this%next
336               ibergs_to_send_s = ibergs_to_send_s + 1
337               IF( nn_verbose_level >= 4 ) THEN
338                  WRITE(numicb,*) 'bergstep ',nktberg,' packing berg ',tmpberg%number(:),' for transfer to south'
339                  CALL flush( numicb )
340               ENDIF
341               CALL icb_pack_into_buffer( tmpberg, obuffer_s, ibergs_to_send_s)
342               CALL icb_utl_delete(first_berg, tmpberg)
343            ELSE
344               this => this%next
345            ENDIF
346         END DO
347      ENDIF
348      if( nn_verbose_level >= 3) then
349         write(numicb,*) 'bergstep ',nktberg,' send ns: ', ibergs_to_send_n, ibergs_to_send_s
350         call flush(numicb)
351      endif
352
353      ! send bergs north
354      ! and receive bergs from south (ie ones sent north)
355     
356      IF( mpinei(jpso) >= 0  )   znsbergs(1) = ibergs_to_send_s
357      IF( mpinei(jpno) >= 0  )   zsnbergs(1) = ibergs_to_send_n
358      IF( mpinei(jpso) >= 0  )   CALL mppsend( 15, znsbergs(1), 1, ipe_S, iml_req2)
359      IF( mpinei(jpno) >= 0  )   CALL mppsend( 16, zsnbergs(1), 1, ipe_N, iml_req3)
360      IF( mpinei(jpno) >= 0  )   CALL mpprecv( 15, znsbergs(2), 1, ipe_N )
361      IF( mpinei(jpso) >= 0  )   CALL mpprecv( 16, zsnbergs(2), 1, ipe_S )
362      IF( mpinei(jpso) >= 0  )   CALL mpi_wait( iml_req2, iml_stat, iml_err )
363      IF( mpinei(jpno) >= 0  )   CALL mpi_wait( iml_req3, iml_stat, iml_err )
364      IF( mpinei(jpno) >= 0  )   ibergs_rcvd_from_n = INT( znsbergs(2) )
365      IF( mpinei(jpso) >= 0  )   ibergs_rcvd_from_s = INT( zsnbergs(2) )
366     
367      IF( nn_verbose_level >= 3) THEN
368         WRITE(numicb,*) 'bergstep ',nktberg,' recv ns: ', ibergs_rcvd_from_s, ibergs_rcvd_from_n
369         CALL FLUSH(numicb)
370      ENDIF
371
372      IF( ibergs_to_send_s > 0 ) CALL mppsend( 17, obuffer_s%data, ibergs_to_send_s*jp_buffer_width, ipe_S, iml_req2 )
373      IF( ibergs_to_send_n > 0 ) CALL mppsend( 18, obuffer_n%data, ibergs_to_send_n*jp_buffer_width, ipe_N, iml_req3 )
374      IF( ibergs_rcvd_from_n > 0 ) THEN
375         CALL icb_increase_ibuffer(ibuffer_n, ibergs_rcvd_from_n)
376         CALL mpprecv( 17, ibuffer_n%data, ibergs_rcvd_from_n*jp_buffer_width )
377      ENDIF
378      IF( ibergs_rcvd_from_s > 0 ) THEN
379         CALL icb_increase_ibuffer(ibuffer_s, ibergs_rcvd_from_s)
380         CALL mpprecv( 18, ibuffer_s%data, ibergs_rcvd_from_s*jp_buffer_width )
381      ENDIF
382      IF( ibergs_to_send_s > 0 ) CALL mpi_wait( iml_req2, iml_stat, iml_err )
383      IF( ibergs_to_send_n > 0 ) CALL mpi_wait( iml_req3, iml_stat, iml_err )
384      DO i = 1, ibergs_rcvd_from_n
385         IF( nn_verbose_level >= 4 ) THEN
386            WRITE(numicb,*) 'bergstep ',nktberg,' unpacking berg ',INT(ibuffer_n%data(16,i)),' from north'
387            CALL FLUSH( numicb )
388         ENDIF
389         CALL icb_unpack_from_buffer(first_berg, ibuffer_n, i)
390      END DO
391      DO i = 1, ibergs_rcvd_from_s
392         IF( nn_verbose_level >= 4 ) THEN
393            WRITE(numicb,*) 'bergstep ',nktberg,' unpacking berg ',INT(ibuffer_s%data(16,i)),' from south'
394            CALL FLUSH( numicb )
395         ENDIF
396         CALL icb_unpack_from_buffer(first_berg, ibuffer_s, i)
397      END DO
398     
399      IF( nn_verbose_level > 0 ) THEN
400         ! compare the number of icebergs on this processor from the start to the end
401         ibergs_end = icb_utl_count()
402         i = ( ibergs_rcvd_from_n + ibergs_rcvd_from_s + ibergs_rcvd_from_e + ibergs_rcvd_from_w ) - &
403             ( ibergs_to_send_n + ibergs_to_send_s + ibergs_to_send_e + ibergs_to_send_w )
404         IF( ibergs_end-(ibergs_start+i) .NE. 0 ) THEN
405            WRITE( numicb,*   ) 'send_bergs_to_other_pes: net change in number of icebergs'
406            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_end=', &
407                                ibergs_end,' on PE',narea
408            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_start=', &
409                                ibergs_start,' on PE',narea
410            WRITE( numicb,1000) 'send_bergs_to_other_pes: delta=', &
411                                i,' on PE',narea
412            WRITE( numicb,1000) 'send_bergs_to_other_pes: error=', &
413                                ibergs_end-(ibergs_start+i),' on PE',narea
414            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_to_send_n=', &
415                                ibergs_to_send_n,' on PE',narea
416            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_to_send_s=', &
417                                ibergs_to_send_s,' on PE',narea
418            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_to_send_e=', &
419                                ibergs_to_send_e,' on PE',narea
420            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_to_send_w=', &
421                                ibergs_to_send_w,' on PE',narea
422            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_rcvd_from_n=', &
423                                ibergs_rcvd_from_n,' on PE',narea
424            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_rcvd_from_s=', &
425                                ibergs_rcvd_from_s,' on PE',narea
426            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_rcvd_from_e=', &
427                                ibergs_rcvd_from_e,' on PE',narea
428            WRITE( numicb,1000) 'send_bergs_to_other_pes: ibergs_rcvd_from_w=', &
429                                ibergs_rcvd_from_w,' on PE',narea
430  1000      FORMAT(a,i5,a,i4)
431            CALL ctl_stop('send_bergs_to_other_pes: lost or gained an iceberg or two')
432         ENDIF
433      ENDIF
434
435      ! deal with north fold if we necessary when there is more than one top row processor
436      ! note that for jpni=1 north fold has been dealt with above in call to icb_lbc
437      IF( l_IdoNFold .AND. jpni > 1 ) CALL icb_lbc_mpp_nfld( )
438
439      IF( nn_verbose_level > 0 ) THEN
440         i = 0
441         this => first_berg
442         DO WHILE (ASSOCIATED(this))
443            pt => this%current_point
444            IF( pt%xi < REAL(mig(nicbdi),wp) - 0.5_wp .OR. &
445                pt%xi > REAL(mig(nicbei),wp) + 0.5_wp .OR. &
446                pt%yj < REAL(mjg(nicbdj),wp) - 0.5_wp .OR. &
447                pt%yj > REAL(mjg(nicbej),wp) + 0.5_wp ) THEN
448               i = i + 1
449               WRITE(numicb,*) 'berg lost in halo: ', this%number(:)
450               WRITE(numicb,*) '                   ', nimpp, njmpp
451               WRITE(numicb,*) '                   ', nicbdi, nicbei, nicbdj, nicbej
452               CALL flush( numicb )
453            ENDIF
454            this => this%next
455         ENDDO ! WHILE
456         CALL mpp_sum('icblbc', i)
457         IF( i .GT. 0 ) THEN
458            WRITE( numicb,'(a,i4)') 'send_bergs_to_other_pes: # of bergs outside computational domain = ',i
459            CALL ctl_stop('send_bergs_to_other_pes:  there are bergs still in halos!')
460         ENDIF ! root_pe
461      ENDIF ! debug
462      !
463      CALL mppsync()
464      !
465   END SUBROUTINE icb_lbc_mpp
466
467
468   SUBROUTINE icb_lbc_mpp_nfld()
469      !!----------------------------------------------------------------------
470      !!                 ***  SUBROUTINE icb_lbc_mpp_nfld  ***
471      !!
472      !! ** Purpose :   north fold treatment in multi processor exchange
473      !!
474      !! ** Method  :   
475      !!----------------------------------------------------------------------
476      TYPE(iceberg)         , POINTER     :: tmpberg, this
477      TYPE(point)           , POINTER     :: pt
478      INTEGER                             :: ibergs_to_send
479      INTEGER                             :: ibergs_to_rcv
480      INTEGER                             :: iiglo, ijglo, jk, jn
481      INTEGER                             :: ifldproc, iproc, ipts
482      INTEGER                             :: iine, ijne
483      INTEGER                             :: jjn
484      REAL(wp), DIMENSION(0:3)            :: zsbergs, znbergs
485      INTEGER                             :: iml_req1, iml_req2, iml_err
486      INTEGER, DIMENSION(MPI_STATUS_SIZE) :: iml_stat
487
488      ! set up indices of neighbouring processors
489
490      ! nicbfldproc is a list of unique processor numbers that this processor
491      ! exchanges with (including itself), so we loop over this array; since
492      ! its of fixed size, the first -1 marks end of list of processors
493      !
494      nicbfldnsend(:) = 0
495      nicbfldexpect(:) = 0
496      nicbfldreq(:) = 0
497      !
498      ! Since each processor may be communicating with more than one northern
499      ! neighbour, cycle through the sends so that the receive order can be
500      ! controlled.
501      !
502      ! First compute how many icebergs each active neighbour should expect
503      DO jn = 1, jpni
504         IF( nicbfldproc(jn) /= -1 ) THEN
505            ifldproc = nicbfldproc(jn)
506            nicbfldnsend(jn) = 0
507
508            ! Find number of bergs that need to be exchanged
509            ! Pick out exchanges with processor ifldproc
510            ! if ifldproc is this processor then don't send
511            !
512            IF( ASSOCIATED(first_berg) ) THEN
513               this => first_berg
514               DO WHILE (ASSOCIATED(this))
515                  pt => this%current_point
516                  iine = INT( pt%xi + 0.5 )
517                  iproc = nicbflddest(mi1(iine))
518                  IF( pt%yj > REAL(mjg(nicbej),wp) + 0.5_wp ) THEN
519                     IF( iproc == ifldproc ) THEN
520                        !
521                        IF( iproc /= narea ) THEN
522                           tmpberg => this
523                           nicbfldnsend(jn) = nicbfldnsend(jn) + 1
524                        ENDIF
525                        !
526                     ENDIF
527                  ENDIF
528                  this => this%next
529               END DO
530            ENDIF
531            !
532         ENDIF
533         !
534      END DO
535      !
536      ! Now tell each active neighbour how many icebergs to expect
537      DO jn = 1, jpni
538         IF( nicbfldproc(jn) /= -1 ) THEN
539            ifldproc = nicbfldproc(jn)
540            IF( ifldproc == narea ) CYCLE
541   
542            zsbergs(0) = narea
543            zsbergs(1) = nicbfldnsend(jn)
544            !IF ( nicbfldnsend(jn) .GT. 0 .AND. nn_verbose_level > 0 ) write(numicb,*) 'ICB sending ',nicbfldnsend(jn),' to ', ifldproc
545            CALL mppsend( 21, zsbergs(0:1), 2, ifldproc-1, nicbfldreq(jn))
546         ENDIF
547         !
548      END DO
549      !
550      ! and receive the heads-up from active neighbours preparing to send
551      DO jn = 1, jpni
552         IF( nicbfldproc(jn) /= -1 ) THEN
553            ifldproc = nicbfldproc(jn)
554            IF( ifldproc == narea ) CYCLE
555
556            CALL mpprecv( 21, znbergs(1:2), 2 )
557            DO jjn = 1,jpni
558             IF( nicbfldproc(jjn) .eq. INT(znbergs(1)) ) EXIT
559            END DO
560            IF( jjn .GT. jpni .AND. nn_verbose_level > 0 ) write(numicb,*) 'ICB ERROR'
561            nicbfldexpect(jjn) = INT( znbergs(2) )
562            !IF ( nicbfldexpect(jjn) .GT. 0 .AND. nn_verbose_level > 0 ) write(numicb,*) 'ICB expecting ',nicbfldexpect(jjn),' from ', nicbfldproc(jjn)
563            !IF (nn_verbose_level > 0) CALL FLUSH(numicb)
564         ENDIF
565         !
566      END DO
567      !
568      ! post the mpi waits if using immediate send protocol
569      DO jn = 1, jpni
570         IF( nicbfldproc(jn) /= -1 ) THEN
571            ifldproc = nicbfldproc(jn)
572            IF( ifldproc == narea ) CYCLE
573            CALL mpi_wait( nicbfldreq(jn), iml_stat, iml_err )
574         ENDIF
575         !
576      END DO
577   
578         !
579         ! Cycle through the icebergs again, this time packing and sending any
580         ! going through the north fold. They will be expected.
581      DO jn = 1, jpni
582         IF( nicbfldproc(jn) /= -1 ) THEN
583            ifldproc = nicbfldproc(jn)
584            ibergs_to_send = 0
585   
586            ! Find number of bergs that need to be exchanged
587            ! Pick out exchanges with processor ifldproc
588            ! if ifldproc is this processor then don't send
589            !
590            IF( ASSOCIATED(first_berg) ) THEN
591               this => first_berg
592               DO WHILE (ASSOCIATED(this))
593                  pt => this%current_point
594                  iine = INT( pt%xi + 0.5 )
595                  ijne = INT( pt%yj + 0.5 )
596                  ipts  = nicbfldpts (mi1(iine))
597                  iproc = nicbflddest(mi1(iine))
598                  IF( pt%yj > REAL(mjg(nicbej),wp) + 0.5_wp ) THEN
599                     IF( iproc == ifldproc ) THEN
600                        !
601                        ! moving across the cut line means both position and
602                        ! velocity must change
603                        ijglo = INT( ipts/nicbpack )
604                        iiglo = ipts - nicbpack*ijglo
605                        pt%xi = iiglo - ( pt%xi - REAL(iine,wp) )
606                        pt%yj = ijglo - ( pt%yj - REAL(ijne,wp) )
607                        pt%uvel = -1._wp * pt%uvel
608                        pt%vvel = -1._wp * pt%vvel
609                        !
610                        ! now remove berg from list and pack it into a buffer
611                        IF( iproc /= narea ) THEN
612                           tmpberg => this
613                           ibergs_to_send = ibergs_to_send + 1
614                           IF( nn_verbose_level >= 4 ) THEN
615                              WRITE(numicb,*) 'bergstep ',nktberg,' packing berg ',tmpberg%number(:),' for north fold'
616                              CALL flush( numicb )
617                           ENDIF
618                           CALL icb_pack_into_buffer( tmpberg, obuffer_f, ibergs_to_send)
619                           CALL icb_utl_delete(first_berg, tmpberg)
620                        ENDIF
621                        !
622                     ENDIF
623                  ENDIF
624                  this => this%next
625               END DO
626            ENDIF
627            if( nn_verbose_level >= 3) then
628               write(numicb,*) 'bergstep ',nktberg,' send nfld: ', ibergs_to_send
629               call flush(numicb)
630            endif
631            !
632            ! if we're in this processor, then we've done everything we need to
633            ! so go on to next element of loop
634            IF( ifldproc == narea ) CYCLE
635   
636            ! send bergs
637   
638            IF( ibergs_to_send > 0 )  &
639                CALL mppsend( 12, obuffer_f%data, ibergs_to_send*jp_buffer_width, ifldproc-1, nicbfldreq(jn) )
640            !
641         ENDIF
642         !
643      END DO
644      !
645      ! Now receive the expected number of bergs from the active neighbours
646      DO jn = 1, jpni
647         IF( nicbfldproc(jn) /= -1 ) THEN
648            ifldproc = nicbfldproc(jn)
649            IF( ifldproc == narea ) CYCLE
650            ibergs_to_rcv = nicbfldexpect(jn)
651
652            IF( ibergs_to_rcv  > 0 ) THEN
653               CALL icb_increase_ibuffer(ibuffer_f, ibergs_to_rcv)
654               CALL mpprecv( 12, ibuffer_f%data, ibergs_to_rcv*jp_buffer_width, ifldproc-1 )
655            ENDIF
656            !
657            DO jk = 1, ibergs_to_rcv
658               IF( nn_verbose_level >= 4 ) THEN
659                  WRITE(numicb,*) 'bergstep ',nktberg,' unpacking berg ',INT(ibuffer_f%data(16,jk)),' from north fold'
660                  CALL flush( numicb )
661               ENDIF
662               CALL icb_unpack_from_buffer(first_berg, ibuffer_f, jk )
663            END DO
664         ENDIF
665         !
666      END DO
667      !
668      ! Finally post the mpi waits if using immediate send protocol
669      DO jn = 1, jpni
670         IF( nicbfldproc(jn) /= -1 ) THEN
671            ifldproc = nicbfldproc(jn)
672            IF( ifldproc == narea ) CYCLE
673            CALL mpi_wait( nicbfldreq(jn), iml_stat, iml_err )
674         ENDIF
675         !
676      END DO
677      !
678   END SUBROUTINE icb_lbc_mpp_nfld
679
680
681   SUBROUTINE icb_pack_into_buffer( berg, pbuff, kb )
682      !!----------------------------------------------------------------------
683      !!----------------------------------------------------------------------
684      TYPE(iceberg), POINTER :: berg
685      TYPE(buffer) , POINTER :: pbuff
686      INTEGER               , INTENT(in) :: kb
687      !
688      INTEGER ::   k   ! local integer
689      !!----------------------------------------------------------------------
690      !
691      IF( .NOT. ASSOCIATED(pbuff) ) CALL icb_increase_buffer( pbuff, jp_delta_buf )
692      IF( kb .GT. pbuff%size ) CALL icb_increase_buffer( pbuff, jp_delta_buf )
693
694      !! pack points into buffer
695
696      pbuff%data( 1,kb) = berg%current_point%lon
697      pbuff%data( 2,kb) = berg%current_point%lat
698      pbuff%data( 3,kb) = berg%current_point%uvel
699      pbuff%data( 4,kb) = berg%current_point%vvel
700      pbuff%data( 5,kb) = berg%current_point%xi
701      pbuff%data( 6,kb) = berg%current_point%yj
702      pbuff%data( 7,kb) = float(berg%current_point%year)
703      pbuff%data( 8,kb) = berg%current_point%day
704      pbuff%data( 9,kb) = berg%current_point%mass
705      pbuff%data(10,kb) = berg%current_point%thickness
706      pbuff%data(11,kb) = berg%current_point%width
707      pbuff%data(12,kb) = berg%current_point%length
708      pbuff%data(13,kb) = berg%current_point%mass_of_bits
709      pbuff%data(14,kb) = berg%current_point%heat_density
710
711      pbuff%data(15,kb) = berg%mass_scaling
712      DO k=1,nkounts
713         pbuff%data(15+k,kb) = REAL( berg%number(k), wp )
714      END DO
715      !
716   END SUBROUTINE icb_pack_into_buffer
717
718
719   SUBROUTINE icb_unpack_from_buffer(first, pbuff, kb)
720      !!----------------------------------------------------------------------
721      !!----------------------------------------------------------------------
722      TYPE(iceberg),             POINTER :: first
723      TYPE(buffer) ,             POINTER :: pbuff
724      INTEGER      , INTENT(in)          :: kb
725      !
726      TYPE(iceberg)                      :: currentberg
727      TYPE(point)                        :: pt
728      INTEGER                            :: ik
729      !!----------------------------------------------------------------------
730      !
731      pt%lon            =      pbuff%data( 1,kb)
732      pt%lat            =      pbuff%data( 2,kb)
733      pt%uvel           =      pbuff%data( 3,kb)
734      pt%vvel           =      pbuff%data( 4,kb)
735      pt%xi             =      pbuff%data( 5,kb)
736      pt%yj             =      pbuff%data( 6,kb)
737      pt%year           = INT( pbuff%data( 7,kb) )
738      pt%day            =      pbuff%data( 8,kb)
739      pt%mass           =      pbuff%data( 9,kb)
740      pt%thickness      =      pbuff%data(10,kb)
741      pt%width          =      pbuff%data(11,kb)
742      pt%length         =      pbuff%data(12,kb)
743      pt%mass_of_bits   =      pbuff%data(13,kb)
744      pt%heat_density   =      pbuff%data(14,kb)
745
746      currentberg%mass_scaling =      pbuff%data(15,kb)
747      DO ik = 1, nkounts
748         currentberg%number(ik) = INT( pbuff%data(15+ik,kb) )
749      END DO
750      !
751      CALL icb_utl_add(currentberg, pt )
752      !
753   END SUBROUTINE icb_unpack_from_buffer
754
755
756   SUBROUTINE icb_increase_buffer(old,kdelta)
757      !!----------------------------------------------------------------------
758      TYPE(buffer), POINTER    :: old
759      INTEGER     , INTENT(in) :: kdelta
760      !
761      TYPE(buffer), POINTER ::   new
762      INTEGER ::   inew_size
763      !!----------------------------------------------------------------------
764      !
765      IF( .NOT. ASSOCIATED(old) ) THEN   ;   inew_size = kdelta
766      ELSE                               ;   inew_size = old%size + kdelta
767      ENDIF
768      ALLOCATE( new )
769      ALLOCATE( new%data( jp_buffer_width, inew_size) )
770      new%size = inew_size
771      IF( ASSOCIATED(old) ) THEN
772         new%data(:,1:old%size) = old%data(:,1:old%size)
773         DEALLOCATE(old%data)
774         DEALLOCATE(old)
775      ENDIF
776      old => new
777      !
778   END SUBROUTINE icb_increase_buffer
779
780
781   SUBROUTINE icb_increase_ibuffer(old,kdelta)
782      !!----------------------------------------------------------------------
783      !!----------------------------------------------------------------------
784      TYPE(buffer),            POINTER :: old
785      INTEGER     , INTENT(in)         :: kdelta
786      !
787      TYPE(buffer),            POINTER :: new
788      INTEGER                          :: inew_size, iold_size
789      !!----------------------------------------------------------------------
790
791      IF( .NOT. ASSOCIATED(old) ) THEN
792         inew_size = kdelta + jp_delta_buf
793         iold_size = 0
794      ELSE
795         iold_size = old%size
796         IF( kdelta .LT. old%size ) THEN
797            inew_size = old%size + kdelta
798         ELSE
799            inew_size = kdelta + jp_delta_buf
800         ENDIF
801      ENDIF
802
803      IF( iold_size .NE. inew_size ) THEN
804         ALLOCATE( new )
805         ALLOCATE( new%data( jp_buffer_width, inew_size) )
806         new%size = inew_size
807         IF( ASSOCIATED(old) ) THEN
808            new%data(:,1:old%size) = old%data(:,1:old%size)
809            DEALLOCATE(old%data)
810            DEALLOCATE(old)
811         ENDIF
812         old => new
813         !IF (nn_verbose_level > 0) WRITE( numicb,*) 'icb_increase_ibuffer',narea,' increased to',inew_size
814      ENDIF
815      !
816   END SUBROUTINE icb_increase_ibuffer
817
818#else
819   !!----------------------------------------------------------------------
820   !!   Default case:            Dummy module        share memory computing
821   !!----------------------------------------------------------------------
822   SUBROUTINE icb_lbc_mpp()
823      WRITE(numout,*) 'icb_lbc_mpp: You should not have seen this message!!'
824   END SUBROUTINE icb_lbc_mpp
825#endif
826
827   !!======================================================================
828END MODULE icblbc
Note: See TracBrowser for help on using the repository browser.