1! 2! Copyright (C) 2003, Northwestern University and Argonne National Laboratory 3! See COPYRIGHT notice in top-level directory. 4! 5! $Id: pnf_test.f 2224 2015-12-16 06:10:36Z wkliao $ 6! 7!============================================================================= 8! 9! CODE DEVELOPER 10! John Tannahill, LLNL 11! jrt@llnl.gov 12! Note that this code was adapted from csnap.c, which was written by: 13! Woo-Sun Yang and Chris Ding 14! NERSC, Lawrence Berkeley National Laboratory 15! 16! FILE 17! pnf_test.F 18! 19! ROUTINES 20! Pnf_Test (Program) 21! 22! The purpose of this program is to test the Fortran interface to the 23! parallel netCDF library being developed at Argonne National Lab and 24! Northwestern Univ. 25! 26! This code writes the array, tt(k)(j)(i), into the file 'pnf_test.nc'. It 27! then reads the array from the file, and compares it with the original 28! values. 29! 30! i=longitude, j=latitude, k=level 31! 32!============================================================================= 33 34 program Pnf_Test 35 36 implicit none 37 38 include "mpif.h" 39 include "pnetcdf.inc" 40 41 42! ----------------------- 43! Parameter declarations. 44! ----------------------- 45 46 integer NREADS, NWRITES 47 parameter (NREADS = 5, NWRITES = 5 ) 48 ! number of read samples 49 ! number of write samples 50 51 integer*8 TOTSIZ_3D(3) ! global sizes of 3D field 52 53 54! ---------------------- 55! Variable declarations. 56! ---------------------- 57 58 logical reorder 59 60 logical isperiodic(3) 61 62 integer comm_cart ! Cartesian communicator 63 integer ierr 64 integer*8 istart, jstart, kstart ! offsets of 3D field 65 integer*8 locsiz 66 integer mype ! rank in comm_cart 67 integer totpes ! total number of PEs 68 69 integer*8 locsiz_3d(3) ! local sizes of 3D fields 70 integer pe_coords(3) ! Cartesian PE coords 71 72 integer numpes(3) ! number of PEs along axes; 73 ! determined by MPI where a 74 ! zero is specified 75 76 double precision filsiz 77 78 real*4 rdt_g(2) 79 real*4 rdt_l(2) 80 real*4 wrt_g(2) 81 real*4 wrt_l(2) 82 83 real*4 rrates_g(2) 84 real*4 rrates_l(2) 85 real*4 wrates_g(2) 86 real*4 wrates_l(2) 87 88 89 data reorder / .false. / 90 data isperiodic / .false., .false., .false. / 91 data numpes / 1, 1, 0 / 92 data TOTSIZ_3D / 256, 256, 256 / 93 94! ---------------- 95! Begin execution. 96! ---------------- 97 98 call MPI_Init (ierr) 99 100 call MPI_Comm_Size (MPI_COMM_WORLD, totpes, ierr) 101 102 call MPI_Dims_Create (totpes, 3, numpes, ierr) 103 104 call MPI_Cart_Create 105 & (MPI_COMM_WORLD, 3, numpes, isperiodic, reorder, 106 & comm_cart, ierr) 107 108 call MPI_Comm_Rank (comm_cart, mype, ierr) 109 110 call MPI_Cart_Coords (comm_cart, mype, 3, pe_coords, ierr) 111 112 113 rdt_l(1) = 1.0e38 114 rdt_l(2) = 1.0e38 115 wrt_l(1) = 1.0e38 116 wrt_l(2) = 1.0e38 117! rdt_l(:) = Huge (rdt_l) ! initialize for timing 118! wrt_l(:) = Huge (wrt_l) 119 120 121! ---------------------------------------- 122! Determine local size for tt (locsiz_3d). 123! ---------------------------------------- 124 125! =============== 126 call Find_Locnx 127 & (TOTSIZ_3D(1), pe_coords(1), numpes(1), locsiz_3d(1), istart) 128 call Find_Locnx 129 & (TOTSIZ_3D(2), pe_coords(2), numpes(2), locsiz_3d(2), jstart) 130 call Find_Locnx 131 & (TOTSIZ_3D(3), pe_coords(3), numpes(3), locsiz_3d(3), kstart) 132! =============== 133 134 135! ------------------------------- 136! Compute file size in 1d6 bytes. 137! ------------------------------- 138 139 filsiz = (TOTSIZ_3D(1) * TOTSIZ_3D(2) * TOTSIZ_3D(3)) * 140 & 1.0d-6 * 4.0d0 141 142 143! ------------------------------------- 144! Print data decomposition information. 145! ------------------------------------- 146 147 if (mype .EQ. 0) Write (6,900) 148 149 call MPI_Barrier (comm_cart, ierr) 150 151 Write (6, 902) 152 & mype, pe_coords(1), pe_coords(2), pe_coords(3), 153 & TOTSIZ_3D(1), TOTSIZ_3D(2), TOTSIZ_3D(3), 154 & locsiz_3d(1), locsiz_3d(2), locsiz_3d(3), 155 & kstart, jstart, istart 156 157 900 format ("mype pe_coords totsiz_3d locsiz_3d ", 158 & "kstart,jstart,istart") 159 902 format (i3,3x,i2,1x,i2,1x,i2,2x,i4,1x,i4,1x,i4,4x,i4,1x,i4,1x,i4, 160 & 3x,i6,1x,i6,1x,i6) 161 162 163! ------------------------- 164! Write and then read back. 165! ------------------------- 166 167 locsiz = locsiz_3d(1) * locsiz_3d(2) * locsiz_3d(3) 168 169! =============== 170 call Write_File 171! =============== 172 & ("pnf_test.nc", NWRITES, mype, comm_cart, istart, jstart, 173 & kstart, locsiz, locsiz_3d, TOTSIZ_3D, wrt_l) 174!!! Write (6,*) wrt_l(1), wrt_l(2) 175 176! ============== 177 call Read_File 178! ============== 179 & ("pnf_test.nc", NREADS, mype, comm_cart, istart, jstart, 180 & kstart, locsiz, locsiz_3d, TOTSIZ_3D, rdt_l) 181 182 183! ---------------------------- 184! Compute and print I/O rates. 185! ---------------------------- 186 187 wrates_l(1) = real(filsiz / wrt_l(2)) ! write rate 188 wrates_l(2) = real(filsiz / (wrt_l(1) + wrt_l(2))) ! effective write rate 189 190 rrates_l(1) = real(filsiz / rdt_l(2)) ! read rate 191 rrates_l(2) = real(filsiz / (rdt_l(1) + rdt_l(2))) ! effective read rate 192 193 194 call MPI_Allreduce 195 & (wrates_l, wrates_g, 2, MPI_REAL, MPI_MIN, comm_cart, ierr) 196 call MPI_Allreduce 197 & (rrates_l, rrates_g, 2, MPI_REAL, MPI_MIN, comm_cart, ierr) 198 199 call MPI_Allreduce 200 & (wrt_l, wrt_g, 2, MPI_REAL, MPI_MAX, comm_cart, ierr) 201 call MPI_Allreduce 202 & (rdt_l, rdt_g, 2, MPI_REAL, MPI_MAX, comm_cart, ierr) 203 204 205 if (mype .EQ. 0) then 206 Write (6,905) filsiz 207 Write (6,910) wrates_g(1), wrates_g(2) 208 Write (6,915) rrates_g(1), rrates_g(2) 209 Write (6,920) totpes 210 Write (6,922) wrt_g(1), wrt_g(2), wrates_g(2), 211 & rdt_g(1), rdt_g(2), rrates_g(2) 212 end if 213 214 905 format ("File size: ", e10.3, " MB") 215 910 format (" Write: ", f9.3, " MB/s (eff., ", f9.3, " MB/s)") 216 915 format (" Read : ", f9.3, " MB/s (eff., ", f9.3, " MB/s)") 217 920 format ("Total number PEs: ", i4) 218 922 format (e11.3, e11.3, f9.3, e11.3, e11.3, f9.3) 219 220 221 call MPI_Comm_Free (comm_cart, ierr) 222 223 call MPI_Finalize (ierr) 224 225 226 Stop 227 228 end ! program Pnf_Test 229 230 231! ------------ 232 233 234 subroutine Write_File 235 & (filename, nwrites, mype, comm_cart, istart, jstart, kstart, 236 & locsiz, locsiz_3d, totsiz_3d, wrt_l) 237 238 implicit none 239 240 include "mpif.h" 241 include "pnetcdf.inc" 242 243 244! ---------------------- 245! Argument declarations. 246! ---------------------- 247 248 character*(*) filename 249 integer nwrites 250 integer mype 251 integer comm_cart 252 integer*8 istart, jstart, kstart 253 integer*8 locsiz 254 integer*8 locsiz_3d(3) 255 integer*8 totsiz_3d(3) 256 real*4 wrt_l(2) 257 258 259! ---------------------------- 260! Local variable declarations. 261! ---------------------------- 262 263 integer ierr 264 integer lon_id, lat_id, lev_id 265 integer ncid 266 integer nw 267 integer tt_id 268 269 integer*8 count_3d(3) 270 integer*8 start_3d(3) 271 272 integer dim_id(3) 273 274 double precision t1, t2, t3 275 276 integer max_loc_size 277 parameter( max_loc_size = 20000000 ) 278 real*4 tt(max_loc_size) ! Need tt(locsiz) 279 280 281 if (locsiz .gt. MAX_LOC_SIZE) then 282 print *, 'locsiz = ', locsiz, ' larger than MAX_LOC_SIZE' 283 stop 284 endif 285! ---------------- 286! Begin execution. 287! ---------------- 288 289! start_3d(1:3) = (/ kstart, jstart, istart /) 290! count_3d(:) = locsiz_3d(:) 291 start_3d(1) = istart 292 start_3d(2) = jstart 293 start_3d(3) = kstart 294 count_3d(1) = locsiz_3d(1) 295 count_3d(2) = locsiz_3d(2) 296 count_3d(3) = locsiz_3d(3) 297 298 299 do nw = 1, nwrites 300 301! ============== 302 call Get_Field 303! ============== 304 & (istart, jstart, kstart, locsiz, locsiz_3d, totsiz_3d, tt) 305 306 307 call MPI_Barrier (comm_cart, ierr) 308 t1 = MPI_Wtime ( ) 309 310 311! ================= 312 ierr = Nfmpi_Create 313! ================= 314 & (comm_cart, filename, NF_CLOBBER, MPI_INFO_NULL, ncid) 315 316 317! ================== 318 ierr = Nfmpi_Def_Dim 319 & (ncid, "level", totsiz_3d(1), lon_id) 320 ierr = Nfmpi_Def_Dim 321 & (ncid, "latitude", totsiz_3d(2), lat_id) 322 ierr = Nfmpi_Def_Dim 323 & (ncid, "longitude", totsiz_3d(3), lev_id) 324! ================== 325 326 327 dim_id(1) = lon_id 328 dim_id(2) = lat_id 329 dim_id(3) = lev_id 330 331! ================== 332 ierr = Nfmpi_Def_Var 333! ================== 334 & (ncid, "tt", NF_REAL, 3, dim_id, tt_id) 335 336 337! ================= 338 ierr = Nfmpi_Enddef (ncid) 339! ================= 340 341 342 t2 = MPI_Wtime ( ) 343 344 345! ============================= 346 ierr = Nfmpi_Put_Vara_Real_All 347! ============================= 348 & (ncid, tt_id, start_3d, count_3d, tt) 349 350 351! ================ 352 ierr = Nfmpi_Close (ncid) 353! ================ 354 355 356 call MPI_Barrier (comm_cart, ierr) 357 t3 = MPI_Wtime ( ) 358 359 360 if (t2 - t1 .LT. wrt_l(1)) wrt_l(1) = real(t2 - t1) 361 if (t3 - t2 .LT. wrt_l(2)) wrt_l(2) = real(t3 - t2) 362 363 if (mype .EQ. 0) Write (6,950) nw, t2-t1, t3-t2 364 365 end do 366 367 368 950 format ("write ", i1, ": ", e10.3, 1x, e10.3) 369 370 371 Return 372 373 end 374 375 376! ------------ 377 378 379 subroutine Read_File 380 & (filename, nreads, mype, comm_cart, istart, jstart, kstart, 381 & locsiz, locsiz_3d, totsiz_3d, rdt_l) 382 383 implicit none 384 385 include "mpif.h" 386 include "pnetcdf.inc" 387 388 389! ---------------------- 390! Argument declarations. 391! ---------------------- 392 393 character*(*) filename 394 integer nreads 395 integer mype 396 integer comm_cart 397 integer*8 istart, jstart, kstart 398 integer*8 locsiz 399 integer*8 locsiz_3d(3) 400 integer*8 totsiz_3d(3) 401 real*4 rdt_l(2) 402 403 404! ---------------------------- 405! Local variable declarations. 406! ---------------------------- 407 408 integer ierr 409 integer ncid 410 integer nr 411 integer tt_id 412 integer ii 413 414 integer*8 count_3d(3) 415 integer*8 start_3d(3) 416 417 double precision t1, t2, t3 418 419 integer max_loc_size 420 parameter( max_loc_size = 20000000 ) 421 real*4 buf(max_loc_size) 422 real*4 tt (max_loc_size) 423 424 425! ---------------- 426! Begin execution. 427! ---------------- 428 429! ============== 430 call Get_Field 431! ============== 432 & (istart, jstart, kstart, locsiz, locsiz_3d, totsiz_3d, tt) 433 434 435! start_3d(1:3) = (/ kstart, jstart, istart /) 436! count_3d(:) = locsiz_3d(:) 437 start_3d(1) = istart 438 start_3d(2) = jstart 439 start_3d(3) = kstart 440 count_3d(1) = locsiz_3d(1) 441 count_3d(2) = locsiz_3d(2) 442 count_3d(3) = locsiz_3d(3) 443 444 445 do nr = 1, nreads 446 447 do ii=1, int(locsiz) 448 buf(ii) = -99.99 449 enddo 450 451 452 call MPI_Barrier (comm_cart, ierr) 453 t1 = MPI_Wtime ( ) 454 455 456! =============== 457 ierr = Nfmpi_Open 458! =============== 459 & (comm_cart, filename, NF_NOWRITE, MPI_INFO_NULL, ncid) 460 461 462! ==================== 463 ierr = Nfmpi_Inq_Varid 464! ==================== 465 & (ncid, "tt", tt_id) 466 467 468 t2 = MPI_Wtime ( ) 469 470 471! ============================= 472 ierr = Nfmpi_Get_Vara_Real_All 473! ============================= 474 & (ncid, tt_id, start_3d, count_3d, buf) 475 476 477! ================ 478 ierr = Nfmpi_Close (ncid) 479! ================ 480 481 482 call MPI_Barrier (comm_cart, ierr) 483 t3 = MPI_Wtime ( ) 484 485 486 if (t2 - t1 .LT. rdt_l(1)) rdt_l(1) = real(t2 - t1) 487 if (t3 - t2 .LT. rdt_l(2)) rdt_l(2) = real(t3 - t2) 488 489 if (mype .EQ. 0) Write (6,970) nr, t2-t1, t3-t2 490 491 492 if (nr .EQ. 1) then 493! ================ 494 call Compare_Vec 495! ================ 496 & (mype, comm_cart, locsiz, tt, buf) 497 end if 498 499 end do 500 501 502 970 format (" read ", i1, ": ", e10.3, 1x, e10.3) 503 504 505 Return 506 507 end 508 509 510! ------------ 511 512 513 subroutine Find_Locnx 514 & (nx, mype, totpes, locnx, ibegin) 515 516 implicit none 517 518 include "mpif.h" 519 include "pnetcdf.inc" 520 521! ---------------------- 522! Argument declarations. 523! ---------------------- 524 525 integer*8 nx 526 integer mype 527 integer totpes 528 integer*8 locnx 529 integer*8 ibegin 530 531 532! ---------------------------- 533! Local variable declarations. 534! ---------------------------- 535 536 integer*8 iremain 537 538 539! ---------------- 540! Begin execution. 541! ---------------- 542 543 locnx = nx / totpes 544 545 iremain = nx - (totpes * locnx) 546 547 if (mype .LT. iremain) locnx = locnx + 1 548 549 ibegin = mype * (nx / totpes) + iremain + 1 550 551 if (mype .LT. iremain) ibegin = ibegin + (mype - iremain) 552 553 554 Return 555 556 end 557 558 559! ------------ 560 561 562 subroutine Get_Field 563 & (istart, jstart, kstart, locsiz, locsiz_3d, totsiz_3d, tt) 564 565 implicit none 566 567 include "mpif.h" 568 include "pnetcdf.inc" 569 570 571! ---------------------- 572! Argument declarations. 573! ---------------------- 574 575 integer*8 istart, jstart, kstart 576 integer*8 locsiz 577 integer*8 locsiz_3d(3) 578 integer*8 totsiz_3d(3) 579 real*4 tt(locsiz) 580 581 582! ---------------------------- 583! Local variable declarations. 584! ---------------------------- 585 586 integer ii, jj, kk 587 integer ind 588 589 590! ---------------- 591! Begin execution. 592! ---------------- 593 594 ind = 1 595 596 597 do kk = 1, int(locsiz_3d(3)) 598 do jj = 1, int(locsiz_3d(2)) 599 do ii = 1, int(locsiz_3d(1)) 600 601 tt(ind) = real( 602 & (istart-1 +(ii - 1) + 1 + totsiz_3d(3)*(jstart-1 + 603 & (jj - 1) + totsiz_3d(2)*(kstart-1 + 604 & (kk-1)))) * 1.0d-3) 605 ind = ind + 1 606 607 end do 608 end do 609 end do 610 611 612 Return 613 614 end 615 616 617! ------------ 618 619 620 subroutine Compare_Vec 621 & (mype, comm_cart, locsiz, tt, buf) 622 623 implicit none 624 625 include "mpif.h" 626 627 628! ---------------------- 629! Argument declarations. 630! ---------------------- 631 632 integer mype 633 integer comm_cart 634 integer*8 locsiz 635 real*4 tt (locsiz) 636 real*4 buf(locsiz) 637 638 639! ---------------------------- 640! Local variable declarations. 641! ---------------------------- 642 643 integer ierr 644 integer ii 645 646 real*4 delmax, delmin, delta 647 real*4 diff 648 649 real*4 wr(5) 650 real*4 ws(5) 651 652 653! ---------------- 654! Begin execution. 655! ---------------- 656 657 ws(1) = 0.0d0 ! diff 658 ws(2) = 0.0d0 ! sumsq 659 ws(3) = int(locsiz) ! locsiz 660 ws(4) = 0.0d0 ! delmax 661 ws(5) = 1.0d38 ! Huge (ws) ! delmin 662 663 664 do ii = 1, int(locsiz) 665 delta = (tt(ii) - buf(ii)) * (tt(ii) - buf(ii)) 666 ws(1) = ws(1) + delta 667 ws(2) = ws(2) + tt(ii) * tt(ii) 668 if (delta .GT. ws(4)) ws(4) = delta 669 if (delta .LT. ws(5)) ws(5) = delta 670 end do 671 672 673 call MPI_Allreduce 674 & (ws, wr, 3, MPI_REAL, MPI_SUM, comm_cart, ierr) 675 call MPI_Allreduce 676 & (ws(4), delmax, 1, MPI_REAL, MPI_MAX, comm_cart, ierr) 677 call MPI_Allreduce 678 & (ws(5), delmin, 1, MPI_REAL, MPI_MIN, comm_cart, ierr) 679 680 681 diff = Sqrt (wr(1) / wr(2)) ! normalized error 682 delmax = Sqrt (wr(3) * delmax/wr(2)) ! normalized max difference 683 delmin = Sqrt (wr(3) * delmin/wr(2)) ! normalized min difference 684 685 686 if (mype .EQ. 0) Write (6,990) diff, delmax, delmin 687 688 990 format ("diff, delmax, delmin = ", 689 & e10.3, 1x, e10.3, 1x, e10.3) 690 691 692 693 Return 694 695 end 696 697