- Timestamp:
- 2015-11-06T19:07:02+01:00 (9 years ago)
- Location:
- branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 2 deleted
- 13 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/ASM/asmbkg.F90
r5836 r5868 34 34 USE zdfmxl ! Mixed layer depth 35 35 USE dom_oce, ONLY : ndastp 36 USE sol_oce, ONLY : gcx ! Solver variables defined in memory37 36 USE in_out_manager ! I/O manager 38 37 USE iom ! I/O module … … 114 113 CALL iom_rstput( kt, nitbkg_r, inum, 'en' , en ) 115 114 #endif 116 CALL iom_rstput( kt, nitbkg_r, inum, 'gcx' , gcx )117 115 ! 118 116 CALL iom_close( inum ) -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90
r4689 r5868 35 35 PRIVATE 36 36 37 PUBLIC bdy_dyn ! routine called in dynspg_flt (if lk_dynspg_flt) or 38 ! dyn_nxt (if lk_dynspg_ts or lk_dynspg_exp) 37 PUBLIC bdy_dyn ! routine called in dyn_nxt 39 38 40 39 # include "domzgr_substitute.h90" -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90
r5836 r5868 10 10 !! 3.4 ! 2011 (D. Storkey) rewrite in preparation for OBC-BDY merge 11 11 !!---------------------------------------------------------------------- 12 #if defined key_bdy && defined key_dynspg_flt12 #if defined key_bdy 13 13 !!---------------------------------------------------------------------- 14 !! 'key_bdy' AND unstructured open boundary conditions 15 !! 'key_dynspg_flt' filtered free surface 14 !! 'key_bdy' unstructured open boundary conditions 16 15 !!---------------------------------------------------------------------- 17 16 USE oce ! ocean dynamics and tracers … … 30 29 PRIVATE 31 30 32 PUBLIC bdy_vol ! routine called by dynspg_flt.h9031 PUBLIC bdy_vol 33 32 34 33 !! * Substitutions … … 45 44 !! *** ROUTINE bdyvol *** 46 45 !! 47 !! ** Purpose : This routine is called in dynspg_flt to control48 !! the volume of the system.A correction velocity is calculated46 !! ** Purpose : This routine controls the volume of the system. 47 !! A correction velocity is calculated 49 48 !! to correct the total transport through the unstructured OBC. 50 49 !! The total depth used is constant (H0) to be consistent with the -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r5836 r5868 30 30 USE zdf_oce ! ocean vertical physics 31 31 USE ldftra ! lateral physics: eddy diffusivity coef. 32 USE sol_oce ! solver variables33 32 USE sbc_oce ! Surface boundary condition: ocean fields 34 33 USE sbc_ice ! Surface boundary condition: ice fields -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90
r5836 r5868 35 35 USE dtauvd ! data: U & V current (dta_uvd routine) 36 36 USE domvvl ! varying vertical mesh 37 USE dynspg_oce ! pressure gradient schemes38 USE dynspg_flt ! filtered free surface39 USE sol_oce ! ocean solver variables40 37 ! 41 38 USE in_out_manager ! I/O manager … … 133 130 ENDIF 134 131 ! 135 IF( lk_agrif ) THEN ! read free surface arrays in restart file136 IF( ln_rstart ) THEN137 IF( lk_dynspg_flt ) THEN ! read or initialize the following fields138 ! ! gcx, gcxb for agrif_opa_init139 IF( sol_oce_alloc() > 0 ) CALL ctl_stop('agrif sol_oce_alloc: allocation of arrays failed')140 CALL flt_rst( nit000, 'READ' )141 ENDIF142 ENDIF ! explicit case not coded yet with AGRIF143 ENDIF144 !145 132 ! 146 133 ! Initialize "now" and "before" barotropic velocities: … … 445 432 !! p=integral [ rau*g dz ] 446 433 !!---------------------------------------------------------------------- 447 USE dynspg ! surface pressure gradient (dyn_spg routine)448 434 USE divhor ! hor. divergence (div_hor routine) 449 435 USE lbclnk ! ocean lateral boundary condition (or mpp link) 450 436 ! 451 437 INTEGER :: ji, jj, jk ! dummy loop indices 452 INTEGER :: indic ! ???453 438 REAL(wp) :: zmsv, zphv, zmsu, zphu, zalfg ! temporary scalars 454 439 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprn … … 517 502 vb(:,:,:) = vn(:,:,:) 518 503 519 ! WARNING !!!!!520 ! after initializing u and v, we need to calculate the initial streamfunction bsf.521 ! Otherwise, only the trend will be computed and the model will blow up (inconsistency).522 ! to do that, we call dyn_spg with a special trick:523 ! we fill ua and va with the velocities divided by dt, and the streamfunction will be brought to the524 ! right value assuming the velocities have been set up in one time step.525 ! we then set bsfd to zero (first guess for next step is d(psi)/dt = 0.)526 ! sets up s false trend to calculate the barotropic streamfunction.527 528 ua(:,:,:) = ub(:,:,:) / rdt529 va(:,:,:) = vb(:,:,:) / rdt530 531 ! calls dyn_spg. we assume euler time step, starting from rest.532 indic = 0533 CALL dyn_spg( nit000, indic ) ! surface pressure gradient534 !535 ! the new velocity is ua*rdt536 !537 CALL lbc_lnk( ua, 'U', -1. )538 CALL lbc_lnk( va, 'V', -1. )539 540 ub(:,:,:) = ua(:,:,:) * rdt541 vb(:,:,:) = va(:,:,:) * rdt542 ua(:,:,:) = 0.e0543 va(:,:,:) = 0.e0544 un(:,:,:) = ub(:,:,:)545 vn(:,:,:) = vb(:,:,:)546 504 ! 547 505 !!gm Check here call to div_hor should not be necessary -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r5643 r5868 78 78 !! (lk_vvl=T), the leap-frog is applied on thickness weighted 79 79 !! velocity. 80 !! Note also that in filtered free surface (lk_dynspg_flt=T),81 !! the time stepping has already been done in dynspg module82 80 !! 83 81 !! * Apply lateral boundary conditions on after velocity … … 101 99 INTEGER :: ji, jj, jk ! dummy loop indices 102 100 INTEGER :: iku, ikv ! local integers 103 #if ! defined key_dynspg_flt104 101 REAL(wp) :: z2dt ! temporary scalar 105 #endif106 102 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zec ! local scalars 107 103 REAL(wp) :: zve3a, zve3n, zve3b, zvf, z1_2dt ! - - … … 120 116 IF(lwp) WRITE(numout,*) '~~~~~~~' 121 117 ENDIF 122 123 #if defined key_dynspg_flt124 !125 ! Next velocity : Leap-frog time stepping already done in dynspg_flt.F routine126 ! -------------127 128 ! Update after velocity on domain lateral boundaries (only local domain required)129 ! --------------------------------------------------130 CALL lbc_lnk( ua, 'U', -1. ) ! local domain boundaries131 CALL lbc_lnk( va, 'V', -1. )132 !133 #else134 118 135 119 # if defined key_dynspg_exp … … 201 185 CALL Agrif_dyn( kt ) !* AGRIF zoom boundaries 202 186 # endif 203 #endif204 187 205 188 IF( l_trddyn ) THEN ! prepare the atf trend computation + some diagnostics -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90
r5836 r5868 21 21 USE dynspg_exp ! surface pressure gradient (dyn_spg_exp routine) 22 22 USE dynspg_ts ! surface pressure gradient (dyn_spg_ts routine) 23 USE dynspg_flt ! surface pressure gradient (dyn_spg_flt routine)24 23 USE dynadv ! dynamics: vector invariant versus flux form 25 24 USE dynhpg, ONLY: ln_dynhpg_imp … … 32 31 USE in_out_manager ! I/O manager 33 32 USE lib_mpp ! MPP library 34 USE solver ! solver initialization35 33 USE wrk_nemo ! Memory Allocation 36 34 USE timing ! Timing … … 55 53 CONTAINS 56 54 57 SUBROUTINE dyn_spg( kt , kindic)55 SUBROUTINE dyn_spg( kt ) 58 56 !!---------------------------------------------------------------------- 59 57 !! *** ROUTINE dyn_spg *** … … 80 78 ! 81 79 INTEGER, INTENT(in ) :: kt ! ocean time-step index 82 INTEGER, INTENT( out) :: kindic ! solver flag83 80 ! 84 81 INTEGER :: ji, jj, jk ! dummy loop indices … … 171 168 CASE ( 0 ) ; CALL dyn_spg_exp( kt ) ! explicit 172 169 CASE ( 1 ) ; CALL dyn_spg_ts ( kt ) ! time-splitting 173 CASE ( 2 ) ; CALL dyn_spg_flt( kt, kindic ) ! filtered174 170 ! 175 171 END SELECT 176 172 ! 177 173 IF( l_trddyn ) THEN ! save the surface pressure gradient trends for further diagnostics 178 SELECT CASE ( nspg ) 179 CASE ( 0, 1 ) 180 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 181 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 182 CASE( 2 ) 183 z2dt = 2. * rdt 184 IF( neuler == 0 .AND. kt == nit000 ) z2dt = rdt 185 ztrdu(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) / z2dt - ztrdu(:,:,:) 186 ztrdv(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) / z2dt - ztrdv(:,:,:) 187 END SELECT 174 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 175 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 188 176 CALL trd_dyn( ztrdu, ztrdv, jpdyn_spg, kt ) 189 177 ! … … 217 205 WRITE(numout,*) ' Explicit free surface lk_dynspg_exp = ', lk_dynspg_exp 218 206 WRITE(numout,*) ' Free surface with time splitting lk_dynspg_ts = ', lk_dynspg_ts 219 WRITE(numout,*) ' Filtered free surface cst volume lk_dynspg_flt = ', lk_dynspg_flt220 207 ENDIF 221 208 … … 233 220 IF(lk_dynspg_exp) ioptio = ioptio + 1 234 221 IF(lk_dynspg_ts ) ioptio = ioptio + 1 235 IF(lk_dynspg_flt) ioptio = ioptio + 1236 222 ! 237 223 IF( ioptio > 1 .OR. ( ioptio == 0 .AND. .NOT. lk_c1d ) ) & … … 242 228 IF( lk_dynspg_exp) nspg = 0 243 229 IF( lk_dynspg_ts ) nspg = 1 244 IF( lk_dynspg_flt) nspg = 2 230 245 231 ! 246 232 IF(lwp) THEN … … 248 234 IF( nspg == 0 ) WRITE(numout,*) ' explicit free surface' 249 235 IF( nspg == 1 ) WRITE(numout,*) ' free surface with time splitting scheme' 250 IF( nspg == 2 ) WRITE(numout,*) ' filtered free surface' 251 ENDIF 252 253 #if defined key_dynspg_flt 254 CALL solver_init( nit000 ) ! Elliptic solver initialisation 255 #endif 236 ENDIF 256 237 ! ! Control of hydrostatic pressure choice 257 238 IF( lk_dynspg_ts .AND. ln_dynhpg_imp ) CALL ctl_stop( 'Semi-implicit hpg not compatible with time splitting' ) -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_oce.F90
r5836 r5868 26 26 #else 27 27 LOGICAL, PUBLIC, PARAMETER :: lk_dynspg_ts = .FALSE. !: Free surface with time splitting flag 28 #endif29 #if defined key_dynspg_flt30 LOGICAL, PUBLIC, PARAMETER :: lk_dynspg_flt = .TRUE. !: Filtered free surface cst volume flag31 #else32 LOGICAL, PUBLIC, PARAMETER :: lk_dynspg_flt = .FALSE. !: Filtered free surface cst volume flag33 28 #endif 34 29 ! !!! Time splitting scheme (key_dynspg_ts) -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trd_oce.F90
r5215 r5868 71 71 INTEGER, PUBLIC, PARAMETER :: jpdyn_bfri = 12 !: implicit bottom friction (ln_bfrimp=.TRUE.) 72 72 INTEGER, PUBLIC, PARAMETER :: jpdyn_ken = 13 !: use for calculation of KE 73 INTEGER, PUBLIC, PARAMETER :: jpdyn_spgflt = 14 !: filter contribution to surface pressure gradient (spg_flt)74 INTEGER, PUBLIC, PARAMETER :: jpdyn_spgexp = 15 !: explicit contribution to surface pressure gradient (spg_flt)75 73 ! 76 74 !!---------------------------------------------------------------------- -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trddyn.F90
r5215 r5868 113 113 CASE( jpdyn_spg ) ; CALL iom_put( "utrd_spg", putrd ) ! surface pressure gradient 114 114 CALL iom_put( "vtrd_spg", pvtrd ) 115 CASE( jpdyn_spgexp ); CALL iom_put( "utrd_spgexp", putrd ) ! surface pressure gradient (explicit)116 CALL iom_put( "vtrd_spgexp", pvtrd )117 CASE( jpdyn_spgflt ); CALL iom_put( "utrd_spgflt", putrd ) ! surface pressure gradient (filtered)118 CALL iom_put( "vtrd_spgflt", pvtrd )119 115 CASE( jpdyn_pvo ) ; CALL iom_put( "utrd_pvo", putrd ) ! planetary vorticity 120 116 CALL iom_put( "vtrd_pvo", pvtrd ) -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90
r5836 r5868 120 120 CASE( jpdyn_hpg ) ; CALL iom_put( "ketrd_hpg", zke ) ! hydrostatic pressure gradient 121 121 CASE( jpdyn_spg ) ; CALL iom_put( "ketrd_spg", zke ) ! surface pressure gradient 122 CASE( jpdyn_spgexp ); CALL iom_put( "ketrd_spgexp", zke ) ! surface pressure gradient (explicit)123 CASE( jpdyn_spgflt ); CALL iom_put( "ketrd_spgflt", zke ) ! surface pressure gradient (filter)124 122 CASE( jpdyn_pvo ) ; CALL iom_put( "ketrd_pvo", zke ) ! planetary vorticity 125 123 CASE( jpdyn_rvo ) ; CALL iom_put( "ketrd_rvo", zke ) ! relative vorticity (or metric term) -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/step.F90
r5836 r5868 214 214 #endif 215 215 CALL dyn_hpg( kstp ) ! horizontal gradient of Hydrostatic pressure 216 CALL dyn_spg( kstp , indic )! surface pressure gradient216 CALL dyn_spg( kstp ) ! surface pressure gradient 217 217 218 218 ua_sv(:,:,:) = ua(:,:,:) ! Save trends (barotropic trend has been fully updated at this stage) … … 329 329 CALL dyn_bfr( kstp ) ! bottom friction 330 330 CALL dyn_zdf( kstp ) ! vertical diffusion 331 CALL dyn_spg( kstp , indic )! surface pressure gradient331 CALL dyn_spg( kstp ) ! surface pressure gradient 332 332 ENDIF 333 333 CALL dyn_nxt( kstp ) ! lateral velocity at next time step -
branches/2015/dev_r5847_MERCATOR9_solveur_simplification/NEMOGCM/NEMO/OPA_SRC/stpctl.F90
r3294 r5868 16 16 USE oce ! ocean dynamics and tracers variables 17 17 USE dom_oce ! ocean space and time domain variables 18 USE sol_oce ! ocean space and time domain variables19 18 USE in_out_manager ! I/O manager 20 19 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 21 20 USE lib_mpp ! distributed memory computing 22 USE dynspg_oce ! pressure gradient schemes23 21 USE c1d ! 1D vertical configuration 24 22 … … 43 41 !! ** Method : - Save the time step in numstp 44 42 !! - Print it each 50 time steps 45 !! - Print solver statistics in numsol 46 !! - Stop the run IF problem for the solver ( indec < 0 ) 43 !! - Stop the run IF problem ( indic < 0 ) 47 44 !! 48 45 !! ** Actions : 'time.step' file containing the last ocean time-step … … 50 47 !!---------------------------------------------------------------------- 51 48 INTEGER, INTENT( in ) :: kt ! ocean time-step index 52 INTEGER, INTENT( inout ) :: kindic ! indicator of solver convergence49 INTEGER, INTENT( inout ) :: kindic ! error indicator 53 50 !! 54 51 INTEGER :: ji, jj, jk ! dummy loop indices … … 143 140 IF( lk_c1d ) RETURN ! No log file in case of 1D vertical configuration 144 141 145 ! log file (solver or ssh statistics) 146 ! -------- 147 IF( lk_dynspg_flt ) THEN ! elliptic solver statistics (if required) 148 ! 149 IF(lwp) WRITE(numsol,9200) kt, niter, res, SQRT(epsr)/eps ! Solver 150 ! 151 IF( kindic < 0 .AND. zsmin > 0.e0 .AND. zumax <= 20.e0 ) THEN ! create a abort file if problem found 152 IF(lwp) THEN 153 WRITE(numout,*) ' stpctl: the elliptic solver DO not converge or explode' 154 WRITE(numout,*) ' ====== ' 155 WRITE(numout,9200) kt, niter, res, sqrt(epsr)/eps 156 WRITE(numout,*) 157 WRITE(numout,*) ' stpctl: output of last fields' 158 WRITE(numout,*) ' ====== ' 159 ENDIF 160 ENDIF 161 ! 162 ELSE !* ssh statistics (and others...) 163 IF( kt == nit000 .AND. lwp ) THEN ! open ssh statistics file (put in solver.stat file) 164 CALL ctl_opn( numsol, 'solver.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 165 ENDIF 166 ! 167 zssh2 = SUM( sshn(:,:) * sshn(:,:) * tmask_i(:,:) ) 168 IF( lk_mpp ) CALL mpp_sum( zssh2 ) ! sum over the global domain 169 ! 170 IF(lwp) WRITE(numsol,9300) kt, zssh2, zumax, zsmin ! ssh statistics 171 ! 142 ! log file (ssh statistics) 143 ! -------- !* ssh statistics (and others...) 144 IF( kt == nit000 .AND. lwp ) THEN ! open ssh statistics file (put in solver.stat file) 145 CALL ctl_opn( numsol, 'solver.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 172 146 ENDIF 147 ! 148 zssh2 = SUM( sshn(:,:) * sshn(:,:) * tmask_i(:,:) ) 149 IF( lk_mpp ) CALL mpp_sum( zssh2 ) ! sum over the global domain 150 ! 151 IF(lwp) WRITE(numsol,9300) kt, zssh2, zumax, zsmin ! ssh statistics 152 ! 173 153 174 154 9200 FORMAT('it:', i8, ' iter:', i4, ' r: ',e16.10, ' b: ',e16.10 )
Note: See TracChangeset
for help on using the changeset viewer.