1!* This file is part of MED. 2!* 3!* COPYRIGHT (C) 1999 - 2019 EDF R&D, CEA/DEN 4!* MED is free software: you can redistribute it and/or modify 5!* it under the terms of the GNU Lesser General Public License as published by 6!* the Free Software Foundation, either version 3 of the License, or 7!* (at your option) any later version. 8!* 9!* MED is distributed in the hope that it will be useful, 10!* but WITHOUT ANY WARRANTY; without even the implied warranty of 11!* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 12!* GNU Lesser General Public License for more details. 13!* 14!* You should have received a copy of the GNU Lesser General Public License 15!* along with MED. If not, see <http://www.gnu.org/licenses/>. 16!* 17 18 19! ****************************************************************************** 20! * - Nom du fichier : Parallel_test1.f90 21! * 22! * - Description : lecture de champs de resultats MED en parallele 23! * 24! ***************************************************************************** 25 26 27program parallel_test1 28 29 implicit none 30 include 'med.hf90' 31 include 'mpif.h' 32 33 integer ret, fid 34 integer USER_INTERLACE,USER_MODE 35 integer*4 com,ioe,rank,nprocs 36 integer info,com4_8 37 integer nent 38 integer nvent 39 integer ncent 40 integer start, stride, count, bsize, lbsize, resd 41 character*64 :: pflname 42 integer*8 flt(1) 43 real*8, allocatable,dimension(:) :: val 44 integer i,j,k 45 46 com4_8=MPI_COMM_WORLD 47 info=MPI_INFO_NULL 48 49 call MPI_INIT(ioe) 50 call MPI_COMM_SIZE(MPI_COMM_WORLD,nprocs,ioe) 51 call MPI_COMM_RANK(MPI_COMM_WORLD,rank,ioe) 52 53 ! ** ouverture du fichier ** 54 call mpfope(fid, 'NENT-942_NVAL-008_NCST-007.med', MED_ACC_RDONLY,com4_8, info, ret) 55 56 if (ret .ne. 0) then 57 print *,"Erreur à l'ouverture du fichier" 58 print *,"Process n° ",rank,"/",nprocs," ret :",ret 59 call efexit(ret) 60 endif 61 62 nent = 942 63 nvent = 008 64 ncent = 007 65 pflname = "" 66 bsize = nent/nprocs 67! Etant donné que l'on affecte qu'un bloc par processus lbsize vaut toujours 0 68 lbsize = 0 69 start = rank*(bsize)+1 70 count = 1 71 stride = bsize 72 resd = 0 73 if (rank.eq.(nprocs-1) ) then 74 resd = nent-(nprocs*bsize) 75 bsize = bsize + resd 76 endif 77 print *,"myrank :",rank," resd", resd," bsize ",bsize," lbsize",lbsize 78 79 call mfrall(1,flt,ret) 80 if (ret .ne. 0) then 81 print *,"Erreur à l'allocation du filtre" 82 print *,"Process n° ",rank,"/",nprocs," ret :",ret 83 call efexit(ret) 84 endif 85 86 call mfrblc (fid, nent, nvent, ncent, & 87 & MED_ALL_CONSTITUENT, MED_FULL_INTERLACE,MED_COMPACT_STMODE ,MED_ALLENTITIES_PROFILE, & 88 & start, stride, count, bsize, lbsize, flt, ret) 89 90 if (ret .ne. 0) then 91 print *,"Erreur à la définition du filtre" 92 print *,"Process n° ",rank,"/",nprocs," ret :",ret 93 call efexit(ret) 94 endif 95 96 allocate(val(bsize*nvent*ncent),STAT=ret) 97 val(:)=-1.1 98 99 call mfdrar ( fid, "NENT-942_NVAL-008_NCST-007_NBL-001",& 100 & MED_NO_DT, MED_NO_IT, MED_CELL, MED_TRIA6,& 101 & flt(1), val, ret ) 102 if (ret .ne. 0) then 103 print *,"Erreur à la lecture du champ résultat" 104 print *,"Process n° ",rank,"/",nprocs," ret :",ret 105 call efexit(ret) 106 endif 107 108 open(40+rank) 109 do i=0,bsize-1 110 do j=0,nvent-1 111 do k=0,ncent-1 112 write(40+rank,'(1X,F10.3,1X)',ADVANCE='NO') val(i*(ncent*nvent)+j*ncent+k+1) 113 enddo 114 write(40+rank,'(A)') "/" 115 enddo 116 write(40+rank,'(A)') "//" 117 enddo 118 close(40+rank) 119 120 deallocate(val) 121 122 call mfrdea(1,flt,ret) 123 if (ret .ne. 0) then 124 print *,"Erreur à la desallocation du filtre" 125 print *,"Process n° ",rank,"/",nprocs," ret :",ret 126 call efexit(ret) 127 endif 128 129 print *,"Process n° ",rank,"/",nprocs," ret :",ret 130 131! call MPI_BARRIER(com,ioe) 132 133 call mficlo(fid,ret) 134 135 call MPI_FINALIZE(ioe) 136 137end program parallel_test1 138