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.
sbcflx.F90 in NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/SBC – NEMO

source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/SBC/sbcflx.F90 @ 12340

Last change on this file since 12340 was 12340, checked in by acc, 4 years ago

Branch 2019/dev_r11943_MERGE_2019. This commit introduces basic do loop macro
substitution to the 2019 option 1, merge branch. These changes have been SETTE
tested. The only addition is the do_loop_substitute.h90 file in the OCE directory but
the macros defined therein are used throughout the code to replace identifiable, 2D-
and 3D- nested loop opening and closing statements with single-line alternatives. Code
indents are also adjusted accordingly.

The following explanation is taken from comments in the new header file:

This header file contains preprocessor definitions and macros used in the do-loop
substitutions introduced between version 4.0 and 4.2. The primary aim of these macros
is to assist in future applications of tiling to improve performance. This is expected
to be achieved by alternative versions of these macros in selected locations. The
initial introduction of these macros simply replaces all identifiable nested 2D- and
3D-loops with single line statements (and adjusts indenting accordingly). Do loops
are identifiable if they comform to either:

DO jk = ....

DO jj = .... DO jj = ...

DO ji = .... DO ji = ...
. OR .
. .

END DO END DO

END DO END DO

END DO

and white-space variants thereof.

Additionally, only loops with recognised jj and ji loops limits are treated; these are:
Lower limits of 1, 2 or fs_2
Upper limits of jpi, jpim1 or fs_jpim1 (for ji) or jpj, jpjm1 or fs_jpjm1 (for jj)

The macro naming convention takes the form: DO_2D_BT_LR where:

B is the Bottom offset from the PE's inner domain;
T is the Top offset from the PE's inner domain;
L is the Left offset from the PE's inner domain;
R is the Right offset from the PE's inner domain

So, given an inner domain of 2,jpim1 and 2,jpjm1, a typical example would replace:

DO jj = 2, jpj

DO ji = 1, jpim1
.
.

END DO

END DO

with:

DO_2D_01_10
.
.
END_2D

similar conventions apply to the 3D loops macros. jk loop limits are retained
through macro arguments and are not restricted. This includes the possibility of
strides for which an extra set of DO_3DS macros are defined.

In the example definition below the inner PE domain is defined by start indices of
(kIs, kJs) and end indices of (kIe, KJe)

#define DO_2D_00_00 DO jj = kJs, kJe ; DO ji = kIs, kIe
#define END_2D END DO ; END DO

TO DO:

Only conventional nested loops have been identified and replaced by this step. There are constructs such as:

DO jk = 2, jpkm1

z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk)

END DO

which may need to be considered.
  • Property svn:keywords set to Id
File size: 9.0 KB
Line 
1MODULE sbcflx
2   !!======================================================================
3   !!                       ***  MODULE  sbcflx  ***
4   !! Ocean forcing:  momentum, heat and freshwater flux formulation
5   !!=====================================================================
6   !! History :  1.0  !  2006-06  (G. Madec)  Original code
7   !!            3.3  !  2010-10  (S. Masson)  add diurnal cycle
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   namflx   : flux formulation namlist
12   !!   sbc_flx  : flux formulation as ocean surface boundary condition (forced mode, fluxes read in NetCDF files)
13   !!----------------------------------------------------------------------
14   USE oce             ! ocean dynamics and tracers
15   USE dom_oce         ! ocean space and time domain
16   USE sbc_oce         ! surface boundary condition: ocean fields
17   USE sbcdcy          ! surface boundary condition: diurnal cycle on qsr
18   USE phycst          ! physical constants
19   !
20   USE fldread         ! read input fields
21   USE iom             ! IOM library
22   USE in_out_manager  ! I/O manager
23   USE lib_mpp         ! distribued memory computing library
24   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC sbc_flx       ! routine called by step.F90
30
31   INTEGER , PARAMETER ::   jpfld   = 5   ! maximum number of files to read
32   INTEGER , PARAMETER ::   jp_utau = 1   ! index of wind stress (i-component) file
33   INTEGER , PARAMETER ::   jp_vtau = 2   ! index of wind stress (j-component) file
34   INTEGER , PARAMETER ::   jp_qtot = 3   ! index of total (non solar+solar) heat file
35   INTEGER , PARAMETER ::   jp_qsr  = 4   ! index of solar heat file
36   INTEGER , PARAMETER ::   jp_emp  = 5   ! index of evaporation-precipation file
37   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf    ! structure of input fields (file informations, fields read)
38
39   !! * Substitutions
40#  include "vectopt_loop_substitute.h90"
41#  include "do_loop_substitute.h90"
42   !!----------------------------------------------------------------------
43   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
44   !! $Id$
45   !! Software governed by the CeCILL license (see ./LICENSE)
46   !!----------------------------------------------------------------------
47CONTAINS
48
49   SUBROUTINE sbc_flx( kt )
50      !!---------------------------------------------------------------------
51      !!                    ***  ROUTINE sbc_flx  ***
52      !!                   
53      !! ** Purpose :   provide at each time step the surface ocean fluxes
54      !!                (momentum, heat, freshwater and runoff)
55      !!
56      !! ** Method  : - READ each fluxes in NetCDF files:
57      !!                   i-component of the stress              utau  (N/m2)
58      !!                   j-component of the stress              vtau  (N/m2)
59      !!                   net downward heat flux                 qtot  (watt/m2)
60      !!                   net downward radiative flux            qsr   (watt/m2)
61      !!                   net upward freshwater (evapo - precip) emp   (kg/m2/s)
62      !!
63      !!      CAUTION :  - never mask the surface stress fields
64      !!                 - the stress is assumed to be in the (i,j) mesh referential
65      !!
66      !! ** Action  :   update at each time-step
67      !!              - utau, vtau  i- and j-component of the wind stress
68      !!              - taum        wind stress module at T-point
69      !!              - wndm        10m wind module at T-point
70      !!              - qns         non solar heat flux including heat flux due to emp
71      !!              - qsr         solar heat flux
72      !!              - emp         upward mass flux (evap. - precip.)
73      !!              - sfx         salt flux; set to zero at nit000 but possibly non-zero
74      !!                            if ice is present
75      !!----------------------------------------------------------------------
76      INTEGER, INTENT(in) ::   kt   ! ocean time step
77      !!
78      INTEGER  ::   ji, jj, jf            ! dummy indices
79      INTEGER  ::   ierror                ! return error code
80      INTEGER  ::   ios                   ! Local integer output status for namelist read
81      REAL(wp) ::   zfact                 ! temporary scalar
82      REAL(wp) ::   zrhoa  = 1.22         ! Air density kg/m3
83      REAL(wp) ::   zcdrag = 1.5e-3       ! drag coefficient
84      REAL(wp) ::   ztx, zty, zmod, zcoef ! temporary variables
85      !!
86      CHARACTER(len=100) ::  cn_dir                               ! Root directory for location of flx files
87      TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i                    ! array of namelist information structures
88      TYPE(FLD_N) ::   sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp  ! informations about the fields to be read
89      NAMELIST/namsbc_flx/ cn_dir, sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp
90      !!---------------------------------------------------------------------
91      !
92      IF( kt == nit000 ) THEN                ! First call kt=nit000 
93         ! set file information
94         READ  ( numnam_ref, namsbc_flx, IOSTAT = ios, ERR = 901)
95901      IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc_flx in reference namelist' )
96
97         READ  ( numnam_cfg, namsbc_flx, IOSTAT = ios, ERR = 902 )
98902      IF( ios >  0 )   CALL ctl_nam ( ios , 'namsbc_flx in configuration namelist' )
99         IF(lwm) WRITE ( numond, namsbc_flx ) 
100         !
101         !                                         ! check: do we plan to use ln_dm2dc with non-daily forcing?
102         IF( ln_dm2dc .AND. sn_qsr%freqh /= 24. )   &
103            &   CALL ctl_stop( 'sbc_blk_core: ln_dm2dc can be activated only with daily short-wave forcing' ) 
104         !
105         !                                         ! store namelist information in an array
106         slf_i(jp_utau) = sn_utau   ;   slf_i(jp_vtau) = sn_vtau
107         slf_i(jp_qtot) = sn_qtot   ;   slf_i(jp_qsr ) = sn_qsr 
108         slf_i(jp_emp ) = sn_emp
109         !
110         ALLOCATE( sf(jpfld), STAT=ierror )        ! set sf structure
111         IF( ierror > 0 ) THEN   
112            CALL ctl_stop( 'sbc_flx: unable to allocate sf structure' )   ;   RETURN 
113         ENDIF
114         DO ji= 1, jpfld
115            ALLOCATE( sf(ji)%fnow(jpi,jpj,1) )
116            IF( slf_i(ji)%ln_tint ) ALLOCATE( sf(ji)%fdta(jpi,jpj,1,2) )
117         END DO
118         !                                         ! fill sf with slf_i and control print
119         CALL fld_fill( sf, slf_i, cn_dir, 'sbc_flx', 'flux formulation for ocean surface boundary condition', 'namsbc_flx' )
120         !
121         sfx(:,:) = 0.0_wp                         ! salt flux due to freezing/melting (non-zero only if ice is present)
122         !
123      ENDIF
124
125      CALL fld_read( kt, nn_fsbc, sf )                            ! input fields provided at the current time-step
126     
127      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN                        ! update ocean fluxes at each SBC frequency
128
129         IF( ln_dm2dc ) THEN   ;   qsr(:,:) = sbc_dcy( sf(jp_qsr)%fnow(:,:,1) )   ! modify now Qsr to include the diurnal cycle
130         ELSE                  ;   qsr(:,:) =          sf(jp_qsr)%fnow(:,:,1)
131         ENDIF
132         DO_2D_11_11
133            utau(ji,jj) = sf(jp_utau)%fnow(ji,jj,1)
134            vtau(ji,jj) = sf(jp_vtau)%fnow(ji,jj,1)
135            qns (ji,jj) = sf(jp_qtot)%fnow(ji,jj,1) - sf(jp_qsr)%fnow(ji,jj,1)
136            emp (ji,jj) = sf(jp_emp )%fnow(ji,jj,1)
137         END_2D
138         !                                                        ! add to qns the heat due to e-p
139         qns(:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp        ! mass flux is at SST
140         !
141         qns(:,:) = qns(:,:) * tmask(:,:,1)
142         emp(:,:) = emp(:,:) * tmask(:,:,1)
143         !
144         !                                                        ! module of wind stress and wind speed at T-point
145         zcoef = 1. / ( zrhoa * zcdrag )
146         DO_2D_00_00
147            ztx = utau(ji-1,jj  ) + utau(ji,jj) 
148            zty = vtau(ji  ,jj-1) + vtau(ji,jj) 
149            zmod = 0.5 * SQRT( ztx * ztx + zty * zty )
150            taum(ji,jj) = zmod
151            wndm(ji,jj) = SQRT( zmod * zcoef )
152         END_2D
153         taum(:,:) = taum(:,:) * tmask(:,:,1) ; wndm(:,:) = wndm(:,:) * tmask(:,:,1)
154         CALL lbc_lnk( 'sbcflx', taum(:,:), 'T', 1. )   ;   CALL lbc_lnk( 'sbcflx', wndm(:,:), 'T', 1. )
155
156         IF( nitend-nit000 <= 100 .AND. lwp ) THEN                ! control print (if less than 100 time-step asked)
157            WRITE(numout,*) 
158            WRITE(numout,*) '        read daily momentum, heat and freshwater fluxes OK'
159            DO jf = 1, jpfld
160               IF( jf == jp_utau .OR. jf == jp_vtau )   zfact =     1.
161               IF( jf == jp_qtot .OR. jf == jp_qsr  )   zfact =     0.1
162               IF( jf == jp_emp                     )   zfact = 86400.
163               WRITE(numout,*) 
164               WRITE(numout,*) ' day: ', ndastp , TRIM(sf(jf)%clvar), ' * ', zfact
165            END DO
166         ENDIF
167         !
168      ENDIF
169      !
170   END SUBROUTINE sbc_flx
171
172   !!======================================================================
173END MODULE sbcflx
Note: See TracBrowser for help on using the repository browser.