Changeset 6140 for trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90
- Timestamp:
- 2015-12-21T12:35:23+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90
r5930 r6140 6 6 !! History : 1.0 ! 2005-12 (C. Talandier, G. Madec, V. Garnier) Original code 7 7 !! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option 8 !! 3.7 ! 2015-11 (J. Chanut) Suppression of filtered free surface9 8 !!---------------------------------------------------------------------- 10 9 … … 21 20 USE dynspg_exp ! surface pressure gradient (dyn_spg_exp routine) 22 21 USE dynspg_ts ! surface pressure gradient (dyn_spg_ts routine) 23 USE sbctide 24 USE updtide 22 USE sbctide ! 23 USE updtide ! 25 24 USE trd_oce ! trends: ocean variables 26 25 USE trddyn ! trend manager: dynamics … … 32 31 USE timing ! Timing 33 32 34 35 33 IMPLICIT NONE 36 34 PRIVATE … … 39 37 PUBLIC dyn_spg_init ! routine called by opa module 40 38 41 INTEGER :: nspg = 0 ! type of surface pressure gradient scheme defined from ln_dynspg_... 39 INTEGER :: nspg = 0 ! type of surface pressure gradient scheme defined from lk_dynspg_... 40 41 ! ! Parameter to control the surface pressure gradient scheme 42 INTEGER, PARAMETER :: np_TS = 1 ! split-explicit time stepping (Time-Splitting) 43 INTEGER, PARAMETER :: np_EXP = 0 ! explicit time stepping 44 INTEGER, PARAMETER :: np_NO =-1 ! no surface pressure gradient, no scheme 42 45 43 46 !! * Substitutions 44 # include "domzgr_substitute.h90"45 47 # include "vectopt_loop_substitute.h90" 46 48 !!---------------------------------------------------------------------- … … 55 57 !! *** ROUTINE dyn_spg *** 56 58 !! 57 !! ** Purpose : achieve the momentum time stepping by computing the 58 !! last trend, the surface pressure gradient including the 59 !! atmospheric pressure forcing (ln_apr_dyn=T), and performing 60 !! the Leap-Frog integration. 61 !!gm In the current version only the filtered solution provide 62 !!gm the after velocity, in the 2 other (ua,va) are still the trends 59 !! ** Purpose : compute surface pressure gradient including the 60 !! atmospheric pressure forcing (ln_apr_dyn=T). 63 61 !! 64 62 !! ** Method : Two schemes: 65 !! - explicit computation: the spg is evaluated at now66 !! - split-explicit computation: a time splitting technique is used63 !! - explicit : the spg is evaluated at now 64 !! - split-explicit : a time splitting technique is used 67 65 !! 68 66 !! ln_apr_dyn=T : the atmospheric pressure forcing is applied … … 73 71 !! period is used to prevent the divergence of odd and even time step. 74 72 !!---------------------------------------------------------------------- 75 !76 73 INTEGER, INTENT(in ) :: kt ! ocean time-step index 77 74 ! … … 84 81 IF( nn_timing == 1 ) CALL timing_start('dyn_spg') 85 82 ! 86 87 !!gm NOTA BENE : the dynspg_exp and dynspg_ts should be modified so that88 !!gm they return the after velocity, not the trends (as in trazdf_imp...)89 !!gm In this case, change/simplify dynnxt90 91 92 83 IF( l_trddyn ) THEN ! temporary save of ta and sa trends 93 CALL wrk_alloc( jpi, jpj, jpk,ztrdu, ztrdv )84 CALL wrk_alloc( jpi,jpj,jpk, ztrdu, ztrdv ) 94 85 ztrdu(:,:,:) = ua(:,:,:) 95 86 ztrdv(:,:,:) = va(:,:,:) 96 87 ENDIF 97 88 ! 98 89 IF( ln_apr_dyn & ! atmos. pressure 99 90 .OR. ( .NOT.ln_dynspg_ts .AND. (ln_tide_pot .AND. lk_tide) ) & ! tide potential (no time slitting) … … 107 98 END DO 108 99 ! 109 IF( ln_apr_dyn .AND. (.NOT. ln_dynspg_ts) ) THEN!== Atmospheric pressure gradient (added later in time-split case) ==!100 IF( ln_apr_dyn .AND. .NOT.ln_dynspg_ts ) THEN !== Atmospheric pressure gradient (added later in time-split case) ==! 110 101 zg_2 = grav * 0.5 111 102 DO jj = 2, jpjm1 ! gradient of Patm using inverse barometer ssh … … 133 124 ! 134 125 IF( nn_ice_embd == 2 ) THEN !== embedded sea ice: Pressure gradient due to snow-ice mass ==! 135 CALL wrk_alloc( jpi, jpj,zpice )126 CALL wrk_alloc( jpi,jpj, zpice ) 136 127 ! 137 128 zintp = REAL( MOD( kt-1, nn_fsbc ) ) / REAL( nn_fsbc ) … … 145 136 END DO 146 137 ! 147 CALL wrk_dealloc( jpi, jpj,zpice )138 CALL wrk_dealloc( jpi,jpj, zpice ) 148 139 ENDIF 149 140 ! 150 DO jk = 1, jpkm1 141 DO jk = 1, jpkm1 !== Add all terms to the general trend 151 142 DO jj = 2, jpjm1 152 143 DO ji = fs_2, fs_jpim1 ! vector opt. … … 156 147 END DO 157 148 END DO 158 149 ! 159 150 !!gm add here a call to dyn_trd for ice pressure gradient, the surf pressure trends ???? 160 161 ENDIF 162 163 SELECT CASE ( nspg ) ! compute surf. pressure gradient trend and add it to the general trend 164 ! 165 CASE ( 0 ) ; CALL dyn_spg_exp( kt ) ! explicit 166 CASE ( 1 ) ; CALL dyn_spg_ts ( kt ) ! time-splitting 167 ! 151 ! 152 ENDIF 153 ! 154 SELECT CASE ( nspg ) !== surface pressure gradient computed and add to the general trend ==! 155 CASE ( np_EXP ) ; CALL dyn_spg_exp( kt ) ! explicit 156 CASE ( np_TS ) ; CALL dyn_spg_ts ( kt ) ! time-splitting 168 157 END SELECT 169 158 ! 170 IF( l_trddyn ) THEN 159 IF( l_trddyn ) THEN ! save the surface pressure gradient trends for further diagnostics 171 160 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 172 161 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 173 162 CALL trd_dyn( ztrdu, ztrdv, jpdyn_spg, kt ) 174 ! 175 CALL wrk_dealloc( jpi, jpj, jpk, ztrdu, ztrdv ) 176 ENDIF 177 ! ! print mean trends (used for debugging) 163 CALL wrk_dealloc( jpi,jpj,jpk, ztrdu, ztrdv ) 164 ENDIF 165 ! ! print mean trends (used for debugging) 178 166 IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' spg - Ua: ', mask1=umask, & 179 167 & tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) … … 191 179 !! surface pressure gradient schemes 192 180 !!---------------------------------------------------------------------- 193 INTEGER :: ioptio, ios 194 ! 195 NAMELIST/namdyn_spg/ ln_dynspg_exp , ln_dynspg_ts,&196 & ln_bt_fw, ln_bt_av , ln_bt_auto,&197 & nn_baro , rn_bt_cmax, nn_bt_flt181 INTEGER :: ioptio, ios ! local integers 182 ! 183 NAMELIST/namdyn_spg/ ln_dynspg_exp , ln_dynspg_ts, & 184 & ln_bt_fw, ln_bt_av , ln_bt_auto , & 185 & nn_baro , rn_bt_cmax, nn_bt_flt 198 186 !!---------------------------------------------------------------------- 199 187 ! … … 202 190 REWIND( numnam_ref ) ! Namelist namdyn_spg in reference namelist : Free surface 203 191 READ ( numnam_ref, namdyn_spg, IOSTAT = ios, ERR = 901) 204 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_spg in reference namelist', lwp )205 192 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_spg in reference namelist', lwp ) 193 ! 206 194 REWIND( numnam_cfg ) ! Namelist namdyn_spg in configuration namelist : Free surface 207 195 READ ( numnam_cfg, namdyn_spg, IOSTAT = ios, ERR = 902 ) 208 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_spg in configuration namelist', lwp )196 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_spg in configuration namelist', lwp ) 209 197 IF(lwm) WRITE ( numond, namdyn_spg ) 210 198 ! … … 216 204 WRITE(numout,*) ' Free surface with time splitting ln_dynspg_ts = ', ln_dynspg_ts 217 205 ENDIF 218 219 IF( ln_dynspg_ts ) THEN 220 CALL dyn_spg_ts_init( nit000 ) ! do it first, to set nn_baro, used to allocate some arrays later on 221 ! ! allocate dyn_spg arrays 222 IF( dyn_spg_ts_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_spg_init: failed to allocate dynspg_ts arrays') 223 IF ((neuler/=0).AND.(ln_bt_fw)) CALL ts_rst( nit000, 'READ' ) 224 ENDIF 225 226 ! ! Control of surface pressure gradient scheme options 227 ioptio = 0 228 IF(ln_dynspg_exp) ioptio = ioptio + 1 229 IF(ln_dynspg_ts ) ioptio = ioptio + 1 230 ! 231 IF( ioptio > 1 .OR. ( ioptio == 0 .AND. .NOT. lk_c1d ) ) & 232 & CALL ctl_stop( ' Choose only one surface pressure gradient scheme' ) 233 IF( ln_dynspg_ts .AND. ln_isfcav ) & 234 & CALL ctl_stop( ' dynspg_ts not tested with ice shelf cavity ' ) 235 ! 236 IF( ln_dynspg_exp) nspg = 0 237 IF( ln_dynspg_ts ) nspg = 1 206 ! ! Control of surface pressure gradient scheme options 207 ; nspg = np_NO ; ioptio = 0 208 IF( ln_dynspg_exp ) THEN ; nspg = np_EXP ; ioptio = ioptio + 1 ; ENDIF 209 IF( ln_dynspg_ts ) THEN ; nspg = np_TS ; ioptio = ioptio + 1 ; ENDIF 210 ! 211 IF( ioptio > 1 ) CALL ctl_stop( 'Choose only one surface pressure gradient scheme' ) 212 IF( ioptio == 0 ) CALL ctl_warn( 'NO surface pressure gradient trend in momentum Eqs.' ) 213 IF( ln_dynspg_exp .AND. ln_isfcav ) & 214 & CALL ctl_stop( ' dynspg_exp not tested with ice shelf cavity ' ) 238 215 ! 239 216 IF(lwp) THEN 240 217 WRITE(numout,*) 241 IF( nspg == 0 ) WRITE(numout,*) ' explicit free surface' 242 IF( nspg == 1 ) WRITE(numout,*) ' free surface with time splitting scheme' 218 IF( nspg == np_EXP ) WRITE(numout,*) ' explicit free surface' 219 IF( nspg == np_TS ) WRITE(numout,*) ' free surface with time splitting scheme' 220 IF( nspg == np_NO ) WRITE(numout,*) ' No surface surface pressure gradient trend in momentum Eqs.' 221 ENDIF 222 ! 223 IF( nspg == np_TS ) THEN ! split-explicit scheme initialisation 224 CALL dyn_spg_ts_init ! do it first: set nn_baro used to allocate some arrays later on 225 IF( dyn_spg_ts_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_spg_init: failed to allocate dynspg_ts arrays' ) 226 IF( neuler/=0 .AND. ln_bt_fw ) CALL ts_rst( nit000, 'READ' ) 243 227 ENDIF 244 228 !
Note: See TracChangeset
for help on using the changeset viewer.