Index: branches/2016/dev_r6519_HPC_4/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90
===================================================================
--- branches/2016/dev_r6519_HPC_4/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90	(revision 6520)
+++ branches/2016/dev_r6519_HPC_4/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90	(revision 7037)
@@ -128,5 +128,9 @@
       !------------------------------------------!
 
+!$OMP PARALLEL
+!$OMP WORKSHARE
       zqnsoce(:,:) = qns(:,:)
+!$OMP END WORKSHARE NOWAIT
+!$OMP DO schedule(static) private(jj,ji,zinda,ifvt,i1mfr,idfr,iflt,ial,iadv,ifral,ifrdv,zqsr,zqns,zqhc,zemp,zemp_snw,zfmm,zfsalt,zcd) 
       DO jj = 1, jpj
          DO ji = 1, jpi
@@ -229,8 +233,11 @@
          END DO
       END DO
+!$OMP END DO NOWAIT
+!$OMP END PARALLEL
       !                                !------------------------------------------!
       !                                !    mass of snow and ice per unit area    !
       !                                !------------------------------------------!
       IF( nn_ice_embd /= 0 ) THEN      ! embedded sea-ice (mass required)
+!$OMP PARALLEL WORKSHARE
          snwice_mass_b(:,:) = snwice_mass(:,:)                  ! save mass from the previous ice time step
          !                                                      ! new mass per unit area
@@ -238,4 +245,5 @@
          !                                                      ! time evolution of snow+ice mass
          snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) / rdt_ice
+!$OMP END PARALLEL WORKSHARE
       ENDIF
 
@@ -253,11 +261,15 @@
 
       IF( ln_cpl) THEN
+!$OMP PARALLEL WORKSHARE
          tn_ice(:,:,1) = sist(:,:)          ! sea-ice surface temperature       
          ht_i(:,:,1) = hicif(:,:)
          ht_s(:,:,1) = hsnif(:,:)
          a_i(:,:,1) = fr_i(:,:)
+!$OMP END PARALLEL WORKSHARE
          !                                  ! Computation of snow/ice and ocean albedo
          CALL albedo_ice( tn_ice, ht_i, ht_s, zalbp, zalb )
+!$OMP PARALLEL WORKSHARE
          alb_ice(:,:,1) =  0.5 * ( zalbp(:,:,1) + zalb (:,:,1) )   ! Ice albedo (mean clear and overcast skys)
+!$OMP END PARALLEL WORKSHARE
          IF( iom_use('icealb_cea' ) )   CALL iom_put( 'icealb_cea', alb_ice(:,:,1) * fr_i(:,:) )  ! ice albedo
       ENDIF
@@ -320,4 +332,6 @@
          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==! (i.e. surface module time-step)
             !
+!$OMP PARALLEL
+!$OMP DO schedule(static) private(jj,ji,zu_i,zv_i) 
             DO jj = 1, jpj                               !* modulus of ice-ocean relative velocity at I-point
                DO ji = 1, jpi
@@ -327,4 +341,5 @@
                END DO
             END DO
+!$OMP DO schedule(static) private(jj,ji,zumt) 
             DO jj = 1, jpjm1                             !* update the modulus of stress at ocean surface (T-point)
                DO ji = 1, jpim1   ! NO vector opt.
@@ -336,8 +351,12 @@
                END DO
             END DO
+!$OMP END DO NOWAIT
+!$OMP END PARALLEL
             CALL lbc_lnk( taum, 'T', 1. )
             !
+!$OMP PARALLEL WORKSHARE
             utau_oce(:,:) = utau(:,:)                    !* save the air-ocean stresses at ice time-step
             vtau_oce(:,:) = vtau(:,:)
+!$OMP END PARALLEL WORKSHARE
             !
          ENDIF
@@ -346,4 +365,6 @@
          !
          !                                               !* ice/ocean stress WITH a ice-ocean rotation angle at I-point
+!$OMP PARALLEL
+!$OMP DO schedule(static) private(jj,ji,zsang,zu_i,zv_i,zmodi) 
          DO jj = 2, jpj
             zsang  = SIGN( 1._wp, gphif(1,jj) ) * sangvg          ! change the cosine angle sign in the SH 
@@ -359,4 +380,5 @@
          END DO
          !                                               !* surface ocean stresses at u- and v-points
+!$OMP DO schedule(static) private(jj,ji,zutau_ice,zvtau_ice,zfrldu,zfrldv) 
          DO jj = 2, jpjm1
             DO ji = 2, jpim1   ! NO vector opt.
@@ -372,4 +394,6 @@
             END DO
          END DO
+!$OMP END DO NOWAIT
+!$OMP END PARALLEL
          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )     ! lateral boundary condition
          !
@@ -381,4 +405,5 @@
          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==! (i.e. surface module time-step)
             !
+!$OMP PARALLEL DO schedule(static) private(jj,ji,zu_t,zv_t,zmodt) 
             DO jj = 2, jpjm1                          !* modulus of the ice-ocean velocity at T-point
                DO ji = fs_2, fs_jpim1
@@ -393,6 +418,8 @@
             CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. )
             !
+!$OMP PARALLEL WORKSHARE
             utau_oce(:,:) = utau(:,:)                 !* save the air-ocean stresses at ice time-step
             vtau_oce(:,:) = vtau(:,:)
+!$OMP END PARALLEL WORKSHARE
             !
          ENDIF
@@ -400,4 +427,5 @@
          !                                        !==  at each ocean time-step  ==!
          !
+!$OMP PARALLEL DO schedule(static) private(jj,ji,zutau_ice,zvtau_ice,zfrldu,zfrldv) 
          DO jj = 2, jpjm1                             !* ice stress over ocean WITHOUT a ice-ocean rotation angle
             DO ji = fs_2, fs_jpim1
@@ -435,5 +463,5 @@
       !! ** input   : Namelist namicedia
       !!-------------------------------------------------------------------
-      INTEGER ::   jk   ! local integer
+      INTEGER ::   jk,jj,ji   ! local integer
       !!-------------------------------------------------------------------
       !
@@ -447,6 +475,8 @@
       r1_rdtice = 1._wp / rdt_ice
       !
+!$OMP PARALLEL WORKSHARE
       soce_0(:,:) = soce                     ! constant SSS and ice salinity used in levitating sea-ice case
       sice_0(:,:) = sice
+!$OMP END PARALLEL WORKSHARE
       !
       IF( cp_cfg == "orca" ) THEN            ! decrease ocean & ice reference salinities in the Baltic sea 
@@ -459,24 +489,35 @@
       !                                      ! embedded sea ice
       IF( nn_ice_embd /= 0 ) THEN            ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass
+!$OMP PARALLEL WORKSHARE
          snwice_mass  (:,:) = tms(:,:) * ( rhosn * hsnif(:,:) + rhoic * hicif(:,:)  ) * ( 1.0 - frld(:,:) )
          snwice_mass_b(:,:) = snwice_mass(:,:)
+!$OMP END PARALLEL WORKSHARE
       ELSE
+!$OMP PARALLEL WORKSHARE
          snwice_mass  (:,:) = 0.e0           ! no mass exchanges
          snwice_mass_b(:,:) = 0.e0           ! no mass exchanges
          snwice_fmass (:,:) = 0.e0           ! no mass exchanges
+!$OMP END PARALLEL WORKSHARE
       ENDIF
       IF( nn_ice_embd == 2 .AND.          &  ! full embedment (case 2) & no restart : 
          &   .NOT.ln_rstart ) THEN           ! deplete the initial ssh below sea-ice area
+!$OMP PARALLEL WORKSHARE
          sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0
          sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0
+!$OMP END PARALLEL WORKSHARE
 !!gm I really don't like this staff here...  Find a way to put that elsewhere or differently
 !!gm
          IF( .NOT.ln_linssh ) THEN
 
+!$OMP PARALLEL
+!$OMP DO schedule(static) private(jk)
             do jk = 1,jpkm1                     ! adjust initial vertical scale factors
                e3t_n(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshn(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) )
                e3t_b(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshb(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) )
             end do
+!$OMP WORKSHARE
             e3t_a(:,:,:) = e3t_b(:,:,:)
+!$OMP END WORKSHARE NOWAIT
+!$OMP END PARALLEL
             ! Reconstruction of all vertical scale factors at now and before time steps
             !        ! Horizontal scale factor interpolations
@@ -493,12 +534,22 @@
             CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' )
             !        ! t- and w- points depth
+!$OMP PARALLEL
+!$OMP WORKSHARE
             gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1)
             gdepw_n(:,:,1) = 0.0_wp
             gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:)
+!$OMP END WORKSHARE
             DO jk = 2, jpk
-               gdept_n(:,:,jk) = gdept_n(:,:,jk-1) + e3w_n(:,:,jk)
-               gdepw_n(:,:,jk) = gdepw_n(:,:,jk-1) + e3t_n(:,:,jk-1)
-               gde3w_n(:,:,jk) = gdept_n(:,:,jk  ) - sshn   (:,:)
-            END DO
+!$OMP DO schedule(static) private(jj,ji)
+               DO jj = 1, jpj
+                  DO ji = 1, jpi
+                     gdept_n(ji,jj,jk) = gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk)
+                     gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
+                     gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk  ) - sshn (ji,jj)
+                  END DO
+               END DO
+!$OMP END DO NOWAIT
+            END DO
+!$OMP END PARALLEL
          ENDIF
 !!gm end
