Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
MODULE stpctl
!!======================================================================
!! *** MODULE stpctl ***
!! Ocean run control : gross check of the ocean time stepping
!! version for STATION_ASF test-case
!!======================================================================
!! History : OPA ! 1991-03 (G. Madec) Original code
!! 6.0 ! 1992-06 (M. Imbard)
!! 8.0 ! 1997-06 (A.M. Treguier)
!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
!! 2.0 ! 2009-07 (G. Madec) Add statistic for time-spliting
!! 3.7 ! 2016-09 (G. Madec) Remove solver
!! 4.0 ! 2017-04 (G. Madec) regroup global communications
!!----------------------------------------------------------------------
!!----------------------------------------------------------------------
!! stp_ctl : Control the run
!!----------------------------------------------------------------------
USE dom_oce ! ocean space and time domain variables
USE sbc_oce ! surface fluxes and stuff
!
USE diawri ! Standard run outputs (dia_wri_state routine)
USE in_out_manager ! I/O manager
USE lbclnk ! ocean lateral boundary conditions (or mpp link)
USE lib_mpp ! distributed memory computing
!
USE netcdf ! NetCDF library
IMPLICIT NONE
PRIVATE
PUBLIC stp_ctl ! routine called by step.F90
INTEGER, PARAMETER :: jpvar = 3
INTEGER :: nrunid ! netcdf file id
INTEGER, DIMENSION(jpvar) :: nvarid ! netcdf variable id
!!----------------------------------------------------------------------
!! NEMO/SAS 4.0 , NEMO Consortium (2018)
!! $Id: stpctl.F90 10603 2019-01-29 11:18:09Z acc $
!! Software governed by the CeCILL license (see ./LICENSE)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE stp_ctl( kt, Kmm )
!!----------------------------------------------------------------------
!! *** ROUTINE stp_ctl ***
!!
!! ** Purpose : Control the run
!!
!! ** Method : - Save the time step in numstp
!! - Stop the run IF problem encountered by setting nstop > 0
!! Problems checked: wind stress module max larger than 5 N/m^2
!! non-solar heat flux max larger than 2000 W/m^2
!! Evaporation-Precip max larger than 1.E-3 kg/m^2/s
!!
!! ** Actions : "time.step" file = last ocean time-step
!! "run.stat" file = run statistics
!! nstop indicator sheared among all local domain
!!----------------------------------------------------------------------
INTEGER, INTENT(in ) :: kt ! ocean time-step index
INTEGER, INTENT(in ) :: Kmm ! ocean time level index
!!
INTEGER, PARAMETER :: jptst = 3
INTEGER :: ji ! dummy loop indices
INTEGER :: idtime, istatus
INTEGER , DIMENSION(jptst) :: iareasum, iareamin, iareamax
INTEGER , DIMENSION(3,jptst) :: iloc ! min/max loc indices
REAL(wp) :: zzz ! local real
REAL(wp), DIMENSION(jpvar+1) :: zmax
REAL(wp), DIMENSION(jptst) :: zmaxlocal
LOGICAL :: ll_wrtstp, ll_colruns, ll_wrtruns, ll_0oce
LOGICAL, DIMENSION(jpi,jpj) :: llmsk
CHARACTER(len=20) :: clname
!!----------------------------------------------------------------------
IF( nstop > 0 .AND. ngrdstop > -1 ) RETURN ! stpctl was already called by a child grid
!
ll_wrtstp = ( MOD( kt-nit000, sn_cfctl%ptimincr ) == 0 ) .OR. ( kt == nitend )
ll_colruns = sn_cfctl%l_runstat .AND. ll_wrtstp .AND. jpnij > 1
ll_wrtruns = sn_cfctl%l_runstat .AND. ll_wrtstp .AND. lwm
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
!
IF( kt == nit000 ) THEN
!
IF( lwp ) THEN
WRITE(numout,*)
WRITE(numout,*) 'stp_ctl : time-stepping control'
WRITE(numout,*) '~~~~~~~'
ENDIF
! ! open time.step ascii file, done only by 1st subdomain
IF( lwm ) CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
!
IF( ll_wrtruns ) THEN
! ! open run.stat ascii file, done only by 1st subdomain
CALL ctl_opn( numrun, 'run.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
! ! open run.stat.nc netcdf file, done only by 1st subdomain
clname = 'run.stat.nc'
IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//"_"//TRIM(clname)
istatus = NF90_CREATE( TRIM(clname), NF90_CLOBBER, nrunid )
istatus = NF90_DEF_DIM( nrunid, 'time', NF90_UNLIMITED, idtime )
istatus = NF90_DEF_VAR( nrunid, 'tau_max', NF90_DOUBLE, (/ idtime /), nvarid(1) )
istatus = NF90_DEF_VAR( nrunid, 'qns_max', NF90_DOUBLE, (/ idtime /), nvarid(2) )
istatus = NF90_DEF_VAR( nrunid, 'emp_max', NF90_DOUBLE, (/ idtime /), nvarid(3) )
istatus = NF90_ENDDEF(nrunid)
ENDIF
!
ENDIF
!
! !== write current time step ==!
! !== done only by 1st subdomain at writting timestep ==!
IF( lwm .AND. ll_wrtstp ) THEN
WRITE ( numstp, '(1x, i8)' ) kt
REWIND( numstp )
ENDIF
! !== test of local extrema ==!
! !== done by all processes at every time step ==!
!
llmsk( 1:nn_hls,:) = .FALSE. ! exclude halos from the checked region
llmsk(Nie0+1: jpi,:) = .FALSE.
llmsk(:, 1:nn_hls) = .FALSE.
llmsk(:,Nje0+1: jpj) = .FALSE.
!
llmsk(Nis0:Nie0,Njs0:Nje0) = tmask(Nis0:Nie0,Njs0:Nje0,1) == 1._wp ! test only the inner domain
!
ll_0oce = .NOT. ANY( llmsk(:,:) ) ! no ocean point in the inner domain?
!
zmax(1) = MAXVAL( taum(:,:) , mask = llmsk(A2D(0)) ) ! max wind stress module
zmax(2) = MAXVAL( ABS( qns(:,:) ), mask = llmsk(A2D(0)) ) ! max non-solar heat flux
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
zmax(3) = MAXVAL( ABS( emp(:,:) ), mask = llmsk ) ! max E-P
zmax(jpvar+1) = REAL( nstop, wp ) ! stop indicator
!
! !== get global extrema ==!
! !== done by all processes if writting run.stat ==!
IF( ll_colruns ) THEN
zmaxlocal(:) = zmax(1:jptst)
CALL mpp_max( "stpctl", zmax ) ! max over the global domain: ok even of ll_0oce = .true.
nstop = NINT( zmax(jpvar+1) ) ! update nstop indicator (now sheared among all local domains)
ELSE
! if no ocean point: MAXVAL returns -HUGE => we must overwrite this value to avoid error handling bellow.
IF( ll_0oce ) zmax(1:jptst) = 0._wp ! default "valid" values...
ENDIF
! !== write "run.stat" files ==!
! !== done only by 1st subdomain at writting timestep ==!
IF( ll_wrtruns ) THEN
WRITE(numrun,9500) kt, zmax(1:jptst)
DO ji = 1, jpvar
istatus = NF90_PUT_VAR( nrunid, nvarid(ji), (/zmax(ji)/), (/kt/), (/1/) )
END DO
IF( kt == nitend ) istatus = NF90_CLOSE(nrunid)
END IF
! !== error handling ==!
! !== done by all processes at every time step ==!
!
IF( zmax(1) > 5._wp .OR. & ! too large wind stress ( > 5 N/m^2 )
& zmax(2) > 2000._wp .OR. & ! too large non-solar heat flux ( > 2000 W/m^2 )
& zmax(3) > 1.E-3_wp .OR. & ! too large net freshwater flux ( > 1.E-3 kg/m^2/s )
& ISNAN( SUM(zmax(1:jptst)) ) .OR. & ! NaN encounter in the tests
& ABS( SUM(zmax(1:jptst)) ) > HUGE(1._wp) ) THEN ! Infinity encounter in the tests
!
iloc(:,:) = 0
IF( ll_colruns ) THEN ! zmax is global, so it is the same on all subdomains -> no dead lock with mpp_maxloc
! first: close the netcdf file, so we can read it
IF( lwm .AND. kt /= nitend ) istatus = NF90_CLOSE(nrunid)
! get global loc on the min/max
CALL mpp_maxloc( 'stpctl', taum(:,:) , llmsk(A2D(0)), zzz, iloc(1:2,1) ) ! mpp_maxloc ok if mask = F
CALL mpp_maxloc( 'stpctl',ABS( qns(:,:) ), llmsk(A2D(0)), zzz, iloc(1:2,2) )
CALL mpp_minloc( 'stpctl',ABS( emp(:,:) ), llmsk , zzz, iloc(1:2,3) )
! find which subdomain has the max.
iareamin(:) = jpnij+1 ; iareamax(:) = 0 ; iareasum(:) = 0
DO ji = 1, jptst
IF( zmaxlocal(ji) == zmax(ji) ) THEN
iareamin(ji) = narea ; iareamax(ji) = narea ; iareasum(ji) = 1
ENDIF
END DO
CALL mpp_min( "stpctl", iareamin ) ! min over the global domain
CALL mpp_max( "stpctl", iareamax ) ! max over the global domain
CALL mpp_sum( "stpctl", iareasum ) ! sum over the global domain
ELSE ! find local min and max locations:
! if we are here, this means that the subdomain contains some oce points -> no need to test the mask used in maxloc
iloc(1:2,1) = MAXLOC( taum(:,:) , mask = llmsk(A2D(0)) )
iloc(1:2,2) = MAXLOC( ABS( qns(:,:) ), mask = llmsk(A2D(0)) )
iloc(1:2,3) = MINLOC( ABS( emp(:,:) ), mask = llmsk )
DO ji = 1, jptst ! local domain indices ==> global domain indices, excluding halos
iloc(1:2,ji) = (/ mig(iloc(1,ji),0), mjg(iloc(2,ji),0) /)
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
END DO
iareamin(:) = narea ; iareamax(:) = narea ; iareasum(:) = 1 ! this is local information
ENDIF
!
WRITE(ctmp1,*) ' stp_ctl: |tau_mod| > 5 N/m2 or |qns| > 2000 W/m2 or |emp| > 1.E-3 or NaN encounter in the tests'
CALL wrt_line( ctmp2, kt, '|tau| max', zmax(1), iloc(:,1), iareasum(1), iareamin(1), iareamax(1) )
CALL wrt_line( ctmp3, kt, '|qns| max', zmax(2), iloc(:,2), iareasum(2), iareamin(2), iareamax(2) )
CALL wrt_line( ctmp4, kt, 'emp max', zmax(3), iloc(:,3), iareasum(3), iareamin(3), iareamax(3) )
IF( Agrif_Root() ) THEN
WRITE(ctmp6,*) ' ===> output of last computed fields in output.abort* files'
ELSE
WRITE(ctmp6,*) ' ===> output of last computed fields in '//TRIM(Agrif_CFixed())//'_output.abort* files'
ENDIF
!
CALL dia_wri_state( Kmm, 'output.abort' ) ! create an output.abort file
!
IF( ll_colruns .OR. jpnij == 1 ) THEN ! all processes synchronized -> use lwp to print in opened ocean.output files
IF(lwp) THEN ; CALL ctl_stop( ctmp1, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ' ', ctmp6 )
ELSE ; nstop = MAX(1, nstop) ! make sure nstop > 0 (automatically done when calling ctl_stop)
ENDIF
ELSE ! only mpi subdomains with errors are here -> STOP now
CALL ctl_stop( 'STOP', ctmp1, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ' ', ctmp6 )
ENDIF
!
ENDIF
!
IF( nstop > 0 ) THEN ! an error was detected and we did not abort yet...
ngrdstop = Agrif_Fixed() ! store which grid got this error
IF( .NOT. ll_colruns .AND. jpnij > 1 ) CALL ctl_stop( 'STOP' ) ! we must abort here to avoid MPI deadlock
ENDIF
!
9500 FORMAT(' it :', i8, ' tau_max: ', D23.16, ' |qns|_max: ', D23.16,' |emp|_max: ', D23.16)
!
END SUBROUTINE stp_ctl
SUBROUTINE wrt_line( cdline, kt, cdprefix, pval, kloc, ksum, kmin, kmax )
!!----------------------------------------------------------------------
!! *** ROUTINE wrt_line ***
!!
!! ** Purpose : write information line
!!
!!----------------------------------------------------------------------
CHARACTER(len=*), INTENT( out) :: cdline
CHARACTER(len=*), INTENT(in ) :: cdprefix
REAL(wp), INTENT(in ) :: pval
INTEGER, DIMENSION(3), INTENT(in ) :: kloc
INTEGER, INTENT(in ) :: kt, ksum, kmin, kmax
!
CHARACTER(len=80) :: clsuff
CHARACTER(len=9 ) :: clkt, clsum, clmin, clmax
CHARACTER(len=9 ) :: cli, clj, clk
CHARACTER(len=1 ) :: clfmt
CHARACTER(len=4 ) :: cl4 ! needed to be able to compile with Agrif, I don't know why
INTEGER :: ifmtk
!!----------------------------------------------------------------------
WRITE(clkt , '(i9)') kt
WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpnij ,wp))) + 1 ! how many digits to we need to write ? (we decide max = 9)
!!! WRITE(clsum, '(i'//clfmt//')') ksum ! this is creating a compilation error with AGRIF
cl4 = '(i'//clfmt//')' ; WRITE(clsum, cl4) ksum
WRITE(clfmt, '(i1)') INT(LOG10(REAL(MAX(1,jpnij-1),wp))) + 1 ! how many digits to we need to write ? (we decide max = 9)
cl4 = '(i'//clfmt//')' ; WRITE(clmin, cl4) kmin-1
WRITE(clmax, cl4) kmax-1
!
WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpiglo,wp))) + 1 ! how many digits to we need to write jpiglo? (we decide max = 9)
cl4 = '(i'//clfmt//')' ; WRITE(cli, cl4) kloc(1) ! this is ok with AGRIF
WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpjglo,wp))) + 1 ! how many digits to we need to write jpjglo? (we decide max = 9)
cl4 = '(i'//clfmt//')' ; WRITE(clj, cl4) kloc(2) ! this is ok with AGRIF
!
IF( ksum == 1 ) THEN ; WRITE(clsuff,9100) TRIM(clmin)
ELSE ; WRITE(clsuff,9200) TRIM(clsum), TRIM(clmin), TRIM(clmax)
ENDIF
IF(kloc(3) == 0) THEN
ifmtk = INT(LOG10(REAL(jpk,wp))) + 1 ! how many digits to we need to write jpk? (we decide max = 9)
clk = REPEAT(' ', ifmtk) ! create the equivalent in blank string
WRITE(cdline,9300) TRIM(ADJUSTL(clkt)), TRIM(ADJUSTL(cdprefix)), pval, TRIM(cli), TRIM(clj), clk(1:ifmtk), TRIM(clsuff)
ELSE
WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpk,wp))) + 1 ! how many digits to we need to write jpk? (we decide max = 9)
!!! WRITE(clk, '(i'//clfmt//')') kloc(3) ! this is creating a compilation error with AGRIF
cl4 = '(i'//clfmt//')' ; WRITE(clk, cl4) kloc(3) ! this is ok with AGRIF
WRITE(cdline,9400) TRIM(ADJUSTL(clkt)), TRIM(ADJUSTL(cdprefix)), pval, TRIM(cli), TRIM(clj), TRIM(clk), TRIM(clsuff)
ENDIF
!
9100 FORMAT('MPI rank ', a)
9200 FORMAT('found in ', a, ' MPI tasks, spread out among ranks ', a, ' to ', a)
9300 FORMAT('kt ', a, ' ', a, ' ', 1pg11.4, ' at i j ', a, ' ', a, ' ', a, ' ', a)
9400 FORMAT('kt ', a, ' ', a, ' ', 1pg11.4, ' at i j k ', a, ' ', a, ' ', a, ' ', a)
!
END SUBROUTINE wrt_line
!!======================================================================
END MODULE stpctl