Changeset 71


Ignore:
Timestamp:
06/15/16 17:13:33 (8 years ago)
Author:
dumas
Message:

First parallelization instructions with OpenMP tested in debug : exactly the same results | isostasy called only every 50 years

Location:
trunk/SOURCES
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • trunk/SOURCES/Fichiers-parametres/Makefile.tof-lsce3130.inc

    r65 r71  
    4242                                                                      # (normalement reproductible) 
    4343    ifeq ($(debug), 1) 
    44                 ARITHM    = $(ARITHMi) -CB -g -traceback -warn all 
     44                ARITHM    = $(ARITHMi) -CB -g -p -traceback -warn all -openmp 
    4545        else 
    4646                ARITHM    = $(ARITHMi) -diag-disable warn 
  • trunk/SOURCES/New-remplimat/diagno-L2_mod.f90

    r65 r71  
    11!module diagno_L2_mod               ! Nouvelle version, compatible remplimat 2008 Cat 
    22module diagno_mod                   ! nom pendant les tests 
    3  
     3 !$ USE OMP_LIB 
    44use module3D_phy 
    55use module_choix 
     
    238238! Pour les noeuds posés mais ayant un voisin stream ou flottant, on calcule 
    239239! la viscosité avec stream/shelves 
    240  
    241240! le calcul se fait sur les noeuds majeurs 
     241 
     242!$  integer :: rang ,nb_taches 
     243!$  logical :: paral 
     244 
     245  integer :: t1,t2,ir 
     246  real                           :: temps, t_cpu_0, t_cpu_1, t_cpu, norme   
     247 
    242248if (itracebug.eq.1)  call tracebug(' Calc pvi') 
     249 
     250!$OMP PARALLEL PRIVATE(rang,iplus1,jplus1,sf3,sf1,BT2,phiphi,ttau,d02,discr) 
     251!$ paral = OMP_IN_PARALLEL() 
     252!$ rang=OMP_GET_THREAD_NUM() 
     253!$ nb_taches=OMP_GET_NUM_THREADS() 
     254 
     255!$OMP WORKSHARE 
    243256pvi(:,:)  = pvimin 
    244257Abar(:,:) = 0. 
    245  
    246  
     258!$OMP END WORKSHARE 
     259 
     260!$OMP DO 
    247261do j=1,ny 
    248262   do i=1,nx 
     
    322336   end do 
    323337end do 
     338!$OMP END DO 
    324339 
    325340! cas des noeuds fictifs, si l'épaisseur est très petite  
    326341! pvimin est très petit 
    327  
     342!$OMP WORKSHARE 
    328343where (H(:,:).le.1.) 
    329344   pvi(:,:) = pvimin 
     
    334349end where 
    335350 
     351 
    336352debug_3D(:,:,27)=pvi(:,:) 
    337  
     353!$OMP END WORKSHARE 
    338354! attention run 35 
    339355!--------------------  
     
    345361 
    346362!  calcul de la viscosite integree au milieu des mailles (pvm) 
    347  
     363!$OMP DO 
    348364do i=2,nx 
    349365   do j=2,ny 
     
    355371   end do 
    356372end do 
     373!$OMP END DO 
     374!$OMP END PARALLEL 
     375 
    357376end subroutine calc_pvi 
    358377!------------------------------------------------------------------ 
     
    379398!    -41               -1  Sud                -12  
    380399 
    381  
    382  
     400!$OMP PARALLEL 
     401!$OMP WORKSHARE 
    383402imx_diag(:,:)=0 
    384403imy_diag(:,:)=0 
     
    432451imx_diag(1,:)=0     ! hors domaine a cause des mailles alternees 
    433452imy_diag(:,1)=0     ! hors domaine a cause des mailles alternees 
     453!$OMP END WORKSHARE 
     454!$OMP END PARALLEL 
    434455  
    435456end subroutine imx_imy_nx_ny 
  • trunk/SOURCES/Temperature-routines/Qprod_icetemp.f90

    r70 r71  
    1717 
    1818  use Icetemp_declar 
     19!$ USE OMP_LIB 
     20 
    1921 
    2022  implicit none 
     
    4143  real, dimension(Nx,Ny)      :: Tox           !< Contraintes Sur Maille Mx 
    4244  real, dimension(Nx,Ny)      :: Toy           !< Contraintes Sur Maille Mx 
     45  real      :: Chalk_1       !< Utilise Pour Le Calcul De Chalk : Glace Posée 
     46  real      :: Chalk_2       !< Utilise Pour Le Calcul De Chalk : Ice Streams Et Ice 
    4347 
    4448 
     
    5256     ! Q_prod_demi(T_m%T,T_m%Tpmp,Deform_m%Tobmx,Deform_m%Tobmy,Geom_m%H,Geom_m%Hmx,Geom_m%Hmy,mask_flot_m%Flot,mask_flot_m%Flotmx,mask_flot_m%Flotmy,mask_flot_m%Flgzmx ,mask_flot_m%Flgzmy,Geom_m%Slop_x,Geom_m%Slop_y,Ice_flow_m%Ddx,Ice_flow_m%Ddy,& 
    5357     !      Ice_flow_m%Ddbx,Ice_flow_m%Ddby,Deform_m%Epsxx,Deform_m%Epsyy,Deform_m%Epsxy,mask_flot_m%Gzmx,mask_flot_m%Gzmy,Deform_m%Btt,Ice_flow_m%Uxbar,Ice_flow_m%Uybar,therm_var_m%Phid,Deform_m%Glen,Deform_m%Visc) 
    54  
    55      do K=2,Nz 
    56         do J=2,Ny-1 
    57            do I=2,Nx-1 
     58      
     59 
     60 
     61    !$OMP PARALLEL PRIVATE(chalk_1,chalk_2) 
     62    !$OMP DO COLLAPSE(2) 
     63    do K=2,Nz 
     64    do J=2,Ny-1 
     65            do I=2,Nx-1 
    5866 
    5967              ! Calcul De La Chaleur De Deformation Selon Xx Yy Zz Et Xy 
     
    8088           end do 
    8189        end do 
    82      end do 
     90!        !$ i_min=min(i_min,j) 
     91!        !$ i_max=max(i_max,j) 
     92     end do 
     93        !$OMP END DO NOWAIT 
     94!       !$ print*,'nb_taches = ', nb_taches 
     95!   !$ print*,"Rang : ",rang," i_min : ",i_min," i_max : ",i_max 
    8396 
    8497     !  Partie Sia   Calcul De La Chaleur Produite Sur Chaque Demi Maille 
    85      do L=1,size(Btt,4) !N1poly,N2poly 
    86         do J=2,Ny 
     98                !$OMP DO COLLAPSE(2) 
     99                do L=1,size(Btt,4) !N1poly,N2poly 
     100      do J=2,Ny 
    87101           do I=2,Nx 
    88102              !          Ffx A 3 Dimensions ! 
     
    92106        end do 
    93107     end do 
    94  
    95      do L=1,size(Btt,4)!N1poly,N2poly 
    96         do K=2,Nz 
    97            do J=2,Ny 
    98               do I=2,Nx 
    99                  if ((.not.Flotmx(I,J)).and.(.not.fleuvemx(I,J))) then ! grounded et slowssa 
    100                     Chalx(I,J,K,L)=(Btt(I-1,J,K,L)+Btt(I,J,K,L))*Ffx(I,J,L) !& 
     108     !$OMP END DO 
     109      
     110     !$OMP DO COLLAPSE(2) 
     111          do K=2,Nz 
     112       do J=2,Ny 
     113              do I=2,Nx 
     114                     do L=1,size(Btt,4)!N1poly,N2poly 
     115                    if ((.not.Flotmx(I,J)).and.(.not.fleuvemx(I,J))) then ! grounded et slowssa 
     116                        Chalx(I,J,K,L)=(Btt(I-1,J,K,L)+Btt(I,J,K,L))*Ffx(I,J,L) !& 
    101117                    !                      *Ro*G*Ee(K)**(Glen(L)+1)/Cp(I,J,K) 
    102118 
    103                  else if (fleuvemx(I,J)) then ! Ice Streams 
    104                     Chalx(I,J,K,L)=0. 
    105                  else                     ! Ice Shelves   
    106                     Chalx(I,J,K,L)=0. 
    107                  endif 
    108                  if ((.not.Flotmy(I,J)).and.(.not.fleuvemy(I,J))) then ! grounded et slowssa 
    109                     Chaly(I,J,K,L)=(Btt(I,J-1,K,L)+Btt(I,J,K,L))*Ffy(I,J,L) !& 
     119                    else if (fleuvemx(I,J)) then ! Ice Streams 
     120                        Chalx(I,J,K,L)=0. 
     121                    else                     ! Ice Shelves   
     122                        Chalx(I,J,K,L)=0. 
     123                    endif 
     124                    if ((.not.Flotmy(I,J)).and.(.not.fleuvemy(I,J))) then ! grounded et slowssa 
     125                        Chaly(I,J,K,L)=(Btt(I,J-1,K,L)+Btt(I,J,K,L))*Ffy(I,J,L) !& 
    110126                    !                      *Ro*G*Ee(K)**(Glen(L)+1)/Cp(I,J,K) 
    111                  else if (fleuvemy(I,J)) then  ! Ice Streams 
    112                     Chaly(I,J,K,L)=0. 
    113                  else                      ! Ice Shelves 
    114                     Chaly(I,J,K,L)=0. 
    115                  endif 
    116  
    117               end do 
    118            end do 
    119         end do 
    120      end do 
     127                    else if (fleuvemy(I,J)) then  ! Ice Streams 
     128                        Chaly(I,J,K,L)=0. 
     129                    else                      ! Ice Shelves 
     130                        Chaly(I,J,K,L)=0. 
     131                    endif 
     132 
     133                   end do 
     134            end do 
     135      end do 
     136     end do 
     137     !$OMP END DO 
     138      
    121139 
    122140     !         Nouvelle Formulation De Chaldef_maj(I,J,K), Le 4 Vient Des Moyennes 
     
    125143     !         Ancienne Formulation Chal=(Ro*G*H(I,J))**4*(Sx2+Sy2)*(Sx*Sx+Sy*Sy)  
    126144 
    127      do K=2,Nz 
    128         do J=1,Ny-1 
    129            do I=1,Nx-1 
    130  
     145    !$OMP DO COLLAPSE(2) 
     146    do K=2,Nz 
     147      do J=1,Ny-1 
     148            do I=1,Nx-1 
    131149              ! Modif Christophe Mars 2000 : Chalx Et Chaly Sont A 4 Dim 
    132150              Chaldef_maj(I,J,K)= 0. 
     
    152170        end do 
    153171     end do 
     172     !$OMP END DO NOWAIT 
     173           
    154174     ! Chaleur Produite A La Base Par Le Glissement 
     175     !$OMP DO 
    155176     do J=2,Ny 
    156177        do I=2,Nx 
     
    168189           endif 
    169190 
    170         end do 
    171      end do 
    172  
     191        end do    
     192     end do 
     193     !$OMP END DO 
     194           
    173195     !  Boundary Condition Ice-Rock Interface  
    174196     K=Nz 
    175197 
     198     !$OMP DO 
    176199     do J=2,Ny-1 
    177200        do I=2,Nx-1 
     
    227250        end do 
    228251     end do 
    229  
     252     !$OMP END DO 
     253     !$OMP END PARALLEL 
     254      
    230255 
    231256  case (2)             ! Q_prod_centre : Calcul Avec La Somme Des Carres 
     
    904929 
    905930  end select 
     931 
    906932end subroutine Qprod_ice 
  • trunk/SOURCES/Temperature-routines/icetemp_declar_mod.f90

    r29 r71  
    3232  ! a la Base De L'Ice Shelf 
    3333  Integer :: Ifail1             !< Permet De Detecter Les Erreurs  
    34   Real :: Chalk_1               !< Utilise Pour Le Calcul De Chalk : Glace Posée 
    35   Real :: Chalk_2               !< Utilise Pour Le Calcul De Chalk : Ice Streams Et Ice Shelves 
    3634  Real :: Chalbed_1             !< Utilise Pour Le Calcul De Chalbed 
    3735  Real :: Coefadv               !< Pour Limiter Le Flux De Chaleur Horiz. 
  • trunk/SOURCES/deformation_mod_2lois.f90

    r4 r71  
    179179subroutine flow_general 
    180180 
     181  !$ USE OMP_LIB 
    181182  implicit none 
    182183 
     184  !$OMP PARALLEL 
     185  !$OMP WORKSHARE 
    183186  teta(:,:,:) = t(:,:,1:nz)-tpmp(:,:,:) 
    184  
    185  
     187  !$OMP END WORKSHARE 
     188 
     189  !$OMP DO COLLAPSE(2) 
    186190  do k=nz-1,1,-1 
    187191     do i=2,nx-1 
     
    192196     end do 
    193197  end do 
     198  !$OMP END DO 
     199  !$OMP END PARALLEL 
    194200 
    195201end subroutine flow_general 
     
    197203subroutine flowlaw (iiglen) 
    198204 
     205  !$ USE OMP_LIB 
    199206 
    200207  implicit none 
     
    216223!  real,dimension(nz) :: e          ! vertical coordinate in ice, scaled to h zeta 
    217224  real :: de= 1./(nz-1) 
     225! variables openmp :   
     226  !$  integer :: rang 
     227  !$ integer, dimension(:), allocatable :: tab_sync 
     228  !$ integer :: nb_taches 
     229 
    218230 
    219231  e(1)=0. 
    220232  e(nz)=1. 
     233   
     234  !$OMP PARALLEL PRIVATE(bat,q,aat,ssec,zetat) 
     235  !$ nb_taches = OMP_GET_NUM_THREADS() 
     236  !$OMP DO  
    221237  do k=1,nz 
    222238     if ((k.ne.1).and.(k.ne.nz)) e(k)=(k-1.)/(nz-1.) 
    223239  end do 
     240  !$OMP END DO NOWAIT 
    224241 
    225242  ! exposant de la loi de deformation utilisee : glen(iiglen) 
     
    230247  !    boucle pour btt a k=1 
    231248  ti1=273.15*273.15 
     249  !$OMP DO 
    232250  do j=2,ny-1 
    233251     do i=2,nx-1 
     
    239257     end do 
    240258  end do 
     259  !$OMP END DO 
    241260 
    242261  ! boucle pour tous les autres btt 
     262  !$OMP DO COLLAPSE(2) 
    243263  do k=nz-1,1,-1 
    244      do j=2,ny-1 
     264        do j=2,ny-1   
    245265        do i=2,nx-1   
    246266           if (h(i,j).gt.1.) then 
     
    361381     end do 
    362382  end do 
    363  
     383  !$OMP END DO NOWAIT 
     384   
    364385  !     integration de sa et s2a 
    365  
     386  !$OMP DO 
    366387  do j=1,ny 
    367388     do i=1,nx   
     
    371392     end do 
    372393  end do 
    373  
     394  !$OMP END DO 
     395  !$OMP END PARALLEL 
     396   
     397 
     398  ! Allocation et initialisation du tableau tab_sync 
     399  ! qui gere la synchronisation entre les differents threads 
     400  !$ allocate(tab_sync(0:nb_taches-1)) 
     401  !$ tab_sync(0:nb_taches-1) = 1 
     402  !$OMP PARALLEL private(i,j,k,rang) shared(tab_sync) 
     403  !$ rang = OMP_GET_THREAD_NUM()   
    374404  do k=nz-1,1,-1 
     405       ! Synchronisation des threads 
     406     !$     if (rang /= 0) then 
     407     !$        do 
     408     !$OMP FLUSH(tab_sync)  
     409     !$           if (tab_sync(rang-1)>=tab_sync(rang)+1) exit 
     410     !$        enddo 
     411     !$OMP FLUSH(sa) 
     412     !$OMP FLUSH(s2a) 
     413     !$     endif 
     414     !$OMP DO SCHEDULE(STATIC) 
    375415     do j=1,ny 
    376416        do i=1,nx 
     
    386426        end do 
    387427     end do 
    388   end do 
    389  
     428     !$OMP END DO NOWAIT 
     429!     ! Mis a jour du tableau tab_sync 
     430     !$     tab_sync(rang) = tab_sync(rang) + 1 
     431     !$OMP FLUSH(tab_sync,sa,s2a) 
     432  end do 
     433 
     434  !$OMP WORKSHARE 
    390435  !     cas particulier des bords 
    391436  sa(:,1,:,iiglen)=sa(:,2,:,iiglen) 
     
    402447  s2a(nx,:,:,iiglen)=s2a(nx-1,:,:,iiglen) 
    403448  btt(nx,:,:,iiglen)=btt(nx-1,:,:,iiglen) 
    404  
     449  !$OMP END WORKSHARE 
     450  !$OMP END PARALLEL 
     451   
    405452  ! calcul des moyennes sur les mailles staggered 
    406  
    407453  do k=1,nz 
    408454     tab2d=sa(:,:,k,iiglen) 
     
    417463     s2a_my(:,:,k,iglen)=tab_my 
    418464  end do 
    419  
     465  
    420466! attribution de debug_3d pour les sorties s2a 
    421467! loi 1 -> 190 dans description variable -> 90 dans debug_3d 
    422 debug_3d(:,:,90+iiglen-1) = s2a(:,:,1,iiglen) 
     468!debug_3d(:,:,90+iiglen-1) = s2a(:,:,1,iiglen) 
     469 
     470 
     471!~   ! Temps elapsed final 
     472!~   call system_clock(count=t2, count_rate=ir) 
     473!~   temps=real(t2 - t1,kind=4)/real(ir,kind=4) 
     474!~   ! Temps CPU de calcul final 
     475!~   call cpu_time(t_cpu_1) 
     476!~   t_cpu = t_cpu_1 - t_cpu_0 
     477  
     478!~   ! Impression du resultat. 
     479!~   print '(//,3X,"Valeurs de nx et ny : ",I5,I5/,              & 
     480!~            & 3X,"Temps elapsed       : ",1PE10.3," sec.",/, & 
     481!~            & 3X,"Temps CPU           : ",1PE10.3," sec.",/, & 
     482!~            & 3X,"Norme (PB si /= 0)  : ",1PE10.3,//)', & 
     483!~            nx,ny,temps,t_cpu,norme 
    423484 
    424485end subroutine flowlaw 
  • trunk/SOURCES/flottab2-0.7.f90

    r58 r71  
    1616module flottab_mod 
    1717   
     18  !$ USE OMP_LIB 
    1819USE module3D_phy 
    1920use module_choix 
     
    104105 
    105106 
     107!~   integer :: t1,t2,ir 
     108!~   real                           :: temps, t_cpu_0, t_cpu_1, t_cpu, norme   
     109   
     110!~    ! Temps CPU de calcul initial. 
     111!~    call cpu_time(t_cpu_0) 
     112!~    ! Temps elapsed de reference. 
     113!~    call system_clock(count=t1, count_rate=ir) 
     114 
    106115 
    107116if (itracebug.eq.1)  call tracebug(' Entree dans routine flottab') 
     
    112121! cas particulier des runs paleo ou on impose un masque grounded 
    113122 
     123!$OMP PARALLEL PRIVATE(archim,surnet) 
    114124 if (igrdline.eq.2) then 
     125    !$OMP WORKSHARE 
    115126    where ( mk_init(:,:).eq.1)                 ! pose 
    116127       flot(:,:) = .False. 
     
    120131       flot(:,:) = .True. 
    121132    end where 
     133    !$OMP END WORKSHARE 
    122134 end if 
    123135 
     
    129141 
    130142       appel_new_flot=.false. 
     143       !$OMP DO 
    131144       do j=1,ny 
    132145          do i=1,nx 
     
    136149          enddo 
    137150       enddo 
    138  
     151       !$OMP END DO 
     152        
    139153! ICE(:,:)=(H(:,:).gt.1) ! ice=.true. si epaisseur > 1m 
    140154 
     155      !$OMP WORKSHARE 
    141156      ICE(:,:)=0 
    142157      front(:,:)=0 
     
    148163      cotemy(:,:)=.false. 
    149164      boost=.false. 
     165      !$OMP END WORKSHARE 
    150166 
    151167! fin de l'initialisation 
     
    155171! ------------------------------------- 
    156172 
    157  
     173 !$OMP DO 
    158174 do j=1,ny 
    159175   do i=1,nx 
     
    223239    end do 
    224240 end do 
    225  
     241 !$OMP END DO 
    226242 
    227243!!$ do i=1,nx 
     
    240256 
    241257!----------------------------------------------------------------------- 
    242  
     258 !$OMP DO 
    243259domain_x: do j=1,ny        
    244260         do i=2,nx 
     
    280296          end do 
    281297       end do domain_x 
    282  
    283 if (itracebug.eq.1)  call tracebug('  routine flottab apres domain_x') 
     298 !$OMP END DO 
     299!if (itracebug.eq.1)  call tracebug('  routine flottab apres domain_x') 
    284300 
    285301!     3_y    B- NOUVELLE DEFINITION DE FLOTMY 
    286302!         -------------------------------- 
     303 !$OMP DO 
    287304domain_y: do j=2,ny        
    288305         do i=1,nx 
     
    318335          end do 
    319336       end do domain_y 
    320  
     337 !$OMP END DO 
    321338 
    322339 
     
    351368!      4- determination des iles 
    352369!      ------------------------- 
    353   
     370          !$OMP WORKSHARE 
    354371          ilemx(:,:)=.false. 
    355372          ilemy(:,:)=.false. 
    356  
    357 !       selon x  
     373          !$OMP END WORKSHARE 
     374 
     375!       selon x 
     376  !$OMP DO 
    358377ilesx:  do j=2,ny-1 
    359378           do i=3,nx-2 
     
    404423         end do 
    405424      end do ilesx 
    406  
    407 !       selon y  
     425   !$OMP END DO 
     426 
     427!       selon y 
     428   !$OMP DO 
    408429ilesy: do j=3,ny-2 
    409430           do i=2,nx-1 
     
    452473         end do 
    453474      end do ilesy 
     475    !$OMP END DO 
     476    !$OMP END PARALLEL 
    454477! fin des iles 
    455478 
     
    479502 
    480503!   6- calcule les vitesses des points qui sont devenus gzm 
    481           
     504!$OMP PARALLEL 
     505!$OMP DO          
    482506do j=1,ny 
    483507   do i=2,nx-1 
     
    492516   end do 
    493517end do 
    494  
     518!$OMP END DO 
     519 
     520!$OMP DO  
    495521do j=2,ny-1 
    496522   do i=1,nx 
     
    504530    end do 
    505531end do 
    506  
     532!$OMP END DO 
    507533 
    508534 
     
    511537 
    512538if (nt.ge.2) then   ! pour ne pas faire ce calcul lors du premier passage 
     539   !$OMP WORKSHARE 
    513540   uxbar(:,:)=uxs1(:,:) 
    514    uybar(:,:)=uys1(:,:)  
     541   uybar(:,:)=uys1(:,:) 
     542   !$OMP END WORKSHARE 
    515543endif 
    516544 
     545!$OMP WORKSHARE 
    517546flgzmx(:,:)=(marine.and.(flotmx(:,:).or.gzmx(:,:).or.ilemx(:,:)))   & 
    518547     .or.(.not.marine.and.flotmx(:,:)) 
    519548flgzmy(:,:)=(marine.and.(flotmy(:,:).or.gzmy(:,:).or.ilemy(:,:)))   & 
    520549     .or.(.not.marine.and.flotmy(:,:)) 
    521  
     550!$OMP END WORKSHARE 
    522551 
    523552 
     
    526555!       fbm est vrai si le point est flottant mais un des voisins est pose  
    527556!_________________________________________________________________________ 
     557!$OMP DO 
    528558do j=2,ny-1 
    529559  do i=2,nx-1 
     
    536566  end do 
    537567end do 
    538  
     568!$OMP END DO 
    539569 
    540570 
     
    552582!!$end do 
    553583 
     584!$OMP WORKSHARE 
    554585where (flot(:,:)) 
    555586   where (H(:,:).gt.(1.1))  
     
    565596   end where 
    566597end where 
     598!$OMP END WORKSHARE 
     599!$OMP END PARALLEL 
    567600 
    568601call DETERMIN_TACHE  
     
    581614!----------------------------------------------! 
    582615!On determine les differents ice strean/shelf  ! 
    583       call DETERMIN_TACHE                      ! 
     616!      call DETERMIN_TACHE                      ! 
    584617!----------------------------------------------! 
    585618 
     
    597630 
    598631!On compte comme englacé uniquement les calottes dont une partie est posée       
    599  
     632!$OMP PARALLEL PRIVATE(smax_,smax_coord,smax_i,smax_j) 
     633!$OMP DO 
    600634do i=3,nx-2 
    601635   do j=3,ny-2 
     
    666700   end do 
    667701end do 
    668  
     702!$OMP END DO 
     703!$OMP END PARALLEL 
    669704 
    670705!---------------------------------------------- 
     
    691726!print*, 'front',front(50,30),ice(50,30),flotmx(i,j),uxbar(i,j) 
    692727 
    693  
     728!~   ! Temps elapsed final 
     729!~   call system_clock(count=t2, count_rate=ir) 
     730!~   temps=real(t2 - t1,kind=4)/real(ir,kind=4) 
     731!~   ! Temps CPU de calcul final 
     732!~   call cpu_time(t_cpu_1) 
     733!~   t_cpu = t_cpu_1 - t_cpu_0 
     734  
     735!~   ! Impression du resultat. 
     736!~   print '(//,3X,"Valeurs de nx et ny : ",I5,I5/,              & 
     737!~            & 3X,"Temps elapsed       : ",1PE10.3," sec.",/, & 
     738!~            & 3X,"Temps CPU           : ",1PE10.3," sec.",/, & 
     739!~            & 3X,"Norme (PB si /= 0)  : ",1PE10.3,//)', & 
     740!~            nx,ny,temps,t_cpu,norme 
    694741 
    695742end subroutine flottab 
     
    700747!> 
    701748subroutine determin_tache 
     749 
     750!!$ USE OMP_LIB 
    702751 
    703752implicit none 
     
    720769enddo 
    721770!      table_in  = .false. 
    722  
     771!!$OMP PARALLEL 
     772!!$OMP WORKSHARE 
    723773table_out(:,:) = 0 
    724774iceberg(:)  = .true. 
    725775icetrim (:) = .true. 
    726776nb_pts_tache(:) = 0 
    727  
     777!!$OMP END WORKSHARE 
     778!!$OMP END PARALLEL 
    728779!    open(unit=100,file="tache.data",status='replace') 
    729780 
    730781! 2-reperage des taches 
    731782!---------------------- 
     783!!$OMP PARALLEL PRIVATE(mask,label,indice) 
     784!!$OMP DO 
    732785do i=2,nx-1 
    733786  do j=2,ny-1 
     
    736789 
    737790     IF (ice(i,j).ge.1) THEN ! on est sur la glace-----------------------------! 
    738  
     791                         
    739792        if ((ice(i-1,j).ge.1).or.(ice(i,j-1).ge.1)) then   !masque de 2 cases adjacentes 
    740793           !      un des voisins est deja en glace 
     
    747800              if (mask(indice).gt.0) label=min(label,mask(indice)) 
    748801           enddo 
     802!cdc       label=min(label,minval(mask(:), mask=mask > 0)) 
    749803 
    750804           !on fixe la valeur de la tache voisine minimun au point etudie (via label) 
     
    798852  enddo 
    799853enddo 
    800  
     854!!$OMP END DO 
     855!!$OMP END PARALLEL 
    801856 
    802857 
    803858! On reorganise compt en ecrivant le numero de la tache fondamentale 
    804859! i.e. du plus petit numero present sur la tache (Sans utiliser de recursivité) 
    805 ! On indique aussi le nb de point que contient chaque taches (nb_pts_tache)         
     860! On indique aussi le nb de point que contient chaque taches (nb_pts_tache) 
     861 
    806862do indice=1,label_max 
    807863   vartemp = compt(indice) 
     
    813869enddo 
    814870 
     871!!$OMP PARALLEL 
     872!!$OMP DO REDUCTION(+:nb_pts_tache) 
    815873do i=1,nx 
    816874   do j=1,ny 
     
    821879   enddo 
    822880enddo 
    823  
    824  
     881!!$OMP END DO 
     882!!$OMP END PARALLEL 
    825883 
    826884 
     
    853911!> 
    854912subroutine determin_front 
    855  
     913!!$ USE OMP_LIB 
    856914integer :: i_moins1,i_plus1,i_plus2 
    857915integer :: j_moins1,j_plus1,j_plus2 
    858916       
    859       do i=3,nx-2 
    860        do j=3,ny-2 
     917      !!$OMP PARALLEL 
     918      !!$OMP DO 
     919      do j=3,ny-2 
     920        do i=3,nx-2 
    861921 
    862922 surice:if  (ice(i,j).eq.0) then 
     
    910970       end do 
    911971      end do 
    912  
     972      !!$OMP END DO 
    913973 
    914974!!$call detect_assym(nx,ny,0,41,1,0,1,0,H,itestf) 
     
    923983 
    924984!     print*,'dans remplissage baies',time 
     985        
    925986baies: do k=1,2 
     987!!$OMP DO PRIVATE(i_moins1,j_moins1,i_plus1,j_plus1,i_plus2,j_plus2) 
    926988do j=1,ny 
    927989   do i=1,nx 
     
    9581020   end do 
    9591021end do 
     1022!!$OMP END DO 
    9601023end do baies 
    9611024 
     
    9691032!!$end if 
    9701033 
    971  
    972 do i=2,nx-1 
    973    do j=2,ny-1 
     1034!!$OMP DO 
     1035do j=2,ny-1 
     1036   do i=2,nx-1 
    9741037 
    9751038      if (ice(i,j).eq.1) then     !         test si ice=1 
     
    9831046   end do 
    9841047end do 
     1048!!$OMP END DO 
    9851049 
    9861050! traitement des bords. On considere que l'exterieur n'a pas de glace 
    9871051! attention ce n'est vrai que sur la grande grille 
    9881052 
    989  
     1053!!$OMP DO PRIVATE(i) 
    9901054do j=2,ny-1 
    9911055   i=1 
     
    9941058   front(i,j)=(ice(i-1,j)+ice(i,j+1)+ice(i,j-1)) 
    9951059end do 
    996  
     1060!!$OMP END DO  
     1061 
     1062!!$OMP DO PRIVATE(j) 
    9971063do i=2,nx-1 
    9981064   j=1  
     
    10011067   front(i,j)=(ice(i-1,j)+ice(i+1,j)+ice(i,j-1)) 
    10021068end do 
    1003  
     1069!!$OMP END DO  
    10041070! traitement des coins 
    10051071 
     
    10231089!   si ces deux taches sont posées (ou une des deux), il n'y a pas assez de conditions aux limites 
    10241090 
    1025  
     1091!!$OMP DO  
    10261092do j=1,ny 
    10271093   do i=1,nx-1 
     
    10361102   end do 
    10371103end do 
    1038  
     1104!!$OMP END DO 
     1105 
     1106!!$OMP DO  
    10391107do j=1,ny-1 
    10401108   do i=1,nx 
     
    10491117   end do 
    10501118end do 
     1119!!$OMP END DO 
    10511120 
    10521121!isolx signifie pas de voisins en x 
     
    10561125 
    10571126! calcul de frontfacex et isolx 
     1127!!$OMP DO 
    10581128do j=1,ny 
    10591129   do i=2,nx-1 
     
    10751145   end do 
    10761146end do 
     1147!!$OMP END DO 
    10771148 
    10781149! calcul de frontfacey et isoly 
     1150!!$OMP DO 
    10791151do j=2,ny-1 
    10801152   do i=1,nx 
     
    10961168   end do 
    10971169end do 
    1098  
     1170!!$OMP END DO 
    10991171 
    11001172 
     
    11021174! attention ce n'est vrai que sur la grande grille 
    11031175 
    1104  
     1176!!$OMP DO PRIVATE(i) 
    11051177do j=2,ny-1 
    11061178   i=1 
     
    11211193   end if 
    11221194end do 
    1123  
     1195!!$OMP END DO 
     1196 
     1197!!$OMP DO PRIVATE(j) 
    11241198do i=2,nx-1 
    11251199   j=1  
     
    11401214   end if 
    11411215end do 
    1142  
    1143  
     1216!!$OMP END DO 
     1217!!OMP END PARALLEL 
    11441218 
    11451219return 
  • trunk/SOURCES/mix-SIA-L1_mod.f90

    r4 r71  
    5353subroutine mix_SIA_L1 
    5454 
     55!$ USE OMP_LIB 
     56! variables parallelisation openMP 
     57!!$  integer :: rang ,nb_taches 
    5558 
    5659! subroutine qui associe SIA et L1 
     
    5962 
    6063! use module_choix 
    61     if (itracebug.eq.1)  call tracebug(' routine mix_SIA_L1') 
     64if (itracebug.eq.1)  call tracebug(' routine mix_SIA_L1') 
    6265 
    6366debug_3D(:,:,9)=0. 
     
    7376   ! que celles du L1. On garde le L1 comme vitesse basale. 
    7477   ! version utilisee pour les runs MIS11 de Cairns 
    75  
     78        !$OMP PARALLEL 
     79   !$OMP DO 
    7680   do j=2,ny 
    7781      do i=2,nx 
     
    9599      end do 
    96100   end do 
     101   !$OMP END DO 
     102        !$OMP END PARALLEL 
    97103 
    98104else if (i_resolmeca.eq.2) then 
    99105 
    100106! addition systematique 
    101  
     107!       !$OMP PARALLEL PRIVATE(rang,nb_taches) 
     108!       !$ rang=OMP_GET_THREAD_NUM() 
     109!       !$ nb_taches=OMP_GET_NUM_THREADS() 
     110!       !$OMP DO SCHEDULE(STATIC,NY/nb_taches) 
     111        !$OMP PARALLEL 
     112   !$OMP DO 
    102113   do j=2,ny 
    103114      do i=2,nx 
     115! test sur tout le tableau : 
    104116         if (flgzmx(i,j)) then  
    105             debug_3D(i,j,11)=uxflgz(i,j) 
     117!            debug_3D(i,j,11)=uxflgz(i,j) 
    106118            uxbar(i,j)=uxflgz(i,j)+uxdef(i,j)                       ! on ajoute la deform. SIA 
    107119            ubx(i,j) = uxflgz(i,j) 
    108             do k=1,nz 
     120!            do k=1,nz 
     121!               ux(i,j,k)=ux(i,j,k)-ux(i,j,nz)+uxflgz(i,j)           ! remplace le glissement par uxflgz 
     122!            end do 
     123 
     124!            debug_3D(i,j,9)=uxdef(i,j)                    
     125         endif 
     126 
     127         if (flgzmy(i,j)) then  
     128!            debug_3D(i,j,12)=uyflgz(i,j) 
     129            uybar(i,j)=uyflgz(i,j)+uydef(i,j)                       ! on ajoute la deform. SIA 
     130            uby(i,j) = uyflgz(i,j) 
     131!            do k=1,nz 
     132!               uy(i,j,k)=uy(i,j,k)-uy(i,j,nz)+uyflgz(i,j)           ! remplace le glissement par uyflgz 
     133!            end do 
     134 
     135!            debug_3D(i,j,10)=uydef(i,j)               
     136         endif 
     137      end do 
     138   end do 
     139        !$OMP END DO 
     140        !$OMP END PARALLEL 
     141 
     142! on ne recalcul ux que lorsque uxflgz est modifié soit tous les dtt 
     143        if (isynchro.eq.1) then 
     144!       !$OMP PARALLEL PRIVATE(rang,nb_taches) 
     145!       !$ rang=OMP_GET_THREAD_NUM() 
     146!       !$ nb_taches=OMP_GET_NUM_THREADS() 
     147!       !$OMP DO SCHEDULE(STATIC,NY/nb_taches) 
     148        !$OMP PARALLEL 
     149   !$OMP DO 
     150        do j=2,ny 
     151                do i=2,nx 
     152        if (flgzmx(i,j)) then  
     153                        do k=1,nz 
    109154               ux(i,j,k)=ux(i,j,k)-ux(i,j,nz)+uxflgz(i,j)           ! remplace le glissement par uxflgz 
    110155            end do 
    111  
    112             debug_3D(i,j,9)=uxdef(i,j)                    
    113          end if 
    114  
     156         endif 
    115157         if (flgzmy(i,j)) then  
    116             debug_3D(i,j,12)=uyflgz(i,j) 
    117             uybar(i,j)=uyflgz(i,j)+uydef(i,j)                       ! on ajoute la deform. SIA 
    118             uby(i,j) = uyflgz(i,j) 
    119             do k=1,nz 
     158                        do k=1,nz 
    120159               uy(i,j,k)=uy(i,j,k)-uy(i,j,nz)+uyflgz(i,j)           ! remplace le glissement par uyflgz 
    121160            end do 
    122  
    123             debug_3D(i,j,10)=uydef(i,j)               
    124          end if 
    125       end do 
    126    end do 
    127  
    128  
     161         endif 
     162       enddo 
     163   enddo 
     164        !$OMP END DO 
     165        !$OMP END PARALLEL  
     166endif 
    129167 
    130168else ! SIA et L1 en 2 zones separe : pas d'action (uxbar est deja remis a uxflgz dans diffusiv) 
  • trunk/SOURCES/moy_mxmy.f90

    r4 r71  
    1818subroutine moy_mxmy(n1,n2,X2D,X_mx,X_my) 
    1919! fait la moyenne ponderee d'un tableau X2D sur les mailles staggered 
    20  
     20 !$ USE OMP_LIB 
    2121implicit none 
    2222integer, intent(in) :: n1,n2   !< dimension des tableaux 
     
    2929 
    3030 
     31  !$OMP PARALLEL PRIVATE(i_moins1,j_moins1,i_plus1,j_plus1) 
     32  !$OMP DO 
     33        do j=1,n2 
     34                do i=1,n1 
     35                        i_moins1=max(1,i-1) 
     36                        j_moins1=max(1,j-1) 
     37                        i_plus1=min(n1,i+1) 
     38                        j_plus1=min(n2,j+1) 
    3139 
    32 do j=1,n2 
    33    do i=1,n1 
    34       i_moins1=max(1,i-1) 
    35       j_moins1=max(1,j-1) 
    36       i_plus1=min(n1,i+1) 
    37       j_plus1=min(n2,j+1) 
    38  
    39       X_mx(i,j)=0.25*(X2D(i,j)+X2D(i_moins1,j)) & 
     40                        X_mx(i,j)=0.25*(X2D(i,j)+X2D(i_moins1,j)) & 
    4041           + 0.125*((X2D(i_moins1,j_plus1)+X2D(i,j_plus1))   & 
    4142           +       (X2D(i_moins1,j_moins1)+X2D(i,j_moins1))) 
    4243 
    43       X_my(i,j)=0.25*(X2D(i,j)+X2D(i,j_moins1)) & 
     44                        X_my(i,j)=0.25*(X2D(i,j)+X2D(i,j_moins1)) & 
    4445           + 0.125*((X2D(i_plus1,j_moins1)+X2D(i_plus1,j))   & 
    4546           +       (X2D(i_moins1,j_moins1)+X2D(i_moins1,j))) 
    46  
    47  
    48    end do 
    49 end do 
     47                end do 
     48        end do 
     49  !$OMP END DO 
     50  !$OMP END PARALLEL 
    5051 
    5152return 
  • trunk/SOURCES/steps_time_loop.f90

    r50 r71  
    340340        if (itracebug.eq.1)  call tracebug(' Dans spinup_3_bed') 
    341341 
    342         if ((nbed.eq.1).and.nt.ne.1.and.isynchro.eq.1) then            
    343            call bedrock()                                            !  bedrock adjustment 
     342!        if ((nbed.eq.1).and.nt.ne.1.and.isynchro.eq.1) then 
     343         if ((nbed.eq.1).and.nt.ne.1.and.isynchro.eq.1.and.(mod(abs(TIME),50.).lt.dtmin)) then 
     344           call bedrock                                              !  bedrock adjustment 
    344345        endif 
    345346 
Note: See TracChangeset for help on using the changeset viewer.