New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7698 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90 – NEMO

Ignore:
Timestamp:
2017-02-18T10:02:03+01:00 (7 years ago)
Author:
mocavero
Message:

update trunk with OpenMP parallelization

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r7646 r7698  
    112112      ! --- case we bypass ice thermodynamics --- ! 
    113113      IF( .NOT. ln_limthd ) THEN   ! we suppose ice is impermeable => ocean is isolated from atmosphere 
    114          hfx_in   (:,:)   = pfrld(:,:) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:) 
    115          hfx_out  (:,:)   = pfrld(:,:) *   qns_oce(:,:)                  + qemp_oce(:,:) 
    116          ftr_ice  (:,:,:) = 0._wp 
    117          emp_ice  (:,:)   = 0._wp 
    118          qemp_ice (:,:)   = 0._wp 
    119          qevap_ice(:,:,:) = 0._wp 
     114!$OMP PARALLEL 
     115!$OMP DO schedule(static) private(jj,ji) 
     116         DO jj = 1, jpj 
     117            DO ji = 1, jpi 
     118               hfx_in   (ji,jj)   = pfrld(ji,jj) * ( qns_oce(ji,jj) + qsr_oce(ji,jj) ) + qemp_oce(ji,jj) 
     119               hfx_out  (ji,jj)   = pfrld(ji,jj) *   qns_oce(ji,jj)                  + qemp_oce(ji,jj) 
     120               emp_ice  (ji,jj)   = 0._wp 
     121               qemp_ice (ji,jj)   = 0._wp 
     122            END DO 
     123         END DO 
     124         DO jl = 1, jpl 
     125!$OMP DO schedule(static) private(jj,ji) 
     126            DO jj = 1, jpj 
     127               DO ji = 1, jpi 
     128                  ftr_ice  (ji,jj,jl) = 0._wp 
     129                  qevap_ice(ji,jj,jl) = 0._wp 
     130               END DO 
     131            END DO 
     132         END DO 
     133!$OMP END PARALLEL 
    120134      ENDIF 
    121135       
     
    123137      CALL wrk_alloc( jpi,jpj, zalb )     
    124138 
    125       zalb(:,:) = 0._wp 
    126       WHERE     ( at_i_b <= epsi06 )  ;  zalb(:,:) = 0.066_wp 
    127       ELSEWHERE                       ;  zalb(:,:) = SUM( alb_ice * a_i_b, dim=3 ) / at_i_b 
    128       END WHERE 
     139!$OMP PARALLEL 
     140!$OMP DO schedule(static) private(jj,ji) 
     141      DO jj = 1, jpj 
     142         DO ji = 1, jpi 
     143            zalb(ji,jj) = 0._wp 
     144         END DO 
     145      END DO 
     146!$OMP DO schedule(static) private(jj,ji,jl) 
     147      DO jj = 1, jpj 
     148         DO ji = 1, jpi 
     149            IF ( at_i_b(ji,jj) <= epsi06 ) THEN 
     150               zalb(ji,jj) = 0.066_wp 
     151            ELSE   
     152               DO jl = 1, jpl 
     153                  zalb(ji,jj) = zalb(ji,jj) + ( alb_ice(ji,jj,jl) * a_i_b(ji,jj,jl) ) / at_i_b(ji,jj) 
     154               END DO 
     155            END IF 
     156          END DO 
     157      END DO 
     158!$OMP END PARALLEL 
    129159      IF( iom_use('alb_ice' ) )  CALL iom_put( "alb_ice"  , zalb(:,:) )           ! ice albedo output 
    130160 
    131       zalb(:,:) = SUM( alb_ice * a_i_b, dim=3 ) + 0.066_wp * ( 1._wp - at_i_b )       
     161!$OMP PARALLEL 
     162!$OMP DO schedule(static) private(jj,ji) 
     163      DO jj = 1, jpj 
     164         DO ji = 1, jpi 
     165            zalb(ji,jj) = 0._wp 
     166         END DO 
     167      END DO 
     168      DO jl = 1, jpl 
     169!$OMP DO schedule(static) private(jj,ji) 
     170         DO jj = 1, jpj 
     171            DO ji = 1, jpi 
     172               zalb(ji,jj) = zalb(ji,jj) + ( alb_ice(ji,jj,jl) * a_i_b(ji,jj,jl) ) + 0.066_wp * ( 1._wp - at_i_b(ji,jj) )       
     173            END DO 
     174         END DO 
     175      END DO 
     176!$OMP END PARALLEL 
    132177      IF( iom_use('albedo'  ) )  CALL iom_put( "albedo"  , zalb(:,:) )           ! ice albedo output 
    133178 
    134179      CALL wrk_dealloc( jpi,jpj, zalb )     
    135180 
     181!$OMP PARALLEL 
     182!$OMP DO schedule(static) private(jj,ji,jl,zqsr,zqmass) 
    136183      DO jj = 1, jpj 
    137184         DO ji = 1, jpi 
     
    186233      !      salt flux at the ocean surface      ! 
    187234      !------------------------------------------! 
    188       sfx(:,:) = sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:)   & 
    189          &     + sfx_res(:,:) + sfx_dyn(:,:) + sfx_bri(:,:) + sfx_sub(:,:) + sfx_lam(:,:) 
     235!$OMP DO schedule(static) private(jj,ji) 
     236      DO jj = 1, jpj 
     237         DO ji = 1, jpi 
     238            sfx(ji,jj) = sfx_bog(ji,jj) + sfx_bom(ji,jj) + sfx_sum(ji,jj) + sfx_sni(ji,jj) + sfx_opw(ji,jj)   & 
     239               &     + sfx_res(ji,jj) + sfx_dyn(ji,jj) + sfx_bri(ji,jj) + sfx_sub(ji,jj) + sfx_lam(ji,jj) 
     240         END DO 
     241      END DO 
     242!$OMP END PARALLEL 
    190243 
    191244      !-------------------------------------------------------------! 
     
    193246      !-------------------------------------------------------------! 
    194247      IF( nn_ice_embd /= 0 ) THEN 
    195          ! save mass from the previous ice time step 
    196          snwice_mass_b(:,:) = snwice_mass(:,:)                   
    197          ! new mass per unit area 
    198          snwice_mass  (:,:) = tmask(:,:,1) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  )  
    199          ! time evolution of snow+ice mass 
    200          snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_rdtice 
     248!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     249         DO jj = 1, jpj 
     250            DO ji = 1, jpi 
     251               ! save mass from the previous ice time step 
     252               snwice_mass_b(ji,jj) = snwice_mass(ji,jj)                   
     253               ! new mass per unit area 
     254               snwice_mass  (ji,jj) = tmask(ji,jj,1) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj)  )  
     255               ! time evolution of snow+ice mass 
     256               snwice_fmass (ji,jj) = ( snwice_mass(ji,jj) - snwice_mass_b(ji,jj) ) * r1_rdtice 
     257            END DO 
     258         END DO 
    201259      ENDIF 
    202260 
     
    204262      !   Storing the transmitted variables           ! 
    205263      !-----------------------------------------------! 
    206       fr_i  (:,:)   = at_i(:,:)             ! Sea-ice fraction             
    207       tn_ice(:,:,:) = t_su(:,:,:)           ! Ice surface temperature                       
     264!$OMP PARALLEL 
     265!$OMP DO schedule(static) private(jj,ji) 
     266      DO jj = 1, jpj 
     267         DO ji = 1, jpi 
     268            fr_i  (ji,jj)   = at_i(ji,jj)             ! Sea-ice fraction             
     269         END DO 
     270      END DO 
     271      DO jl = 1, jpl 
     272!$OMP DO schedule(static) private(jj,ji) 
     273         DO jj = 1, jpj 
     274            DO ji = 1, jpi 
     275               tn_ice(ji,jj,jl) = t_su(ji,jj,jl)           ! Ice surface temperature                       
     276            END DO 
     277         END DO 
     278      END DO 
     279!$OMP END PARALLEL 
    208280 
    209281      !------------------------------------------------------------------------! 
     
    212284      CALL wrk_alloc( jpi,jpj,jpl,   zalb_cs, zalb_os )     
    213285      CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos 
    214       alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     286      DO jl = 1, jpl 
     287!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     288         DO jj = 1, jpj 
     289            DO ji = 1, jpi 
     290               alb_ice(ji,jj,jl) = ( 1. - cldf_ice ) * zalb_cs(ji,jj,jl) + cldf_ice * zalb_os(ji,jj,jl) 
     291            END DO 
     292         END DO 
     293      END DO 
    215294      CALL wrk_dealloc( jpi,jpj,jpl,   zalb_cs, zalb_os ) 
    216295 
     
    260339      ! 
    261340      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==!   (i.e. surface module time-step) 
     341!$OMP PARALLEL DO schedule(static) private(jj,ji,zu_t,zv_t,zmodt) 
    262342         DO jj = 2, jpjm1                             !* update the modulus of stress at ocean surface (T-point) 
    263343            DO ji = fs_2, fs_jpim1 
     
    274354         CALL lbc_lnk_multi( taum, 'T', 1., tmod_io, 'T', 1. ) 
    275355         ! 
    276          utau_oce(:,:) = utau(:,:)                    !* save the air-ocean stresses at ice time-step 
    277          vtau_oce(:,:) = vtau(:,:) 
     356!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     357         DO jj = 1, jpj 
     358            DO ji = 1, jpi 
     359               utau_oce(ji,jj) = utau(ji,jj)                    !* save the air-ocean stresses at ice time-step 
     360               vtau_oce(ji,jj) = vtau(ji,jj) 
     361            END DO 
     362         END DO 
    278363         ! 
    279364      ENDIF 
     
    281366      !                                      !==  every ocean time-step  ==! 
    282367      ! 
     368!$OMP PARALLEL DO schedule(static) private(jj,ji,zat_u,zat_v,zutau_ice,zvtau_ice) 
    283369      DO jj = 2, jpjm1                                !* update the stress WITHOUT a ice-ocean rotation angle 
    284370         DO ji = fs_2, fs_jpim1   ! Vect. Opt. 
     
    319405      IF( lim_sbc_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'lim_sbc_init : unable to allocate standard arrays' ) 
    320406      ! 
    321       soce_0(:,:) = soce                     ! constant SSS and ice salinity used in levitating sea-ice case 
    322       sice_0(:,:) = sice 
     407!$OMP PARALLEL 
     408!$OMP DO schedule(static) private(jj,ji) 
     409      DO jj = 1, jpj 
     410         DO ji = 1, jpi 
     411            soce_0(ji,jj) = soce                     ! constant SSS and ice salinity used in levitating sea-ice case 
     412            sice_0(ji,jj) = sice 
     413         END DO 
     414      END DO 
    323415      !                                      ! decrease ocean & ice reference salinities in the Baltic Sea area 
    324       WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.   & 
    325          &   54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp         )  
    326          soce_0(:,:) = 4._wp 
    327          sice_0(:,:) = 2._wp 
    328       END WHERE 
     416!$OMP DO schedule(static) private(jj,ji) 
     417      DO jj = 1, jpj 
     418         DO ji = 1, jpi 
     419            IF ( 14._wp <= glamt(ji,jj) .AND. glamt(ji,jj) <= 32._wp .AND.   & 
     420               &   54._wp <= gphit(ji,jj) .AND. gphit(ji,jj) <= 66._wp         ) THEN 
     421               soce_0(ji,jj) = 4._wp 
     422               sice_0(ji,jj) = 2._wp 
     423            END IF 
     424         END DO 
     425      END DO 
     426!$OMP END PARALLEL 
    329427      ! 
    330428      IF( .NOT. ln_rstart ) THEN 
    331429         !                                      ! embedded sea ice 
    332430         IF( nn_ice_embd /= 0 ) THEN            ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass 
    333             snwice_mass  (:,:) = tmask(:,:,1) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  ) 
    334             snwice_mass_b(:,:) = snwice_mass(:,:) 
     431!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     432            DO jj = 1, jpj 
     433               DO ji = 1, jpi 
     434                  snwice_mass  (ji,jj) = tmask(ji,jj,1) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj)  ) 
     435                  snwice_mass_b(ji,jj) = snwice_mass(ji,jj) 
     436               END DO 
     437            END DO 
    335438         ELSE 
    336             snwice_mass  (:,:) = 0._wp          ! no mass exchanges 
    337             snwice_mass_b(:,:) = 0._wp          ! no mass exchanges 
     439!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     440            DO jj = 1, jpj 
     441               DO ji = 1, jpi 
     442                  snwice_mass  (ji,jj) = 0._wp          ! no mass exchanges 
     443                  snwice_mass_b(ji,jj) = 0._wp          ! no mass exchanges 
     444               END DO 
     445            END DO 
    338446         ENDIF 
    339447         IF( nn_ice_embd == 2 ) THEN            ! full embedment (case 2) deplete the initial ssh below sea-ice area 
    340             sshn(:,:) = sshn(:,:) - snwice_mass(:,:) * r1_rau0 
    341             sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0 
     448!$OMP PARALLEL DO schedule(static) private(jj,ji) 
     449            DO jj = 1, jpj 
     450               DO ji = 1, jpi 
     451                  sshn(ji,jj) = sshn(ji,jj) - snwice_mass(ji,jj) * r1_rau0 
     452                  sshb(ji,jj) = sshb(ji,jj) - snwice_mass(ji,jj) * r1_rau0 
     453               END DO 
     454            END DO 
    342455 
    343456!!gm I really don't like this stuff here...  Find a way to put that elsewhere or differently 
    344457!!gm 
    345458            IF( .NOT.ln_linssh ) THEN 
     459!$OMP PARALLEL 
     460!$OMP DO schedule(static) private(jj,ji) 
    346461               DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
    347                   e3t_n(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshn(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
    348                   e3t_b(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshb(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
    349                END DO 
    350                e3t_a(:,:,:) = e3t_b(:,:,:) 
     462                  DO jj = 1, jpj 
     463                     DO ji = 1, jpi 
     464                        e3t_n(ji,jj,jk) = e3t_0(ji,jj,jk)*( 1._wp + sshn(ji,jj)*tmask(ji,jj,1)/(ht_0(ji,jj) + 1.0 - tmask(ji,jj,1)) ) 
     465                        e3t_b(ji,jj,jk) = e3t_0(ji,jj,jk)*( 1._wp + sshb(ji,jj)*tmask(ji,jj,1)/(ht_0(ji,jj) + 1.0 - tmask(ji,jj,1)) ) 
     466                     END DO 
     467                  END DO 
     468               END DO 
     469!$OMP DO schedule(static) private(jj,ji) 
     470               DO jk = 1,jpk 
     471                  DO jj = 1, jpj 
     472                     DO ji = 1, jpi 
     473                        e3t_a(ji,jj,jk) = e3t_b(ji,jj,jk) 
     474                     END DO 
     475                  END DO 
     476               END DO 
     477!$OMP END PARALLEL 
    351478               ! Reconstruction of all vertical scale factors at now and before time-steps 
    352479               ! ========================================================================= 
     
    368495               ! ---------------------- 
    369496!!gm not sure of that.... 
    370                gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1) 
    371                gdepw_n(:,:,1) = 0.0_wp 
    372                gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:) 
     497!$OMP PARALLEL 
     498!$OMP DO schedule(static) private(jj,ji) 
     499               DO jj = 1, jpj 
     500                  DO ji = 1, jpi 
     501                     gdept_n(ji,jj,1) = 0.5_wp * e3w_n(ji,jj,1) 
     502                     gdepw_n(ji,jj,1) = 0.0_wp 
     503                     gde3w_n(ji,jj,1) = gdept_n(ji,jj,1) - sshn(ji,jj) 
     504                  END DO 
     505               END DO 
    373506               DO jk = 2, jpk 
    374                   gdept_n(:,:,jk) = gdept_n(:,:,jk-1) + e3w_n(:,:,jk) 
    375                   gdepw_n(:,:,jk) = gdepw_n(:,:,jk-1) + e3t_n(:,:,jk-1) 
    376                   gde3w_n(:,:,jk) = gdept_n(:,:,jk  ) - sshn   (:,:) 
    377                END DO 
     507!$OMP DO schedule(static) private(jj,ji) 
     508                  DO jj = 1, jpj 
     509                     DO ji = 1, jpi 
     510                        gdept_n(ji,jj,jk) = gdept_n(ji,jj,jk-1) + e3w_n(ji,jj,jk) 
     511                        gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
     512                        gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk  ) - sshn   (ji,jj) 
     513                     END DO 
     514                  END DO 
     515               END DO 
     516!$OMP END PARALLEL 
    378517            ENDIF 
    379518         ENDIF 
Note: See TracChangeset for help on using the changeset viewer.