00001 #define _SINGLEMSG 1
00002
00003
00004
00005
00006 module ice_gather_scatter
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 use ice_kinds_mod
00027 use ice_communicate
00028 use ice_constants
00029 use ice_blocks
00030 use ice_distribution
00031 use ice_domain
00032 use ice_domain_size
00033 use ice_exit
00034
00035 implicit none
00036 private
00037 save
00038
00039
00040
00041 public :: gather_global, &
00042 scatter_global, &
00043 gatherArray, &
00044 scatter_global_stress
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054 interface gather_global
00055 module procedure gather_global_dbl, &
00056 gather_global_real, &
00057 gather_global_int
00058 end interface
00059
00060 interface scatter_global
00061 module procedure scatter_global_dbl, &
00062 scatter_global_real, &
00063 scatter_global_int
00064 end interface
00065
00066 interface gatherArray
00067 module procedure gatherArray_dbl
00068 end interface
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079 contains
00080
00081
00082 subroutine gatherArray_dbl(array_g,array,length,root)
00083
00084 include 'mpif.h'
00085
00086 real(dbl_kind) :: array_g(:)
00087 real(dbl_kind) :: array(:)
00088 integer(int_kind) :: length
00089 integer(int_kind) :: root
00090
00091 integer(int_kind) :: ierr
00092
00093 call MPI_Gather(array,length,MPI_REAL8,array_g, length,MPI_REAL8,root,MPI_COMM_ICE,ierr)
00094
00095 end subroutine gatherArray_dbl
00096
00097
00098
00099
00100
00101
00102
00103 subroutine gather_global_dbl(ARRAY_G, ARRAY, dst_task, src_dist)
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122 include 'mpif.h'
00123
00124
00125
00126 integer (int_kind), intent(in) ::
00127 dst_task
00128
00129 type (distrb), intent(in) ::
00130 src_dist
00131
00132 real (dbl_kind), dimension(:,:,:), intent(in) ::
00133 ARRAY
00134
00135
00136
00137 real (dbl_kind), dimension(:,:), intent(inout) ::
00138 ARRAY_G
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148 integer (int_kind) ::
00149 i,j,n ,
00150 nsends ,
00151 src_block ,
00152 src_task ,
00153 ierr
00154
00155 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
00156 status
00157
00158 integer (int_kind), dimension(:), allocatable ::
00159 snd_request
00160
00161 integer (int_kind), dimension(:,:), allocatable ::
00162 snd_status
00163
00164 #ifdef _SINGLEMSG
00165 real (dbl_kind), dimension(:,:,:), allocatable ::
00166 msg_buffer
00167 #else
00168 real (dbl_kind), dimension(:,:), allocatable ::
00169 msg_buffer
00170 #endif
00171
00172 type (block) ::
00173 this_block
00174
00175 integer (int_kind) :: ib, ig, itask, nprocs, maxBlocks
00176 integer (int_kind) :: iig,ijg,it
00177 integer (int_kind) :: msgLen, msgTag
00178
00179
00180
00181
00182
00183
00184
00185
00186 nprocs = get_num_procs()
00187
00188 if (my_task == dst_task) then
00189
00190 do n=1,nblocks_tot
00191
00192
00193
00194 if (src_dist%blockLocation(n) == my_task+1) then
00195
00196 this_block = get_block(n,n)
00197
00198 do j=this_block%jlo,this_block%jhi
00199 do i=this_block%ilo,this_block%ihi
00200 ARRAY_G(this_block%i_glob(i), &
00201 this_block%j_glob(j)) = &
00202 ARRAY(i,j,src_dist%blockLocalID(n))
00203 end do
00204 end do
00205
00206
00207
00208 else if (src_dist%blockLocation(n) == 0) then
00209
00210 this_block = get_block(n,n)
00211
00212 do j=this_block%jlo,this_block%jhi
00213 do i=this_block%ilo,this_block%ihi
00214 ARRAY_G(this_block%i_glob(i), &
00215 this_block%j_glob(j)) = c0
00216 end do
00217 end do
00218 endif
00219
00220 end do
00221
00222
00223 #ifdef _SINGLEMSG
00224 allocate (msg_buffer(nx_block,ny_block,max_blocks))
00225 do itask = 0,nprocs-1
00226 if(itask /= dst_task) then
00227 maxBlocks = src_dist%BlockCnt(itask+1)
00228 msgLen = nx_block*ny_block*maxBlocks
00229 msgTag = 3*mpitag_gs+itask
00230 if(maxBlocks>0) then
00231
00232 call MPI_RECV(msg_buffer(:,:,1:maxBlocks),msgLen, &
00233 mpiR8, itask, msgTag, MPI_COMM_ICE, status, ierr)
00234
00235 do ib=1,maxBlocks
00236 ig = src_dist%blockIndex(itask+1,ib)
00237 this_block = get_block(ig,ig)
00238 do j=this_block%jlo,this_block%jhi
00239 do i=this_block%ilo,this_block%ihi
00240 iig = this_block%i_glob(i)
00241 ijg = this_block%j_glob(j)
00242 ARRAY_G(iig,ijg) = msg_buffer(i,j,ib)
00243 end do
00244 end do
00245 enddo
00246 endif
00247 endif
00248 enddo
00249
00250 deallocate(msg_buffer)
00251 #else
00252 allocate (msg_buffer(nx_block,ny_block))
00253
00254 do n=1,nblocks_tot
00255 if (src_dist%blockLocation(n) > 0 .and. &
00256 src_dist%blockLocation(n) /= my_task+1) then
00257
00258 this_block = get_block(n,n)
00259
00260 call MPI_RECV(msg_buffer, size(msg_buffer), &
00261 mpiR8, src_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
00262 MPI_COMM_ICE, status, ierr)
00263
00264 do j=this_block%jlo,this_block%jhi
00265 do i=this_block%ilo,this_block%ihi
00266 ARRAY_G(this_block%i_glob(i), &
00267 this_block%j_glob(j)) = msg_buffer(i,j)
00268 end do
00269 end do
00270 endif
00271 end do
00272
00273 deallocate(msg_buffer)
00274 #endif
00275
00276
00277
00278
00279
00280
00281
00282 else
00283
00284 #ifdef _SINGLEMSG
00285
00286 if(nblocks>0) then
00287 msgLen = nx_block*ny_block*nblocks
00288 msgTag = 3*mpitag_gs+my_task
00289 call MPI_SEND(ARRAY(:,:,1:nblocks),msgLen,mpiR8,dst_task,msgTag,MPI_COMM_ICE,ierr)
00290 endif
00291
00292 #else
00293 allocate(snd_request(nblocks_tot), &
00294 snd_status (MPI_STATUS_SIZE, nblocks_tot))
00295
00296 nsends = 0
00297 do n=1,nblocks_tot
00298 if (src_dist%blockLocation(n) == my_task+1) then
00299
00300 nsends = nsends + 1
00301 src_block = src_dist%blockLocalID(n)
00302 call MPI_ISEND(ARRAY(1,1,src_block), nx_block*ny_block, &
00303 mpiR8, dst_task, 3*mpitag_gs+n, &
00304 MPI_COMM_ICE, snd_request(nsends), ierr)
00305 endif
00306 end do
00307
00308 if (nsends > 0) &
00309 call MPI_WAITALL(nsends, snd_request, snd_status, ierr)
00310 deallocate(snd_request, snd_status)
00311 #endif
00312
00313 endif
00314
00315
00316
00317 end subroutine gather_global_dbl
00318
00319
00320
00321 subroutine gather_global_real(ARRAY_G, ARRAY, dst_task, src_dist)
00322
00323
00324
00325
00326
00327
00328
00329
00330 include 'mpif.h'
00331
00332
00333
00334
00335
00336
00337
00338 integer (int_kind), intent(in) ::
00339 dst_task
00340
00341 type (distrb), intent(in) ::
00342 src_dist
00343
00344 real (real_kind), dimension(:,:,:), intent(in) ::
00345 ARRAY
00346
00347
00348
00349
00350
00351
00352
00353 real (real_kind), dimension(:,:), intent(inout) ::
00354 ARRAY_G
00355
00356
00357
00358
00359
00360
00361
00362 integer (int_kind) ::
00363 i,j,n ,
00364 nsends ,
00365 src_block ,
00366 ierr
00367
00368 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
00369 status
00370
00371 integer (int_kind), dimension(:), allocatable ::
00372 snd_request
00373
00374 integer (int_kind), dimension(:,:), allocatable ::
00375 snd_status
00376
00377 #ifdef _SINGLEMSG
00378 real (real_kind), dimension(:,:,:), allocatable ::
00379 msg_buffer
00380 #else
00381 real (real_kind), dimension(:,:), allocatable ::
00382 msg_buffer
00383 #endif
00384
00385 type (block) ::
00386 this_block
00387
00388 integer (int_kind) :: ib, ig, itask, nprocs, maxBlocks
00389 integer (int_kind) :: iig,ijg
00390 integer (int_kind) :: msgLen, msgTag
00391
00392
00393
00394
00395
00396
00397
00398
00399 nprocs = get_num_procs()
00400
00401 if (my_task == dst_task) then
00402
00403 do n=1,nblocks_tot
00404
00405
00406
00407 if (src_dist%blockLocation(n) == my_task+1) then
00408
00409 this_block = get_block(n,n)
00410
00411 do j=this_block%jlo,this_block%jhi
00412 do i=this_block%ilo,this_block%ihi
00413 ARRAY_G(this_block%i_glob(i), &
00414 this_block%j_glob(j)) = &
00415 ARRAY(i,j,src_dist%blockLocalID(n))
00416 end do
00417 end do
00418
00419
00420
00421 else if (src_dist%blockLocation(n) == 0) then
00422
00423 this_block = get_block(n,n)
00424
00425 do j=this_block%jlo,this_block%jhi
00426 do i=this_block%ilo,this_block%ihi
00427 ARRAY_G(this_block%i_glob(i), &
00428 this_block%j_glob(j)) = spval
00429 end do
00430 end do
00431 endif
00432
00433 end do
00434
00435
00436
00437 #ifdef _SINGLEMSG
00438
00439 allocate (msg_buffer(nx_block,ny_block,max_blocks))
00440 do itask = 0,nprocs-1
00441 if( itask /= dst_task) then
00442 maxBlocks = src_dist%BlockCnt(itask+1)
00443 msgLen = nx_block*ny_block*maxBlocks
00444 msgTag = 3*mpitag_gs+itask
00445 if( maxBlocks>0) then
00446
00447 call MPI_RECV(msg_buffer(:,:,1:maxBlocks), msgLen, &
00448 mpiR4, itask, msgTag ,MPI_COMM_ICE, status, ierr)
00449
00450 do ib=1,maxBlocks
00451 ig = src_dist%blockIndex(itask+1,ib)
00452 this_block = get_block(ig,ig)
00453 do j=this_block%jlo,this_block%jhi
00454 do i=this_block%ilo,this_block%ihi
00455 iig = this_block%i_glob(i)
00456 ijg = this_block%j_glob(j)
00457 ARRAY_G(iig,ijg) = msg_buffer(i,j,ib)
00458 end do
00459 end do
00460 enddo
00461 endif
00462 endif
00463 enddo
00464
00465 deallocate(msg_buffer)
00466
00467 #else
00468 allocate (msg_buffer(nx_block,ny_block))
00469
00470 do n=1,nblocks_tot
00471 if (src_dist%blockLocation(n) > 0 .and. &
00472 src_dist%blockLocation(n) /= my_task+1) then
00473
00474 this_block = get_block(n,n)
00475
00476 call MPI_RECV(msg_buffer, size(msg_buffer), &
00477 mpiR4, src_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
00478 MPI_COMM_ICE, status, ierr)
00479
00480 do j=this_block%jlo,this_block%jhi
00481 do i=this_block%ilo,this_block%ihi
00482 ARRAY_G(this_block%i_glob(i), &
00483 this_block%j_glob(j)) = msg_buffer(i,j)
00484 end do
00485 end do
00486 endif
00487 end do
00488 deallocate(msg_buffer)
00489 #endif
00490
00491
00492
00493
00494
00495
00496
00497
00498 else
00499 #ifdef _SINGLEMSG
00500
00501 if(nblocks>0) then
00502 msgLen = nx_block*ny_block*nblocks
00503 msgTag = 3*mpitag_gs+my_task
00504 call MPI_SEND(ARRAY(:,:,1:nblocks),msgLen,mpiR4,dst_task,msgTag,MPI_COMM_ICE,ierr)
00505 endif
00506
00507 #else
00508
00509 allocate(snd_request(nblocks_tot), &
00510 snd_status (MPI_STATUS_SIZE, nblocks_tot))
00511
00512 nsends = 0
00513 do n=1,nblocks_tot
00514 if (src_dist%blockLocation(n) == my_task+1) then
00515
00516 nsends = nsends + 1
00517 src_block = src_dist%blockLocalID(n)
00518 call MPI_ISEND(ARRAY(1,1,src_block), nx_block*ny_block, &
00519 mpiR4, dst_task, 3*mpitag_gs+n, &
00520 MPI_COMM_ICE, snd_request(nsends), ierr)
00521 endif
00522 end do
00523
00524 if (nsends > 0) &
00525 call MPI_WAITALL(nsends, snd_request, snd_status, ierr)
00526 deallocate(snd_request, snd_status)
00527 #endif
00528
00529 endif
00530
00531
00532
00533 end subroutine gather_global_real
00534
00535
00536
00537 subroutine gather_global_int(ARRAY_G, ARRAY, dst_task, src_dist)
00538
00539
00540
00541
00542
00543
00544
00545
00546 include 'mpif.h'
00547
00548
00549
00550
00551
00552
00553
00554 integer (int_kind), intent(in) ::
00555 dst_task
00556
00557 type (distrb), intent(in) ::
00558 src_dist
00559
00560 integer (int_kind), dimension(:,:,:), intent(in) ::
00561 ARRAY
00562
00563
00564
00565
00566
00567
00568
00569 integer (int_kind), dimension(:,:), intent(inout) ::
00570 ARRAY_G
00571
00572
00573
00574
00575
00576
00577
00578 integer (int_kind) ::
00579 i,j,n ,
00580 nsends ,
00581 src_block ,
00582 ierr
00583
00584 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
00585 status
00586
00587 integer (int_kind), dimension(:), allocatable ::
00588 snd_request
00589
00590 integer (int_kind), dimension(:,:), allocatable ::
00591 snd_status
00592
00593 #ifdef _SINGLEMSG
00594 integer (int_kind), dimension(:,:,:), allocatable ::
00595 msg_buffer
00596 #else
00597 integer (int_kind), dimension(:,:), allocatable ::
00598 msg_buffer
00599 #endif
00600
00601 type (block) ::
00602 this_block
00603
00604 integer (int_kind) :: ib, ig, len, itask, nprocs, maxBlocks
00605 integer (int_kind) :: iig,ijg
00606 integer (int_kind) :: msgLen, msgTag
00607
00608
00609
00610
00611
00612
00613
00614
00615 nprocs = get_num_procs()
00616
00617 if (my_task == dst_task) then
00618
00619 do n=1,nblocks_tot
00620
00621
00622
00623 if (src_dist%blockLocation(n) == my_task+1) then
00624
00625 this_block = get_block(n,n)
00626
00627 do j=this_block%jlo,this_block%jhi
00628 do i=this_block%ilo,this_block%ihi
00629 ARRAY_G(this_block%i_glob(i), &
00630 this_block%j_glob(j)) = &
00631 ARRAY(i,j,src_dist%blockLocalID(n))
00632 end do
00633 end do
00634
00635
00636
00637 else if (src_dist%blockLocation(n) == 0) then
00638
00639 this_block = get_block(n,n)
00640
00641 do j=this_block%jlo,this_block%jhi
00642 do i=this_block%ilo,this_block%ihi
00643 ARRAY_G(this_block%i_glob(i), &
00644 this_block%j_glob(j)) = 0
00645 end do
00646 end do
00647 endif
00648
00649 end do
00650
00651
00652
00653 #ifdef _SINGLEMSG
00654
00655 allocate (msg_buffer(nx_block,ny_block,max_blocks))
00656 do itask = 0,nprocs-1
00657 if( itask /= dst_task) then
00658 maxBlocks = src_dist%BlockCnt(itask+1)
00659 msgLen = nx_block*ny_block*max_blocks
00660 msgTag = 3*mpitag_gs+itask
00661 if(maxBLocks>0) then
00662 call MPI_RECV(msg_buffer(:,:,1:maxBlocks),msgLen, &
00663 mpi_integer, itask, msgTag, MPI_COMM_ICE, status, ierr)
00664 do ib=1,maxBlocks
00665 ig = src_dist%blockIndex(itask+1,ib)
00666 this_block = get_block(ig,ig)
00667 do j=this_block%jlo,this_block%jhi
00668 do i=this_block%ilo,this_block%ihi
00669 iig = this_block%i_glob(i)
00670 ijg = this_block%j_glob(j)
00671 ARRAY_G(iig,ijg) = msg_buffer(i,j,ib)
00672 end do
00673 end do
00674 enddo
00675 endif
00676 endif
00677 enddo
00678 deallocate(msg_buffer)
00679
00680 #else
00681 allocate (msg_buffer(nx_block,ny_block))
00682
00683 do n=1,nblocks_tot
00684 if (src_dist%blockLocation(n) > 0 .and. &
00685 src_dist%blockLocation(n) /= my_task+1) then
00686
00687 this_block = get_block(n,n)
00688
00689 call MPI_RECV(msg_buffer, size(msg_buffer), &
00690 mpi_integer, src_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
00691 MPI_COMM_ICE, status, ierr)
00692
00693 do j=this_block%jlo,this_block%jhi
00694 do i=this_block%ilo,this_block%ihi
00695 ARRAY_G(this_block%i_glob(i), &
00696 this_block%j_glob(j)) = msg_buffer(i,j)
00697 end do
00698 end do
00699 endif
00700 end do
00701 deallocate(msg_buffer)
00702 #endif
00703
00704
00705
00706
00707
00708
00709
00710
00711 else
00712
00713 #ifdef _SINGLEMSG
00714
00715 if(nblocks>0) then
00716 msgLen = nx_block*ny_block*nblocks
00717 msgTag = 3*mpitag_gs+my_task
00718 call MPI_SEND(ARRAY(:,:,1:nblocks),msgLen, &
00719 mpi_integer, dst_task, msgTag, MPI_COMM_ICE, ierr)
00720 endif
00721
00722 #else
00723
00724 allocate(snd_request(nblocks_tot), &
00725 snd_status (MPI_STATUS_SIZE, nblocks_tot))
00726
00727 nsends = 0
00728 do n=1,nblocks_tot
00729 if (src_dist%blockLocation(n) == my_task+1) then
00730
00731 nsends = nsends + 1
00732 src_block = src_dist%blockLocalID(n)
00733 call MPI_ISEND(ARRAY(1,1,src_block), nx_block*ny_block, &
00734 mpi_integer, dst_task, 3*mpitag_gs+n, &
00735 MPI_COMM_ICE, snd_request(nsends), ierr)
00736 endif
00737 end do
00738
00739 if (nsends > 0) &
00740 call MPI_WAITALL(nsends, snd_request, snd_status, ierr)
00741 deallocate(snd_request, snd_status)
00742 #endif
00743
00744 endif
00745
00746
00747
00748 end subroutine gather_global_int
00749
00750
00751
00752
00753
00754
00755
00756 subroutine scatter_global_dbl(ARRAY, ARRAY_G, src_task, dst_dist, &
00757 field_loc, field_type)
00758
00759
00760
00761
00762
00763
00764
00765
00766
00767
00768
00769
00770
00771 include 'mpif.h'
00772
00773
00774
00775 integer (int_kind), intent(in) ::
00776 src_task
00777
00778 type (distrb), intent(in) ::
00779 dst_dist
00780
00781 real (dbl_kind), dimension(:,:), intent(in) ::
00782 ARRAY_G
00783
00784 integer (int_kind), intent(in) ::
00785 field_type,
00786 field_loc
00787
00788
00789
00790
00791 real (dbl_kind), dimension(:,:,:), intent(inout) ::
00792 ARRAY
00793
00794
00795
00796
00797
00798
00799
00800
00801
00802 integer (int_kind) ::
00803 i,j,n,bid,
00804 nrecvs,
00805 isrc, jsrc,
00806 dst_block,
00807 xoffset, yoffset,
00808 yoffset2,
00809 isign,
00810 ierr
00811
00812 type (block) ::
00813 this_block
00814
00815 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
00816 status
00817
00818 integer (int_kind), dimension(:), allocatable ::
00819 rcv_request
00820
00821 integer (int_kind), dimension(:,:), allocatable ::
00822 rcv_status
00823
00824 real (dbl_kind), dimension(:,:), allocatable ::
00825 msg_buffer
00826
00827
00828
00829
00830
00831
00832
00833 ARRAY = c0
00834
00835 this_block = get_block(1,1)
00836 if (this_block%tripoleTFlag) then
00837 select case (field_loc)
00838 case (field_loc_center)
00839 xoffset = 2
00840 yoffset = 0
00841 case (field_loc_NEcorner)
00842 xoffset = 1
00843 yoffset = -1
00844 case (field_loc_Eface)
00845 xoffset = 1
00846 yoffset = 0
00847 case (field_loc_Nface)
00848 xoffset = 2
00849 yoffset = -1
00850 case (field_loc_noupdate)
00851 xoffset = 1
00852 yoffset = 1
00853 end select
00854 else
00855 select case (field_loc)
00856 case (field_loc_center)
00857 xoffset = 1
00858 yoffset = 1
00859 case (field_loc_NEcorner)
00860 xoffset = 0
00861 yoffset = 0
00862 case (field_loc_Eface)
00863 xoffset = 0
00864 yoffset = 1
00865 case (field_loc_Nface)
00866 xoffset = 1
00867 yoffset = 0
00868 case (field_loc_noupdate)
00869 xoffset = 1
00870 yoffset = 1
00871 end select
00872 endif
00873
00874 select case (field_type)
00875 case (field_type_scalar)
00876 isign = 1
00877 case (field_type_vector)
00878 isign = -1
00879 case (field_type_angle)
00880 isign = -1
00881 case (field_type_noupdate)
00882 isign = 1
00883 case default
00884 call abort_ice('Unknown field type in scatter')
00885 end select
00886
00887
00888
00889
00890
00891
00892
00893
00894 if (my_task == src_task) then
00895
00896
00897
00898 allocate (msg_buffer(nx_block,ny_block))
00899
00900 do n=1,nblocks_tot
00901 if (dst_dist%blockLocation(n) > 0 .and. &
00902 dst_dist%blockLocation(n)-1 /= my_task) then
00903
00904 msg_buffer = c0
00905 this_block = get_block(n,n)
00906
00907
00908
00909
00910 if (this_block%iblock > 1 .and. &
00911 this_block%iblock < nblocks_x .and. &
00912 this_block%jblock > 1 .and. &
00913 this_block%jblock < nblocks_y) then
00914
00915 do j=1,ny_block
00916 do i=1,nx_block
00917 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
00918 this_block%j_glob(j))
00919 end do
00920 end do
00921
00922
00923
00924
00925
00926 else if (this_block%jblock /= nblocks_y) then
00927
00928 do j=1,ny_block
00929 if (this_block%j_glob(j) /= 0) then
00930 do i=1,nx_block
00931 if (this_block%i_glob(i) /= 0) then
00932 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
00933 this_block%j_glob(j))
00934 endif
00935 end do
00936 endif
00937 end do
00938
00939
00940
00941
00942 else
00943
00944 do j=1,ny_block
00945 if (this_block%j_glob(j) > 0) then
00946
00947 do i=1,nx_block
00948 if (this_block%i_glob(i) /= 0) then
00949 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
00950 this_block%j_glob(j))
00951 endif
00952 end do
00953
00954 else if (this_block%j_glob(j) < 0) then
00955
00956
00957
00958 do yoffset2=0,max(yoffset,0)-yoffset
00959 jsrc = ny_global + yoffset + yoffset2 + &
00960 (this_block%j_glob(j) + ny_global)
00961 do i=1,nx_block
00962 if (this_block%i_glob(i) /= 0) then
00963 isrc = nx_global + xoffset - this_block%i_glob(i)
00964 if (isrc < 1) isrc = isrc + nx_global
00965 if (isrc > nx_global) isrc = isrc - nx_global
00966 msg_buffer(i,j-yoffset2) = isign * ARRAY_G(isrc,jsrc)
00967 endif
00968 end do
00969 end do
00970
00971 endif
00972 end do
00973
00974 endif
00975
00976 call MPI_SEND(msg_buffer, nx_block*ny_block, &
00977 mpiR8, dst_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
00978 MPI_COMM_ICE, status, ierr)
00979
00980 endif
00981 end do
00982
00983 deallocate(msg_buffer)
00984
00985
00986
00987 do n=1,nblocks_tot
00988 if (dst_dist%blockLocation(n) == my_task+1) then
00989 dst_block = dst_dist%blockLocalID(n)
00990 this_block = get_block(n,n)
00991
00992
00993
00994
00995 if (this_block%iblock > 1 .and. &
00996 this_block%iblock < nblocks_x .and. &
00997 this_block%jblock > 1 .and. &
00998 this_block%jblock < nblocks_y) then
00999
01000 do j=1,ny_block
01001 do i=1,nx_block
01002 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01003 this_block%j_glob(j))
01004 end do
01005 end do
01006
01007
01008
01009
01010
01011 else if (this_block%jblock /= nblocks_y) then
01012
01013 do j=1,ny_block
01014 if (this_block%j_glob(j) /= 0) then
01015 do i=1,nx_block
01016 if (this_block%i_glob(i) /= 0) then
01017 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01018 this_block%j_glob(j))
01019 endif
01020 end do
01021 endif
01022 end do
01023
01024
01025
01026
01027 else
01028
01029 do j=1,ny_block
01030 if (this_block%j_glob(j) > 0) then
01031
01032 do i=1,nx_block
01033 if (this_block%i_glob(i) /= 0) then
01034 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01035 this_block%j_glob(j))
01036 endif
01037 end do
01038
01039 else if (this_block%j_glob(j) < 0) then
01040
01041
01042
01043 do yoffset2=0,max(yoffset,0)-yoffset
01044 jsrc = ny_global + yoffset + yoffset2 + &
01045 (this_block%j_glob(j) + ny_global)
01046 do i=1,nx_block
01047 if (this_block%i_glob(i) /= 0) then
01048 isrc = nx_global + xoffset - this_block%i_glob(i)
01049 if (isrc < 1) isrc = isrc + nx_global
01050 if (isrc > nx_global) isrc = isrc - nx_global
01051 ARRAY(i,j-yoffset2,dst_block) &
01052 = isign * ARRAY_G(isrc,jsrc)
01053 endif
01054 end do
01055 end do
01056
01057 endif
01058 end do
01059
01060 endif
01061 endif
01062 end do
01063
01064
01065
01066
01067
01068
01069
01070 else
01071
01072 allocate (rcv_request(nblocks_tot), &
01073 rcv_status(MPI_STATUS_SIZE, nblocks_tot))
01074
01075 rcv_request = 0
01076 rcv_status = 0
01077
01078 nrecvs = 0
01079 do n=1,nblocks_tot
01080 if (dst_dist%blockLocation(n) == my_task+1) then
01081 nrecvs = nrecvs + 1
01082 dst_block = dst_dist%blockLocalID(n)
01083 call MPI_IRECV(ARRAY(1,1,dst_block), nx_block*ny_block, &
01084 mpiR8, src_task, 3*mpitag_gs+n, &
01085 MPI_COMM_ICE, rcv_request(nrecvs), ierr)
01086 endif
01087 end do
01088
01089 if (nrecvs > 0) &
01090 call MPI_WAITALL(nrecvs, rcv_request, rcv_status, ierr)
01091
01092 deallocate(rcv_request, rcv_status)
01093 endif
01094
01095
01096
01097
01098
01099 if (field_loc == field_loc_noupdate) then
01100 do n=1,nblocks_tot
01101 dst_block = dst_dist%blockLocalID(n)
01102 this_block = get_block(n,n)
01103
01104 if (dst_block > 0) then
01105
01106
01107 do j = this_block%jhi+1,ny_block
01108 do i = 1, nx_block
01109 ARRAY (i,j,dst_block) = c0
01110 enddo
01111 enddo
01112
01113 do j = 1, ny_block
01114 do i = this_block%ihi+1,nx_block
01115 ARRAY (i,j,dst_block) = c0
01116 enddo
01117 enddo
01118
01119 do j = 1, this_block%jlo-1
01120 do i = 1, nx_block
01121 ARRAY (i,j,dst_block) = c0
01122 enddo
01123 enddo
01124
01125 do j = 1, ny_block
01126 do i = 1, this_block%ilo-1
01127 ARRAY (i,j,dst_block) = c0
01128 enddo
01129 enddo
01130
01131 endif
01132 enddo
01133 endif
01134
01135
01136
01137 end subroutine scatter_global_dbl
01138
01139
01140
01141 subroutine scatter_global_real(ARRAY, ARRAY_G, src_task, dst_dist, &
01142 field_loc, field_type)
01143
01144
01145
01146
01147
01148
01149
01150 include 'mpif.h'
01151
01152
01153
01154
01155
01156
01157
01158 integer (int_kind), intent(in) ::
01159 src_task
01160
01161 type (distrb), intent(in) ::
01162 dst_dist
01163
01164 real (real_kind), dimension(:,:), intent(in) ::
01165 ARRAY_G
01166
01167 integer (int_kind), intent(in) ::
01168 field_type,
01169 field_loc
01170
01171
01172
01173
01174
01175
01176
01177
01178 real (real_kind), dimension(:,:,:), intent(inout) ::
01179 ARRAY
01180
01181
01182
01183
01184
01185
01186
01187 integer (int_kind) ::
01188 i,j,n,bid,
01189 nrecvs,
01190 isrc, jsrc,
01191 dst_block,
01192 xoffset, yoffset,
01193 yoffset2,
01194 isign,
01195 ierr
01196
01197 type (block) ::
01198 this_block
01199
01200 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
01201 status
01202
01203 integer (int_kind), dimension(:), allocatable ::
01204 rcv_request
01205
01206 integer (int_kind), dimension(:,:), allocatable ::
01207 rcv_status
01208
01209 real (real_kind), dimension(:,:), allocatable ::
01210 msg_buffer
01211
01212
01213
01214
01215
01216
01217
01218 ARRAY = 0._real_kind
01219
01220 this_block = get_block(1,1)
01221 if (this_block%tripoleTFlag) then
01222 select case (field_loc)
01223 case (field_loc_center)
01224 xoffset = 2
01225 yoffset = 0
01226 case (field_loc_NEcorner)
01227 xoffset = 1
01228 yoffset = 1
01229 case (field_loc_Eface)
01230 xoffset = 1
01231 yoffset = 0
01232 case (field_loc_Nface)
01233 xoffset = 2
01234 yoffset = 1
01235 case (field_loc_noupdate)
01236 xoffset = 1
01237 yoffset = 1
01238 end select
01239 else
01240 select case (field_loc)
01241 case (field_loc_center)
01242 xoffset = 1
01243 yoffset = 1
01244 case (field_loc_NEcorner)
01245 xoffset = 0
01246 yoffset = 0
01247 case (field_loc_Eface)
01248 xoffset = 0
01249 yoffset = 1
01250 case (field_loc_Nface)
01251 xoffset = 1
01252 yoffset = 0
01253 case (field_loc_noupdate)
01254 xoffset = 1
01255 yoffset = 1
01256 end select
01257 endif
01258
01259 select case (field_type)
01260 case (field_type_scalar)
01261 isign = 1
01262 case (field_type_vector)
01263 isign = -1
01264 case (field_type_angle)
01265 isign = -1
01266 case (field_type_noupdate)
01267 isign = 1
01268 case default
01269 call abort_ice('Unknown field type in scatter')
01270 end select
01271
01272
01273
01274
01275
01276
01277
01278
01279 if (my_task == src_task) then
01280
01281
01282
01283 allocate (msg_buffer(nx_block,ny_block))
01284
01285 do n=1,nblocks_tot
01286 if (dst_dist%blockLocation(n) > 0 .and. &
01287 dst_dist%blockLocation(n)-1 /= my_task) then
01288
01289 msg_buffer = 0._real_kind
01290 this_block = get_block(n,n)
01291
01292
01293
01294
01295 if (this_block%iblock > 1 .and. &
01296 this_block%iblock < nblocks_x .and. &
01297 this_block%jblock > 1 .and. &
01298 this_block%jblock < nblocks_y) then
01299
01300 do j=1,ny_block
01301 do i=1,nx_block
01302 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01303 this_block%j_glob(j))
01304 end do
01305 end do
01306
01307
01308
01309
01310
01311 else if (this_block%jblock /= nblocks_y) then
01312
01313 do j=1,ny_block
01314 if (this_block%j_glob(j) /= 0) then
01315 do i=1,nx_block
01316 if (this_block%i_glob(i) /= 0) then
01317 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01318 this_block%j_glob(j))
01319 endif
01320 end do
01321 endif
01322 end do
01323
01324
01325
01326
01327 else
01328
01329 do j=1,ny_block
01330 if (this_block%j_glob(j) > 0) then
01331
01332 do i=1,nx_block
01333 if (this_block%i_glob(i) /= 0) then
01334 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01335 this_block%j_glob(j))
01336 endif
01337 end do
01338
01339 else if (this_block%j_glob(j) < 0) then
01340
01341
01342
01343 do yoffset2=0,max(yoffset,0)-yoffset
01344 jsrc = ny_global + yoffset + yoffset2 + &
01345 (this_block%j_glob(j) + ny_global)
01346 do i=1,nx_block
01347 if (this_block%i_glob(i) /= 0) then
01348 isrc = nx_global + xoffset - this_block%i_glob(i)
01349 if (isrc < 1) isrc = isrc + nx_global
01350 if (isrc > nx_global) isrc = isrc - nx_global
01351 msg_buffer(i,j-yoffset2) = isign * ARRAY_G(isrc,jsrc)
01352 endif
01353 end do
01354 end do
01355
01356 endif
01357 end do
01358
01359 endif
01360
01361 call MPI_SEND(msg_buffer, nx_block*ny_block, &
01362 mpiR4, dst_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
01363 MPI_COMM_ICE, status, ierr)
01364
01365 endif
01366 end do
01367
01368 deallocate(msg_buffer)
01369
01370
01371
01372 do n=1,nblocks_tot
01373 if (dst_dist%blockLocation(n) == my_task+1) then
01374 dst_block = dst_dist%blockLocalID(n)
01375 this_block = get_block(n,n)
01376
01377
01378
01379
01380 if (this_block%iblock > 1 .and. &
01381 this_block%iblock < nblocks_x .and. &
01382 this_block%jblock > 1 .and. &
01383 this_block%jblock < nblocks_y) then
01384
01385 do j=1,ny_block
01386 do i=1,nx_block
01387 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01388 this_block%j_glob(j))
01389 end do
01390 end do
01391
01392
01393
01394
01395
01396 else if (this_block%jblock /= nblocks_y) then
01397
01398 do j=1,ny_block
01399 if (this_block%j_glob(j) /= 0) then
01400 do i=1,nx_block
01401 if (this_block%i_glob(i) /= 0) then
01402 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01403 this_block%j_glob(j))
01404 endif
01405 end do
01406 endif
01407 end do
01408
01409
01410
01411
01412 else
01413
01414 do j=1,ny_block
01415 if (this_block%j_glob(j) > 0) then
01416
01417 do i=1,nx_block
01418 if (this_block%i_glob(i) /= 0) then
01419 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01420 this_block%j_glob(j))
01421 endif
01422 end do
01423
01424 else if (this_block%j_glob(j) < 0) then
01425
01426
01427
01428 do yoffset2=0,max(yoffset,0)-yoffset
01429 jsrc = ny_global + yoffset + yoffset2 + &
01430 (this_block%j_glob(j) + ny_global)
01431 do i=1,nx_block
01432 if (this_block%i_glob(i) /= 0) then
01433 isrc = nx_global + xoffset - this_block%i_glob(i)
01434 if (isrc < 1) isrc = isrc + nx_global
01435 if (isrc > nx_global) isrc = isrc - nx_global
01436 ARRAY(i,j-yoffset2,dst_block) &
01437 = isign * ARRAY_G(isrc,jsrc)
01438 endif
01439 end do
01440 end do
01441
01442 endif
01443 end do
01444
01445 endif
01446 endif
01447 end do
01448
01449
01450
01451
01452
01453
01454
01455 else
01456
01457 allocate (rcv_request(nblocks_tot), &
01458 rcv_status(MPI_STATUS_SIZE, nblocks_tot))
01459
01460 rcv_request = 0
01461 rcv_status = 0
01462
01463 nrecvs = 0
01464 do n=1,nblocks_tot
01465 if (dst_dist%blockLocation(n) == my_task+1) then
01466 nrecvs = nrecvs + 1
01467 dst_block = dst_dist%blockLocalID(n)
01468 call MPI_IRECV(ARRAY(1,1,dst_block), nx_block*ny_block, &
01469 mpiR4, src_task, 3*mpitag_gs+n, &
01470 MPI_COMM_ICE, rcv_request(nrecvs), ierr)
01471 endif
01472 end do
01473
01474 if (nrecvs > 0) &
01475 call MPI_WAITALL(nrecvs, rcv_request, rcv_status, ierr)
01476
01477 deallocate(rcv_request, rcv_status)
01478 endif
01479
01480
01481
01482
01483
01484 if (field_loc == field_loc_noupdate) then
01485 do n=1,nblocks_tot
01486 dst_block = dst_dist%blockLocalID(n)
01487 this_block = get_block(n,n)
01488
01489 if (dst_block > 0) then
01490
01491
01492 do j = this_block%jhi+1,ny_block
01493 do i = 1, nx_block
01494 ARRAY (i,j,dst_block) = 0._real_kind
01495 enddo
01496 enddo
01497
01498 do j = 1, ny_block
01499 do i = this_block%ihi+1,nx_block
01500 ARRAY (i,j,dst_block) = 0._real_kind
01501 enddo
01502 enddo
01503
01504 do j = 1, this_block%jlo-1
01505 do i = 1, nx_block
01506 ARRAY (i,j,dst_block) = 0._real_kind
01507 enddo
01508 enddo
01509
01510 do j = 1, ny_block
01511 do i = 1, this_block%ilo-1
01512 ARRAY (i,j,dst_block) = 0._real_kind
01513 enddo
01514 enddo
01515
01516 endif
01517 enddo
01518 endif
01519
01520
01521
01522 end subroutine scatter_global_real
01523
01524
01525
01526 subroutine scatter_global_int(ARRAY, ARRAY_G, src_task, dst_dist, &
01527 field_loc, field_type)
01528
01529
01530
01531
01532
01533
01534
01535 include 'mpif.h'
01536
01537
01538
01539
01540
01541
01542
01543 integer (int_kind), intent(in) ::
01544 src_task
01545
01546 integer (int_kind), intent(in) ::
01547 field_type,
01548 field_loc
01549
01550
01551 type (distrb), intent(in) ::
01552 dst_dist
01553
01554 integer (int_kind), dimension(:,:), intent(in) ::
01555 ARRAY_G
01556
01557
01558
01559
01560
01561
01562
01563 integer (int_kind), dimension(:,:,:), intent(inout) ::
01564 ARRAY
01565
01566
01567
01568
01569
01570
01571
01572 integer (int_kind) ::
01573 i,j,n,bid,
01574 nrecvs,
01575 isrc, jsrc,
01576 dst_block,
01577 xoffset, yoffset,
01578 yoffset2,
01579 isign,
01580 ierr
01581
01582 type (block) ::
01583 this_block
01584
01585 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
01586 status
01587
01588 integer (int_kind), dimension(:), allocatable ::
01589 rcv_request
01590
01591 integer (int_kind), dimension(:,:), allocatable ::
01592 rcv_status
01593
01594 integer (int_kind), dimension(:,:), allocatable ::
01595 msg_buffer
01596
01597
01598
01599
01600
01601
01602
01603 ARRAY = 0
01604
01605 this_block = get_block(1,1)
01606 if (this_block%tripoleTFlag) then
01607 select case (field_loc)
01608 case (field_loc_center)
01609 xoffset = 2
01610 yoffset = 0
01611 case (field_loc_NEcorner)
01612 xoffset = 1
01613 yoffset = 1
01614 case (field_loc_Eface)
01615 xoffset = 1
01616 yoffset = 0
01617 case (field_loc_Nface)
01618 xoffset = 2
01619 yoffset = 1
01620 case (field_loc_noupdate)
01621 xoffset = 1
01622 yoffset = 1
01623 end select
01624 else
01625 select case (field_loc)
01626 case (field_loc_center)
01627 xoffset = 1
01628 yoffset = 1
01629 case (field_loc_NEcorner)
01630 xoffset = 0
01631 yoffset = 0
01632 case (field_loc_Eface)
01633 xoffset = 0
01634 yoffset = 1
01635 case (field_loc_Nface)
01636 xoffset = 1
01637 yoffset = 0
01638 case (field_loc_noupdate)
01639 xoffset = 1
01640 yoffset = 1
01641 end select
01642 endif
01643
01644 select case (field_type)
01645 case (field_type_scalar)
01646 isign = 1
01647 case (field_type_vector)
01648 isign = -1
01649 case (field_type_angle)
01650 isign = -1
01651 case (field_type_noupdate)
01652 isign = 1
01653 case default
01654 call abort_ice('Unknown field type in scatter')
01655 end select
01656
01657
01658
01659
01660
01661
01662
01663
01664 if (my_task == src_task) then
01665
01666
01667
01668 allocate (msg_buffer(nx_block,ny_block))
01669
01670 do n=1,nblocks_tot
01671 if (dst_dist%blockLocation(n) > 0 .and. &
01672 dst_dist%blockLocation(n)-1 /= my_task) then
01673
01674 msg_buffer = 0
01675 this_block = get_block(n,n)
01676
01677
01678
01679
01680 if (this_block%iblock > 1 .and. &
01681 this_block%iblock < nblocks_x .and. &
01682 this_block%jblock > 1 .and. &
01683 this_block%jblock < nblocks_y) then
01684
01685 do j=1,ny_block
01686 do i=1,nx_block
01687 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01688 this_block%j_glob(j))
01689 end do
01690 end do
01691
01692
01693
01694
01695
01696 else if (this_block%jblock /= nblocks_y) then
01697
01698 do j=1,ny_block
01699 if (this_block%j_glob(j) /= 0) then
01700 do i=1,nx_block
01701 if (this_block%i_glob(i) /= 0) then
01702 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01703 this_block%j_glob(j))
01704 endif
01705 end do
01706 endif
01707 end do
01708
01709
01710
01711
01712 else
01713
01714 do j=1,ny_block
01715 if (this_block%j_glob(j) > 0) then
01716
01717 do i=1,nx_block
01718 if (this_block%i_glob(i) /= 0) then
01719 msg_buffer(i,j) = ARRAY_G(this_block%i_glob(i),&
01720 this_block%j_glob(j))
01721 endif
01722 end do
01723
01724 else if (this_block%j_glob(j) < 0) then
01725
01726
01727
01728 do yoffset2=0,max(yoffset,0)-yoffset
01729 jsrc = ny_global + yoffset + yoffset2 + &
01730 (this_block%j_glob(j) + ny_global)
01731 do i=1,nx_block
01732 if (this_block%i_glob(i) /= 0) then
01733 isrc = nx_global + xoffset - this_block%i_glob(i)
01734 if (isrc < 1) isrc = isrc + nx_global
01735 if (isrc > nx_global) isrc = isrc - nx_global
01736 msg_buffer(i,j-yoffset2) = isign * ARRAY_G(isrc,jsrc)
01737 endif
01738 end do
01739 end do
01740
01741 endif
01742 end do
01743
01744 endif
01745
01746 call MPI_SEND(msg_buffer, nx_block*ny_block, &
01747 mpi_integer, dst_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
01748 MPI_COMM_ICE, status, ierr)
01749
01750 endif
01751 end do
01752
01753 deallocate(msg_buffer)
01754
01755
01756
01757 do n=1,nblocks_tot
01758 if (dst_dist%blockLocation(n) == my_task+1) then
01759 dst_block = dst_dist%blockLocalID(n)
01760 this_block = get_block(n,n)
01761
01762
01763
01764
01765 if (this_block%iblock > 1 .and. &
01766 this_block%iblock < nblocks_x .and. &
01767 this_block%jblock > 1 .and. &
01768 this_block%jblock < nblocks_y) then
01769
01770 do j=1,ny_block
01771 do i=1,nx_block
01772 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01773 this_block%j_glob(j))
01774 end do
01775 end do
01776
01777
01778
01779
01780
01781 else if (this_block%jblock /= nblocks_y) then
01782
01783 do j=1,ny_block
01784 if (this_block%j_glob(j) /= 0) then
01785 do i=1,nx_block
01786 if (this_block%i_glob(i) /= 0) then
01787 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01788 this_block%j_glob(j))
01789 endif
01790 end do
01791 endif
01792 end do
01793
01794
01795
01796
01797 else
01798
01799 do j=1,ny_block
01800 if (this_block%j_glob(j) > 0) then
01801
01802 do i=1,nx_block
01803 if (this_block%i_glob(i) /= 0) then
01804 ARRAY(i,j,dst_block) = ARRAY_G(this_block%i_glob(i),&
01805 this_block%j_glob(j))
01806 endif
01807 end do
01808
01809 else if (this_block%j_glob(j) < 0) then
01810
01811
01812
01813 do yoffset2=0,max(yoffset,0)-yoffset
01814 jsrc = ny_global + yoffset + yoffset2 + &
01815 (this_block%j_glob(j) + ny_global)
01816 do i=1,nx_block
01817 if (this_block%i_glob(i) /= 0) then
01818 isrc = nx_global + xoffset - this_block%i_glob(i)
01819 if (isrc < 1) isrc = isrc + nx_global
01820 if (isrc > nx_global) isrc = isrc - nx_global
01821 ARRAY(i,j-yoffset2,dst_block) &
01822 = isign * ARRAY_G(isrc,jsrc)
01823 endif
01824 end do
01825 end do
01826
01827 endif
01828 end do
01829
01830 endif
01831 endif
01832 end do
01833
01834
01835
01836
01837
01838
01839
01840 else
01841
01842 allocate (rcv_request(nblocks_tot), &
01843 rcv_status(MPI_STATUS_SIZE, nblocks_tot))
01844
01845 rcv_request = 0
01846 rcv_status = 0
01847
01848 nrecvs = 0
01849 do n=1,nblocks_tot
01850 if (dst_dist%blockLocation(n) == my_task+1) then
01851 nrecvs = nrecvs + 1
01852 dst_block = dst_dist%blockLocalID(n)
01853 call MPI_IRECV(ARRAY(1,1,dst_block), nx_block*ny_block, &
01854 mpi_integer, src_task, 3*mpitag_gs+n, &
01855 MPI_COMM_ICE, rcv_request(nrecvs), ierr)
01856 endif
01857 end do
01858
01859 if (nrecvs > 0) &
01860 call MPI_WAITALL(nrecvs, rcv_request, rcv_status, ierr)
01861
01862 deallocate(rcv_request, rcv_status)
01863 endif
01864
01865
01866
01867
01868
01869 if (field_loc == field_loc_noupdate) then
01870 do n=1,nblocks_tot
01871 dst_block = dst_dist%blockLocalID(n)
01872 this_block = get_block(n,n)
01873
01874 if (dst_block > 0) then
01875
01876
01877 do j = this_block%jhi+1,ny_block
01878 do i = 1, nx_block
01879 ARRAY (i,j,dst_block) = 0
01880 enddo
01881 enddo
01882
01883 do j = 1, ny_block
01884 do i = this_block%ihi+1,nx_block
01885 ARRAY (i,j,dst_block) = 0
01886 enddo
01887 enddo
01888
01889 do j = 1, this_block%jlo-1
01890 do i = 1, nx_block
01891 ARRAY (i,j,dst_block) = 0
01892 enddo
01893 enddo
01894
01895 do j = 1, ny_block
01896 do i = 1, this_block%ilo-1
01897 ARRAY (i,j,dst_block) = 0
01898 enddo
01899 enddo
01900
01901 endif
01902 enddo
01903 endif
01904
01905
01906
01907 end subroutine scatter_global_int
01908
01909
01910
01911
01912
01913
01914
01915 subroutine scatter_global_stress(ARRAY, ARRAY_G1, ARRAY_G2, &
01916 src_task, dst_dist)
01917
01918
01919
01920
01921
01922
01923
01924
01925
01926
01927
01928
01929
01930 include 'mpif.h'
01931
01932
01933
01934 integer (int_kind), intent(in) ::
01935 src_task
01936
01937 type (distrb), intent(in) ::
01938 dst_dist
01939
01940 real (dbl_kind), dimension(:,:), intent(in) ::
01941 ARRAY_G1,
01942 ARRAY_G2
01943
01944
01945
01946 real (dbl_kind), dimension(:,:,:), intent(inout) ::
01947 ARRAY
01948
01949
01950
01951
01952
01953
01954
01955
01956
01957 integer (int_kind) ::
01958 i,j,n,bid,
01959 nrecvs,
01960 isrc, jsrc,
01961 dst_block,
01962 xoffset, yoffset,
01963 yoffset2,
01964 isign,
01965 ierr
01966
01967 type (block) ::
01968 this_block
01969
01970 integer (int_kind), dimension(MPI_STATUS_SIZE) ::
01971 status
01972
01973 integer (int_kind), dimension(:), allocatable ::
01974 rcv_request
01975
01976 integer (int_kind), dimension(:,:), allocatable ::
01977 rcv_status
01978
01979 real (dbl_kind), dimension(:,:), allocatable ::
01980 msg_buffer
01981
01982
01983
01984
01985
01986
01987
01988 ARRAY = c0
01989
01990 this_block = get_block(1,1)
01991 if (this_block%tripoleTFlag) then
01992 xoffset = 2
01993 yoffset = 0
01994 else
01995 xoffset = 1
01996 yoffset = 1
01997 endif
01998 isign = 1
01999
02000
02001
02002
02003
02004
02005
02006
02007 if (my_task == src_task) then
02008
02009
02010
02011 allocate (msg_buffer(nx_block,ny_block))
02012
02013 do n=1,nblocks_tot
02014 if (dst_dist%blockLocation(n) > 0 .and. &
02015 dst_dist%blockLocation(n)-1 /= my_task) then
02016
02017 msg_buffer = c0
02018 this_block = get_block(n,n)
02019
02020
02021
02022
02023 if (this_block%iblock > 1 .and. &
02024 this_block%iblock < nblocks_x .and. &
02025 this_block%jblock > 1 .and. &
02026 this_block%jblock < nblocks_y) then
02027
02028 do j=1,ny_block
02029 do i=1,nx_block
02030 msg_buffer(i,j) = ARRAY_G1(this_block%i_glob(i),&
02031 this_block%j_glob(j))
02032 end do
02033 end do
02034
02035
02036
02037
02038
02039 else if (this_block%jblock /= nblocks_y) then
02040
02041 do j=1,ny_block
02042 if (this_block%j_glob(j) /= 0) then
02043 do i=1,nx_block
02044 if (this_block%i_glob(i) /= 0) then
02045 msg_buffer(i,j) = ARRAY_G1(this_block%i_glob(i),&
02046 this_block%j_glob(j))
02047 endif
02048 end do
02049 endif
02050 end do
02051
02052
02053
02054
02055 else
02056
02057 do j=1,ny_block
02058 if (this_block%j_glob(j) > 0) then
02059
02060 do i=1,nx_block
02061 if (this_block%i_glob(i) /= 0) then
02062 msg_buffer(i,j) = ARRAY_G1(this_block%i_glob(i),&
02063 this_block%j_glob(j))
02064 endif
02065 end do
02066
02067 else if (this_block%j_glob(j) < 0) then
02068
02069 jsrc = ny_global + yoffset + &
02070 (this_block%j_glob(j) + ny_global)
02071 do i=1,nx_block
02072 if (this_block%i_glob(i) /= 0) then
02073 isrc = nx_global + xoffset - this_block%i_glob(i)
02074 if (isrc < 1) isrc = isrc + nx_global
02075 if (isrc > nx_global) isrc = isrc - nx_global
02076 msg_buffer(i,j) = isign * ARRAY_G2(isrc,jsrc)
02077 endif
02078 end do
02079
02080 endif
02081 end do
02082
02083 endif
02084
02085 call MPI_SEND(msg_buffer, nx_block*ny_block, &
02086 mpiR8, dst_dist%blockLocation(n)-1, 3*mpitag_gs+n, &
02087 MPI_COMM_ICE, status, ierr)
02088
02089 endif
02090 end do
02091
02092 deallocate(msg_buffer)
02093
02094
02095
02096 do n=1,nblocks_tot
02097 if (dst_dist%blockLocation(n) == my_task+1) then
02098 dst_block = dst_dist%blockLocalID(n)
02099 this_block = get_block(n,n)
02100
02101
02102
02103
02104 if (this_block%iblock > 1 .and. &
02105 this_block%iblock < nblocks_x .and. &
02106 this_block%jblock > 1 .and. &
02107 this_block%jblock < nblocks_y) then
02108
02109 do j=1,ny_block
02110 do i=1,nx_block
02111 ARRAY(i,j,dst_block) = ARRAY_G1(this_block%i_glob(i),&
02112 this_block%j_glob(j))
02113 end do
02114 end do
02115
02116
02117
02118
02119
02120 else if (this_block%jblock /= nblocks_y) then
02121
02122 do j=1,ny_block
02123 if (this_block%j_glob(j) /= 0) then
02124 do i=1,nx_block
02125 if (this_block%i_glob(i) /= 0) then
02126 ARRAY(i,j,dst_block) = ARRAY_G1(this_block%i_glob(i),&
02127 this_block%j_glob(j))
02128 endif
02129 end do
02130 endif
02131 end do
02132
02133
02134
02135
02136 else
02137
02138 do j=1,ny_block
02139 if (this_block%j_glob(j) > 0) then
02140
02141 do i=1,nx_block
02142 if (this_block%i_glob(i) /= 0) then
02143 ARRAY(i,j,dst_block) = ARRAY_G1(this_block%i_glob(i),&
02144 this_block%j_glob(j))
02145 endif
02146 end do
02147
02148 else if (this_block%j_glob(j) < 0) then
02149
02150
02151
02152 do yoffset2=0,max(yoffset,0)-yoffset
02153 jsrc = ny_global + yoffset + yoffset2 + &
02154 (this_block%j_glob(j) + ny_global)
02155 do i=1,nx_block
02156 if (this_block%i_glob(i) /= 0) then
02157 isrc = nx_global + xoffset - this_block%i_glob(i)
02158 if (isrc < 1) isrc = isrc + nx_global
02159 if (isrc > nx_global) isrc = isrc - nx_global
02160 ARRAY(i,j-yoffset2,dst_block) &
02161 = isign * ARRAY_G2(isrc,jsrc)
02162 endif
02163 end do
02164 end do
02165
02166 endif
02167 end do
02168
02169 endif
02170 endif
02171 end do
02172
02173
02174
02175
02176
02177
02178
02179 else
02180
02181 allocate (rcv_request(nblocks_tot), &
02182 rcv_status(MPI_STATUS_SIZE, nblocks_tot))
02183
02184 rcv_request = 0
02185 rcv_status = 0
02186
02187 nrecvs = 0
02188 do n=1,nblocks_tot
02189 if (dst_dist%blockLocation(n) == my_task+1) then
02190 nrecvs = nrecvs + 1
02191 dst_block = dst_dist%blockLocalID(n)
02192 call MPI_IRECV(ARRAY(1,1,dst_block), nx_block*ny_block, &
02193 mpiR8, src_task, 3*mpitag_gs+n, &
02194 MPI_COMM_ICE, rcv_request(nrecvs), ierr)
02195 endif
02196 end do
02197
02198 if (nrecvs > 0) &
02199 call MPI_WAITALL(nrecvs, rcv_request, rcv_status, ierr)
02200
02201 deallocate(rcv_request, rcv_status)
02202 endif
02203
02204
02205
02206 end subroutine scatter_global_stress
02207
02208
02209
02210 end module ice_gather_scatter
02211
02212