[26] | 1 | MODULE transfert_mpi_mod |
---|
| 2 | USE genmod |
---|
| 3 | |
---|
| 4 | TYPE array |
---|
| 5 | INTEGER,POINTER :: value(:) |
---|
| 6 | INTEGER :: domain |
---|
| 7 | INTEGER :: rank |
---|
| 8 | INTEGER :: size |
---|
| 9 | INTEGER,POINTER :: buffer(:) |
---|
| 10 | REAL,POINTER :: buffer_r2(:) |
---|
| 11 | REAL,POINTER :: buffer_r3(:,:) |
---|
| 12 | REAL,POINTER :: buffer_r4(:,:,:) |
---|
| 13 | END TYPE array |
---|
| 14 | |
---|
| 15 | TYPE t_request |
---|
| 16 | INTEGER :: type_field |
---|
| 17 | INTEGER :: max_size |
---|
| 18 | INTEGER :: size |
---|
| 19 | INTEGER,POINTER :: src_domain(:) |
---|
| 20 | INTEGER,POINTER :: src_i(:) |
---|
| 21 | INTEGER,POINTER :: src_j(:) |
---|
| 22 | INTEGER,POINTER :: src_ind(:) |
---|
| 23 | INTEGER,POINTER :: target_ind(:) |
---|
| 24 | INTEGER,POINTER :: target_i(:) |
---|
| 25 | INTEGER,POINTER :: target_j(:) |
---|
| 26 | INTEGER :: nrecv |
---|
| 27 | TYPE(ARRAY),POINTER :: recv(:) |
---|
| 28 | INTEGER :: nsend |
---|
| 29 | TYPE(ARRAY),POINTER :: send(:) |
---|
| 30 | END TYPE t_request |
---|
| 31 | |
---|
| 32 | TYPE(t_request),POINTER :: req_i1(:) |
---|
| 33 | TYPE(t_request),POINTER :: req_e1(:) |
---|
| 34 | |
---|
| 35 | |
---|
| 36 | CONTAINS |
---|
| 37 | |
---|
| 38 | SUBROUTINE init_transfert |
---|
| 39 | USE domain_mod |
---|
| 40 | USE dimensions |
---|
| 41 | USE field_mod |
---|
| 42 | USE metric |
---|
| 43 | USE mpipara |
---|
| 44 | IMPLICIT NONE |
---|
| 45 | INTEGER :: ind,i,j |
---|
| 46 | |
---|
| 47 | CALL create_request(field_t,req_i1) |
---|
| 48 | |
---|
| 49 | DO ind=1,ndomain |
---|
| 50 | CALL swap_dimensions(ind) |
---|
| 51 | DO i=ii_begin,ii_end+1 |
---|
| 52 | CALL request_add_point(ind,i,jj_begin-1,req_i1) |
---|
| 53 | ENDDO |
---|
| 54 | |
---|
| 55 | DO j=jj_begin,jj_end |
---|
| 56 | CALL request_add_point(ind,ii_end+1,j,req_i1) |
---|
| 57 | ENDDO |
---|
| 58 | DO i=ii_begin,ii_end |
---|
| 59 | CALL request_add_point(ind,i,jj_end+1,req_i1) |
---|
| 60 | ENDDO |
---|
| 61 | |
---|
| 62 | DO j=jj_begin,jj_end+1 |
---|
| 63 | CALL request_add_point(ind,ii_begin-1,j,req_i1) |
---|
| 64 | ENDDO |
---|
| 65 | |
---|
| 66 | DO i=ii_begin,ii_end |
---|
| 67 | CALL request_add_point(ind,i,jj_begin,req_i1) |
---|
| 68 | ENDDO |
---|
| 69 | |
---|
| 70 | DO j=jj_begin,jj_end |
---|
| 71 | CALL request_add_point(ind,ii_end,j,req_i1) |
---|
| 72 | ENDDO |
---|
| 73 | |
---|
| 74 | DO i=ii_begin,ii_end |
---|
| 75 | CALL request_add_point(ind,i,jj_end,req_i1) |
---|
| 76 | ENDDO |
---|
| 77 | |
---|
| 78 | DO j=jj_begin,jj_end |
---|
| 79 | CALL request_add_point(ind,ii_begin,j,req_i1) |
---|
| 80 | ENDDO |
---|
| 81 | |
---|
| 82 | ENDDO |
---|
| 83 | |
---|
| 84 | CALL finalize_request(req_i1) |
---|
| 85 | |
---|
| 86 | CALL create_request(field_u,req_e1) |
---|
| 87 | DO ind=1,ndomain |
---|
| 88 | CALL swap_dimensions(ind) |
---|
| 89 | DO i=ii_begin,ii_end |
---|
| 90 | CALL request_add_point(ind,i,jj_begin-1,req_e1,rup) |
---|
| 91 | CALL request_add_point(ind,i+1,jj_begin-1,req_e1,lup) |
---|
| 92 | ENDDO |
---|
| 93 | |
---|
| 94 | DO j=jj_begin,jj_end |
---|
| 95 | CALL request_add_point(ind,ii_end+1,j,req_e1,left) |
---|
| 96 | CALL request_add_point(ind,ii_end+1,j-1,req_e1,lup) |
---|
| 97 | ENDDO |
---|
| 98 | |
---|
| 99 | DO i=ii_begin,ii_end |
---|
| 100 | CALL request_add_point(ind,i,jj_end+1,req_e1,ldown) |
---|
| 101 | CALL request_add_point(ind,i-1,jj_end+1,req_e1,rdown) |
---|
| 102 | ENDDO |
---|
| 103 | |
---|
| 104 | DO j=jj_begin,jj_end |
---|
| 105 | CALL request_add_point(ind,ii_begin-1,j,req_e1,right) |
---|
| 106 | CALL request_add_point(ind,ii_begin-1,j+1,req_e1,rdown) |
---|
| 107 | ENDDO |
---|
| 108 | |
---|
| 109 | DO i=ii_begin+1,ii_end-1 |
---|
| 110 | CALL request_add_point(ind,i,jj_begin,req_e1,right) |
---|
| 111 | CALL request_add_point(ind,i,jj_end,req_e1,right) |
---|
| 112 | ENDDO |
---|
| 113 | |
---|
| 114 | DO j=jj_begin+1,jj_end-1 |
---|
| 115 | CALL request_add_point(ind,ii_begin,j,req_e1,rup) |
---|
| 116 | CALL request_add_point(ind,ii_end,j,req_e1,rup) |
---|
| 117 | ENDDO |
---|
| 118 | |
---|
| 119 | CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1,left) |
---|
| 120 | CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1,ldown) |
---|
| 121 | CALL request_add_point(ind,ii_begin+1,jj_end,req_e1,left) |
---|
| 122 | CALL request_add_point(ind,ii_end,jj_begin+1,req_e1,ldown) |
---|
| 123 | |
---|
| 124 | CALL finalize_request(req_e1) |
---|
| 125 | |
---|
| 126 | ENDDO |
---|
| 127 | |
---|
| 128 | END SUBROUTINE init_transfert |
---|
| 129 | |
---|
| 130 | SUBROUTINE create_request(type_field,request) |
---|
| 131 | USE domain_mod |
---|
| 132 | USE field_mod |
---|
| 133 | IMPLICIT NONE |
---|
| 134 | INTEGER :: type_field |
---|
| 135 | TYPE(t_request),POINTER :: request(:) |
---|
| 136 | TYPE(t_request),POINTER :: req |
---|
| 137 | TYPE(t_domain),POINTER :: d |
---|
| 138 | INTEGER :: ind |
---|
| 139 | INTEGER :: max_size |
---|
| 140 | |
---|
| 141 | ALLOCATE(request(ndomain)) |
---|
| 142 | |
---|
| 143 | DO ind=1,ndomain |
---|
| 144 | req=>request(ind) |
---|
| 145 | d=>domain(ind) |
---|
| 146 | IF (type_field==field_t) THEN |
---|
| 147 | Max_size=2*(d%iim+2)+2*(d%jjm+2) |
---|
| 148 | ELSE IF (type_field==field_u) THEN |
---|
| 149 | Max_size=3*(2*(d%iim+2)+2*(d%jjm+2)) |
---|
| 150 | ELSE IF (type_field==field_z) THEN |
---|
| 151 | Max_size=2*(2*(d%iim+2)+2*(d%jjm+2)) |
---|
| 152 | ENDIF |
---|
| 153 | |
---|
| 154 | req%type_field=type_field |
---|
| 155 | req%max_size=max_size*2 |
---|
| 156 | req%size=0 |
---|
| 157 | ALLOCATE(req%src_domain(req%max_size)) |
---|
| 158 | ALLOCATE(req%src_ind(req%max_size)) |
---|
| 159 | ALLOCATE(req%target_ind(req%max_size)) |
---|
| 160 | ALLOCATE(req%src_i(req%max_size)) |
---|
| 161 | ALLOCATE(req%src_j(req%max_size)) |
---|
| 162 | ALLOCATE(req%target_i(req%max_size)) |
---|
| 163 | ALLOCATE(req%target_j(req%max_size)) |
---|
| 164 | ENDDO |
---|
| 165 | |
---|
| 166 | END SUBROUTINE create_request |
---|
| 167 | |
---|
| 168 | SUBROUTINE reallocate_request(req) |
---|
| 169 | IMPLICIT NONE |
---|
| 170 | TYPE(t_request),POINTER :: req |
---|
| 171 | |
---|
| 172 | INTEGER,POINTER :: src_domain(:) |
---|
| 173 | INTEGER,POINTER :: src_ind(:) |
---|
| 174 | INTEGER,POINTER :: target_ind(:) |
---|
| 175 | INTEGER,POINTER :: src_i(:) |
---|
| 176 | INTEGER,POINTER :: src_j(:) |
---|
| 177 | INTEGER,POINTER :: target_i(:) |
---|
| 178 | INTEGER,POINTER :: target_j(:) |
---|
| 179 | |
---|
| 180 | PRINT *,"REALLOCATE_REQUEST" |
---|
| 181 | src_domain=>req%src_domain |
---|
| 182 | src_ind=>req%src_ind |
---|
| 183 | target_ind=>req%target_ind |
---|
| 184 | src_i=>req%src_i |
---|
| 185 | src_j=>req%src_j |
---|
| 186 | target_i=>req%target_i |
---|
| 187 | target_j=>req%target_j |
---|
| 188 | ! req%max_size=req%max_size*2 |
---|
| 189 | ALLOCATE(req%src_domain(req%max_size*2)) |
---|
| 190 | ALLOCATE(req%src_ind(req%max_size*2)) |
---|
| 191 | ALLOCATE(req%target_ind(req%max_size*2)) |
---|
| 192 | ALLOCATE(req%src_i(req%max_size*2)) |
---|
| 193 | ALLOCATE(req%src_j(req%max_size*2)) |
---|
| 194 | ALLOCATE(req%target_i(req%max_size*2)) |
---|
| 195 | ALLOCATE(req%target_j(req%max_size*2)) |
---|
| 196 | |
---|
| 197 | req%src_domain(1:req%max_size)=src_domain(:) |
---|
| 198 | req%src_ind(1:req%max_size)=src_ind(:) |
---|
| 199 | req%target_ind(1:req%max_size)=target_ind(:) |
---|
| 200 | req%src_i(1:req%max_size)=src_i(:) |
---|
| 201 | req%src_j(1:req%max_size)=src_j(:) |
---|
| 202 | req%target_i(1:req%max_size)=target_i(:) |
---|
| 203 | req%target_j(1:req%max_size)=target_j(:) |
---|
| 204 | |
---|
| 205 | req%max_size=req%max_size*2 |
---|
| 206 | |
---|
| 207 | DEALLOCATE(src_domain) |
---|
| 208 | DEALLOCATE(src_ind) |
---|
| 209 | DEALLOCATE(target_ind) |
---|
| 210 | DEALLOCATE(src_i) |
---|
| 211 | DEALLOCATE(src_j) |
---|
| 212 | DEALLOCATE(target_i) |
---|
| 213 | DEALLOCATE(target_j) |
---|
| 214 | |
---|
| 215 | END SUBROUTINE reallocate_request |
---|
| 216 | |
---|
| 217 | |
---|
| 218 | SUBROUTINE request_add_point(ind,i,j,request,pos) |
---|
| 219 | USE domain_mod |
---|
| 220 | USE field_mod |
---|
| 221 | IMPLICIT NONE |
---|
| 222 | INTEGER,INTENT(IN) :: ind |
---|
| 223 | INTEGER,INTENT(IN) :: i |
---|
| 224 | INTEGER,INTENT(IN) :: j |
---|
| 225 | TYPE(t_request),POINTER :: request(:) |
---|
| 226 | INTEGER,INTENT(IN),OPTIONAL :: pos |
---|
| 227 | |
---|
| 228 | INTEGER :: src_domain |
---|
| 229 | INTEGER :: src_iim,src_i,src_j,src_n,src_pos,src_delta |
---|
| 230 | TYPE(t_request),POINTER :: req |
---|
| 231 | TYPE(t_domain),POINTER :: d |
---|
| 232 | |
---|
| 233 | req=>request(ind) |
---|
| 234 | d=>domain(ind) |
---|
| 235 | |
---|
| 236 | IF (req%max_size==req%size) CALL reallocate_request(req) |
---|
| 237 | req%size=req%size+1 |
---|
| 238 | IF (req%type_field==field_t) THEN |
---|
| 239 | src_domain=domain(ind)%assign_domain(i,j) |
---|
| 240 | src_iim=domain_glo(src_domain)%iim |
---|
| 241 | src_i=domain(ind)%assign_i(i,j) |
---|
| 242 | src_j=domain(ind)%assign_j(i,j) |
---|
| 243 | |
---|
| 244 | req%target_ind(req%size)=(j-1)*d%iim+i |
---|
| 245 | req%src_domain(req%size)=src_domain |
---|
| 246 | req%src_ind(req%size)=(src_j-1)*src_iim+src_i |
---|
| 247 | ELSE IF (req%type_field==field_u) THEN |
---|
| 248 | IF (.NOT. PRESENT(pos)) STOP 'argument request_add_point non conforme' |
---|
| 249 | |
---|
| 250 | src_domain=domain(ind)%edge_assign_domain(pos-1,i,j) |
---|
| 251 | src_iim=domain_glo(src_domain)%iim |
---|
| 252 | src_i=domain(ind)%edge_assign_i(pos-1,i,j) |
---|
| 253 | src_j=domain(ind)%edge_assign_j(pos-1,i,j) |
---|
| 254 | src_n=(src_j-1)*src_iim+src_i |
---|
| 255 | src_delta=domain(ind)%delta(i,j) |
---|
| 256 | |
---|
| 257 | ! src_pos=MOD(pos-1+src_delta+6,6)+1 |
---|
| 258 | src_pos=domain(ind)%edge_assign_pos(pos-1,i,j)+1 |
---|
| 259 | |
---|
| 260 | req%target_ind(req%size)=(j-1)*d%iim+i+d%u_pos(pos) |
---|
| 261 | req%src_domain(req%size)=src_domain |
---|
| 262 | req%src_ind(req%size)=src_n+domain_glo(src_domain)%u_pos(src_pos) |
---|
| 263 | |
---|
| 264 | ! req%target_i(req%size)=i |
---|
| 265 | ! req%target_j(req%size)=j |
---|
| 266 | ! req%src_i(req%size)=domain(ind)%assign_i(i,j) |
---|
| 267 | ! req%src_j(req%size)=domain(ind)%assign_j(i,j) |
---|
| 268 | |
---|
| 269 | ! PRINT *,"1--->",ind,i,j,pos |
---|
| 270 | ! PRINT *,"2--->",src_domain,src_i,src_j,src_pos |
---|
| 271 | |
---|
| 272 | ELSE IF (req%type_field==field_z) THEN |
---|
| 273 | IF (.NOT. PRESENT(pos)) STOP 'argument request_add_point non conforme' |
---|
| 274 | |
---|
| 275 | src_domain=domain(ind)%assign_domain(i,j) |
---|
| 276 | src_iim=domain_glo(src_domain)%iim |
---|
| 277 | src_i=domain(ind)%assign_i(i,j) |
---|
| 278 | src_j=domain(ind)%assign_j(i,j) |
---|
| 279 | src_n=(src_j-1)*src_iim+src_i |
---|
| 280 | src_delta=domain(ind)%delta(i,j) |
---|
| 281 | |
---|
| 282 | src_pos=MOD(pos-1+src_delta+6,6)+1 |
---|
| 283 | |
---|
| 284 | req%target_ind(req%size)=(j-1)*d%iim+i+d%z_pos(pos) |
---|
| 285 | req%src_domain(req%size)=src_domain |
---|
| 286 | req%src_ind(req%size)=src_n+domain_glo(src_domain)%z_pos(src_pos) |
---|
| 287 | ENDIF |
---|
| 288 | END SUBROUTINE request_add_point |
---|
| 289 | |
---|
| 290 | |
---|
| 291 | SUBROUTINE Finalize_request(request) |
---|
| 292 | USE mpipara |
---|
| 293 | USE domain_mod |
---|
| 294 | USE mpi_mod |
---|
| 295 | IMPLICIT NONE |
---|
| 296 | TYPE(t_request),POINTER :: request(:) |
---|
| 297 | TYPE(t_request),POINTER :: req |
---|
| 298 | INTEGER :: nb_domain_recv(0:mpi_size-1) |
---|
| 299 | INTEGER :: nb_domain_send(0:mpi_size-1) |
---|
| 300 | INTEGER :: nb_data_domain_recv(ndomain_glo) |
---|
| 301 | INTEGER :: list_domain_recv(ndomain_glo) |
---|
| 302 | INTEGER,ALLOCATABLE :: list_domain_send(:) |
---|
| 303 | INTEGER :: list_domain(ndomain) |
---|
| 304 | |
---|
| 305 | INTEGER :: rank,i,j |
---|
| 306 | INTEGER :: size,ind_glo,ind_loc |
---|
| 307 | INTEGER :: isend, irecv, ireq, nreq |
---|
| 308 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 309 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 310 | |
---|
| 311 | IF (.NOT. using_mpi) RETURN |
---|
| 312 | |
---|
| 313 | DO ind_loc=1,ndomain |
---|
| 314 | req=>request(ind_loc) |
---|
| 315 | |
---|
| 316 | nb_data_domain_recv(:) = 0 |
---|
| 317 | nb_domain_recv(:) = 0 |
---|
| 318 | |
---|
| 319 | DO i=1,req%size |
---|
| 320 | ind_glo=req%src_domain(i) |
---|
| 321 | nb_data_domain_recv(ind_glo)=nb_data_domain_recv(ind_glo)+1 |
---|
| 322 | ENDDO |
---|
| 323 | |
---|
| 324 | DO ind_glo=1,ndomain_glo |
---|
| 325 | IF ( nb_data_domain_recv(ind_glo) > 0 ) nb_domain_recv(domglo_rank(ind_glo))=nb_domain_recv(domglo_rank(ind_glo))+1 |
---|
| 326 | ENDDO |
---|
| 327 | |
---|
| 328 | req%nrecv=sum(nb_domain_recv(:)) |
---|
| 329 | ALLOCATE(req%recv(req%nrecv)) |
---|
| 330 | |
---|
| 331 | irecv=0 |
---|
| 332 | DO ind_glo=1,ndomain_glo |
---|
| 333 | IF (nb_data_domain_recv(ind_glo)>0) THEN |
---|
| 334 | irecv=irecv+1 |
---|
| 335 | list_domain_recv(ind_glo)=irecv |
---|
| 336 | req%recv(irecv)%rank=domglo_rank(ind_glo) |
---|
| 337 | req%recv(irecv)%size=nb_data_domain_recv(ind_glo) |
---|
| 338 | req%recv(irecv)%domain=domglo_loc_ind(ind_glo) |
---|
| 339 | ALLOCATE(req%recv(irecv)%value(req%recv(irecv)%size)) |
---|
| 340 | ALLOCATE(req%recv(irecv)%buffer(req%recv(irecv)%size)) |
---|
| 341 | ENDIF |
---|
| 342 | ENDDO |
---|
| 343 | |
---|
| 344 | req%recv(:)%size=0 |
---|
| 345 | irecv=0 |
---|
| 346 | DO i=1,req%size |
---|
| 347 | irecv=list_domain_recv(req%src_domain(i)) |
---|
| 348 | req%recv(irecv)%size=req%recv(irecv)%size+1 |
---|
| 349 | size=req%recv(irecv)%size |
---|
| 350 | req%recv(irecv)%value(size)=req%src_ind(i) |
---|
| 351 | req%recv(irecv)%buffer(size)=req%target_ind(i) |
---|
| 352 | ENDDO |
---|
| 353 | ENDDO |
---|
| 354 | |
---|
| 355 | nb_domain_recv(:) = 0 |
---|
| 356 | DO ind_loc=1,ndomain |
---|
| 357 | req=>request(ind_loc) |
---|
| 358 | |
---|
| 359 | DO irecv=1,req%nrecv |
---|
| 360 | rank=req%recv(irecv)%rank |
---|
| 361 | nb_domain_recv(rank)=nb_domain_recv(rank)+1 |
---|
| 362 | ENDDO |
---|
| 363 | ENDDO |
---|
| 364 | |
---|
| 365 | |
---|
| 366 | CALL MPI_Alltoall(nb_domain_recv,1,MPI_INTEGER,nb_domain_send,1,MPI_INTEGER,comm_icosa,ierr) |
---|
| 367 | |
---|
| 368 | |
---|
| 369 | ALLOCATE(list_domain_send(sum(nb_domain_send))) |
---|
| 370 | |
---|
| 371 | nreq=sum(nb_domain_recv(:))+sum(nb_domain_send(:)) |
---|
| 372 | ALLOCATE(mpi_req(nreq)) |
---|
| 373 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 374 | |
---|
| 375 | ireq=0 |
---|
| 376 | DO ind_loc=1,ndomain |
---|
| 377 | req=>request(ind_loc) |
---|
| 378 | DO irecv=1,req%nrecv |
---|
| 379 | ireq=ireq+1 |
---|
| 380 | CALL MPI_ISEND(req%recv(irecv)%domain,1,MPI_INTEGER,req%recv(irecv)%rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 381 | ENDDO |
---|
| 382 | ENDDO |
---|
| 383 | |
---|
| 384 | j=0 |
---|
| 385 | DO rank=0,mpi_size-1 |
---|
| 386 | DO i=1,nb_domain_send(rank) |
---|
| 387 | j=j+1 |
---|
| 388 | ireq=ireq+1 |
---|
| 389 | CALL MPI_IRECV(list_domain_send(j),1,MPI_INTEGER,rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 390 | ENDDO |
---|
| 391 | ENDDO |
---|
| 392 | |
---|
| 393 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 394 | |
---|
| 395 | list_domain(:)=0 |
---|
| 396 | DO i=1,sum(nb_domain_send) |
---|
| 397 | ind_loc=list_domain_send(i) |
---|
| 398 | list_domain(ind_loc)=list_domain(ind_loc)+1 |
---|
| 399 | ENDDO |
---|
| 400 | |
---|
| 401 | DO ind_loc=1,ndomain |
---|
| 402 | req=>request(ind_loc) |
---|
| 403 | req%nsend=list_domain(ind_loc) |
---|
| 404 | ALLOCATE(req%send(req%nsend)) |
---|
| 405 | ENDDO |
---|
| 406 | |
---|
| 407 | ireq=0 |
---|
| 408 | DO ind_loc=1,ndomain |
---|
| 409 | req=>request(ind_loc) |
---|
| 410 | |
---|
| 411 | DO irecv=1,req%nrecv |
---|
| 412 | ireq=ireq+1 |
---|
| 413 | CALL MPI_ISEND(mpi_rank,1,MPI_INTEGER,req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 414 | ENDDO |
---|
| 415 | |
---|
| 416 | DO isend=1,req%nsend |
---|
| 417 | ireq=ireq+1 |
---|
| 418 | CALL MPI_IRECV(req%send(isend)%rank,1,MPI_INTEGER,MPI_ANY_SOURCE,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
| 419 | ENDDO |
---|
| 420 | ENDDO |
---|
| 421 | |
---|
| 422 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 423 | CALL MPI_BARRIER(comm_icosa,ierr) |
---|
| 424 | |
---|
| 425 | ireq=0 |
---|
| 426 | DO ind_loc=1,ndomain |
---|
| 427 | req=>request(ind_loc) |
---|
| 428 | |
---|
| 429 | DO irecv=1,req%nrecv |
---|
| 430 | ireq=ireq+1 |
---|
| 431 | CALL MPI_ISEND(req%recv(irecv)%size,1,MPI_INTEGER,req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 432 | ENDDO |
---|
| 433 | |
---|
| 434 | DO isend=1,req%nsend |
---|
| 435 | ireq=ireq+1 |
---|
| 436 | CALL MPI_IRECV(req%send(isend)%size,1,MPI_INTEGER,req%send(isend)%rank,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
| 437 | ENDDO |
---|
| 438 | ENDDO |
---|
| 439 | |
---|
| 440 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 441 | |
---|
| 442 | ireq=0 |
---|
| 443 | DO ind_loc=1,ndomain |
---|
| 444 | req=>request(ind_loc) |
---|
| 445 | |
---|
| 446 | DO irecv=1,req%nrecv |
---|
| 447 | ireq=ireq+1 |
---|
[44] | 448 | CALL MPI_ISEND(req%recv(irecv)%value,req%recv(irecv)%size,MPI_INTEGER,& |
---|
| 449 | req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
[26] | 450 | ENDDO |
---|
| 451 | |
---|
| 452 | DO isend=1,req%nsend |
---|
| 453 | ireq=ireq+1 |
---|
| 454 | ALLOCATE(req%send(isend)%value(req%send(isend)%size)) |
---|
[44] | 455 | CALL MPI_IRECV(req%send(isend)%value,req%send(isend)%size,MPI_INTEGER,& |
---|
| 456 | req%send(isend)%rank,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
[26] | 457 | ENDDO |
---|
| 458 | ENDDO |
---|
| 459 | |
---|
| 460 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 461 | |
---|
| 462 | DO ind_loc=1,ndomain |
---|
| 463 | req=>request(ind_loc) |
---|
| 464 | |
---|
| 465 | DO irecv=1,req%nrecv |
---|
| 466 | req%recv(irecv)%value(:)=req%recv(irecv)%buffer(:) |
---|
| 467 | DEALLOCATE(req%recv(irecv)%buffer) |
---|
| 468 | ENDDO |
---|
| 469 | ENDDO |
---|
| 470 | |
---|
| 471 | |
---|
| 472 | ! nb_domain_recv(:)=0 |
---|
| 473 | ! nb_data_domain_recv(:)=0 |
---|
| 474 | ! |
---|
| 475 | ! DO ind_loc=1,ndomain |
---|
| 476 | ! |
---|
| 477 | ! DO i=1,req%size |
---|
| 478 | ! ind_glo=req%src_domain(i) |
---|
| 479 | ! nb_data_domain_recv(ind_glo)=nb_data_domain_recv(ind_glo)+1 |
---|
| 480 | ! ENDDO |
---|
| 481 | ! |
---|
| 482 | ! DO ind_glo=1,ndomain_glo |
---|
| 483 | ! IF ( nb_data_domain_recv(ind_glo) > 0 ) nb_domain_recv(domglo_rank(ind_glo))=nb_domain_recv(domglo_rank(ind_glo))+1 |
---|
| 484 | ! ENDDO |
---|
| 485 | ! |
---|
| 486 | ! CALL MPI_Alltoall(nb_domain_recv,1,MPI_INTEGER,nb_domain_send,1,MPI_INTEGER,comm_icosa,ierr) |
---|
| 487 | ! ENDDO |
---|
| 488 | ! |
---|
| 489 | ! DO |
---|
| 490 | ! recv=sum(nb_domain_recv(:)) |
---|
| 491 | ! send=sum(nb_domain_send(:)) |
---|
| 492 | |
---|
| 493 | ! ALLOCATE(req%recv(recv)) |
---|
| 494 | ! ALLOCATE(req%send(send)) |
---|
| 495 | |
---|
| 496 | ! ALLOCATE(mpi_req(2*(send+recv))) |
---|
| 497 | ! ALLOCATE(status(MPI_STATUS_SIZE,2*(send+recv))) |
---|
| 498 | ! |
---|
| 499 | ! recv=0 |
---|
| 500 | ! ireq=0 |
---|
| 501 | ! DO ind_glo=1,ndomain_glo |
---|
| 502 | ! IF (nb_data_domain_recv(ind_glo)>0) THEN |
---|
| 503 | ! recv=recv+1 |
---|
| 504 | ! list_domain_recv(ind_glo)=recv |
---|
| 505 | ! req%recv(recv)%rank=domglo_rank(ind_glo) |
---|
| 506 | ! req%recv(recv)%size=nb_data_domain_recv(ind_glo) |
---|
| 507 | ! req%recv(recv)%domain=domglo_loc_ind(ind_glo) |
---|
| 508 | ! ALLOCATE(req%recv(recv)%value(req%recv(recv)%size)) |
---|
| 509 | ! ireq=ireq+1 |
---|
| 510 | ! CALL MPI_ISEND(req%recv(recv)%domain,1,MPI_INTEGER,req%recv(recv)%rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 511 | ! ireq=ireq+1 |
---|
| 512 | ! CALL MPI_ISEND(req%recv(recv)%size,1,MPI_INTEGER,req%recv(recv)%rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 513 | ! ENDIF |
---|
| 514 | ! ENDDO |
---|
| 515 | ! |
---|
| 516 | ! |
---|
| 517 | ! send=0 |
---|
| 518 | ! DO rank=0,mpi_size-1 |
---|
| 519 | ! DO j=1,nb_domain_send(rank) |
---|
| 520 | ! send=send+1 |
---|
| 521 | ! req%send(send)%rank=rank |
---|
| 522 | ! ireq=ireq+1 |
---|
| 523 | ! CALL MPI_IRECV(req%send(send)%domain,1,MPI_INTEGER,rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 524 | ! ireq=ireq+1 |
---|
| 525 | ! CALL MPI_IRECV(req%send(send)%size,1,MPI_INTEGER,rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 526 | ! ENDDO |
---|
| 527 | ! ENDDO |
---|
| 528 | ! |
---|
| 529 | ! CALL MPI_WAITALL(2*(send+recv),mpi_req,status,ierr) |
---|
| 530 | |
---|
| 531 | ! req%recv(:)%size=0 |
---|
| 532 | ! |
---|
| 533 | ! DO i=1,req%size |
---|
| 534 | ! j=list_domain_recv(req%src_domain(i)) |
---|
| 535 | ! req%recv(j)%size=req%recv(j)%size+1 |
---|
| 536 | ! size=req%recv(j)%size |
---|
| 537 | ! req%recv(j)%value(size)=req%src_ind(i) |
---|
| 538 | ! ENDDO |
---|
| 539 | ! |
---|
| 540 | ! ireq=0 |
---|
| 541 | ! DO i=1,recv |
---|
| 542 | ! ireq=ireq+1 |
---|
| 543 | ! CALL MPI_ISEND(req%recv(i)%value,req%recv(i)%size,MPI_INTEGER,req%recv(i)%rank,req%recv(i)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 544 | ! ENDDO |
---|
| 545 | |
---|
| 546 | ! DO i=1,send |
---|
| 547 | ! ireq=ireq+1 |
---|
| 548 | ! ALLOCATE(req%send(i)%value(req%send(i)%size)) |
---|
| 549 | ! CALL MPI_IRECV(req%send(i)%value,req%send(i)%size,MPI_INTEGER,req%send(i)%rank,req%send(i)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 550 | ! ENDDO |
---|
| 551 | ! |
---|
| 552 | ! CALL MPI_WAITALL(send+recv,mpi_req,status,ierr) |
---|
| 553 | |
---|
| 554 | |
---|
| 555 | END SUBROUTINE Finalize_request |
---|
| 556 | |
---|
| 557 | |
---|
| 558 | SUBROUTINE transfert_request_mpi(field,request) |
---|
| 559 | USE field_mod |
---|
| 560 | USE domain_mod |
---|
| 561 | USE mpi_mod |
---|
| 562 | USE mpipara |
---|
| 563 | IMPLICIT NONE |
---|
| 564 | TYPE(t_field),POINTER :: field(:) |
---|
| 565 | TYPE(t_request),POINTER :: request(:) |
---|
| 566 | REAL(rstd),POINTER :: rval2d(:) |
---|
| 567 | REAL(rstd),POINTER :: rval3d(:,:) |
---|
| 568 | REAL(rstd),POINTER :: rval4d(:,:,:) |
---|
| 569 | REAL(rstd),POINTER :: buffer_r2(:) |
---|
| 570 | REAL(rstd),POINTER :: buffer_r3(:,:) |
---|
| 571 | REAL(rstd),POINTER :: buffer_r4(:,:,:) |
---|
| 572 | INTEGER,POINTER :: value(:) |
---|
| 573 | TYPE(ARRAY),POINTER :: recv,send |
---|
| 574 | TYPE(t_request),POINTER :: req |
---|
| 575 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 576 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 577 | INTEGER :: irecv,isend |
---|
| 578 | INTEGER :: ireq,nreq |
---|
| 579 | INTEGER :: ind,n |
---|
| 580 | INTEGER :: dim3 |
---|
| 581 | |
---|
| 582 | IF (field(1)%data_type==type_real) THEN |
---|
| 583 | IF (field(1)%ndim==2) THEN |
---|
| 584 | |
---|
| 585 | nreq=sum(request(:)%nsend)+sum(request(:)%nrecv) |
---|
| 586 | ALLOCATE(mpi_req(nreq)) |
---|
| 587 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 588 | |
---|
| 589 | ireq=0 |
---|
| 590 | DO ind=1,ndomain |
---|
| 591 | rval2d=>field(ind)%rval2d |
---|
| 592 | |
---|
| 593 | req=>request(ind) |
---|
| 594 | DO isend=1,req%nsend |
---|
| 595 | send=>req%send(isend) |
---|
| 596 | |
---|
| 597 | ALLOCATE(send%buffer_r2(send%size)) |
---|
| 598 | buffer_r2=>send%buffer_r2 |
---|
| 599 | value=>send%value |
---|
| 600 | DO n=1,send%size |
---|
| 601 | buffer_r2(n)=rval2d(value(n)) |
---|
| 602 | ENDDO |
---|
| 603 | |
---|
| 604 | ireq=ireq+1 |
---|
| 605 | CALL MPI_ISEND(buffer_r2,send%size,MPI_REAL8,send%rank,ind,comm_icosa, mpi_req(ireq),ierr) |
---|
| 606 | ENDDO |
---|
| 607 | |
---|
| 608 | DO irecv=1,req%nrecv |
---|
| 609 | recv=>req%recv(irecv) |
---|
| 610 | ALLOCATE(recv%buffer_r2(recv%size)) |
---|
| 611 | |
---|
| 612 | ireq=ireq+1 |
---|
| 613 | CALL MPI_IRECV(recv%buffer_r2,recv%size,MPI_REAL8,recv%rank,recv%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 614 | ENDDO |
---|
| 615 | |
---|
| 616 | ENDDO |
---|
| 617 | |
---|
| 618 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 619 | ! DO ind=1,ndomain |
---|
| 620 | ! field(ind)%rval2d(:)=0. |
---|
| 621 | ! ENDDO |
---|
| 622 | |
---|
| 623 | DO ind=1,ndomain |
---|
| 624 | rval2d=>field(ind)%rval2d |
---|
| 625 | |
---|
| 626 | req=>request(ind) |
---|
| 627 | DO isend=1,req%nsend |
---|
| 628 | send=>req%send(isend) |
---|
| 629 | DEALLOCATE(send%buffer_r2) |
---|
| 630 | ENDDO |
---|
| 631 | |
---|
| 632 | DO irecv=1,req%nrecv |
---|
| 633 | recv=>req%recv(irecv) |
---|
| 634 | buffer_r2=>recv%buffer_r2 |
---|
| 635 | value=>recv%value |
---|
| 636 | DO n=1,recv%size |
---|
| 637 | rval2d(value(n))=buffer_r2(n) |
---|
| 638 | ENDDO |
---|
| 639 | DEALLOCATE(recv%buffer_r2) |
---|
| 640 | ENDDO |
---|
| 641 | |
---|
| 642 | ENDDO |
---|
| 643 | |
---|
| 644 | |
---|
| 645 | ELSE IF (field(1)%ndim==3) THEN |
---|
| 646 | |
---|
| 647 | nreq=sum(request(:)%nsend)+sum(request(:)%nrecv) |
---|
| 648 | ALLOCATE(mpi_req(nreq)) |
---|
| 649 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 650 | |
---|
| 651 | ireq=0 |
---|
| 652 | DO ind=1,ndomain |
---|
| 653 | dim3=size(field(ind)%rval3d,2) |
---|
| 654 | rval3d=>field(ind)%rval3d |
---|
| 655 | |
---|
| 656 | req=>request(ind) |
---|
| 657 | DO isend=1,req%nsend |
---|
| 658 | send=>req%send(isend) |
---|
| 659 | |
---|
| 660 | ALLOCATE(send%buffer_r3(send%size,dim3)) |
---|
| 661 | buffer_r3=>send%buffer_r3 |
---|
| 662 | value=>send%value |
---|
| 663 | DO n=1,send%size |
---|
| 664 | buffer_r3(n,:)=rval3d(value(n),:) |
---|
| 665 | ENDDO |
---|
| 666 | |
---|
| 667 | ireq=ireq+1 |
---|
| 668 | CALL MPI_ISEND(buffer_r3,send%size*dim3,MPI_REAL8,send%rank,ind,comm_icosa, mpi_req(ireq),ierr) |
---|
| 669 | ENDDO |
---|
| 670 | |
---|
| 671 | DO irecv=1,req%nrecv |
---|
| 672 | recv=>req%recv(irecv) |
---|
| 673 | ALLOCATE(recv%buffer_r3(recv%size,dim3)) |
---|
| 674 | |
---|
| 675 | ireq=ireq+1 |
---|
| 676 | CALL MPI_IRECV(recv%buffer_r3,recv%size*dim3,MPI_REAL8,recv%rank,recv%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 677 | ENDDO |
---|
| 678 | |
---|
| 679 | ENDDO |
---|
| 680 | |
---|
| 681 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 682 | ! DO ind=1,ndomain |
---|
| 683 | ! field(ind)%rval2d(:)=0. |
---|
| 684 | ! ENDDO |
---|
| 685 | |
---|
| 686 | DO ind=1,ndomain |
---|
| 687 | rval3d=>field(ind)%rval3d |
---|
| 688 | |
---|
| 689 | req=>request(ind) |
---|
| 690 | DO isend=1,req%nsend |
---|
| 691 | send=>req%send(isend) |
---|
| 692 | DEALLOCATE(send%buffer_r3) |
---|
| 693 | ENDDO |
---|
| 694 | |
---|
| 695 | DO irecv=1,req%nrecv |
---|
| 696 | recv=>req%recv(irecv) |
---|
| 697 | buffer_r3=>recv%buffer_r3 |
---|
| 698 | value=>recv%value |
---|
| 699 | DO n=1,recv%size |
---|
| 700 | rval3d(value(n),:)=buffer_r3(n,:) |
---|
| 701 | ENDDO |
---|
| 702 | DEALLOCATE(recv%buffer_r3) |
---|
| 703 | ENDDO |
---|
| 704 | |
---|
| 705 | ENDDO |
---|
| 706 | |
---|
| 707 | ENDIF |
---|
| 708 | |
---|
| 709 | |
---|
| 710 | ENDIF |
---|
| 711 | |
---|
| 712 | END SUBROUTINE transfert_request_mpi |
---|
| 713 | |
---|
| 714 | SUBROUTINE transfert_request(field,request) |
---|
| 715 | USE field_mod |
---|
| 716 | USE domain_mod |
---|
| 717 | IMPLICIT NONE |
---|
| 718 | TYPE(t_field),POINTER :: field(:) |
---|
| 719 | TYPE(t_request),POINTER :: request(:) |
---|
| 720 | REAL(rstd),POINTER :: rval2d(:) |
---|
| 721 | REAL(rstd),POINTER :: rval3d(:,:) |
---|
| 722 | REAL(rstd),POINTER :: rval4d(:,:,:) |
---|
| 723 | INTEGER :: ind |
---|
| 724 | TYPE(t_request),POINTER :: req |
---|
| 725 | INTEGER :: n |
---|
| 726 | REAL(rstd) :: var1,var2 |
---|
| 727 | |
---|
| 728 | DO ind=1,ndomain |
---|
| 729 | req=>request(ind) |
---|
| 730 | rval2d=>field(ind)%rval2d |
---|
| 731 | rval3d=>field(ind)%rval3d |
---|
| 732 | rval4d=>field(ind)%rval4d |
---|
| 733 | |
---|
| 734 | IF (field(ind)%data_type==type_real) THEN |
---|
| 735 | IF (field(ind)%ndim==2) THEN |
---|
| 736 | DO n=1,req%size |
---|
| 737 | rval2d(req%target_ind(n))=field(req%src_domain(n))%rval2d(req%src_ind(n)) |
---|
| 738 | ENDDO |
---|
| 739 | ELSE IF (field(ind)%ndim==3) THEN |
---|
| 740 | DO n=1,req%size |
---|
| 741 | rval3d(req%target_ind(n),:)=field(req%src_domain(n))%rval3d(req%src_ind(n),:) |
---|
| 742 | ENDDO |
---|
| 743 | ELSE IF (field(ind)%ndim==4) THEN |
---|
| 744 | DO n=1,req%size |
---|
| 745 | rval4d(req%target_ind(n),:,:)=field(req%src_domain(n))%rval4d(req%src_ind(n),:,:) |
---|
| 746 | ENDDO |
---|
| 747 | ENDIF |
---|
| 748 | ENDIF |
---|
| 749 | |
---|
| 750 | ENDDO |
---|
| 751 | |
---|
| 752 | END SUBROUTINE transfert_request |
---|
| 753 | |
---|
| 754 | |
---|
| 755 | SUBROUTINE gather_field(field_loc,field_glo) |
---|
| 756 | USE field_mod |
---|
| 757 | USE domain_mod |
---|
| 758 | USE mpi_mod |
---|
| 759 | USE mpipara |
---|
| 760 | IMPLICIT NONE |
---|
| 761 | TYPE(t_field),POINTER :: field_loc(:) |
---|
| 762 | TYPE(t_field),POINTER :: field_glo(:) |
---|
| 763 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 764 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 765 | INTEGER :: ireq,nreq |
---|
| 766 | INTEGER :: ind_glo,ind_loc |
---|
| 767 | |
---|
| 768 | IF (.NOT. using_mpi) THEN |
---|
| 769 | |
---|
| 770 | DO ind_loc=1,ndomain |
---|
| 771 | IF (field_loc(ind_loc)%ndim==2) field_glo(ind_loc)%rval2d=field_loc(ind_loc)%rval2d |
---|
| 772 | IF (field_loc(ind_loc)%ndim==3) field_glo(ind_loc)%rval3d=field_loc(ind_loc)%rval3d |
---|
| 773 | IF (field_loc(ind_loc)%ndim==4) field_glo(ind_loc)%rval4d=field_loc(ind_loc)%rval4d |
---|
| 774 | ENDDO |
---|
| 775 | |
---|
| 776 | ELSE |
---|
| 777 | |
---|
| 778 | nreq=ndomain |
---|
| 779 | IF (mpi_rank==0) nreq=nreq+ndomain_glo |
---|
| 780 | ALLOCATE(mpi_req(nreq)) |
---|
| 781 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 782 | |
---|
| 783 | |
---|
| 784 | ireq=0 |
---|
| 785 | IF (mpi_rank==0) THEN |
---|
| 786 | DO ind_glo=1,ndomain_glo |
---|
| 787 | ireq=ireq+1 |
---|
| 788 | |
---|
| 789 | IF (field_glo(ind_glo)%ndim==2) THEN |
---|
| 790 | CALL MPI_IRECV(field_glo(ind_glo)%rval2d,size(field_glo(ind_glo)%rval2d) , MPI_REAL8 , & |
---|
| 791 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 792 | |
---|
| 793 | ELSE IF (field_glo(ind_glo)%ndim==3) THEN |
---|
| 794 | CALL MPI_IRECV(field_glo(ind_glo)%rval3d,size(field_glo(ind_glo)%rval3d) , MPI_REAL8 , & |
---|
| 795 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 796 | |
---|
| 797 | ELSE IF (field_glo(ind_glo)%ndim==4) THEN |
---|
[31] | 798 | CALL MPI_IRECV(field_glo(ind_glo)%rval4d,size(field_glo(ind_glo)%rval4d) , MPI_REAL8 , & |
---|
[26] | 799 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 800 | ENDIF |
---|
| 801 | |
---|
| 802 | ENDDO |
---|
| 803 | ENDIF |
---|
| 804 | |
---|
| 805 | DO ind_loc=1,ndomain |
---|
| 806 | ireq=ireq+1 |
---|
| 807 | |
---|
| 808 | IF (field_loc(ind_loc)%ndim==2) THEN |
---|
| 809 | CALL MPI_ISEND(field_loc(ind_loc)%rval2d,size(field_loc(ind_loc)%rval2d) , MPI_REAL8 , & |
---|
| 810 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 811 | ELSE IF (field_loc(ind_loc)%ndim==3) THEN |
---|
| 812 | CALL MPI_ISEND(field_loc(ind_loc)%rval3d,size(field_loc(ind_loc)%rval3d) , MPI_REAL8 , & |
---|
| 813 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 814 | ELSE IF (field_loc(ind_loc)%ndim==4) THEN |
---|
| 815 | CALL MPI_ISEND(field_loc(ind_loc)%rval4d,size(field_loc(ind_loc)%rval4d) , MPI_REAL8 , & |
---|
| 816 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 817 | ENDIF |
---|
| 818 | |
---|
| 819 | ENDDO |
---|
| 820 | |
---|
| 821 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 822 | |
---|
| 823 | ENDIF |
---|
| 824 | |
---|
| 825 | END SUBROUTINE gather_field |
---|
| 826 | |
---|
| 827 | |
---|
| 828 | END MODULE transfert_mpi_mod |
---|
| 829 | |
---|
| 830 | |
---|
| 831 | |
---|
| 832 | |
---|
| 833 | |
---|