Changeset 295 for codes/icosagcm/trunk


Ignore:
Timestamp:
10/31/14 14:52:01 (10 years ago)
Author:
ymipsl
Message:

Merging OpenMP parallisme mode : by subdomain and on vertical level.
This feature is actually experimental but may be retro-compatible with the last method based only on subdomain

YM

Location:
codes/icosagcm/trunk/src
Files:
1 added
7 deleted
27 edited

Legend:

Unmodified
Added
Removed

  • codes/icosagcm/trunk/src/advect_tracer.f90

    r252 r295  
    2626  SUBROUTINE init_advect_tracer 
    2727    USE advect_mod 
     28    USE omp_para 
    2829    REAL(rstd),POINTER :: tangent(:,:) 
    2930    REAL(rstd),POINTER :: normal(:,:) 
    … …  
    4849       tangent=f_tangent(ind) 
    4950       sqrt_leng=f_sqrt_leng(ind) 
    50        CALL init_advect(normal,tangent,sqrt_leng) 
     51       IF (is_omp_level_master) CALL init_advect(normal,tangent,sqrt_leng) 
    5152    END DO 
    5253 
    … …  
    238239 
    239240!--> flush dzqw, adzqw 
    240 !!$OMP BARRIER 
     241!$OMP BARRIER 
    241242 
    242243    ! minmod-limited slope of q 
    … …  
    258259 
    259260    ! 0 slope in top and bottom layers 
    260     IF (omp_first) THEN 
     261    IF (is_omp_first_level) THEN 
    261262      DO ij=ijb,ije 
    262263           dzq(ij,1)=0. 
    … …  
    264265    ENDIF 
    265266       
    266     IF (omp_last) THEN 
     267    IF (is_omp_last_level) THEN 
    267268      DO ij=ijb,ije 
    268269          dzq(ij,llm)=0. 
    … …  
    271272 
    272273!---> flush dzq 
    273 !!$OMP BARRIER   
     274!$OMP BARRIER   
    274275 
    275276    ! sigw = fraction of mass that leaves level l/l+1 
    … …  
    290291    END DO 
    291292    ! wq = 0 at top and bottom 
    292     IF (omp_first) THEN 
     293    IF (is_omp_first_level) THEN 
    293294       DO ij=ijb,ije 
    294295            wq(ij,1)=0. 
    … …  
    296297    ENDIF 
    297298     
    298     IF (omp_last) THEN 
     299    IF (is_omp_last_level) THEN 
    299300      DO ij=ijb,ije 
    300301            wq(ij,llm+1)=0. 
    … …  
    303304 
    304305! --> flush wq 
    305 !!$OMP BARRIER 
     306!$OMP BARRIER 
    306307 
    307308 

  • codes/icosagcm/trunk/src/caldyn_gcm.f90

    r286 r295  
    3333    USE exner_mod 
    3434    USE mpipara 
     35    USE omp_para 
    3536    IMPLICIT NONE 
    3637    CHARACTER(len=255) :: def 
    … …  
    5051       STOP 
    5152    END SELECT 
    52     IF (is_mpi_root) PRINT *, 'caldyn_conserv=',def 
     53    IF (is_master) PRINT *, 'caldyn_conserv=',def 
    5354 
    5455    CALL allocate_caldyn 
    … …  
    7273    CALL allocate_field(f_qv,field_z,type_real,llm)  
    7374   
    74     CALL allocate_field(f_buf_i,   field_t,type_real,llm) 
     75    CALL allocate_field(f_buf_i,   field_t,type_real,llm,name="buffer_i") 
    7576    CALL allocate_field(f_buf_p,   field_t,type_real,llm+1)  
    7677    CALL allocate_field(f_buf_u3d, field_t,type_real,3,llm)  ! 3D vel at cell centers 
    … …  
    102103    INTEGER :: ind,i,j,ij,l 
    103104 
    104     IF (omp_first) THEN 
     105    IF (is_omp_first_level) THEN 
    105106       DO ind=1,ndomain 
    106107          IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    128129    ENDIF 
    129130 
    130 !    !$OMP BARRIER 
     131    !$OMP BARRIER 
    131132  END SUBROUTINE caldyn_BC 
    132133    
    … …  
    143144    USE omp_para 
    144145    USE output_field_mod 
     146    USE checksum_mod 
    145147    IMPLICIT NONE 
    146148    LOGICAL,INTENT(IN)    :: write_out 
    … …  
    291293    END SELECT 
    292294 
    293 !!$OMP BARRIER 
     295!$OMP BARRIER 
    294296    IF (write_out) THEN 
    295297 
    296        IF (is_mpi_root) PRINT *,'CALL write_output_fields' 
     298       IF (is_master) PRINT *,'CALL write_output_fields' 
    297299 
    298300! ---> for openMP test to fix later 
    … …  
    450452 
    451453!!! Compute exner function and geopotential 
    452        DO l = 1,llm 
     454       IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency 
     455         DO l = 1,llm 
    453456!         !$OMP DO SCHEDULE(STATIC)  
    454457          !DIR$ SIMD 
    455           DO ij=ij_begin_ext,ij_end_ext          
    456              p_ik = ptop + mass_ak(l) + mass_bk(l)*ps(ij) ! FIXME : leave ps for the moment ; change ps to Ms later 
    457              !         p_ik = ptop + g*(mass_ak(l)+ mass_bk(l)*ps(i,j)) 
    458              exner_ik = cpp * (p_ik/preff) ** kappa 
    459              pk(ij,l) = exner_ik 
     458            DO ij=ij_begin_ext,ij_end_ext          
     459               p_ik = ptop + mass_ak(l) + mass_bk(l)*ps(ij) ! FIXME : leave ps for the moment ; change ps to Ms later 
     460               !         p_ik = ptop + g*(mass_ak(l)+ mass_bk(l)*ps(i,j)) 
     461               exner_ik = cpp * (p_ik/preff) ** kappa 
     462               pk(ij,l) = exner_ik 
    460463             ! specific volume v = kappa*theta*pi/p = dphi/g/rhodz 
    461              geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik 
     464               geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik 
    462465          ENDDO 
    463        ENDDO 
    464  
     466         ENDDO 
     467       ENDIF 
    465468    ELSE  
    466469       ! We are using a Lagrangian vertical coordinate 
    … …  
    471474       IF(boussinesq) THEN ! compute only geopotential : pressure pk will be computed in compute_caldyn_horiz 
    472475          ! specific volume 1 = dphi/g/rhodz 
    473           DO l = 1,llm 
     476         IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency 
     477           DO l = 1,llm 
    474478!             !$OMP DO SCHEDULE(STATIC)  
    475479             !DIR$ SIMD 
    … …  
    477481                geopot(ij,l+1) = geopot(ij,l) + g*rhodz(ij,l) 
    478482             ENDDO 
    479           ENDDO 
     483           ENDDO 
     484         ENDIF 
    480485       ELSE ! non-Boussinesq, compute geopotential and Exner pressure 
    481486          ! uppermost layer 
    482           !DIR$ SIMD 
    483           DO ij=ij_begin_ext,ij_end_ext          
    484              pk(ij,llm) = ptop + (.5*g)*rhodz(ij,llm) 
    485           END DO 
    486           ! other layers 
    487           DO l = llm-1, 1, -1 
    488  
    489 !          !$OMP DO SCHEDULE(STATIC)  
    490              !DIR$ SIMD 
    491              DO ij=ij_begin_ext,ij_end_ext          
    492                 pk(ij,l) = pk(ij,l+1) + (.5*g)*(rhodz(ij,l)+rhodz(ij,l+1)) 
    493              END DO 
    494           END DO 
     487         IF (is_omp_level_master) THEN  ! no openMP on vertical due to dependency 
     488 
     489           !DIR$ SIMD 
     490           DO ij=ij_begin_ext,ij_end_ext          
     491              pk(ij,llm) = ptop + (.5*g)*rhodz(ij,llm) 
     492           END DO 
     493           ! other layers 
     494           DO l = llm-1, 1, -1 
     495 
     496!           !$OMP DO SCHEDULE(STATIC)  
     497              !DIR$ SIMD 
     498              DO ij=ij_begin_ext,ij_end_ext          
     499                 pk(ij,l) = pk(ij,l+1) + (.5*g)*(rhodz(ij,l)+rhodz(ij,l+1)) 
     500              END DO 
     501           END DO 
    495502          ! surface pressure (for diagnostics) 
    496           DO ij=ij_begin_ext,ij_end_ext          
    497              ps(ij) = pk(ij,1) + (.5*g)*rhodz(ij,1) 
    498           END DO 
     503           DO ij=ij_begin_ext,ij_end_ext          
     504              ps(ij) = pk(ij,1) + (.5*g)*rhodz(ij,1) 
     505           END DO 
    499506 
    500507          ! specific volume v = kappa*theta*pi/p = dphi/g/rhodz 
    501           DO l = 1,llm 
     508           DO l = 1,llm 
    502509 
    503510!             !$OMP DO SCHEDULE(STATIC)  
    504511             !DIR$ SIMD 
    505              DO ij=ij_begin_ext,ij_end_ext          
    506                 p_ik = pk(ij,l) 
    507                 exner_ik = cpp * (p_ik/preff) ** kappa 
    508                 geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik  
    509                 pk(ij,l) = exner_ik 
    510              ENDDO 
    511           ENDDO 
     512              DO ij=ij_begin_ext,ij_end_ext          
     513                 p_ik = pk(ij,l) 
     514                 exner_ik = cpp * (p_ik/preff) ** kappa 
     515                 geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik  
     516                 pk(ij,l) = exner_ik 
     517              ENDDO 
     518            ENDDO 
     519          ENDIF 
    512520       END IF 
    513521 
    514522    END IF 
     523 
     524!ym flush geopot 
     525!$OMP BARRIER 
    515526 
    516527  CALL trace_end("compute_geopot") 
    … …  
    799810  CALL trace_start("compute_caldyn_vert") 
    800811 
    801 !!$OMP BARRIER    
     812!$OMP BARRIER    
    802813!!! cumulate mass flux convergence from top to bottom 
    803   DO  l = llm-1, 1, -1 
     814  IF (is_omp_level_master) THEN 
     815    DO  l = llm-1, 1, -1 
    804816!    IF (caldyn_conserv==energy) CALL test_message(req_qu)  
    805817 
    806818!!$OMP DO SCHEDULE(STATIC)  
    807819!DIR$ SIMD 
    808     DO ij=ij_begin,ij_end          
    809         convm(ij,l) = convm(ij,l) + convm(ij,l+1) 
     820      DO ij=ij_begin,ij_end          
     821          convm(ij,l) = convm(ij,l) + convm(ij,l+1) 
     822      ENDDO 
    810823    ENDDO 
    811   ENDDO 
    812  
    813 ! IMPLICIT FLUSH on convm 
     824  ENDIF 
     825 
     826!$OMP BARRIER 
     827! FLUSH on convm 
    814828!!!!!!!!!!!!!!!!!!!!!!!!! 
    815829 
    816830! compute dps 
    817   IF (omp_first) THEN 
     831  IF (is_omp_first_level) THEN 
    818832!DIR$ SIMD 
    819833     DO ij=ij_begin,ij_end          
    … …  
    834848  ENDDO 
    835849 
     850 !--> flush wflux   
     851 !$OMP BARRIER 
     852 
    836853  DO l=ll_begin,ll_endm1 
    837854!DIR$ SIMD 
    … …  
    847864    ENDDO 
    848865  ENDDO 
     866 
    849867   
    850868! Compute vertical transport 
    … …  
    859877 
    860878 !--> flush wwuu   
    861  ! !$OMP BARRIER 
     879 !$OMP BARRIER 
    862880 
    863881! Add vertical transport to du 
    … …  
    10181036       ps = f_ps(ind) 
    10191037       p  = f_p(ind) 
     1038!$OMP BARRIER 
    10201039       CALL compute_pression(ps,p,0) 
    10211040       pk = f_pk(ind) 
    10221041       pks = f_pks(ind) 
     1042!$OMP BARRIER 
    10231043       CALL compute_exner(ps,p,pks,pk,0) 
     1044!$OMP BARRIER 
    10241045       theta_rhodz = f_theta_rhodz(ind) 
    10251046       theta = f_theta(ind) 

  • codes/icosagcm/trunk/src/check_conserve.f90

    r286 r295  
    99    
    1010  PUBLIC init_check_conserve, check_conserve  
    11     REAL(rstd),SAVE:: mtot0,ztot0,etot0,ang0,stot0,rmsv0 
    12 !$OMP THREADPRIVATE(mtot0,ztot0,etot0,ang0,stot0,rmsv0) 
    13     REAL(rstd),SAVE:: etot,ang,stot,rmsv 
    14 !$OMP THREADPRIVATE(etot,ang,stot,rmsv) 
    15     REAL(rstd),SAVE:: ztot 
    16 !$OMP THREADPRIVATE(ztot)        
     11    REAL(rstd),SAVE:: mtot0,ztot0,etot0,angtot0,stot0,rmsvtot0 
     12!$OMP THREADPRIVATE(mtot0,ztot0,etot0,angtot0,stot0,rmsvtot0) 
    1713      
    1814CONTAINS  
    … …  
    3632  USE caldyn_gcm_mod 
    3733  USE exner_mod  
    38   USE mpipara, ONLY : is_mpi_root, comm_icosa 
     34  USE mpipara, ONLY : is_mpi_master, comm_icosa 
     35  USE omp_para 
    3936  IMPLICIT NONE 
    4037    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    4845    INTEGER::ind,ierr 
    4946    REAL(rstd) :: mtot, rmsdpdt 
    50  
    51     etot=0.0; ang=0.0;stot=0.0;rmsv=0.0  
    52     ztot = 0.0  
    53  
     47    REAL(rstd) :: etot, stot, angtot, rmsvtot, ztot 
     48 
     49    CALL transfert_request(f_ue,req_e1_vect) 
    5450    CALL pression(f_ps,f_p) 
    5551 
    … …  
    6561    CALL vorticity(f_ue,f_vort) 
    6662    CALL check_mass_conserve(f_ps,f_dps,mtot,rmsdpdt) 
    67     CALL check_PV  
     63    CALL check_PV(ztot)  
    6864    CALL exner(f_ps,f_p,f_pks,f_pk) 
    69     CALL check_EN(f_ue,f_theta_rhodz,f_phis)  
    70  
    71     IF (is_mpi_root) THEN  
    72 !$OMP MASTER        
    73        IF ( it == itau0  ) Then  
     65    CALL check_EN(f_ue,f_theta_rhodz,f_phis, etot, stot, angtot, rmsvtot)  
     66 
     67    IF (is_master) THEN  
     68       IF ( it == itau0  ) THEN  
    7469          ztot0 = ztot 
    7570          mtot0 = mtot 
    7671          etot0 = etot  
    77           ang0  = ang  
     72          angtot0  = angtot  
    7873          stot0 = stot  
    7974       END IF 
    8075 
    81        rmsv=SQRT(rmsv/mtot)  
     76       rmsvtot=SQRT(rmsvtot/mtot)  
    8277       ztot=ztot/ztot0-1. ; mtot=mtot/mtot0-1.  
    83        etot=etot/etot0-1. ; ang=ang/ang0-1. ; stot=stot/stot0-1.  
     78       etot=etot/etot0-1. ; angtot=angtot/angtot0-1. ; stot=stot/stot0-1.  
    8479       rmsdpdt= daysec*1.e-2*sqrt(rmsdpdt/ncell_glo)   
    8580        
    8681       OPEN(134,file="checkconsicosa.txt",position='append') 
    87        WRITE(134,4000)mtot,rmsdpdt,etot,ztot,stot,rmsv,ang 
    88        WRITE(134,*)mtot,rmsdpdt,etot,ztot,stot,rmsv,ang  
     82       WRITE(134,4000)mtot,rmsdpdt,etot,ztot,stot,rmsvtot,angtot 
     83       WRITE(134,*)mtot,rmsdpdt,etot,ztot,stot,rmsvtot,angtot  
    8984       WRITE(134,*)"==================================================" 
    90        WRITE(*,4000)mtot,rmsdpdt,etot,ztot,stot,rmsv,ang 
     85       WRITE(*,4000)mtot,rmsdpdt,etot,ztot,stot,rmsvtot,angtot 
    9186        
    92874000   FORMAT(10x,'masse',5x,'rmsdpdt',5x,'energie',5x,'enstrophie'  & 
    9388            ,5x,'entropie',5x,'rmsv',5x,'mt.ang',/,'GLOB  '                & 
    94             ,e10.3,e13.6,5e10.3/)      
     89            ,e10.3,e13.6,5e13.3/)      
    9590       close(134) 
    96 !$OMP END MASTER  
     91       
    9792    END IF 
    9893  END SUBROUTINE check_conserve 
    … …  
    105100  USE transfert_omp_mod 
    106101  USE icosa 
     102  USE omp_para 
    107103  IMPLICIT NONE 
    108104    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    117113    mloc=0.0; rmsloc=0.0 
    118114    DO ind=1,ndomain 
    119       IF (.NOT. assigned_domain(ind)) CYCLE 
     115      IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    120116      CALL swap_dimensions(ind) 
    121117      CALL swap_geometry(ind) 
    … …  
    141137!$OMP END MASTER 
    142138    ELSE 
    143        
     139      mtot=mloc 
     140      rmsdpdt=rmsloc 
    144141    ENDIF 
    145142     
    … …  
    148145!--------------------------------------------------------------------- 
    149146 
    150   SUBROUTINE check_en(f_ue,f_theta_rhodz,f_phis)  
     147  SUBROUTINE check_en(f_ue,f_theta_rhodz,f_phis, etot, stot, angtot, rmsvtot)  
    151148  USE icosa 
    152149  USE pression_mod 
    153150  USE vorticity_mod 
     151  USE mpi_mod 
     152  USE mpipara 
     153  USE transfert_omp_mod 
    154154  IMPLICIT NONE 
    155155    TYPE(t_field), POINTER :: f_ue(:) 
    156156    TYPE(t_field), POINTER :: f_theta_rhodz(:) 
    157157    TYPE(t_field), POINTER :: f_phis(:) 
     158    REAL(rstd),INTENT(OUT) :: etot, stot, angtot, rmsvtot 
    158159    REAL(rstd), POINTER :: ue(:,:) 
    159160    REAL(rstd), POINTER :: p(:,:)  
    … …  
    162163    REAL(rstd), POINTER :: phis(:)  
    163164    REAL(rstd), POINTER :: rhodz(:,:)  
     165    REAL(rstd) :: e, s, ang, rmsv 
     166    REAL(rstd) :: e_mpi, s_mpi, ang_mpi, rmsv_mpi 
    164167    INTEGER :: ind 
    165168 
     169    e = 0 ; s = 0 ; ang = 0 ; rmsv = 0 
    166170 
    167171    DO ind=1,ndomain  
    … …  
    175179      pk=f_pk(ind)  
    176180      phis=f_phis(ind)  
    177       CALL compute_en(ind,ue,p,rhodz,theta_rhodz,pk,phis) 
     181      CALL compute_en(ind,ue,p,rhodz,theta_rhodz,pk,phis, e, s, ang, rmsv) 
    178182    END DO  
    179   
     183 
     184    IF (using_mpi) THEN 
     185      CALL reduce_sum_omp(e, e_mpi) 
     186      CALL reduce_sum_omp(s, s_mpi) 
     187      CALL reduce_sum_omp(ang, ang_mpi) 
     188      CALL reduce_sum_omp(rmsv, rmsv_mpi) 
     189!$OMP MASTER     
     190      CALL MPI_REDUCE(e_mpi, etot, 1, MPI_REAL8, MPI_SUM, 0, comm_icosa, ierr) 
     191      CALL MPI_REDUCE(s_mpi, stot, 1, MPI_REAL8, MPI_SUM, 0, comm_icosa, ierr) 
     192      CALL MPI_REDUCE(ang_mpi, angtot, 1, MPI_REAL8, MPI_SUM, 0, comm_icosa, ierr) 
     193      CALL MPI_REDUCE(rmsv_mpi, rmsvtot, 1, MPI_REAL8, MPI_SUM, 0, comm_icosa, ierr) 
     194!$OMP END MASTER 
     195    ELSE 
     196      etot=e 
     197      stot=s 
     198      angtot=ang 
     199      rmsvtot=rmsv 
     200    ENDIF  
     201 
    180202  END SUBROUTINE check_en  
    181203 
    182   SUBROUTINE compute_en(ind,u,p,rhodz,theta_rhodz,pk,phis) 
     204  SUBROUTINE compute_en(ind,u,p,rhodz,theta_rhodz,pk,phis, e, s, ang, rmsv) 
    183205  USE icosa 
    184206  USE disvert_mod 
    185207  USE wind_mod 
     208  USE omp_para 
    186209  IMPLICIT NONE 
    187210    INTEGER,INTENT(IN)::ind  
    … …  
    192215    REAL(rstd),INTENT(IN) :: pk(iim*jjm,llm) 
    193216    REAL(rstd),INTENT(IN) :: phis(iim*jjm) 
     217    REAL(rstd),INTENT(INOUT) :: e 
     218    REAL(rstd),INTENT(INOUT) :: s 
     219    REAL(rstd),INTENT(INOUT) :: ang 
     220    REAL(rstd),INTENT(INOUT) :: rmsv 
     221     
    194222    REAL(rstd):: theta(iim*jjm,llm) 
    195223    REAL(rstd)::KE(iim*jjm,llm)  
    … …  
    204232 
    205233   
    206     DO l = 1, llm 
     234    DO l = ll_begin, ll_end 
    207235      DO j=jj_begin,jj_end 
    208236        DO i=ii_begin,ii_end 
    … …  
    211239          theta(ij,l) = theta_rhodz(ij,l)/rhodz(ij,l) 
    212240          IF (domain(ind)%own(i,j)) THEN 
    213             stot = stot + Ai(ij)*theta_rhodz(ij,l) 
     241            s = s + Ai(ij)*theta_rhodz(ij,l) 
    214242          END IF  
    215243        END DO  
    … …  
    217245    END DO  
    218246 
    219     DO l = 1, llm 
     247    DO l = ll_begin,ll_end 
    220248      DO j=jj_begin,jj_end 
    221249        DO i=ii_begin,ii_end 
    … …  
    234262    END DO   
    235263 
    236     DO l = 1, llm 
     264    DO l = ll_begin,ll_end 
    237265      DO j=jj_begin,jj_end 
    238266        DO i=ii_begin,ii_end 
    239267          ij=(j-1)*iim+i 
    240268          IF (domain(ind)%own(i,j)) THEN 
    241             etot = etot + masse(ij,l)*(phis(ij)+theta(ij,l)*pk(ij,l)+KE(ij,l)) 
     269            e = e + masse(ij,l)*(phis(ij)+theta(ij,l)*pk(ij,l)+KE(ij,l)) 
    242270          END IF  
    243271        END DO  
    … …  
    251279    radomeg = rad*omega 
    252280 
    253     DO l = 1, llm 
     281    DO l = ll_begin,ll_end 
    254282      DO j=jj_begin,jj_end 
    255283        DO i=ii_begin,ii_end 
    … …  
    267295!--------------------------------------------------------------------- 
    268296 
    269   SUBROUTINE check_PV  
    270   USE icosa 
    271   IMPLICIT NONE 
     297  SUBROUTINE check_PV(ztot) 
     298  USE icosa 
     299  USE mpi_mod 
     300  USE mpipara 
     301  USE transfert_omp_mod 
     302  IMPLICIT NONE 
     303     REAL(rstd),INTENT(OUT) :: ztot 
    272304     REAL(rstd), POINTER :: vort(:,:) 
    273305     REAL(rstd), POINTER :: rhodz(:,:)  
    274306     INTEGER :: ind 
    275  
     307     REAL(rstd) :: z,z_mpi 
     308 
     309     z=0 
    276310     DO ind=1,ndomain  
    277311       IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    280314       vort=f_vort(ind) 
    281315       rhodz=f_rhodz(ind)  
    282        CALL compute_PV(vort,rhodz)  
     316       CALL compute_PV(vort,rhodz,z)  
    283317     ENDDO 
     318 
     319    IF (using_mpi) THEN 
     320      CALL reduce_sum_omp(z, z_mpi) 
     321!$OMP MASTER     
     322      CALL MPI_REDUCE(z_mpi, ztot, 1, MPI_REAL8, MPI_SUM, 0, comm_icosa, ierr) 
     323!$OMP END MASTER 
     324    ELSE 
     325      ztot=z 
     326    ENDIF  
    284327  
    285328  END SUBROUTINE check_PV  
    286329 
    287   SUBROUTINE compute_PV(vort,rhodz)  
     330  SUBROUTINE compute_PV(vort,rhodz,z)  
    288331  USE icosa 
    289332  USE disvert_mod 
     333  USE omp_para 
    290334  IMPLICIT NONE  
    291335    REAL(rstd),INTENT(IN) :: vort(2*iim*jjm,llm) 
    292336    REAL(rstd),INTENT(IN) :: rhodz(iim*jjm,llm)  
     337    REAL(rstd),INTENT(INOUT) :: z 
    293338    REAL(rstd)::qv1,qv2  
    294339    REAL(rstd)::hv1,hv2   
    … …  
    297342    hv1 = 0.0 ; hv2 = 0.0  
    298343    
    299     DO l = 1,llm 
     344    DO l = ll_begin,ll_end 
    300345      DO j=jj_begin+1,jj_end-1 
    301346        DO i=ii_begin+1,ii_end-1 
    … …  
    313358          qv2 =( vort(ij+z_down,l)+fv(ij+z_down) )/hv2 
    314359 
    315           ztot = ztot + (qv1*qv1*hv1 + qv2*qv2*hv2)  
     360          z = z + (qv1*qv1*hv1 + qv2*qv2*hv2)  
    316361 
    317362        ENDDO 
    318363      ENDDO 
     364       
    319365    ENDDO 
    320366 

  • codes/icosagcm/trunk/src/dissip_gcm.f90

    r286 r295  
    7171  USE time_mod 
    7272  USE transfert_omp_mod 
     73  USE omp_para 
    7374  IMPLICIT NONE 
    7475   
    … …  
    100101   CASE('none') 
    101102      rayleigh_friction_type=0 
    102       IF (is_mpi_root) PRINT *, 'No Rayleigh friction' 
     103      IF (is_master) PRINT *, 'No Rayleigh friction' 
    103104   CASE('dcmip2_schaer_noshear') 
    104105      rayleigh_friction_type=1 
    105106      rayleigh_shear=0 
    106       IF (is_mpi_root) PRINT *, 'Rayleigh friction : Schaer-like mountain without shear DCMIP2.1' 
     107      IF (is_master) PRINT *, 'Rayleigh friction : Schaer-like mountain without shear DCMIP2.1' 
    107108   CASE('dcmip2_schaer_shear') 
    108109      rayleigh_shear=1 
    109110      rayleigh_friction_type=2 
    110       IF (is_mpi_root) PRINT *, 'Rayleigh friction : Schaer-like mountain with shear DCMIP2.2' 
     111      IF (is_master) PRINT *, 'Rayleigh friction : Schaer-like mountain with shear DCMIP2.2' 
    111112   CASE DEFAULT 
    112       IF (is_mpi_root) PRINT *, 'Bad selector : rayleigh_friction_type =', TRIM(rayleigh_friction_key), ' in dissip_gcm.f90/init_dissip' 
     113      IF (is_master) PRINT *, 'Bad selector : rayleigh_friction_type =', TRIM(rayleigh_friction_key), ' in dissip_gcm.f90/init_dissip' 
    113114      STOP 
    114115   END SELECT 
    … …  
    119120      rayleigh_tau = rayleigh_tau / scale_factor 
    120121      IF(rayleigh_tau<=0) THEN 
    121          IF (is_mpi_root) PRINT *, 'Forbidden : negative value for rayleigh_friction_tau =',rayleigh_tau 
     122         IF (is_master) PRINT *, 'Forbidden : negative value for rayleigh_friction_tau =',rayleigh_tau 
    122123         STOP 
    123124      END IF 
    … …  
    156157    cdivgrad=-1 
    157158    cgradrot=-1 
    158      
     159 
     160!$OMP BARRIER 
     161!$OMP MASTER 
    159162    DO ind=1,ndomain 
    160       IF (.NOT. assigned_domain(ind)) CYCLE 
    161163      CALL swap_dimensions(ind) 
    162164      CALL swap_geometry(ind) 
    … …  
    175177      ENDDO         
    176178    ENDDO 
    177    
    178  
     179!$OMP END MASTER   
     180!$OMP BARRIER 
    179181 
    180182    DO it=1,20 
    … …  
    184186        CALL transfert_request(f_u,req_e1_vect) 
    185187        DO ind=1,ndomain 
    186           IF (.NOT. assigned_domain(ind)) CYCLE 
     188          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    187189          CALL swap_dimensions(ind) 
    188190          CALL swap_geometry(ind) 
    … …  
    197199       
    198200      DO ind=1,ndomain 
    199         IF (.NOT. assigned_domain(ind)) CYCLE 
     201        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    200202        CALL swap_dimensions(ind) 
    201203        CALL swap_geometry(ind) 
    … …  
    214216 
    215217      IF (using_mpi) THEN  
    216         CALL reduce_sum_omp(dumax,dumax1) 
     218        CALL reduce_max_omp(dumax,dumax1) 
    217219!$OMP MASTER         
    218220        CALL MPI_ALLREDUCE(dumax1,dumax,1,MPI_REAL8,MPI_MAX,comm_icosa,ierr) 
    … …  
    220222        CALL bcast_omp(dumax)   
    221223      ELSE 
    222         CALL allreduce_sum_omp(dumax,dumax1) 
     224        CALL allreduce_max_omp(dumax,dumax1) 
    223225        dumax=dumax1 
    224226      ENDIF   
    225227                         
    226228      DO ind=1,ndomain 
    227         IF (.NOT. assigned_domain(ind)) CYCLE 
     229        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    228230        CALL swap_dimensions(ind) 
    229231        CALL swap_geometry(ind) 
    … …  
    232234        u=du/dumax 
    233235      ENDDO 
    234       IF (is_mpi_root) PRINT *,"gradiv : it :",it ,": dumax",dumax 
     236      IF (is_master) PRINT *,"gradiv : it :",it ,": dumax",dumax 
    235237 
    236238    ENDDO  
    237     IF (is_mpi_root) PRINT *,"gradiv : dumax",dumax 
    238     IF (is_mpi_root) PRINT *, 'mean T-cell edge size (km)', 1.45*radius/iim_glo/1000., & 
     239    IF (is_master) PRINT *,"gradiv : dumax",dumax 
     240    IF (is_master) PRINT *, 'mean T-cell edge size (km)', 1.45*radius/iim_glo/1000., & 
    239241                              'effective T-cell half-edge size (km)', dumax**(-.5/nitergdiv)/1000 
    240     IF (is_mpi_root)  PRINT *, 'Max. time step assuming c=340 m/s and Courant number=2.8 :', & 
     242    IF (is_master)  PRINT *, 'Max. time step assuming c=340 m/s and Courant number=2.8 :', & 
    241243                               2.8/340.*dumax**(-.5/nitergdiv) 
    242244 
    243245    cgraddiv=dumax**(-1./nitergdiv) 
    244     IF (is_mpi_root) PRINT *,"cgraddiv : ",cgraddiv 
    245  
     246    IF (is_master) PRINT *,"cgraddiv : ",cgraddiv 
     247 
     248!$OMP BARRIER 
     249!$OMP MASTER 
    246250    DO ind=1,ndomain 
    247       IF (.NOT. assigned_domain(ind)) CYCLE 
    248251      CALL swap_dimensions(ind) 
    249252      CALL swap_geometry(ind) 
    … …  
    262265      ENDDO         
    263266    ENDDO 
     267!$OMP END MASTER   
     268!$OMP BARRIER 
    264269 
    265270 
    … …  
    270275        CALL transfert_request(f_u,req_e1_vect) 
    271276        DO ind=1,ndomain 
    272           IF (.NOT. assigned_domain(ind)) CYCLE 
     277          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    273278          CALL swap_dimensions(ind) 
    274279          CALL swap_geometry(ind) 
    … …  
    283288       
    284289      DO ind=1,ndomain 
    285         IF (.NOT. assigned_domain(ind)) CYCLE 
     290        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    286291        CALL swap_dimensions(ind) 
    287292        CALL swap_geometry(ind) 
    … …  
    300305 
    301306      IF (using_mpi) THEN  
    302         CALL reduce_sum_omp(dumax,dumax1) 
     307        CALL reduce_max_omp(dumax,dumax1) 
    303308!$OMP MASTER         
    304309        CALL MPI_ALLREDUCE(dumax1,dumax,1,MPI_REAL8,MPI_MAX,comm_icosa,ierr) 
    … …  
    306311        CALL bcast_omp(dumax)   
    307312      ELSE 
    308         CALL allreduce_sum_omp(dumax,dumax1) 
     313        CALL allreduce_max_omp(dumax,dumax1) 
    309314        dumax=dumax1 
    310315      ENDIF   
    … …  
    312317       
    313318      DO ind=1,ndomain 
    314         IF (.NOT. assigned_domain(ind)) CYCLE 
     319        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    315320        CALL swap_dimensions(ind) 
    316321        CALL swap_geometry(ind) 
    … …  
    320325      ENDDO 
    321326       
    322       IF (is_mpi_root) PRINT *,"gradrot : it :",it ,": dumax",dumax 
     327      IF (is_master) PRINT *,"gradrot : it :",it ,": dumax",dumax 
    323328 
    324329    ENDDO  
    325     IF (is_mpi_root) PRINT *,"gradrot : dumax",dumax 
     330    IF (is_master) PRINT *,"gradrot : dumax",dumax 
    326331   
    327332    cgradrot=dumax**(-1./nitergrot)  
    328     IF (is_mpi_root) PRINT *,"cgradrot : ",cgradrot 
     333    IF (is_master) PRINT *,"cgradrot : ",cgradrot 
    329334    
    330335 
    331336 
     337!$OMP BARRIER 
     338!$OMP MASTER 
    332339    DO ind=1,ndomain 
    333       IF (.NOT. assigned_domain(ind)) CYCLE 
    334340      CALL swap_dimensions(ind) 
    335341      CALL swap_geometry(ind) 
    … …  
    344350      ENDDO   
    345351    ENDDO 
     352!$OMP END MASTER 
     353!$OMP BARRIER 
    346354 
    347355    DO it=1,20 
    … …  
    351359        CALL transfert_request(f_theta,req_i1) 
    352360        DO ind=1,ndomain 
    353           IF (.NOT. assigned_domain(ind)) CYCLE 
     361          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    354362          CALL swap_dimensions(ind) 
    355363          CALL swap_geometry(ind) 
    … …  
    364372       
    365373      DO ind=1,ndomain 
    366         IF (.NOT. assigned_domain(ind)) CYCLE 
     374        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    367375        CALL swap_dimensions(ind) 
    368376        CALL swap_geometry(ind) 
    … …  
    379387 
    380388      IF (using_mpi) THEN  
    381         CALL reduce_sum_omp(dthetamax ,dthetamax1) 
     389        CALL reduce_max_omp(dthetamax ,dthetamax1) 
    382390!$OMP MASTER         
    383391        CALL MPI_ALLREDUCE(dthetamax1,dthetamax,1,MPI_REAL8,MPI_MAX,comm_icosa,ierr) 
    … …  
    385393        CALL bcast_omp(dthetamax)   
    386394      ELSE 
    387         CALL allreduce_sum_omp(dthetamax,dthetamax1) 
     395        CALL allreduce_max_omp(dthetamax,dthetamax1) 
    388396        dumax=dumax1 
    389397      ENDIF   
    390398       
    391       IF (is_mpi_root) PRINT *,"divgrad : it :",it ,": dthetamax",dthetamax 
     399      IF (is_master) PRINT *,"divgrad : it :",it ,": dthetamax",dthetamax 
    392400 
    393401      DO ind=1,ndomain 
    394         IF (.NOT. assigned_domain(ind)) CYCLE 
     402        IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    395403        CALL swap_dimensions(ind) 
    396404        CALL swap_geometry(ind) 
    … …  
    401409    ENDDO  
    402410 
    403     IF (is_mpi_root) PRINT *,"divgrad : divgrad",dthetamax 
     411    IF (is_master) PRINT *,"divgrad : divgrad",dthetamax 
    404412   
    405413    cdivgrad=dthetamax**(-1./niterdivgrad)  
    406     IF (is_mpi_root) PRINT *,"cdivgrad : ",cdivgrad 
     414    IF (is_master) PRINT *,"cdivgrad : ",cdivgrad 
    407415 
    408416      
    … …  
    431439       dtdissip=itau_dissip*dt 
    432440    ELSE 
    433        IF (is_mpi_root) PRINT *,"No dissipation time set, setting itau_dissip to 1000000000" 
     441       IF (is_master) PRINT *,"No dissipation time set, setting itau_dissip to 1000000000" 
    434442       itau_dissip=100000000 
    435443    END IF 
    436444    itau_dissip=MAX(1,itau_dissip) 
    437     IF (is_mpi_root) PRINT *,"mintau ",mintau,"itau_dissip",itau_dissip," dtdissip ",dtdissip 
     445    IF (is_master) PRINT *,"mintau ",mintau,"itau_dissip",itau_dissip," dtdissip ",dtdissip 
    438446 
    439447  END SUBROUTINE init_dissip  
    … …  
    713721 
    714722    CALL trace_start("divgrad") 
    715         
     723 
    716724    DO ind=1,ndomain 
    717725      IF (.NOT. assigned_domain(ind)) CYCLE 

  • codes/icosagcm/trunk/src/domain.f90

    r266 r295  
    449449     
    450450    ind=0 
    451     DO rank=0,omp_size-1 
    452       nb_domain=ndomain/omp_size 
    453       IF ( rank < MOD(ndomain,omp_size) ) nb_domain=nb_domain+1 
     451    DO rank=0,omp_domain_size-1 
     452      nb_domain=ndomain/omp_domain_size 
     453      IF ( rank < MOD(ndomain,omp_domain_size) ) nb_domain=nb_domain+1 
    454454     
    455455      DO i=1,nb_domain 
    456456       ind=ind+1 
    457        IF (rank==omp_rank) THEN  
     457       IF (rank==omp_domain_rank) THEN  
    458458         assigned_domain(ind)=.TRUE. 
    459459         PRINT *,"Rank ",mpi_rank," task ",rank," local domain : ",ind," global domain ",domloc_glo_ind(ind) 

  • codes/icosagcm/trunk/src/etat0.f90

    r286 r295  
    2525    USE etat0_dcmip4_mod, ONLY : etat0_dcmip4=>etat0   
    2626    USE etat0_heldsz_mod, ONLY : etat0_heldsz=>etat0   
    27     USE dynetat0_gcm_mod, ONLY : dynetat0_start=>etat0  
    28     USE dynetat0_hz_mod,  ONLY : dynetat0_hz=>etat0  
    2927    USE etat0_start_file_mod, ONLY : etat0_start_file=>etat0   
    3028 
    … …  
    8078    CASE ('dcmip3') 
    8179       CALL etat0_dcmip3(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 
    82      CASE ('dcmip4') 
     80    CASE ('dcmip4') 
    8381        IF(nqtot<2) THEN 
    8482           IF (is_mpi_root)  THEN 
    … …  
    8785           STOP 
    8886        END IF 
    89        CALL etat0_dcmip4(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 
    90      CASE ('readnf_start')  
    91           print*,"readnf_start used"     
    92        CALL dynetat0_start(f_ps,f_phis,f_theta_rhodz,f_u,f_q)  
    93         CASE ('readnf_hz')  
    94           print*,"readnf_hz used" 
    95        CALL dynetat0_hz(f_ps,f_phis,f_theta_rhodz,f_u,f_q)  
     87        CALL etat0_dcmip4(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 
    9688   CASE DEFAULT 
    9789       PRINT*, 'Bad selector for variable etat0 <',etat0_type, & 
    … …  
    114106 
    115107  SUBROUTINE etat0_collocated(f_phis,f_ps,f_mass,f_theta_rhodz,f_u, f_q) 
     108    USE theta2theta_rhodz_mod 
    116109    IMPLICIT NONE 
    117110    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    119112    TYPE(t_field),POINTER :: f_phis(:) 
    120113    TYPE(t_field),POINTER :: f_theta_rhodz(:) 
     114    TYPE(t_field),POINTER :: f_temp(:) 
    121115    TYPE(t_field),POINTER :: f_u(:) 
    122116    TYPE(t_field),POINTER :: f_q(:) 
    … …  
    126120    REAL(rstd),POINTER :: phis(:) 
    127121    REAL(rstd),POINTER :: theta_rhodz(:,:) 
     122    REAL(rstd),POINTER :: temp(:,:) 
    128123    REAL(rstd),POINTER :: u(:,:) 
    129124    REAL(rstd),POINTER :: q(:,:,:) 
    … …  
    138133      phis=f_phis(ind) 
    139134      theta_rhodz=f_theta_rhodz(ind) 
     135      temp=f_temp(ind) 
    140136      u=f_u(ind) 
    141137      q=f_q(ind) 
    142       CALL compute_etat0_collocated(ps,mass, phis, theta_rhodz, u, q) 
     138 
     139      IF( TRIM(etat0_type)=='williamson91.6' ) THEN  
     140       CALL compute_etat0_collocated(ps,mass, phis, theta_rhodz, u, q) 
     141      ELSE 
     142       CALL compute_etat0_collocated(ps,mass, phis, temp, u, q) 
     143      ENDIF 
    143144    ENDDO 
     145     
     146    IF( TRIM(etat0_type)/='williamson91.6' ) CALL temperature2theta_rhodz(f_ps,f_temp,f_theta_rhodz) 
     147     
     148     
    144149  END SUBROUTINE etat0_collocated 
    145150 
    146   SUBROUTINE compute_etat0_collocated(ps,mass, phis, theta_rhodz, u, q) 
    147     USE theta2theta_rhodz_mod 
     151  SUBROUTINE compute_etat0_collocated(ps,mass, phis, temp_i, u, q) 
    148152    USE wind_mod 
    149153    USE etat0_jablonowsky06_mod, ONLY : compute_jablonowsky06 => compute_etat0 
    … …  
    154158    REAL(rstd),INTENT(INOUT) :: mass(iim*jjm,llm) 
    155159    REAL(rstd),INTENT(OUT) :: phis(iim*jjm) 
    156     REAL(rstd),INTENT(OUT) :: theta_rhodz(iim*jjm,llm) 
     160    REAL(rstd),INTENT(OUT) :: temp_i(iim*jjm,llm) 
    157161    REAL(rstd),INTENT(OUT) :: u(3*iim*jjm,llm) 
    158162    REAL(rstd),INTENT(OUT) :: q(iim*jjm,llm,nqtot) 
    159163 
    160     REAL(rstd) :: temp_i(iim*jjm,llm) 
     164    REAL(rstd) :: lon_i(iim*jjm) 
     165    REAL(rstd) :: lat_i(iim*jjm) 
    161166    REAL(rstd) :: ulon_i(iim*jjm,llm) 
    162167    REAL(rstd) :: ulat_i(iim*jjm,llm) 
    163168 
     169    REAL(rstd) :: lon_e(3*iim*jjm) 
     170    REAL(rstd) :: lat_e(3*iim*jjm) 
    164171    REAL(rstd) :: ps_e(3*iim*jjm) 
    165172    REAL(rstd) :: mass_e(3*iim*jjm,llm) 
    … …  
    183190       CALL compute_dcmip5(3*iim*jjm,lon_e,lat_e, phis_e, ps_e, temp_e, ulon_e, ulat_e, q_e) 
    184191    CASE('williamson91.6') 
    185        CALL compute_w91_6(iim*jjm,lon_i,lat_i, phis, mass(:,1), theta_rhodz(:,1), ulon_i(:,1), ulat_i(:,1)) 
     192       CALL compute_w91_6(iim*jjm,lon_i,lat_i, phis, mass(:,1), temp_i(:,1), ulon_i(:,1), ulat_i(:,1)) 
    186193       CALL compute_w91_6(3*iim*jjm,lon_e,lat_e, phis_e, mass_e(:,1), temp_e(:,1), ulon_e(:,1), ulat_e(:,1)) 
    187194    END SELECT 
    188195 
    189     SELECT CASE (TRIM(etat0_type)) 
    190     CASE('williamson91.6') 
    191        ! Do nothing 
    192     CASE DEFAULT 
    193        CALL compute_temperature2theta_rhodz(ps,temp_i,theta_rhodz,0)     
    194     END SELECT 
    195   
    196196    CALL compute_wind_perp_from_lonlat_compound(ulon_e, ulat_e, u) 
    197197 

  • codes/icosagcm/trunk/src/etat0_academic.f90

    r286 r295  
    134134      
    135135 
     136!$OMP BARRIER 
    136137    CALL compute_pression(ps,p,1)      
     138!$OMP BARRIER 
    137139    CALL compute_exner(ps,p,pks,pk,1)   
     140!$OMP BARRIER 
    138141    CALL compute_geopotential(phis,pks,pk,theta,phi,1) 
    139142 

  • codes/icosagcm/trunk/src/etat0_dcmip2.f90

    r286 r295  
    1414  SUBROUTINE etat0(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 
    1515  USE icosa 
     16  USE theta2theta_rhodz_mod 
    1617  IMPLICIT NONE 
    1718    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    2021    TYPE(t_field),POINTER :: f_u(:) 
    2122    TYPE(t_field),POINTER :: f_q(:) 
    22      
     23         
     24    TYPE(t_field),POINTER,SAVE :: f_temp(:) 
    2325    REAL(rstd),POINTER :: ps(:) 
    2426    REAL(rstd),POINTER :: phis(:) 
    2527    REAL(rstd),POINTER :: u(:,:) 
    26     REAL(rstd),POINTER :: theta_rhodz(:,:) 
     28    REAL(rstd),POINTER :: Temp(:,:) 
    2729    REAL(rstd),POINTER :: q(:,:,:) 
    2830 
    … …  
    4547    PRINT *, 'Orographic gravity-wave test case :', TRIM(etat0_type) 
    4648 
     49    CALL allocate_field(f_temp,field_t,type_real,llm,name='temp')  
     50 
    4751    DO ind=1,ndomain 
    4852      IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    5357      phis=f_phis(ind) 
    5458      u=f_u(ind) 
    55       theta_rhodz=f_theta_rhodz(ind) 
    5659      q=f_q(ind) 
    57       CALL compute_etat0_DCMIP2(icase,ps,phis,u,theta_rhodz,q) 
     60      temp=f_temp(ind) 
     61      CALL compute_etat0_DCMIP2(icase,ps,phis,u,Temp,q) 
    5862    ENDDO 
     63     
     64    CALL temperature2theta_rhodz(f_ps,f_temp,f_theta_rhodz) 
    5965             
     66    CALL deallocate_field(f_temp) 
     67    
    6068  END SUBROUTINE etat0 
    6169   
    62   SUBROUTINE compute_etat0_DCMIP2(icase, ps, phis, u, theta_rhodz,q) 
     70  SUBROUTINE compute_etat0_DCMIP2(icase, ps, phis, u, temp,q) 
    6371    USE icosa 
    6472    USE disvert_mod, ONLY : ap,bp 
    … …  
    7179    REAL(rstd), INTENT(OUT) :: phis(iim*jjm) 
    7280    REAL(rstd), INTENT(OUT) :: u(3*iim*jjm,llm) 
    73     REAL(rstd), INTENT(OUT) :: theta_rhodz(iim*jjm,llm) 
     81    REAL(rstd), INTENT(OUT) :: temp(iim*jjm,llm) 
    7482    REAL(rstd), INTENT(OUT) :: q(iim*jjm,llm) 
    75     REAL(rstd) :: ulon(3*iim*jjm,llm), ulat(3*iim*jjm,llm), temp(iim*jjm,llm) 
     83    REAL(rstd) :: ulon(3*iim*jjm,llm), ulat(3*iim*jjm,llm) 
    7684     
    7785    INTEGER :: i,j,l,ij 
    … …  
    9098       END DO 
    9199    END DO 
    92     CALL compute_temperature2theta_rhodz(ps,temp,theta_rhodz,1) 
    93100     
    94101    ! Edges : ulon,ulat 

  • codes/icosagcm/trunk/src/etat0_dcmip3.f90

    r286 r295  
    1818  SUBROUTINE etat0(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 
    1919  USE icosa 
     20  USE theta2theta_rhodz_mod 
    2021  IMPLICIT NONE 
    2122    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    2425    TYPE(t_field),POINTER :: f_u(:) 
    2526    TYPE(t_field),POINTER :: f_q(:) 
     27    TYPE(t_field),POINTER,SAVE :: f_temp(:) 
    2628     
    2729    REAL(rstd),POINTER :: ps(:) 
    2830    REAL(rstd),POINTER :: phis(:) 
    2931    REAL(rstd),POINTER :: u(:,:) 
    30     REAL(rstd),POINTER :: theta_rhodz(:,:) 
     32    REAL(rstd),POINTER :: Temp(:,:) 
    3133    REAL(rstd),POINTER :: q(:,:,:) 
    3234 
    3335    INTEGER :: ind 
     36 
     37    CALL allocate_field(f_temp,field_t,type_real,llm,name='temp')  
    3438 
    3539    DO ind=1,ndomain 
    … …  
    4145      phis=f_phis(ind) 
    4246      u=f_u(ind) 
    43       theta_rhodz=f_theta_rhodz(ind) 
    4447      q=f_q(ind) 
    45       CALL compute_etat0_DCMIP3(ps,phis,u,theta_rhodz,q) 
     48      temp=f_temp(ind) 
     49      CALL compute_etat0_DCMIP3(ps,phis,u,Temp,q) 
    4650    ENDDO 
     51 
     52    CALL temperature2theta_rhodz(f_ps,f_temp,f_theta_rhodz) 
     53    CALL deallocate_field(f_temp) 
    4754             
    4855  END SUBROUTINE etat0 
    4956   
    5057 
    51   SUBROUTINE compute_etat0_DCMIP3(ps, phis, u, theta_rhodz,q) 
     58  SUBROUTINE compute_etat0_DCMIP3(ps, phis, u, temp,q) 
    5259  USE icosa 
    5360  USE pression_mod 
    … …  
    6875  REAL(rstd), INTENT(OUT) :: phis(iim*jjm) 
    6976  REAL(rstd), INTENT(OUT) :: u(3*iim*jjm,llm) 
    70   REAL(rstd), INTENT(OUT) :: theta_rhodz(iim*jjm,llm) 
     77  REAL(rstd), INTENT(OUT) :: Temp(iim*jjm,llm) 
    7178  REAL(rstd), INTENT(OUT) :: q(iim*jjm,llm,nqtot) 
    7279   
    7380  REAL(rstd) :: Ts(iim*jjm) 
    7481  REAL(rstd) :: s(iim*jjm) 
    75   REAL(rstd) :: T(iim*jjm,llm) 
    7682  REAL(rstd) :: p(iim*jjm,llm+1) 
    7783  REAL(rstd) :: theta(iim*jjm,llm) 
    … …  
    125131  END DO 
    126132   
     133!$OMP BARRIER 
    127134  CALL compute_pression(ps,p,0) 
     135!$OMP BARRIER 
    128136   
    129137  DO l=1,llm 
    … …  
    134142           IF(use_dcmip_routine) THEN 
    135143              CALL test3_gravity_wave(lon_i(ij),lat_i(ij),pp,dummy,0, & 
    136                    dummy,dummy,dummy,T(ij,l),dummy,dummy,dummy,dummy) 
     144                   dummy,dummy,dummy,Temp(ij,l),dummy,dummy,dummy,dummy) 
    137145           ELSE 
    138146              pspsk=(pp/ps(ij))**kappa 
    … …  
    142150              thetap = dtheta *sin(2*Pi*zz/Lz) * s(ij)                ! perturbation pot. temp. 
    143151              theta(ij,l) = thetab + thetap 
    144               T(ij,l) = theta(ij,l)* ((pp/peq)**kappa) 
     152              Temp(ij,l) = theta(ij,l)* ((pp/peq)**kappa) 
    145153              ! T(ij,l) = Ts(ij)*pspsk / ( Ts(ij) / GG * ( pspsk-1) +1)  ! background temp. 
    146154           END IF 
    … …  
    149157  ENDDO 
    150158   
    151   IF(use_dcmip_routine) THEN 
    152      CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,0) 
    153   ELSE 
    154      CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,0) 
    155   END IF 
     159!  IF(use_dcmip_routine) THEN 
     160!     CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,0) 
     161!  ELSE 
     162!     CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,0) 
     163!  END IF 
    156164   
    157165  pp=peq 

  • codes/icosagcm/trunk/src/exner.f90

    r186 r295  
    2222    INTEGER :: ind 
    2323 
     24!$OMP BARRIER 
    2425    DO ind=1,ndomain 
    2526      IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    3233      CALL compute_exner(ps, p, pks, pk, 0) 
    3334    ENDDO 
     35!$OMP BARRIER 
    3436   
    3537  END SUBROUTINE exner 
    … …  
    3941  USE disvert_mod 
    4042  USE pression_mod 
     43  USE omp_para 
    4144  IMPLICIT NONE 
    4245    REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
    … …  
    5255       !! Compute Alpha and Beta 
    5356        
     57       IF (is_omp_level_master) THEN 
    5458       ! for llm layer 
    55        DO j=jj_begin-offset,jj_end+offset 
    56           DO i=ii_begin-offset,ii_end+offset 
    57              ij=(j-1)*iim+i 
    58              alpha(ij,llm) = 0. 
    59              beta (ij,llm) = 1./ (1+ 2*kappa) 
    60           ENDDO 
    61        ENDDO 
     59         DO j=jj_begin-offset,jj_end+offset 
     60            DO i=ii_begin-offset,ii_end+offset 
     61               ij=(j-1)*iim+i 
     62               alpha(ij,llm) = 0. 
     63               beta (ij,llm) = 1./ (1+ 2*kappa) 
     64            ENDDO 
     65         ENDDO 
    6266        
    6367       ! for other layer    
    64        DO l = llm-1 , 2 , -1 
     68         DO l = llm-1 , 2 , -1 
     69            DO j=jj_begin-offset,jj_end+offset 
     70               DO i=ii_begin-offset,ii_end+offset 
     71                  ij=(j-1)*iim+i 
     72                  delta = p(ij,l)* (1+2*kappa) + p(ij,l+1)* ( beta(ij,l+1)- (1+2*kappa) ) 
     73                  alpha(ij,l)  = - p(ij,l+1) / delta * alpha(ij,l+1) 
     74                  beta (ij,l)  =   p(ij,l  ) / delta    
     75               ENDDO 
     76            ENDDO 
     77         ENDDO 
     78        
     79         !! Compute pk 
     80        
     81         ! for first layer 
     82         DO j=jj_begin-offset,jj_end+offset 
     83            DO i=ii_begin-offset,ii_end+offset 
     84               ij=(j-1)*iim+i 
     85               pks(ij) = cpp * ( ps(ij)/preff ) ** kappa 
     86               pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )     /    & 
     87                    (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-(1.+2*kappa) )* p(ij,2) ) 
     88            ENDDO 
     89         ENDDO 
     90        
     91       ! for other layers 
     92        
     93         DO l = 2, llm 
     94            DO j=jj_begin-offset,jj_end+offset 
     95               DO i=ii_begin-offset,ii_end+offset 
     96                  ij=(j-1)*iim+i 
     97                  pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1) 
     98               ENDDO 
     99            ENDDO 
     100         ENDDO 
     101 
     102       ENDIF 
     103 
     104    ELSE ! Simple calculation of Exner pressure based on centered average 
     105       ! surface : pks 
     106       IF (is_omp_level_master) THEN 
     107 
     108         DO j=jj_begin-offset,jj_end+offset 
     109            DO i=ii_begin-offset,ii_end+offset 
     110               ij=(j-1)*iim+i 
     111               pks(ij) = cpp * ( ps(ij)/preff ) ** kappa 
     112            ENDDO 
     113         ENDDO 
     114 
     115       ENDIF 
     116 
     117         ! 3D : pk 
     118       DO l = 1, llm 
    65119          DO j=jj_begin-offset,jj_end+offset 
    66120             DO i=ii_begin-offset,ii_end+offset 
    67                 ij=(j-1)*iim+i 
    68                 delta = p(ij,l)* (1+2*kappa) + p(ij,l+1)* ( beta(ij,l+1)- (1+2*kappa) ) 
    69                 alpha(ij,l)  = - p(ij,l+1) / delta * alpha(ij,l+1) 
    70                 beta (ij,l)  =   p(ij,l  ) / delta    
     121               ij=(j-1)*iim+i 
     122               pk(ij,l) = cpp * ((.5/preff)*(p(ij,l)+p(ij,l+1))) ** kappa 
    71123             ENDDO 
    72124          ENDDO 
    73125       ENDDO 
    74         
    75        !! Compute pk 
    76         
    77        ! for first layer 
    78        DO j=jj_begin-offset,jj_end+offset 
    79           DO i=ii_begin-offset,ii_end+offset 
    80              ij=(j-1)*iim+i 
    81              pks(ij) = cpp * ( ps(ij)/preff ) ** kappa 
    82              pk(ij,1) = ( p(ij,1)*pks(ij) - 0.5*alpha(ij,2)*p(ij,2) )     /    & 
    83                   (  p(ij,1)* (1.+kappa) + 0.5*( beta(ij,2)-(1.+2*kappa) )* p(ij,2) ) 
    84           ENDDO 
    85        ENDDO 
    86         
    87        ! for other layers 
    88         
    89        DO l = 2, llm 
    90           DO j=jj_begin-offset,jj_end+offset 
    91              DO i=ii_begin-offset,ii_end+offset 
    92                 ij=(j-1)*iim+i 
    93                 pk(ij,l) = alpha(ij,l) + beta(ij,l) * pk(ij,l-1) 
    94              ENDDO 
    95           ENDDO 
    96        ENDDO 
    97         
    98     ELSE ! Simple calculation of Exner pressure based on centered average 
    99        ! surface : pks 
    100        DO j=jj_begin-offset,jj_end+offset 
    101           DO i=ii_begin-offset,ii_end+offset 
    102              ij=(j-1)*iim+i 
    103              pks(ij) = cpp * ( ps(ij)/preff ) ** kappa 
    104           ENDDO 
    105        ENDDO 
    106126 
    107        ! 3D : pk 
    108        DO l = 1, llm 
    109           DO j=jj_begin-offset,jj_end+offset 
    110              DO i=ii_begin-offset,ii_end+offset 
    111                 ij=(j-1)*iim+i 
    112                 pk(ij,l) = cpp * ((.5/preff)*(p(ij,l)+p(ij,l+1))) ** kappa 
    113              ENDDO 
    114           ENDDO 
    115        ENDDO 
    116127    END IF 
    117128     

  • codes/icosagcm/trunk/src/field.f90

    r266 r295  
    4848  SUBROUTINE allocate_field(field,field_type,data_type,dim1,dim2,name) 
    4949  USE domain_mod 
     50  USE omp_para 
    5051  IMPLICIT NONE 
    5152    TYPE(t_field),POINTER :: field(:) 
    … …  
    6465 
    6566    DO ind=1,ndomain 
    66       IF (.NOT. assigned_domain(ind)) CYCLE 
     67      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE 
    6768 
    6869      IF(PRESENT(name)) THEN 
    … …  
    181182  SUBROUTINE deallocate_field(field) 
    182183  USE domain_mod 
     184  USE omp_para 
    183185  IMPLICIT NONE 
    184186    TYPE(t_field),POINTER :: field(:) 
    … …  
    188190!$OMP BARRIER 
    189191    DO ind=1,ndomain 
    190       IF (.NOT. assigned_domain(ind)) CYCLE 
     192      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE 
    191193 
    192194      data_type=field(ind)%data_type 

  • codes/icosagcm/trunk/src/geometry.f90

    r294 r295  
    166166    USE vector 
    167167    USE transfert_mod 
     168    USE omp_para 
    168169 
    169170    IMPLICIT NONE 
    … …  
    186187     
    187188    DO ind=1,ndomain 
    188       IF (.NOT. assigned_domain(ind)) CYCLE 
     189      IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    189190      CALL swap_dimensions(ind) 
    190191      CALL swap_geometry(ind) 
    … …  
    212213  USE vector 
    213214  USE getin_mod 
     215  USE omp_para 
    214216  IMPLICIT NONE 
    215217    INTEGER :: nb_it=0 
    … …  
    226228     
    227229    DO ind=1,ndomain 
    228       IF (.NOT. assigned_domain(ind)) CYCLE 
     230      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master) CYCLE 
    229231      d=>domain(ind) 
    230232      CALL swap_dimensions(ind) 
    … …  
    241243 
    242244    DO ind=1,ndomain 
    243       IF (.NOT. assigned_domain(ind)) CYCLE 
     245      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master ) CYCLE 
    244246      d=>domain(ind) 
    245247      CALL swap_dimensions(ind) 
    … …  
    271273      sum=0 
    272274      DO ind=1,ndomain 
    273       IF (.NOT. assigned_domain(ind)) CYCLE 
     275      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master ) CYCLE 
    274276        CALL swap_dimensions(ind) 
    275277        CALL swap_geometry(ind) 
    … …  
    311313  USE transfert_mod 
    312314  USE getin_mod 
     315  USE omp_para 
    313316  IMPLICIT NONE 
    314317 
    … …  
    330333  
    331334    DO ind=1,ndomain 
    332       IF (.NOT. assigned_domain(ind)) CYCLE 
     335      IF (.NOT. assigned_domain(ind)  .OR. .NOT. is_omp_level_master ) CYCLE 
    333336      d=>domain(ind) 
    334337      CALL swap_dimensions(ind) 
    … …  
    559562    CALL transfert_request(geom%centroid,req_i1) 
    560563 
    561     CALL surf_triangle(d%xyz(:,ii_begin,jj_begin),d%xyz(:,ii_begin,jj_end),d%xyz(:,ii_end,jj_begin),S) 
     564!    CALL surf_triangle(d%xyz(:,ii_begin,jj_begin),d%xyz(:,ii_begin,jj_end),d%xyz(:,ii_end,jj_begin),S) 
    562565  
    563566  END SUBROUTINE set_geometry 

  • codes/icosagcm/trunk/src/geopotential_mod.f90

    r186 r295  
    3535  SUBROUTINE compute_geopotential(phis,pks,pk,theta,phi,offset) 
    3636  USE icosa 
     37  USE omp_para 
    3738  IMPLICIT NONE 
    3839    REAL(rstd),INTENT(IN) :: phis(iim*jjm) 
    … …  
    5152    
    5253  ! for first layer 
    53    DO j=jj_begin-offset,jj_end+offset 
    54      DO i=ii_begin-offset,ii_end+offset 
    55        ij=(j-1)*iim+i 
    56        phi( ij,1 ) = phis( ij ) + theta(ij,1) * ( pks(ij) - pk(ij,1) ) 
    57      ENDDO 
    58    ENDDO 
    59            
    60   ! for other layers 
    61    DO l = 2, llm 
     54 
     55! flush pks, pk thetha, phis 
     56!$OMP BARRIER 
     57   IF(is_omp_level_master) THEN 
    6258     DO j=jj_begin-offset,jj_end+offset 
    6359       DO i=ii_begin-offset,ii_end+offset 
    6460         ij=(j-1)*iim+i 
    65          phi(ij,l) = phi(ij,l-1) + 0.5 * ( theta(ij,l)  + theta(ij,l-1) )  &  
    66                                        * (  pk(ij,l-1) -  pk(ij,l)    ) 
     61         phi( ij,1 ) = phis( ij ) + theta(ij,1) * ( pks(ij) - pk(ij,1) ) 
    6762       ENDDO 
    6863     ENDDO 
    69    ENDDO 
     64           
     65    ! for other layers 
     66     DO l = 2, llm 
     67       DO j=jj_begin-offset,jj_end+offset 
     68         DO i=ii_begin-offset,ii_end+offset 
     69           ij=(j-1)*iim+i 
     70           phi(ij,l) = phi(ij,l-1) + 0.5 * ( theta(ij,l)  + theta(ij,l-1) )  &  
     71                                         * (  pk(ij,l-1) -  pk(ij,l)    ) 
     72         ENDDO 
     73       ENDDO 
     74     ENDDO 
     75   ENDIF 
     76! flush phi 
     77!$OMP BARRIER 
    7078    
    7179  END SUBROUTINE compute_geopotential 

  • codes/icosagcm/trunk/src/icosa_gcm.f90

    r280 r295  
    1919  CALL xios_init 
    2020  CALL init_earth_const  
    21   CALL init_grid_param(is_mpi_root) 
    22   CALL init_omp_para 
     21  CALL init_grid_param(is_mpi_master) 
     22  CALL init_omp_para(is_mpi_master) 
    2323  CALL compute_metric 
    2424  CALL compute_domain 
    … …  
    2828 
    2929!$OMP PARALLEL   
     30  CALL switch_omp_no_distrib_level 
    3031  CALL compute_geometry 
    3132  CALL check_total_area 
    … …  
    4344   
    4445  CALL timeloop 
    45  
     46  CALL switch_omp_no_distrib_level 
    4647!$OMP END PARALLEL 
    4748 

  • codes/icosagcm/trunk/src/kinetic.f90

    r186 r295  
    3030  SUBROUTINE compute_kinetic(ue, Ki) 
    3131  USE icosa 
     32  USE omp_para 
    3233  IMPLICIT NONE 
    3334    REAL(rstd),INTENT(IN) :: ue(3*iim*jjm,llm) 
    … …  
    3536    INTEGER :: i,j,ij,l 
    3637     
    37     DO l=1,llm 
     38    DO l=ll_begin,ll_end 
    3839      DO j=jj_begin,jj_end 
    3940        DO i=ii_begin,ii_end 

  • codes/icosagcm/trunk/src/omega.f90

    r186 r295  
    2323  END SUBROUTINE W_omega 
    2424 
     25 
     26 
    2527  SUBROUTINE compute_omega(ps,u, w) 
    2628    USE disvert_mod, ONLY : ap,bp 
     29    USE omp_para 
     30    IMPLICIT NONE 
    2731    REAL(rstd),INTENT(IN) :: u(iim*3*jjm,llm), ps(iim*jjm) 
    2832    REAL(rstd),INTENT(OUT):: w(iim*jjm,llm) 
    … …  
    3034    REAL(rstd):: p(iim*jjm,llm+1), rhodz(iim*jjm,llm), Fe(iim*3*jjm,llm) 
    3135    REAL(rstd):: ugradps 
    32     DO    l    = 1, llm+1 
    33        DO j=jj_begin-1,jj_end+1 
    34           DO i=ii_begin-1,ii_end+1 
    35              ij=(j-1)*iim+i 
    36              p(ij,l) = ap(l) + bp(l) * ps(ij) 
    37           ENDDO 
    38        ENDDO 
    39     ENDDO 
     36     
     37    INTEGER :: i,j,l,ij 
     38 
     39!$OMP BARRIER     
     40    IF (is_omp_level_master) THEN 
     41      DO    l    = 1, llm+1 
     42         DO j=jj_begin-1,jj_end+1 
     43           DO i=ii_begin-1,ii_end+1 
     44               ij=(j-1)*iim+i 
     45               p(ij,l) = ap(l) + bp(l) * ps(ij) 
     46            ENDDO 
     47         ENDDO 
     48      ENDDO 
    4049     
    4150!!! Compute mass 
    42     DO l = 1, llm 
    43        DO j=jj_begin-1,jj_end+1 
    44           DO i=ii_begin-1,ii_end+1 
    45              ij=(j-1)*iim+i 
    46              rhodz(ij,l) = ( p(ij,l) - p(ij,l+1) ) / g 
    47           ENDDO 
    48        ENDDO 
    49     ENDDO 
     51      DO l = 1, llm 
     52         DO j=jj_begin-1,jj_end+1 
     53            DO i=ii_begin-1,ii_end+1 
     54               ij=(j-1)*iim+i 
     55               rhodz(ij,l) = ( p(ij,l) - p(ij,l+1) ) / g 
     56            ENDDO 
     57         ENDDO 
     58      ENDDO 
    5059     
    5160!!!  Compute mass flux 
    52     DO l = 1, llm 
    53        DO j=jj_begin-1,jj_end+1 
    54           DO i=ii_begin-1,ii_end+1 
    55              ij=(j-1)*iim+i 
    56              Fe(ij+u_right,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_right,l))*u(ij+u_right,l)*le(ij+u_right) 
    57              Fe(ij+u_lup,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_lup,l))*u(ij+u_lup,l)*le(ij+u_lup) 
    58              Fe(ij+u_ldown,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_ldown,l))*u(ij+u_ldown,l)*le(ij+u_ldown) 
    59           ENDDO 
    60        ENDDO 
    61     ENDDO 
     61      DO l = 1, llm 
     62         DO j=jj_begin-1,jj_end+1 
     63            DO i=ii_begin-1,ii_end+1 
     64               ij=(j-1)*iim+i 
     65               Fe(ij+u_right,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_right,l))*u(ij+u_right,l)*le(ij+u_right) 
     66               Fe(ij+u_lup,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_lup,l))*u(ij+u_lup,l)*le(ij+u_lup) 
     67               Fe(ij+u_ldown,l)=0.5*(rhodz(ij,l)+rhodz(ij+t_ldown,l))*u(ij+u_ldown,l)*le(ij+u_ldown) 
     68            ENDDO 
     69         ENDDO 
     70      ENDDO 
    6271     
    6372!!! mass flux convergence computation   
    6473     
    6574    ! horizontal convergence 
    66     DO l = 1, llm 
    67        DO j=jj_begin,jj_end 
    68           DO i=ii_begin,ii_end 
    69              ij=(j-1)*iim+i 
    70              ! convm = +div(mass flux), sign convention as in Ringler et al. 2012, eq. 21 
    71              convm(ij,l)= 1./Ai(ij)*(ne(ij,right)*Fe(ij+u_right,l)   +  & 
    72                   ne(ij,rup)*Fe(ij+u_rup,l)       +  &   
    73                   ne(ij,lup)*Fe(ij+u_lup,l)       +  &   
    74                   ne(ij,left)*Fe(ij+u_left,l)     +  &   
    75                   ne(ij,ldown)*Fe(ij+u_ldown,l)   +  &   
    76                   ne(ij,rdown)*Fe(ij+u_rdown,l)) 
    77           ENDDO 
    78        ENDDO 
    79     ENDDO 
     75      DO l = 1, llm 
     76         DO j=jj_begin,jj_end 
     77            DO i=ii_begin,ii_end 
     78               ij=(j-1)*iim+i 
     79               ! convm = +div(mass flux), sign convention as in Ringler et al. 2012, eq. 21 
     80               convm(ij,l)= 1./Ai(ij)*(ne(ij,right)*Fe(ij+u_right,l)   +  & 
     81                    ne(ij,rup)*Fe(ij+u_rup,l)       +  &   
     82                    ne(ij,lup)*Fe(ij+u_lup,l)       +  &   
     83                    ne(ij,left)*Fe(ij+u_left,l)     +  &   
     84                    ne(ij,ldown)*Fe(ij+u_ldown,l)   +  &   
     85                    ne(ij,rdown)*Fe(ij+u_rdown,l)) 
     86            ENDDO 
     87         ENDDO 
     88      ENDDO 
    8089     
    81     ! vertical integration from up to down 
    82     DO  l = llm-1, 1, -1 
    83        DO j=jj_begin,jj_end 
    84           DO i=ii_begin,ii_end 
    85              ij=(j-1)*iim+i 
    86              convm(ij,l) = convm(ij,l) + convm(ij,l+1) 
    87           ENDDO 
    88        ENDDO 
    89     ENDDO 
    90     convm(:,llm+1)=0. 
     90      ! vertical integration from up to down 
     91      DO  l = llm-1, 1, -1 
     92         DO j=jj_begin,jj_end 
     93            DO i=ii_begin,ii_end 
     94               ij=(j-1)*iim+i 
     95               convm(ij,l) = convm(ij,l) + convm(ij,l+1) 
     96            ENDDO 
     97         ENDDO 
     98      ENDDO 
     99      convm(:,llm+1)=0. 
    91100 
    92101!!! Compute dps 
    … …  
    125134    ! -grad ps : ( ne(ij,ldown)*ps(ij,l) + ne(ij+t_ldown,rup)*ps(ij+t_ldown,l) ) ) / de(ij+u_ldown) 
    126135 
    127     DO l = 1,llm 
    128        DO j=jj_begin,jj_end 
    129           DO i=ii_begin,ii_end 
    130              toto = 1 
     136       DO l = 1,llm 
     137         DO j=jj_begin,jj_end  
     138           DO i=ii_begin,ii_end 
    131139             ij=(j-1)*iim+i 
    132140             ugradps = & 
    … …  
    140148             w( ij, l) =  ugradps - .5*(convm( ij,l+1)+convm(ij,l))  
    141149          ENDDO 
    142        ENDDO 
    143     ENDDO 
     150        ENDDO 
     151      ENDDO 
     152    ENDIF 
     153!$OMP BARRIER 
    144154     
    145155  END SUBROUTINE compute_omega 

  • codes/icosagcm/trunk/src/omp_para.F90

    r186 r295  
    55!$OMP THREADPRIVATE(omp_rank) 
    66 
    7   LOGICAL,SAVE :: omp_first 
    8   LOGICAL,SAVE :: omp_last 
    9   LOGICAL,SAVE :: omp_master 
    10 !$OMP THREADPRIVATE(omp_first, omp_last,omp_master) 
     7  LOGICAL,SAVE :: is_omp_first_level 
     8  LOGICAL,SAVE :: is_omp_last_level 
     9  LOGICAL,SAVE :: is_omp_master 
     10!$OMP THREADPRIVATE(is_omp_first_level, is_omp_last_level,is_omp_master) 
    1111 
    1212  INTEGER,SAVE :: ll_begin 
    … …  
    1717!$OMP THREADPRIVATE(ll_begin,ll_beginp1,ll_end,ll_endm1,ll_endp1) 
    1818  LOGICAL,SAVE :: using_openmp 
    19  
     19   
     20  INTEGER,SAVE :: omp_domain_size 
     21  INTEGER,SAVE :: omp_domain_rank 
     22  INTEGER,SAVE :: omp_level_size 
     23  INTEGER,SAVE :: omp_level_rank 
     24!$OMP THREADPRIVATE( omp_domain_size, omp_level_size,omp_domain_rank,omp_level_rank) 
     25  LOGICAL,SAVE :: is_omp_domain_master 
     26  LOGICAL,SAVE :: is_omp_level_master 
     27!$OMP THREADPRIVATE(is_omp_domain_master,is_omp_level_master ) 
     28   
    2029  LOGICAL,PARAMETER :: omp_by_domain=.TRUE.  
     30  LOGICAL,SAVE :: is_master 
     31!$OMP THREADPRIVATE(is_master) 
     32 
     33 
     34   LOGICAL,SAVE :: is_omp_first_level_full 
     35   LOGICAL,SAVE :: is_omp_last_level_full 
     36   INTEGER,SAVE :: ll_begin_full 
     37   INTEGER,SAVE :: ll_beginp1_full 
     38   INTEGER,SAVE :: ll_end_full 
     39   INTEGER,SAVE :: ll_endm1_full 
     40   INTEGER,SAVE :: ll_endp1_full 
     41!$OMP THREADPRIVATE(is_omp_first_level_full,is_omp_last_level_full)    
     42!$OMP THREADPRIVATE( ll_begin_full, ll_beginp1_full, ll_end_full, ll_endm1_full, ll_endp1_full) 
     43  PRIVATE :: is_omp_first_level_full,is_omp_last_level_full 
     44  PRIVATE :: ll_begin_full, ll_beginp1_full, ll_end_full, ll_endm1_full, ll_endp1_full  
     45   
     46 
     47   LOGICAL,SAVE :: is_omp_first_level_distrib 
     48   LOGICAL,SAVE :: is_omp_last_level_distrib 
     49   INTEGER,SAVE :: ll_begin_distrib 
     50   INTEGER,SAVE :: ll_beginp1_distrib 
     51   INTEGER,SAVE :: ll_end_distrib 
     52   INTEGER,SAVE :: ll_endm1_distrib 
     53   INTEGER,SAVE :: ll_endp1_distrib 
     54!$OMP THREADPRIVATE(is_omp_first_level_distrib,is_omp_last_level_distrib)    
     55!$OMP THREADPRIVATE( ll_begin_distrib, ll_beginp1_distrib, ll_end_distrib, ll_endm1_distrib, ll_endp1_distrib) 
     56 
     57  PRIVATE :: is_omp_first_level_distrib,is_omp_last_level_distrib 
     58  PRIVATE :: ll_begin_distrib, ll_beginp1_distrib, ll_end_distrib, ll_endm1_distrib, ll_endp1_distrib  
     59 
    2160 
    2261CONTAINS 
    2362 
    2463 
    25   SUBROUTINE init_omp_para 
     64  SUBROUTINE init_omp_para(is_mpi_master) 
    2665  USE grid_param 
     66  USE getin_mod 
    2767#ifdef CPP_USING_OMP 
    2868  USE omp_lib 
    2969#endif 
    3070  IMPLICIT NONE 
     71  LOGICAL, INTENT(IN) :: is_mpi_master 
    3172  INTEGER :: ll_nb,i 
    3273 
    … …  
    4990    omp_rank=OMP_GET_THREAD_NUM() 
    5091#endif 
    51  
    52     IF (omp_by_domain) THEN 
    53       omp_first=.TRUE. 
    54       omp_last=.TRUE. 
    55       IF (omp_rank==0) THEN 
    56         omp_master=.TRUE. 
    57       ELSE 
    58         omp_master=.FALSE. 
    59       ENDIF 
    60       
    61       ll_begin=1 
    62       ll_beginp1=2 
    63       ll_end=llm 
    64       ll_endm1=llm-1 
    65       ll_endp1=llm+1 
    66       
    67     ELSE     
    68      
    69       omp_first=.FALSE. 
    70       omp_last=.FALSE. 
    71       omp_master=.FALSE. 
    72      
    73       IF (omp_rank==0) THEN 
    74         omp_first=.TRUE. 
    75         omp_master=.TRUE. 
    76       ENDIF 
    77      
    78       IF (omp_rank==omp_size-1) omp_last=.TRUE. 
    79      
    80       ll_end=0 
    81       DO i=0,omp_rank 
    82         ll_begin=ll_end+1 
    83         ll_nb=llm/omp_size 
    84         IF (MOD(llm,omp_size)>i) ll_nb=ll_nb+1 
    85         ll_end=ll_begin+ll_nb-1 
    86       ENDDO 
    87      
    88       ll_beginp1=ll_begin 
    89       ll_endp1=ll_end 
    90       ll_endm1=ll_end 
    91  
    92       IF (omp_first) ll_beginp1=ll_begin+1 
    93       IF (omp_last) ll_endp1=ll_endp1+1 
    94       IF (omp_last) ll_endm1=ll_endm1-1 
    95      
     92     
     93    is_omp_master=.FALSE. 
     94    is_master=.FALSE. 
     95       
     96    IF (omp_rank==0) THEN 
     97      is_omp_master=.TRUE. 
     98      IF (is_mpi_master) is_master=.TRUE. 
    9699    ENDIF 
     100 
     101    omp_level_size=1  
     102    CALL getin("omp_level_size",omp_level_size) 
     103    IF (MOD(omp_size,omp_level_size)/=0) THEN 
     104      IF (is_mpi_master) PRINT*,"omp_size /= omp_level_size x omp_domain_size  => disable omp threads on vertical layers" 
     105      omp_level_size=1 
     106    ENDIF 
     107    omp_domain_size=omp_size/omp_level_size 
     108    omp_domain_rank = omp_rank / omp_level_size     
     109    omp_level_rank = MOD(omp_rank, omp_level_size)     
     110     
     111    IF (is_mpi_master) PRINT*,"omp_domain_size",omp_domain_size,"omp_domain_rank", omp_domain_rank 
     112    IF (is_mpi_master) PRINT*,"omp_level_size",omp_level_size,"omp_level_rank", omp_level_rank 
     113     
     114    is_omp_first_level=.FALSE. 
     115    is_omp_last_level= .FALSE. 
     116    is_omp_domain_master=.FALSE. 
     117    is_omp_level_master=.FALSE. 
     118 
     119    IF (omp_domain_rank==0) is_omp_domain_master = .TRUE. 
     120    IF (omp_level_rank==0) is_omp_level_master = .TRUE. 
     121    IF (omp_level_rank==0) is_omp_first_level=.TRUE. 
     122     
     123    IF (omp_level_rank==omp_level_size-1) is_omp_last_level=.TRUE. 
     124     
     125    ll_end=0 
     126    DO i=0,omp_level_rank 
     127      ll_begin=ll_end+1 
     128      ll_nb=llm/omp_level_size 
     129      IF (MOD(llm,omp_level_size)>i) ll_nb=ll_nb+1 
     130      ll_end=ll_begin+ll_nb-1 
     131    ENDDO 
     132     
     133    ll_beginp1=ll_begin 
     134    ll_endp1=ll_end 
     135    ll_endm1=ll_end 
     136 
     137    IF (is_omp_first_level) ll_beginp1=ll_begin+1 
     138    IF (is_omp_last_level) ll_endp1=ll_endp1+1 
     139    IF (is_omp_last_level) ll_endm1=ll_endm1-1 
     140     
     141     
     142     
     143    is_omp_first_level_distrib = is_omp_first_level 
     144    is_omp_last_level_distrib  = is_omp_last_level 
     145    ll_begin_distrib           = ll_begin 
     146    ll_beginp1_distrib         = ll_beginp1 
     147    ll_end_distrib             = ll_end 
     148    ll_endm1_distrib           = ll_endm1 
     149    ll_endp1_distrib           = ll_endp1 
     150     
     151    is_omp_first_level_full = .TRUE. 
     152    is_omp_last_level_full  = .TRUE. 
     153    ll_begin_full           = 1 
     154    ll_beginp1_full         = 2 
     155    ll_end_full             = llm 
     156    ll_endm1_full           = llm-1 
     157    ll_endp1_full           = llm+1     
     158 
    97159!$OMP END PARALLEL 
    98160 
    99161   ELSE 
    100162     omp_size=1 
     163     omp_level_size=1 
     164     omp_domain_size=1 
    101165     omp_rank=0 
    102      omp_first=.TRUE. 
    103      omp_last=.TRUE. 
    104      omp_master=.TRUE. 
     166     omp_level_rank=0 
     167     omp_domain_rank=0 
     168     is_omp_first_level=.TRUE. 
     169     is_omp_last_level=.TRUE. 
     170     is_omp_master=.TRUE. 
     171     is_omp_domain_master=.TRUE. 
     172     is_omp_level_master=.TRUE. 
    105173     ll_begin=1 
    106174     ll_beginp1=2 
    … …  
    108176     ll_endm1=llm-1 
    109177     ll_endp1=llm+1 
     178      
     179     is_omp_first_level_distrib = is_omp_first_level 
     180     is_omp_last_level_distrib  = is_omp_last_level 
     181     ll_begin_distrib           = ll_begin 
     182     ll_beginp1_distrib         = ll_beginp1 
     183     ll_end_distrib             = ll_end 
     184     ll_endm1_distrib           = ll_endm1 
     185     ll_endp1_distrib           = ll_endp1 
     186     
     187     is_omp_first_level_full = .TRUE. 
     188     is_omp_last_level_full  = .TRUE. 
     189     ll_begin_full           = 1 
     190     ll_beginp1_full         = 2 
     191     ll_end_full             = llm 
     192     ll_endm1_full           = llm-1 
     193     ll_endp1_full           = llm+1     
     194 
    110195   ENDIF 
    111196    
    112197  END SUBROUTINE init_omp_para 
    113198 
     199 
     200  SUBROUTINE switch_omp_distrib_level 
     201  IMPLICIT NONE 
     202     is_omp_first_level = is_omp_first_level_distrib  
     203     is_omp_last_level  = is_omp_last_level_distrib  
     204     ll_begin           = ll_begin_distrib  
     205     ll_beginp1         = ll_beginp1_distrib  
     206     ll_end             = ll_end_distrib  
     207     ll_endm1           = ll_endm1_distrib  
     208     ll_endp1           = ll_endp1_distrib  
     209   
     210  END SUBROUTINE switch_omp_distrib_level 
     211   
     212 
     213  SUBROUTINE switch_omp_no_distrib_level 
     214  IMPLICIT NONE 
     215 
     216     is_omp_first_level = is_omp_first_level_full  
     217     is_omp_last_level  = is_omp_last_level_full  
     218     ll_begin           = ll_begin_full  
     219     ll_beginp1         = ll_beginp1_full  
     220     ll_end             = ll_end_full  
     221     ll_endm1           = ll_endm1_full  
     222     ll_endp1           = ll_endp1_full  
     223   
     224  END SUBROUTINE switch_omp_no_distrib_level 
     225   
     226 
    114227  FUNCTION omp_in_parallel() 
    115228#ifdef CPP_USING_OMP 

  • codes/icosagcm/trunk/src/physics.f90

    r281 r295  
    1010  TYPE(t_field),POINTER :: f_extra_physics_2D(:), f_extra_physics_3D(:) 
    1111  TYPE(t_field),POINTER :: f_dulon(:), f_dulat(:) 
     12  TYPE(t_field),POINTER :: f_temp(:) 
    1213 
    1314  CHARACTER(LEN=255) :: physics_type 
    … …  
    3839       CALL allocate_field(f_dulon,field_t,type_real,llm, name='dulon') 
    3940       CALL allocate_field(f_dulat,field_t,type_real,llm, name='dulat') 
     41       CALL allocate_field(f_temp,field_t,type_real,llm, name='temp') 
    4042       CALL init_pack_before ! Compute physics_inout%ngrid and offsets used by pack/unpack 
    4143       CALL init_physics_dcmip 
    … …  
    102104    USE physics_interface_mod 
    103105    USE physics_dcmip_mod, ONLY : full_physics_dcmip => full_physics 
     106    USE theta2theta_rhodz_mod 
    104107    USE mpipara 
    105108    IMPLICIT NONE 
    … …  
    111114    REAL(rstd),POINTER :: phis(:) 
    112115    REAL(rstd),POINTER :: ps(:) 
     116    REAL(rstd),POINTER :: temp(:,:) 
    113117    REAL(rstd),POINTER :: theta_rhodz(:,:) 
    114118    REAL(rstd),POINTER :: ue(:,:) 
    … …  
    118122    INTEGER :: it, ind 
    119123 
     124    CALL theta_rhodz2temperature(f_ps,f_theta_rhodz,f_temp) 
     125 
    120126    DO ind=1,ndomain 
    121127       IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    124130       phis=f_phis(ind) 
    125131       ps=f_ps(ind) 
    126        theta_rhodz=f_theta_rhodz(ind) 
     132       temp=f_temp(ind) 
    127133       ue=f_ue(ind) 
    128134       q=f_q(ind) 
    129        CALL pack_physics(pack_info(ind), phis, ps, theta_rhodz, ue, q) 
     135       CALL pack_physics(pack_info(ind), phis, ps, temp, ue, q) 
    130136    END DO 
    131137 
    … …  
    134140       CALL full_physics_dcmip 
    135141    CASE DEFAULT 
    136        IF(is_mpi_root) PRINT *,'Internal error : illegal value of phys_type', phys_type 
     142       IF(is_mpi_master) PRINT *,'Internal error : illegal value of phys_type', phys_type 
    137143       STOP 
    138144    END SELECT 
    … …  
    143149       CALL swap_geometry(ind) 
    144150       ps=f_ps(ind) 
    145        theta_rhodz=f_theta_rhodz(ind) 
     151       temp=f_temp(ind) 
    146152       q=f_q(ind) 
    147153       dulon=f_dulon(ind) 
    148154       dulat=f_dulat(ind) 
    149        CALL unpack_physics(pack_info(ind), ps, theta_rhodz, q, dulon, dulat) 
     155       CALL unpack_physics(pack_info(ind), ps, temp, q, dulon, dulat) 
    150156    END DO 
     157    CALL temperature2theta_rhodz(f_ps,f_temp,f_theta_rhodz) 
    151158 
    152159    ! Transfer dulon, dulat 
    … …  
    166173  END SUBROUTINE physics_column 
    167174 
    168   SUBROUTINE pack_physics(info, phis, ps, theta_rhodz, ue, q) 
     175  SUBROUTINE pack_physics(info, phis, ps, temp, ue, q) 
    169176    USE icosa 
    170177    USE wind_mod 
    … …  
    176183    REAL(rstd) :: phis(iim*jjm) 
    177184    REAL(rstd) :: ps(iim*jjm) 
    178     REAL(rstd) :: theta_rhodz(iim*jjm,llm) 
     185    REAL(rstd) :: temp(iim*jjm,llm) 
    179186    REAL(rstd) :: ue(3*iim*jjm,llm) 
    180187    REAL(rstd) :: q(iim*jjm,llm,nqtot) 
    181188 
    182189    REAL(rstd) :: p(iim*jjm,llm+1) 
    183     REAL(rstd) :: Temp(iim*jjm,llm) 
    184190    REAL(rstd) :: uc(iim*jjm,3,llm) 
    185191    REAL(rstd) :: ulon(iim*jjm,llm) 
    186192    REAL(rstd) :: ulat(iim*jjm,llm) 
    187193 
     194!$OMP BARRIER 
    188195    CALL compute_pression(ps,p,0) 
    189     CALL compute_theta_rhodz2temperature(ps,theta_rhodz,Temp,0) 
     196!$OMP BARRIER 
    190197    CALL compute_wind_centered(ue,uc) 
    191198    CALL compute_wind_centered_lonlat_compound(uc, ulon, ulat) 
    … …  
    199206  END SUBROUTINE pack_physics 
    200207 
    201   SUBROUTINE unpack_physics(info, ps,theta_rhodz, q, dulon, dulat) 
     208  SUBROUTINE unpack_physics(info, ps,temp, q, dulon, dulat) 
    202209    USE icosa 
    203210    USE physics_interface_mod 
    … …  
    206213    TYPE(t_pack_info) :: info 
    207214    REAL(rstd) :: ps(iim*jjm) 
    208     REAL(rstd) :: theta_rhodz(iim*jjm,llm) 
    209     REAL(rstd) :: Temp(iim*jjm,llm) 
     215    REAL(rstd) :: temp(iim*jjm,llm) 
    210216    REAL(rstd) :: q(iim*jjm,llm,nqtot) 
    211217    REAL(rstd) :: dulon(iim*jjm,llm) 
    … …  
    221227    q = q + physics_inout%dt_phys * dq 
    222228    Temp = Temp + physics_inout%dt_phys * dTemp 
    223     CALL compute_temperature2theta_rhodz(ps,Temp,theta_rhodz,0) 
     229!    CALL compute_temperature2theta_rhodz(ps,Temp,theta_rhodz,0) 
    224230  END SUBROUTINE unpack_physics 
    225231 

  • codes/icosagcm/trunk/src/pression.f90

    r186 r295  
    1313    INTEGER :: ind 
    1414 
     15!$OMP BARRIER 
    1516    DO ind=1,ndomain 
    1617      IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    2122      CALL compute_pression(ps, p,0) 
    2223    ENDDO 
     24!$OMP BARRIER 
    2325   
    2426  END SUBROUTINE pression 
    … …  
    2729  USE icosa 
    2830  USE disvert_mod 
     31  USE omp_para 
    2932  IMPLICIT NONE 
    3033    REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
    … …  
    3437 
    3538    IF(ap_bp_present) THEN 
    36     DO    l    = 1, llm+1 
    37       DO j=jj_begin-offset,jj_end+offset 
    38         DO i=ii_begin-offset,ii_end+offset 
    39           ij=(j-1)*iim+i 
    40           p(ij,l) = ap(l) + bp(l) * ps(ij) 
     39      DO    l    = ll_begin, ll_endp1 
     40!      DO    l    = 1, llm + 1 
     41        DO j=jj_begin-offset,jj_end+offset 
     42          DO i=ii_begin-offset,ii_end+offset 
     43            ij=(j-1)*iim+i 
     44            p(ij,l) = ap(l) + bp(l) * ps(ij) 
     45          ENDDO 
    4146        ENDDO 
    4247      ENDDO 
    43     ENDDO 
    4448    END IF 
     49 
    4550  END SUBROUTINE compute_pression 
    4651 

  • codes/icosagcm/trunk/src/theta_rhodz.f90

    r186 r295  
    11MODULE theta2theta_rhodz_mod 
     2  USE field_mod 
     3   
     4  TYPE(t_field), POINTER, SAVE  :: f_p(:) 
     5  TYPE(t_field), POINTER, SAVE  :: f_pks(:) 
     6  TYPE(t_field), POINTER, SAVE  :: f_pk(:) 
     7 
     8  PRIVATE :: f_p,f_pk,f_pks  
    29 
    310CONTAINS 
    411   
     12  SUBROUTINE init_theta2theta_rhodz 
     13  USE icosa 
     14  USE field_mod 
     15  IMPLICIT NONE 
     16    CALL allocate_field(f_p,field_t,type_real,llm+1,name='p') 
     17    CALL allocate_field(f_pk,field_t,type_real,llm,name='pk') 
     18    CALL allocate_field(f_pks,field_t,type_real,name='pks') 
     19     
     20  END SUBROUTINE init_theta2theta_rhodz 
     21 
     22 
     23 
    524  SUBROUTINE theta_rhodz2theta(f_ps,f_theta_rhodz,f_theta) 
    625  USE icosa 
    … …  
    1534    INTEGER :: ind 
    1635 
     36!$OMP BARRIER 
    1737    DO ind=1,ndomain 
    1838      IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    2444      CALL compute_theta_rhodz2theta(ps, theta_rhodz,theta,0) 
    2545    ENDDO 
     46!$OMP BARRIER 
    2647   
    2748  END SUBROUTINE theta_rhodz2theta 
    … …  
    2950  SUBROUTINE theta_rhodz2temperature(f_ps,f_theta_rhodz,f_temp) 
    3051  USE icosa 
     52  USE pression_mod 
     53  USE exner_mod 
    3154  IMPLICIT NONE 
    3255    TYPE(t_field), POINTER :: f_ps(:) 
    … …  
    3760    REAL(rstd), POINTER :: theta_rhodz(:,:) 
    3861    REAL(rstd), POINTER :: temp(:,:) 
    39     INTEGER :: ind 
    40  
    41     DO ind=1,ndomain 
    42       IF (.NOT. assigned_domain(ind)) CYCLE 
    43       CALL swap_dimensions(ind) 
    44       CALL swap_geometry(ind) 
    45       ps=f_ps(ind) 
     62    REAL(rstd), POINTER :: p(:) 
     63    REAL(rstd), POINTER :: pk(:,:) 
     64    REAL(rstd), POINTER :: pks(:,:) 
     65    INTEGER :: ind 
     66 
     67    DO ind=1,ndomain 
     68      IF (.NOT. assigned_domain(ind)) CYCLE 
     69      CALL swap_dimensions(ind) 
     70      CALL swap_geometry(ind) 
     71      ps=f_ps(ind) 
     72      p=f_p(ind) 
     73      pks=f_pks(ind) 
     74      pk=f_pk(ind) 
    4675      theta_rhodz=f_theta_rhodz(ind) 
    4776      temp=f_temp(ind) 
    48       CALL compute_theta_rhodz2temperature(ps, theta_rhodz,temp,0) 
    49     ENDDO 
     77 
     78!$OMP BARRIER 
     79      CALL compute_pression(ps,p,0) 
     80!$OMP BARRIER 
     81      CALL compute_exner(ps,p,pks,pk,0) 
     82!$OMP BARRIER 
     83      CALL compute_theta_rhodz2temperature(p, pk, theta_rhodz,temp,0) 
     84    ENDDO 
     85!$OMP BARRIER 
    5086   
    5187  END SUBROUTINE theta_rhodz2temperature 
     88  
     89  SUBROUTINE temperature2theta_rhodz(f_ps,f_temp,f_theta_rhodz) 
     90  USE icosa 
     91  USE pression_mod 
     92  USE exner_mod 
     93  IMPLICIT NONE 
     94    TYPE(t_field), POINTER :: f_ps(:) 
     95    TYPE(t_field), POINTER :: f_theta_rhodz(:) 
     96    TYPE(t_field), POINTER :: f_temp(:) 
     97   
     98    REAL(rstd), POINTER :: ps(:) 
     99    REAL(rstd), POINTER :: theta_rhodz(:,:) 
     100    REAL(rstd), POINTER :: temp(:,:) 
     101    REAL(rstd), POINTER :: p(:) 
     102    REAL(rstd), POINTER :: pk(:,:) 
     103    REAL(rstd), POINTER :: pks(:,:) 
     104    INTEGER :: ind 
     105 
     106    DO ind=1,ndomain 
     107      IF (.NOT. assigned_domain(ind)) CYCLE 
     108      CALL swap_dimensions(ind) 
     109      CALL swap_geometry(ind) 
     110      ps=f_ps(ind) 
     111      p=f_p(ind) 
     112      pks=f_pks(ind) 
     113      pk=f_pk(ind) 
     114      theta_rhodz=f_theta_rhodz(ind) 
     115      temp=f_temp(ind) 
     116 
     117!$OMP BARRIER 
     118      CALL compute_pression(ps,p,0) 
     119!$OMP BARRIER 
     120      CALL compute_exner(ps,p,pks,pk,0) 
     121!$OMP BARRIER 
     122      CALL compute_temperature2theta_rhodz(p, pk, temp, theta_rhodz, 0) 
     123    ENDDO 
     124!$OMP BARRIER 
     125   
     126  END SUBROUTINE temperature2theta_rhodz 
     127  
     128  
    52129     
    53130  SUBROUTINE theta2theta_rhodz(f_ps,f_theta,f_theta_rhodz) 
    … …  
    63140    INTEGER :: ind 
    64141 
     142!$OMP BARRIER 
    65143    DO ind=1,ndomain 
    66144      IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    72150      CALL compute_theta2theta_rhodz(ps, theta, theta_rhodz,0) 
    73151    ENDDO 
     152!$OMP BARRIER 
    74153   
    75154  END SUBROUTINE theta2theta_rhodz 
    … …  
    77156  SUBROUTINE compute_theta2theta_rhodz(ps,theta, theta_rhodz,offset) 
    78157  USE icosa 
    79   USE pression_mod 
     158  USE disvert_mod 
     159  USE omp_para 
    80160  IMPLICIT NONE 
    81161    REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
    … …  
    84164    INTEGER,INTENT(IN) :: offset 
    85165    INTEGER :: i,j,ij,l 
    86     REAL(rstd),ALLOCATABLE,SAVE :: p(:,:) 
    87 !$OMP THREADPRIVATE(p) 
    88  
    89     ALLOCATE( p(iim*jjm,llm+1)) 
    90     CALL compute_pression(ps,p,offset) 
    91      
    92     DO    l    = 1, llm 
    93       DO j=jj_begin-offset,jj_end+offset 
    94         DO i=ii_begin-offset,ii_end+offset 
    95           ij=(j-1)*iim+i 
    96           theta_rhodz(ij,l) = theta(ij,l) * (p(ij,l)-p(ij,l+1))/g 
    97         ENDDO 
    98       ENDDO 
    99     ENDDO 
    100  
    101     DEALLOCATE( p) 
     166     
     167!$OMP BARRIER 
     168    DO    l    = ll_begin, ll_end 
     169      DO j=jj_begin-offset,jj_end+offset 
     170        DO i=ii_begin-offset,ii_end+offset 
     171          ij=(j-1)*iim+i 
     172          theta_rhodz(ij,l) = theta(ij,l) * ( (ap(l)-ap(l+1)) + ( bp(l)- bp(l+1))* ps(ij) )/g 
     173        ENDDO 
     174      ENDDO 
     175    ENDDO 
     176!$OMP BARRIER 
     177 
    102178 
    103179  END SUBROUTINE compute_theta2theta_rhodz 
    … …  
    105181  SUBROUTINE compute_theta_rhodz2theta(ps,theta_rhodz,theta,offset) 
    106182  USE icosa 
    107   USE pression_mod 
     183  USE disvert_mod 
     184  USE omp_para 
    108185  IMPLICIT NONE 
    109186    REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
    … …  
    112189    INTEGER,INTENT(IN) :: offset 
    113190    INTEGER :: i,j,ij,l 
    114     REAL(rstd),SAVE,ALLOCATABLE :: p(:,:) 
    115 !$OMP THREADPRIVATE(p) 
    116  
    117     ALLOCATE( p(iim*jjm,llm+1)) 
    118  
    119     CALL compute_pression(ps,p,offset) 
    120      
    121     DO    l    = 1, llm 
    122       DO j=jj_begin-offset,jj_end+offset 
    123         DO i=ii_begin-offset,ii_end+offset 
    124           ij=(j-1)*iim+i 
    125           theta(ij,l) = theta_rhodz(ij,l) / ((p(ij,l)-p(ij,l+1))/g) 
    126         ENDDO 
    127       ENDDO 
    128     ENDDO 
    129  
    130     DEALLOCATE( p) 
     191 
     192!$OMP BARRIER 
     193    DO    l    = ll_begin, ll_end 
     194      DO j=jj_begin-offset,jj_end+offset 
     195        DO i=ii_begin-offset,ii_end+offset 
     196          ij=(j-1)*iim+i 
     197          theta(ij,l) = theta_rhodz(ij,l) / ( (ap(l)-ap(l+1)) + ( bp(l)- bp(l+1))* ps(ij) )/g 
     198        ENDDO 
     199      ENDDO 
     200    ENDDO 
     201!$OMP BARRIER 
     202 
    131203     
    132204  END SUBROUTINE compute_theta_rhodz2theta 
    133205 
    134   SUBROUTINE compute_theta_rhodz2temperature(ps,theta_rhodz,temp,offset) 
    135   USE icosa 
    136   USE pression_mod 
    137   USE exner_mod 
    138   IMPLICIT NONE 
    139     REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
     206 
     207 
     208 
     209 
     210 
     211  SUBROUTINE compute_theta_rhodz2temperature(p,pk,theta_rhodz,temp,offset) 
     212  USE icosa 
     213  USE pression_mod 
     214  USE exner_mod 
     215  USE omp_para 
     216  IMPLICIT NONE 
     217    REAL(rstd),INTENT(IN) :: p(iim*jjm,llm+1) 
     218    REAL(rstd),INTENT(IN) :: pk(iim*jjm,llm) 
    140219    REAL(rstd),INTENT(IN) :: theta_rhodz(iim*jjm,llm) 
    141220    REAL(rstd),INTENT(OUT) :: temp(iim*jjm,llm) 
    142221    INTEGER,INTENT(IN) :: offset 
    143222    INTEGER :: i,j,ij,l 
    144     REAL(rstd) :: p(iim*jjm,llm+1) 
    145     REAL(rstd) :: pk(iim*jjm,llm) 
    146     REAL(rstd) :: pks(iim*jjm) 
    147  
    148     CALL compute_pression(ps,p,offset) 
    149     CALL compute_exner(ps,p,pks,pk,offset) 
    150223         
    151     DO    l    = 1, llm 
     224! flush p 
     225!$OMP BARRIER 
     226    DO    l    = ll_begin, ll_end 
    152227      DO j=jj_begin-offset,jj_end+offset 
    153228        DO i=ii_begin-offset,ii_end+offset 
    … …  
    157232      ENDDO 
    158233    ENDDO 
     234!$OMP BARRIER 
     235     
    159236     
    160237  END SUBROUTINE compute_theta_rhodz2temperature 
    161238 
    162   SUBROUTINE compute_temperature2theta_rhodz(ps,temp,theta_rhodz,offset) 
    163   USE icosa 
    164   USE pression_mod 
    165   USE exner_mod 
    166   IMPLICIT NONE 
    167     REAL(rstd),INTENT(IN) :: ps(iim*jjm) 
     239  SUBROUTINE compute_temperature2theta_rhodz(p,pk,temp,theta_rhodz,offset) 
     240  USE icosa 
     241  USE pression_mod 
     242  USE exner_mod 
     243  USE omp_para 
     244  IMPLICIT NONE 
     245    REAL(rstd),INTENT(IN)  :: p(iim*jjm,llm+1) 
     246    REAL(rstd),INTENT(IN)  :: pk(iim*jjm,llm) 
    168247    REAL(rstd),INTENT(OUT) :: theta_rhodz(iim*jjm,llm) 
    169     REAL(rstd),INTENT(IN) :: temp(iim*jjm,llm) 
    170     INTEGER,INTENT(IN) :: offset 
    171     INTEGER :: i,j,ij,l 
    172     REAL(rstd) :: p(iim*jjm,llm+1) 
    173     REAL(rstd) :: pk(iim*jjm,llm) 
    174     REAL(rstd) :: pks(iim*jjm) 
    175  
    176     CALL compute_pression(ps,p,offset) 
    177     CALL compute_exner(ps,p,pks,pk,offset) 
     248    REAL(rstd),INTENT(IN)  :: temp(iim*jjm,llm) 
     249    INTEGER,INTENT(IN) :: offset 
     250    INTEGER :: i,j,ij,l 
     251 
    178252         
    179     DO    l    = 1, llm 
     253! flush p 
     254!$OMP BARRIER 
     255 
     256    DO    l    = ll_begin, ll_end 
    180257      DO j=jj_begin-offset,jj_end+offset 
    181258        DO i=ii_begin-offset,ii_end+offset 
    … …  
    185262      ENDDO 
    186263    ENDDO 
     264!$OMP BARRIER 
    187265     
    188266  END SUBROUTINE compute_temperature2theta_rhodz 

  • codes/icosagcm/trunk/src/time.f90

    r278 r295  
    1919  INTEGER,SAVE    :: itau_out, itau_adv, itau_dissip, itau_physics, itaumax 
    2020!$OMP THREADPRIVATE(itau_out, itau_adv, itau_dissip, itau_physics, itaumax)   
     21 
     22  INTEGER,SAVE    :: itau_check_conserv 
     23!$OMP THREADPRIVATE(itau_check_conserv)   
    2124   
    2225  INTEGER,SAVE :: day_step,ndays 
    … …  
    3841 
    3942  PUBLIC create_time_counter_header, update_time_counter, close_time_counter, init_time,  & 
    40          dt, write_period, itau_out, itau_adv, itau_dissip, itau_physics, itaumax, &  
    41 day_step,ndays,jD_ref,jH_ref,day_ini,day_end,annee_ref,day_ref,an, mois, jour,heure, & 
    42             calend,time_style,itau0 
     43         dt, write_period, itau_out, itau_adv, itau_dissip, itau_physics, itaumax,  & 
     44         itau_check_conserv,  &  
     45         day_step,ndays,jD_ref,jH_ref,day_ini,day_end,annee_ref,day_ref, & 
     46         an, mois, jour,heure, calend,time_style,itau0 
    4347 
    4448 
    … …  
    8286    itau_adv=1 
    8387    CALL getin('itau_adv',itau_adv) 
    84      
     88 
    8589    itau_dissip=1 
    8690    CALL getin('itau_dissip',itau_dissip) 
    87      
     91 
    8892    itau_physics=1 
    8993    CALL getin('itau_physics',itau_physics) 
    90      
     94 
     95    itau_check_conserv=HUGE(itau_check_conserv) 
     96    CALL getin('itau_check_conserv',itau_check_conserv) 
     97 
    9198    IF (is_mpi_root)  THEN 
    9299       PRINT *, 'itaumax=',itaumax 

  • codes/icosagcm/trunk/src/timeloop_gcm.f90

    r281 r295  
    3939  USE output_field_mod 
    4040  USE write_field 
     41  USE theta2theta_rhodz_mod 
    4142  IMPLICIT NONE 
    4243 
    … …  
    117118    END SELECT 
    118119 
    119  
     120    CALL init_theta2theta_rhodz 
    120121    CALL init_dissip 
    121122    CALL init_caldyn 
    … …  
    124125    CALL init_check_conserve 
    125126    CALL init_physics 
     127     
    126128 
    127129    CALL etat0(f_ps,f_mass,f_phis,f_theta_rhodz,f_u, f_q) 
    … …  
    157159  USE output_field_mod 
    158160  USE write_etat0_mod 
     161  USE checksum_mod 
    159162  IMPLICIT NONE   
    160163    REAL(rstd),POINTER :: q(:,:,:) 
    … …  
    172175    INTEGER :: stop_clock 
    173176    INTEGER :: rate_clock 
     177    INTEGER :: l 
    174178     
    175179     
    … …  
    178182!    CALL write_restart(f_ps,f_mass,f_phis,f_theta_rhodz,f_u,f_q)  
    179183     
     184    CALL switch_omp_distrib_level 
    180185    CALL caldyn_BC(f_phis, f_wflux) ! set constant values in first/last interfaces 
    181  
    182 !!$OMP BARRIER 
     186     
     187!$OMP BARRIER 
    183188  DO ind=1,ndomain 
    184189     IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    189194        CALL compute_rhodz(.TRUE., ps, rhodz) ! save rhodz for transport scheme before dynamics update ps 
    190195     ELSE 
    191         rhodz(:,:)=mass(:,:) 
     196       DO l=ll_begin,ll_end 
     197         rhodz(:,l)=mass(:,l) 
     198       ENDDO 
    192199     END IF 
    193200  END DO 
     201!$OMP BARRIER 
    194202  fluxt_zero=.TRUE. 
    195203 
    … …  
    224232    ENDIF 
    225233 
    226 !$OMP MASTER     
    227     IF (is_mpi_root) PRINT *,"It No :",It,"   t :",dt*It 
    228 !$OMP END MASTER     
     234    IF (is_master) PRINT *,"It No :",It,"   t :",dt*It 
     235 
    229236    IF (mod(it,itau_out)==0 ) THEN 
    230237      CALL update_time_counter(dt*it) 
    231238      CALL output_field("q",f_q) 
    232       CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)   
    233239    ENDIF 
    234240     
    … …  
    256262        
    257263       IF(caldyn_eta==eta_mass) THEN 
     264!ym flush ps 
     265!$OMP BARRIER 
    258266          DO ind=1,ndomain 
    259267             IF (.NOT. assigned_domain(ind)) CYCLE 
    … …  
    267275!       CALL wait_message(req_mass)   
    268276       CALL dissip(f_u,f_du,f_mass,f_phis, f_theta_rhodz,f_dtheta_rhodz) 
     277 
    269278!       CALL send_message(f_mass,req_mass) 
    270279!       CALL wait_message(req_mass)   
    … …  
    291300 
    292301    CALL physics(it,f_phis, f_ps, f_theta_rhodz, f_u,  f_q) 
     302 
     303    IF (MOD(it,itau_check_conserv)==0) THEN 
     304      CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)   
     305    ENDIF 
    293306     
    294307  ENDDO 
    … …  
    324337          IF(caldyn_eta==eta_mass) THEN ! update ps 
    325338             ps=f_ps(ind) ; dps=f_dps(ind) ;               
    326              IF (omp_first) THEN 
     339             IF (is_omp_first_level) THEN 
    327340!$SIMD 
    328341                DO ij=ij_begin,ij_end 
    … …  
    376389      ! if mass coordinate, deal with ps first on one core 
    377390      IF(caldyn_eta==eta_mass) THEN 
    378          IF (omp_first) THEN 
     391         IF (is_omp_first_level) THEN 
    379392 
    380393            DO ind=1,ndomain 

  • codes/icosagcm/trunk/src/transfert_mpi.f90

    r266 r295  
    11191119 
    11201120        DO ind=1,ndomain 
    1121           IF (.NOT. assigned_domain(ind)) CYCLE 
     1121          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    11221122           
    11231123          rval2d=>field(ind)%rval2d 
    … …  
    11521152         
    11531153        DO ind=1,ndomain 
    1154           IF (.NOT. assigned_domain(ind)) CYCLE 
     1154          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    11551155          rval2d=>field(ind)%rval2d 
    11561156          req=>message%request(ind)         
    … …  
    11921192       
    11931193        DO ind=1,ndomain 
    1194           IF (.NOT. assigned_domain(ind)) CYCLE 
     1194          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    11951195 
    11961196          dim3=size(field(ind)%rval3d,2) 
    … …  
    12261226          
    12271227        DO ind=1,ndomain 
    1228           IF (.NOT. assigned_domain(ind)) CYCLE 
     1228          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    12291229          dim3=size(field(ind)%rval3d,2) 
    12301230          rval3d=>field(ind)%rval3d 
    … …  
    12631263     
    12641264        DO ind=1,ndomain 
    1265           IF (.NOT. assigned_domain(ind)) CYCLE 
     1265          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master ) CYCLE 
    12661266 
    12671267          dim3=size(field(ind)%rval4d,2) 
    … …  
    13021302         
    13031303        DO ind=1,ndomain 
    1304           IF (.NOT. assigned_domain(ind)) CYCLE 
     1304          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    13051305           
    13061306          dim3=size(field(ind)%rval4d,2) 
    … …  
    14211421         
    14221422        DO ind=1,ndomain 
    1423           IF (.NOT. assigned_domain(ind)) CYCLE 
     1423          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    14241424           
    14251425          rval2d=>field(ind)%rval2d 
    … …  
    14521452         
    14531453        DO ind=1,ndomain 
    1454           IF (.NOT. assigned_domain(ind)) CYCLE 
     1454          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    14551455 
    14561456          rval3d=>field(ind)%rval3d 
    … …  
    14851485                 
    14861486        DO ind=1,ndomain 
    1487           IF (.NOT. assigned_domain(ind)) CYCLE 
     1487          IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE 
    14881488 
    14891489          rval4d=>field(ind)%rval4d 

  • codes/icosagcm/trunk/src/transfert_omp.f90

    r238 r295  
    3434  END INTERFACE  
    3535 
    36   PUBLIC bcast_omp, reduce_sum_omp, allreduce_sum_omp 
     36 
     37  INTERFACE reduce_max_omp 
     38    MODULE PROCEDURE reduce_max_omp_i,reduce_max_omp_i1,reduce_max_omp_i2,reduce_max_omp_i3,reduce_max_omp_i4, & 
     39                     reduce_max_omp_r,reduce_max_omp_r1,reduce_max_omp_r2,reduce_max_omp_r3,reduce_max_omp_r4 
     40  END INTERFACE  
     41 
     42  INTERFACE allreduce_max_omp 
     43    MODULE PROCEDURE allreduce_max_omp_i,allreduce_max_omp_i1,allreduce_max_omp_i2,allreduce_max_omp_i3,allreduce_max_omp_i4, & 
     44                     allreduce_max_omp_r,allreduce_max_omp_r1,allreduce_max_omp_r2,allreduce_max_omp_r3,allreduce_max_omp_r4 
     45  END INTERFACE  
     46 
     47 
     48  PUBLIC bcast_omp, reduce_sum_omp, allreduce_sum_omp, reduce_max_omp, allreduce_max_omp 
    3749   
    3850CONTAINS 
    … …  
    521533  END SUBROUTINE allreduce_sum_omp_r4   
    522534   
     535 
     536 
     537 
     538 
     539 
     540 
     541 
     542 
     543 
     544 
     545 
     546  SUBROUTINE reduce_max_omp_i(VarIn, VarOut) 
     547    IMPLICIT NONE 
     548   
     549    INTEGER,INTENT(IN)  :: VarIn 
     550    INTEGER,INTENT(OUT) :: VarOut 
     551    INTEGER             :: VarIn_tmp(1) 
     552    INTEGER             :: VarOut_tmp(1) 
     553     
     554    VarIn_tmp(1)=VarIn 
     555    CALL Check_buffer_i(1)    
     556    CALL reduce_max_omp_igen(VarIn_tmp,Varout_tmp,1,buffer_i) 
     557    VarOut=VarOut_tmp(1) 
     558     
     559  END SUBROUTINE reduce_max_omp_i 
     560 
     561  SUBROUTINE reduce_max_omp_i1(VarIn, VarOut) 
     562    IMPLICIT NONE 
     563   
     564    INTEGER,INTENT(IN),DIMENSION(:)  :: VarIn 
     565    INTEGER,INTENT(OUT),DIMENSION(:) :: VarOut 
     566     
     567    CALL Check_buffer_i(size(VarIn))    
     568    CALL reduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     569    
     570  END SUBROUTINE reduce_max_omp_i1 
     571   
     572   
     573  SUBROUTINE reduce_max_omp_i2(VarIn, VarOut) 
     574    IMPLICIT NONE 
     575   
     576    INTEGER,INTENT(IN),DIMENSION(:,:)  :: VarIn 
     577    INTEGER,INTENT(OUT),DIMENSION(:,:) :: VarOut 
     578 
     579    CALL Check_buffer_i(size(VarIn))    
     580    CALL reduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     581   
     582  END SUBROUTINE reduce_max_omp_i2 
     583 
     584 
     585  SUBROUTINE reduce_max_omp_i3(VarIn, VarOut) 
     586    IMPLICIT NONE 
     587   
     588    INTEGER,INTENT(IN),DIMENSION(:,:,:)  :: VarIn 
     589    INTEGER,INTENT(OUT),DIMENSION(:,:,:) :: VarOut 
     590     
     591    CALL Check_buffer_i(size(VarIn))    
     592    CALL reduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     593   
     594  END SUBROUTINE reduce_max_omp_i3 
     595 
     596 
     597  SUBROUTINE reduce_max_omp_i4(VarIn, VarOut) 
     598    IMPLICIT NONE 
     599 
     600    INTEGER,INTENT(IN),DIMENSION(:,:,:,:)  :: VarIn 
     601    INTEGER,INTENT(OUT),DIMENSION(:,:,:,:) :: VarOut 
     602   
     603    CALL Check_buffer_i(size(VarIn))    
     604    CALL reduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     605   
     606  END SUBROUTINE reduce_max_omp_i4 
     607 
     608 
     609  SUBROUTINE reduce_max_omp_r(VarIn, VarOut) 
     610    IMPLICIT NONE 
     611   
     612    REAL,INTENT(IN)  :: VarIn 
     613    REAL,INTENT(OUT) :: VarOut 
     614    REAL             :: VarIn_tmp(1) 
     615    REAL             :: VarOut_tmp(1) 
     616     
     617    VarIn_tmp(1)=VarIn 
     618    CALL Check_buffer_r(1)    
     619    CALL reduce_max_omp_rgen(VarIn_tmp,Varout_tmp,1,buffer_r) 
     620    VarOut=VarOut_tmp(1) 
     621   
     622  END SUBROUTINE reduce_max_omp_r 
     623 
     624  SUBROUTINE reduce_max_omp_r1(VarIn, VarOut) 
     625    IMPLICIT NONE 
     626   
     627    REAL,INTENT(IN),DIMENSION(:)  :: VarIn 
     628    REAL,INTENT(OUT),DIMENSION(:) :: VarOut 
     629     
     630    CALL Check_buffer_r(size(VarIn))    
     631    CALL reduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     632    
     633  END SUBROUTINE reduce_max_omp_r1 
     634   
     635   
     636  SUBROUTINE reduce_max_omp_r2(VarIn, VarOut) 
     637    IMPLICIT NONE 
     638   
     639    REAL,INTENT(IN),DIMENSION(:,:)  :: VarIn 
     640    REAL,INTENT(OUT),DIMENSION(:,:) :: VarOut 
     641     
     642    CALL Check_buffer_r(size(VarIn))    
     643    CALL reduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     644   
     645  END SUBROUTINE reduce_max_omp_r2 
     646 
     647 
     648  SUBROUTINE reduce_max_omp_r3(VarIn, VarOut) 
     649    IMPLICIT NONE 
     650   
     651    REAL,INTENT(IN),DIMENSION(:,:,:)  :: VarIn 
     652    REAL,INTENT(OUT),DIMENSION(:,:,:) :: VarOut 
     653     
     654    CALL Check_buffer_r(size(VarIn))    
     655    CALL reduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     656   
     657  END SUBROUTINE reduce_max_omp_r3 
     658 
     659 
     660  SUBROUTINE reduce_max_omp_r4(VarIn, VarOut) 
     661    IMPLICIT NONE 
     662 
     663    REAL,INTENT(IN),DIMENSION(:,:,:,:)  :: VarIn 
     664    REAL,INTENT(OUT),DIMENSION(:,:,:,:) :: VarOut 
     665   
     666    CALL Check_buffer_r(size(VarIn))    
     667    CALL reduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     668   
     669  END SUBROUTINE reduce_max_omp_r4   
     670   
     671   
     672   
     673   
     674    SUBROUTINE allreduce_max_omp_i(VarIn, VarOut) 
     675    IMPLICIT NONE 
     676   
     677    INTEGER,INTENT(IN)  :: VarIn 
     678    INTEGER,INTENT(OUT) :: VarOut 
     679    INTEGER             :: VarIn_tmp(1) 
     680    INTEGER             :: VarOut_tmp(1) 
     681     
     682    VarIn_tmp(1)=VarIn 
     683    CALL Check_buffer_i(1)    
     684    CALL allreduce_max_omp_igen(VarIn_tmp,Varout_tmp,1,buffer_i) 
     685    VarOut=VarOut_tmp(1) 
     686     
     687  END SUBROUTINE allreduce_max_omp_i 
     688 
     689  SUBROUTINE allreduce_max_omp_i1(VarIn, VarOut) 
     690    IMPLICIT NONE 
     691   
     692    INTEGER,INTENT(IN),DIMENSION(:)  :: VarIn 
     693    INTEGER,INTENT(OUT),DIMENSION(:) :: VarOut 
     694     
     695    CALL Check_buffer_i(size(VarIn))    
     696    CALL allreduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     697    
     698  END SUBROUTINE allreduce_max_omp_i1 
     699   
     700   
     701  SUBROUTINE allreduce_max_omp_i2(VarIn, VarOut) 
     702    IMPLICIT NONE 
     703   
     704    INTEGER,INTENT(IN),DIMENSION(:,:)  :: VarIn 
     705    INTEGER,INTENT(OUT),DIMENSION(:,:) :: VarOut 
     706 
     707    CALL Check_buffer_i(size(VarIn))    
     708    CALL allreduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     709   
     710  END SUBROUTINE allreduce_max_omp_i2 
     711 
     712 
     713  SUBROUTINE allreduce_max_omp_i3(VarIn, VarOut) 
     714    IMPLICIT NONE 
     715   
     716    INTEGER,INTENT(IN),DIMENSION(:,:,:)  :: VarIn 
     717    INTEGER,INTENT(OUT),DIMENSION(:,:,:) :: VarOut 
     718     
     719    CALL Check_buffer_i(size(VarIn))    
     720    CALL allreduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     721   
     722  END SUBROUTINE allreduce_max_omp_i3 
     723 
     724 
     725  SUBROUTINE allreduce_max_omp_i4(VarIn, VarOut) 
     726    IMPLICIT NONE 
     727 
     728    INTEGER,INTENT(IN),DIMENSION(:,:,:,:)  :: VarIn 
     729    INTEGER,INTENT(OUT),DIMENSION(:,:,:,:) :: VarOut 
     730   
     731    CALL Check_buffer_i(size(VarIn))    
     732    CALL allreduce_max_omp_igen(VarIn,Varout,Size(VarIn),buffer_i) 
     733   
     734  END SUBROUTINE allreduce_max_omp_i4 
     735 
     736 
     737  SUBROUTINE allreduce_max_omp_r(VarIn, VarOut) 
     738    IMPLICIT NONE 
     739   
     740    REAL,INTENT(IN)  :: VarIn 
     741    REAL,INTENT(OUT) :: VarOut 
     742    REAL             :: VarIn_tmp(1) 
     743    REAL             :: VarOut_tmp(1) 
     744     
     745    VarIn_tmp(1)=VarIn 
     746    CALL Check_buffer_r(1)    
     747    CALL allreduce_max_omp_rgen(VarIn_tmp,Varout_tmp,1,buffer_r) 
     748    VarOut=VarOut_tmp(1) 
     749   
     750  END SUBROUTINE allreduce_max_omp_r 
     751 
     752  SUBROUTINE allreduce_max_omp_r1(VarIn, VarOut) 
     753    IMPLICIT NONE 
     754   
     755    REAL,INTENT(IN),DIMENSION(:)  :: VarIn 
     756    REAL,INTENT(OUT),DIMENSION(:) :: VarOut 
     757     
     758    CALL Check_buffer_r(size(VarIn))    
     759    CALL allreduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     760    
     761  END SUBROUTINE allreduce_max_omp_r1 
     762   
     763   
     764  SUBROUTINE allreduce_max_omp_r2(VarIn, VarOut) 
     765    IMPLICIT NONE 
     766   
     767    REAL,INTENT(IN),DIMENSION(:,:)  :: VarIn 
     768    REAL,INTENT(OUT),DIMENSION(:,:) :: VarOut 
     769     
     770    CALL Check_buffer_r(size(VarIn))    
     771    CALL allreduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     772   
     773  END SUBROUTINE allreduce_max_omp_r2 
     774 
     775 
     776  SUBROUTINE allreduce_max_omp_r3(VarIn, VarOut) 
     777    IMPLICIT NONE 
     778   
     779    REAL,INTENT(IN),DIMENSION(:,:,:)  :: VarIn 
     780    REAL,INTENT(OUT),DIMENSION(:,:,:) :: VarOut 
     781     
     782    CALL Check_buffer_r(size(VarIn))    
     783    CALL allreduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     784   
     785  END SUBROUTINE allreduce_max_omp_r3 
     786 
     787 
     788  SUBROUTINE allreduce_max_omp_r4(VarIn, VarOut) 
     789    IMPLICIT NONE 
     790 
     791    REAL,INTENT(IN),DIMENSION(:,:,:,:)  :: VarIn 
     792    REAL,INTENT(OUT),DIMENSION(:,:,:,:) :: VarOut 
     793   
     794    CALL Check_buffer_r(size(VarIn))    
     795    CALL allreduce_max_omp_rgen(VarIn,Varout,Size(VarIn),buffer_r) 
     796   
     797  END SUBROUTINE allreduce_max_omp_r4   
     798 
    523799   
    524800!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
    … …  
    7401016   
    7411017  END SUBROUTINE allreduce_sum_omp_rgen 
     1018 
     1019 
     1020 
     1021 
     1022 
     1023 
     1024  SUBROUTINE reduce_max_omp_igen(VarIn,VarOut,dimsize,Buff) 
     1025  IMPLICIT NONE 
     1026 
     1027    INTEGER,INTENT(IN) :: dimsize 
     1028    INTEGER,INTENT(IN),DIMENSION(dimsize) :: VarIn 
     1029    INTEGER,INTENT(OUT),DIMENSION(dimsize) :: VarOut 
     1030    INTEGER,INTENT(INOUT),DIMENSION(dimsize) :: Buff 
     1031 
     1032    INTEGER :: i 
     1033 
     1034  !$OMP MASTER 
     1035    Buff(:)=VarIn(1) 
     1036  !$OMP END MASTER 
     1037  !$OMP BARRIER 
     1038   
     1039  !$OMP CRITICAL      
     1040    DO i=1,dimsize 
     1041      Buff(i)=MAX(Buff(i),VarIn(i)) 
     1042    ENDDO 
     1043  !$OMP END CRITICAL 
     1044  !$OMP BARRIER   
     1045   
     1046  !$OMP MASTER 
     1047    DO i=1,dimsize 
     1048      VarOut(i)=Buff(i) 
     1049    ENDDO 
     1050  !$OMP END MASTER 
     1051  !$OMP BARRIER 
     1052   
     1053  END SUBROUTINE reduce_max_omp_igen 
     1054 
     1055  SUBROUTINE reduce_max_omp_rgen(VarIn,VarOut,dimsize,Buff) 
     1056  IMPLICIT NONE 
     1057 
     1058    INTEGER,INTENT(IN) :: dimsize 
     1059    REAL,INTENT(IN),DIMENSION(dimsize) :: VarIn 
     1060    REAL,INTENT(OUT),DIMENSION(dimsize) :: VarOut 
     1061    REAL,INTENT(INOUT),DIMENSION(dimsize) :: Buff 
     1062 
     1063    INTEGER :: i 
     1064 
     1065  !$OMP MASTER 
     1066    Buff(:)=VarIn(1) 
     1067  !$OMP END MASTER 
     1068  !$OMP BARRIER 
     1069   
     1070  !$OMP CRITICAL      
     1071    DO i=1,dimsize 
     1072      Buff(i)=MAX(Buff(i),VarIn(i)) 
     1073    ENDDO 
     1074  !$OMP END CRITICAL 
     1075  !$OMP BARRIER   
     1076   
     1077    DO i=1,dimsize 
     1078      VarOut(i)=Buff(i) 
     1079    ENDDO 
     1080  !$OMP BARRIER 
     1081   
     1082  END SUBROUTINE reduce_max_omp_rgen 
     1083 
     1084 
     1085 
     1086  SUBROUTINE allreduce_max_omp_igen(VarIn,VarOut,dimsize,Buff) 
     1087  IMPLICIT NONE 
     1088 
     1089    INTEGER,INTENT(IN) :: dimsize 
     1090    INTEGER,INTENT(IN),DIMENSION(dimsize) :: VarIn 
     1091    INTEGER,INTENT(OUT),DIMENSION(dimsize) :: VarOut 
     1092    INTEGER,INTENT(INOUT),DIMENSION(dimsize) :: Buff 
     1093 
     1094    INTEGER :: i 
     1095 
     1096  !$OMP MASTER 
     1097    Buff(:)=VarIn(1) 
     1098  !$OMP END MASTER 
     1099  !$OMP BARRIER 
     1100   
     1101  !$OMP CRITICAL      
     1102    DO i=1,dimsize 
     1103      Buff(i)=MAX(Buff(i),VarIn(i)) 
     1104    ENDDO 
     1105  !$OMP END CRITICAL 
     1106  !$OMP BARRIER   
     1107   
     1108    DO i=1,dimsize 
     1109      VarOut(i)=Buff(i) 
     1110    ENDDO 
     1111  !$OMP BARRIER 
     1112   
     1113  END SUBROUTINE allreduce_max_omp_igen 
     1114 
     1115  SUBROUTINE allreduce_max_omp_rgen(VarIn,VarOut,dimsize,Buff) 
     1116  IMPLICIT NONE 
     1117 
     1118    INTEGER,INTENT(IN) :: dimsize 
     1119    REAL,INTENT(IN),DIMENSION(dimsize) :: VarIn 
     1120    REAL,INTENT(OUT),DIMENSION(dimsize) :: VarOut 
     1121    REAL,INTENT(INOUT),DIMENSION(dimsize) :: Buff 
     1122 
     1123    INTEGER :: i 
     1124 
     1125  !$OMP MASTER 
     1126    Buff(:)=VarIn(1) 
     1127  !$OMP END MASTER 
     1128  !$OMP BARRIER 
     1129   
     1130  !$OMP CRITICAL      
     1131    DO i=1,dimsize 
     1132      Buff(i)=MAX(Buff(i),VarIn(i)) 
     1133    ENDDO 
     1134  !$OMP END CRITICAL 
     1135  !$OMP BARRIER   
     1136   
     1137    DO i=1,dimsize 
     1138      VarOut(i)=Buff(i) 
     1139    ENDDO 
     1140 
     1141  !$OMP BARRIER 
     1142   
     1143  END SUBROUTINE allreduce_max_omp_rgen 
     1144     
     1145 
     1146 
     1147 
     1148 
     1149 
     1150 
     1151 
     1152 
     1153 
     1154 
    7421155     
    7431156END MODULE transfert_omp_mod 

  • codes/icosagcm/trunk/src/vertical_interp.f90

    r186 r295  
    2121  USE icosa 
    2222  USE pression_mod 
     23  USE omp_para 
    2324  IMPLICIT NONE 
    2425    TYPE(t_field),POINTER :: f_ps(:) 
    … …  
    4849 
    4950  SUBROUTINE compute_vertical_interp(p,in,out,pval) 
     51  USE omp_para 
    5052  IMPLICIT NONE 
    5153    REAL(rstd),INTENT(IN) :: p(iim*jjm,llm+1) 
    … …  
    5658    INTEGER :: i,j,ij,l 
    5759         
    58     DO j=jj_begin-1,jj_end+1 
    59       DO i=ii_begin-1,ii_end+1 
    60         ij=(j-1)*iim+i 
    61         l=llm-1 
    62         DO WHILE(0.5*(p(ij,l)+p(ij,l+1))<pval .AND. l>1) 
    63           l=l-1 
     60!$OMP BARRIER     
     61    IF (is_omp_level_master) THEN 
     62     
     63      DO j=jj_begin-1,jj_end+1 
     64        DO i=ii_begin-1,ii_end+1 
     65          ij=(j-1)*iim+i 
     66          l=llm-1 
     67          DO WHILE(0.5*(p(ij,l)+p(ij,l+1))<pval .AND. l>1) 
     68            l=l-1 
     69          ENDDO 
     70          pmid=0.5*(p(ij,l)+p(ij,l+1)) 
     71          pmidp1=0.5*(p(ij,l+1)+p(ij,l+2)) 
     72 
     73          coeff=(pval-pmid)/(pmid-pmidp1) 
     74         
     75          out(ij)=in(ij,l)+coeff*(in(ij,l)-in(ij,l+1)) 
    6476        ENDDO 
    65         pmid=0.5*(p(ij,l)+p(ij,l+1)) 
    66         pmidp1=0.5*(p(ij,l+1)+p(ij,l+2)) 
    67  
    68         coeff=(pval-pmid)/(pmid-pmidp1) 
    69          
    70         out(ij)=in(ij,l)+coeff*(in(ij,l)-in(ij,l+1)) 
    7177      ENDDO 
    72     ENDDO 
     78    
     79    ENDIF 
     80!$OMP BARRIER     
    7381         
    7482  END SUBROUTINE compute_vertical_interp 

  • codes/icosagcm/trunk/src/vorticity.f90

    r186 r295  
    2929  USE icosa 
    3030  USE disvert_mod 
     31  USE omp_para 
    3132  IMPLICIT NONE 
    3233    REAL(rstd),INTENT(IN)  :: ue(3*iim*jjm,llm) 
    … …  
    3435    INTEGER :: i,j,ij,l 
    3536 
    36     DO l = 1,llm 
     37    DO l = ll_begin,ll_end 
    3738      DO j=jj_begin-1,jj_end+1 
    3839        DO i=ii_begin-1,ii_end+1 

  • codes/icosagcm/trunk/src/wind.f90

    r266 r295  
    4747  SUBROUTINE compute_wind_centered(ue,ucenter) 
    4848  USE icosa 
    49    
     49  USE omp_para 
    5050  IMPLICIT NONE 
    5151  REAL(rstd) :: ue(3*iim*jjm,llm) 
    … …  
    5353  INTEGER :: i,j,ij,l     
    5454  
    55     DO l=1,llm 
     55    DO l=ll_begin,ll_end 
    5656      DO j=jj_begin,jj_end 
    5757        DO i=ii_begin,ii_end 
    … …  
    7373 SUBROUTINE compute_wind_on_edge(ue,uedge) 
    7474  USE icosa 
     75  USE omp_para 
    7576     
    7677  IMPLICIT NONE 
    … …  
    8384    CALL compute_tangential_compound(ue,ut) 
    8485   
    85     DO l=1,llm 
     86    DO l=ll_begin,ll_end 
    8687      DO j=jj_begin,jj_end 
    8788        DO i=ii_begin,ii_end 
    … …  
    100101 SUBROUTINE compute_tangential_compound(ue,ut) 
    101102  USE icosa   
     103  USE omp_para 
    102104  IMPLICIT NONE 
    103105  REAL(rstd) :: ue(3*iim*jjm,llm) 
    … …  
    105107  INTEGER :: i,j,l,ij 
    106108     
    107   DO l=1,llm 
     109  DO l=ll_begin,ll_end 
    108110    DO j=jj_begin,jj_end 
    109111      DO i=ii_begin,ii_end 
    … …  
    155157 SUBROUTINE compute_wind_lonlat_compound(u, ulon, ulat) 
    156158  USE icosa   
     159  USE omp_para 
    157160     
    158161  IMPLICIT NONE 
    … …  
    164167     
    165168   
    166     DO l=1,llm 
     169    DO l=ll_begin,ll_end 
    167170      DO j=jj_begin-1,jj_end+1 
    168171        DO i=ii_begin-1,ii_end+1 
    … …  
    184187  SUBROUTINE compute_wind_from_lonlat_compound(ulon, ulat, u) 
    185188  USE icosa   
     189  USE omp_para 
    186190     
    187191  IMPLICIT NONE 
    … …  
    192196  INTEGER :: i,j,ij,l      
    193197   
    194     DO l=1,llm 
     198    DO l=ll_begin,ll_end 
    195199      DO j=jj_begin-1,jj_end+1 
    196200        DO i=ii_begin-1,ii_end+1 
    … …  
    207211  SUBROUTINE compute_wind_centered_from_lonlat_compound(ulon, ulat, u) 
    208212  USE icosa   
     213  USE omp_para 
    209214     
    210215  IMPLICIT NONE 
    … …  
    214219 
    215220  INTEGER :: i,j,ij,l      
    216   DO l=1,llm 
     221  DO l=ll_begin,ll_end 
    217222      DO j=jj_begin-1,jj_end+1 
    218223        DO i=ii_begin-1,ii_end+1 
    … …  
    248253  SUBROUTINE compute_wind_perp_from_lonlat_compound(ulon, ulat, up) 
    249254  USE icosa   
     255  USE omp_para 
    250256     
    251257  IMPLICIT NONE 
    … …  
    259265   CALL compute_wind_from_lonlat_compound(ulon, ulat, u) 
    260266 
    261     DO l=1,llm 
     267    DO l=ll_begin,ll_end 
    262268      DO j=jj_begin-1,jj_end+1 
    263269        DO i=ii_begin-1,ii_end+1 
    … …  
    297303 SUBROUTINE compute_wind_centered_lonlat_compound(uc, ulon, ulat) 
    298304  USE icosa   
     305  USE omp_para 
    299306     
    300307  IMPLICIT NONE 
    … …  
    306313     
    307314   
    308     DO l=1,llm 
     315    DO l=ll_begin,ll_end 
    309316      DO j=jj_begin,jj_end 
    310317        DO i=ii_begin,ii_end 
    … …  
    320327 SUBROUTINE compute_wind_centered_from_wind_lonlat_centered(ulon, ulat,uc) 
    321328  USE icosa   
     329  USE omp_para 
    322330     
    323331  IMPLICIT NONE 
    … …  
    329337     
    330338   
    331     DO l=1,llm 
     339    DO l=ll_begin,ll_end 
    332340      DO j=jj_begin,jj_end 
    333341        DO i=ii_begin,ii_end 
    … …  
    344352 SUBROUTINE compute_wind_perp_from_wind_centered(uc,un) 
    345353  USE icosa   
     354  USE omp_para 
    346355     
    347356  IMPLICIT NONE 
    … …  
    352361     
    353362   
    354     DO l=1,llm 
     363    DO l=ll_begin,ll_end 
    355364      DO j=jj_begin,jj_end 
    356365        DO i=ii_begin,ii_end 
Note: See TracChangeset for help on using the changeset viewer.