- Timestamp:
- 2015-07-21T10:55:28+02:00 (9 years ago)
- Location:
- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM
- Files:
-
- 12 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90
r4162 r5620 72 72 !!---------------------------------------------------------------------- 73 73 INTEGER :: jc ! dummy loop indices 74 INTEGER :: isrow ! local index 74 75 !!---------------------------------------------------------------------- 75 76 … … 91 92 CASE ( 1 ) ! ORCA_R1 configuration 92 93 ! ! ======================= 94 ! This dirty section will be suppressed by simplification process: 95 ! all this will come back in input files 96 ! Currently these hard-wired indices relate to configuration with 97 ! extend grid (jpjglo=332) 98 isrow = 332 - jpjglo 99 ! 93 100 ncsnr(1) = 1 ; ncstt(1) = 0 ! Caspian Sea 94 ncsi1(1) = 332 ; ncsj1(1) = 2 0395 ncsi2(1) = 344 ; ncsj2(1) = 2 35101 ncsi1(1) = 332 ; ncsj1(1) = 243 - isrow 102 ncsi2(1) = 344 ; ncsj2(1) = 275 - isrow 96 103 ncsir(1,1) = 1 ; ncsjr(1,1) = 1 97 104 ! -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90
r5038 r5620 73 73 !!---------------------------------------------------------------------- 74 74 ! 75 ! max number of seconds between each restart 76 IF( REAL( nitend - nit000 + 1 ) * rdt > REAL( HUGE( nsec1jan000 ) ) ) THEN 77 CALL ctl_stop( 'The number of seconds between each restart exceeds the integer 4 max value: 2^31-1. ', & 78 & 'You must do a restart at higher frequency (or remove this stop and recompile the code in I8)' ) 79 ENDIF 75 80 ! all calendar staff is based on the fact that MOD( rday, rdttra(1) ) == 0 76 81 IF( MOD( rday , rdttra(1) ) /= 0. ) CALL ctl_stop( 'the time step must devide the number of second of in a day' ) … … 238 243 nday_year = 1 239 244 nsec_year = ndt05 240 IF( nsec1jan000 >= 2 * (2**30 - nsecd * nyear_len(1) / 2 ) ) THEN ! test integer 4 max value241 CALL ctl_stop( 'The number of seconds between Jan. 1st 00h of nit000 year and Jan. 1st 00h ', &242 & 'of the current year is exceeding the INTEGER 4 max VALUE: 2^31-1 -> 68.09 years in seconds', &243 & 'You must do a restart at higher frequency (or remove this STOP and recompile everything in I8)' )244 ENDIF245 245 nsec1jan000 = nsec1jan000 + nsecd * nyear_len(1) 246 246 IF( nleapy == 1 ) CALL day_mth -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r5038 r5620 162 162 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: gphit, gphiu !: latitude of t-, u-, v- and f-points (degre) 163 163 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: gphiv, gphif !: 164 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1t, e2t !: horizontal scale factorsat t-point (m)165 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1u, e2u !: horizontal scale factorsat u-point (m)166 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1v, e2v !: horizontal scale factorsat v-point (m)167 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1f, e2f !: horizontal scale factorsat f-point (m)164 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1t, e2t, r1_e1t, r1_e2t !: horizontal scale factors and inverse at t-point (m) 165 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1u, e2u, r1_e1u, r1_e2u !: horizontal scale factors and inverse at u-point (m) 166 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1v, e2v, r1_e1v, r1_e2v !: horizontal scale factors and inverse at v-point (m) 167 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: e1f, e2f, r1_e1f, r1_e2f !: horizontal scale factors and inverse at f-point (m) 168 168 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: e1e2t !: surface at t-point (m2) 169 169 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ff !: coriolis factor (2.*omega*sin(yphi) ) (s-1) … … 262 262 263 263 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: tmask, umask, vmask, fmask !: land/ocean mask at T-, U-, V- and F-pts 264 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: wmask, wumask, wvmask !: land/ocean mask at WT-, WU- and WV-pts 264 265 265 266 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: tpol, fpol !: north fold mask (jperio= 3 or 4) … … 332 333 INTEGER FUNCTION dom_oce_alloc() 333 334 !!---------------------------------------------------------------------- 334 INTEGER, DIMENSION(1 1) :: ierr335 INTEGER, DIMENSION(12) :: ierr 335 336 !!---------------------------------------------------------------------- 336 337 ierr(:) = 0 … … 345 346 & tpol(jpiglo) , fpol(jpiglo) , STAT=ierr(2) ) 346 347 ! 347 ALLOCATE( glamt(jpi,jpj) , gphit(jpi,jpj) , e1t(jpi,jpj) , e2t(jpi,jpj) , & 348 & glamu(jpi,jpj) , gphiu(jpi,jpj) , e1u(jpi,jpj) , e2u(jpi,jpj) , & 349 & glamv(jpi,jpj) , gphiv(jpi,jpj) , e1v(jpi,jpj) , e2v(jpi,jpj) , e1e2t(jpi,jpj) , & 350 & glamf(jpi,jpj) , gphif(jpi,jpj) , e1f(jpi,jpj) , e2f(jpi,jpj) , ff (jpi,jpj) , STAT=ierr(3) ) 348 ALLOCATE( glamt(jpi,jpj) , gphit(jpi,jpj) , e1t(jpi,jpj) , e2t(jpi,jpj) , r1_e1t(jpi,jpj) , r1_e2t(jpi,jpj) , & 349 & glamu(jpi,jpj) , gphiu(jpi,jpj) , e1u(jpi,jpj) , e2u(jpi,jpj) , r1_e1u(jpi,jpj) , r1_e2u(jpi,jpj) , & 350 & glamv(jpi,jpj) , gphiv(jpi,jpj) , e1v(jpi,jpj) , e2v(jpi,jpj) , r1_e1v(jpi,jpj) , r1_e2v(jpi,jpj) , & 351 & glamf(jpi,jpj) , gphif(jpi,jpj) , e1f(jpi,jpj) , e2f(jpi,jpj) , r1_e1f(jpi,jpj) , r1_e2f(jpi,jpj) , & 352 & e1e2t(jpi,jpj) , ff (jpi,jpj) , STAT=ierr(3) ) 351 353 ! 352 354 ALLOCATE( gdep3w_0(jpi,jpj,jpk) , e3v_0(jpi,jpj,jpk) , e3f_0 (jpi,jpj,jpk) , & … … 400 402 & vmask(jpi,jpj,jpk) , fmask(jpi,jpj,jpk), STAT=ierr(11) ) 401 403 404 ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(12) ) 405 402 406 #if defined key_noslip_accurate 403 ALLOCATE( npcoa(4,jpk), nicoa(2*(jpi+jpj),4,jpk), njcoa(2*(jpi+jpj),4,jpk), STAT=ierr(1 1) )407 ALLOCATE( npcoa(4,jpk), nicoa(2*(jpi+jpj),4,jpk), njcoa(2*(jpi+jpj),4,jpk), STAT=ierr(12) ) 404 408 #endif 405 409 ! -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r5038 r5620 135 135 !!---------------------------------------------------------------------- 136 136 USE ioipsl 137 NAMELIST/namrun/ nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl, & 137 NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list, & 138 & nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl, & 138 139 & nn_it000, nn_itend , nn_date0 , nn_leapy , nn_istate , nn_stock , & 139 & nn_write, ln_dimgnnn, ln_mskland , ln_c lobber, nn_chunksz, nn_euler140 & nn_write, ln_dimgnnn, ln_mskland , ln_cfmeta , ln_clobber, nn_chunksz, nn_euler 140 141 NAMELIST/namdom/ nn_bathy, rn_bathy , rn_e3zps_min, rn_e3zps_rat, nn_msh, rn_hmin, & 141 142 & nn_acc , rn_atfp , rn_rdt , rn_rdtmin , & … … 169 170 WRITE(numout,*) ' experiment name for output cn_exp = ', cn_exp 170 171 WRITE(numout,*) ' file prefix restart input cn_ocerst_in= ', cn_ocerst_in 172 WRITE(numout,*) ' restart input directory cn_ocerst_indir= ', cn_ocerst_indir 171 173 WRITE(numout,*) ' file prefix restart output cn_ocerst_out= ', cn_ocerst_out 174 WRITE(numout,*) ' restart output directory cn_ocerst_outdir= ', cn_ocerst_outdir 172 175 WRITE(numout,*) ' restart logical ln_rstart = ', ln_rstart 173 176 WRITE(numout,*) ' start with forward time step nn_euler = ', nn_euler … … 178 181 WRITE(numout,*) ' leap year calendar (0/1) nn_leapy = ', nn_leapy 179 182 WRITE(numout,*) ' initial state output nn_istate = ', nn_istate 180 WRITE(numout,*) ' frequency of restart file nn_stock = ', nn_stock 183 IF( ln_rst_list ) THEN 184 WRITE(numout,*) ' list of restart dump times nn_stocklist =', nn_stocklist 185 ELSE 186 WRITE(numout,*) ' frequency of restart file nn_stock = ', nn_stock 187 ENDIF 181 188 WRITE(numout,*) ' frequency of output file nn_write = ', nn_write 182 189 WRITE(numout,*) ' multi file dimgout ln_dimgnnn = ', ln_dimgnnn 183 190 WRITE(numout,*) ' mask land points ln_mskland = ', ln_mskland 191 WRITE(numout,*) ' additional CF standard metadata ln_cfmeta = ', ln_cfmeta 184 192 WRITE(numout,*) ' overwrite an existing file ln_clobber = ', ln_clobber 185 193 WRITE(numout,*) ' NetCDF chunksize (bytes) nn_chunksz = ', nn_chunksz … … 195 203 ninist = nn_istate 196 204 nstock = nn_stock 205 nstocklist = nn_stocklist 197 206 nwrite = nn_write 198 207 neuler = nn_euler 199 IF ( neuler == 1 .AND. .NOT. ln_rstart ) THEN208 IF ( neuler == 1 .AND. .NOT. ln_rstart ) THEN 200 209 WRITE(ctmp1,*) 'ln_rstart =.FALSE., nn_euler is forced to 0 ' 201 210 CALL ctl_warn( ctmp1 ) -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90
r5038 r5620 105 105 REAL(wp) :: zlam1, zcos_alpha, zim1 , zjm1 , ze1, ze1deg 106 106 REAL(wp) :: zphi1, zsin_alpha, zim05, zjm05 107 INTEGER :: isrow ! index for ORCA1 starting row 108 107 109 !!---------------------------------------------------------------------- 108 110 ! … … 159 161 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration 160 162 ! ! ===================== 161 162 ii0 = 281 ; ii1 = 282 ! Gibraltar Strait (e2u = 20 km) 163 ij0 = 200 ; ij1 = 200 ; e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 163 ! This dirty section will be suppressed by simplification process: all this will come back in input files 164 ! Currently these hard-wired indices relate to configuration with 165 ! extend grid (jpjglo=332) 166 ! which had a grid-size of 362x292. 167 ! 168 isrow = 332 - jpjglo 169 ! 170 ii0 = 282 ; ii1 = 283 ! Gibraltar Strait (e2u = 20 km) 171 ij0 = 241 - isrow ; ij1 = 241 - isrow ; e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 164 172 IF(lwp) WRITE(numout,*) 165 173 IF(lwp) WRITE(numout,*) ' orca_r1: Gibraltar : e2u reduced to 20 km' 166 174 167 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u = 10 km)168 ij0 = 2 08 ; ij1 = 208; e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e3175 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u = 10 km) 176 ij0 = 248 - isrow ; ij1 = 248 - isrow ; e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 10.e3 169 177 IF(lwp) WRITE(numout,*) 170 178 IF(lwp) WRITE(numout,*) ' orca_r1: Bhosporus : e2u reduced to 10 km' 171 179 172 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v = 13 km)173 ij0 = 1 24 ; ij1 = 125; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3180 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v = 13 km) 181 ij0 = 164 - isrow ; ij1 = 165 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3 174 182 IF(lwp) WRITE(numout,*) 175 183 IF(lwp) WRITE(numout,*) ' orca_r1: Lombok : e1v reduced to 10 km' 176 184 177 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v = 8 km) [closed from bathy_11 on]178 ij0 = 1 24 ; ij1 = 125; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 8.e3185 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v = 8 km) [closed from bathy_11 on] 186 ij0 = 164 - isrow ; ij1 = 165 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 8.e3 179 187 IF(lwp) WRITE(numout,*) 180 188 IF(lwp) WRITE(numout,*) ' orca_r1: Sumba : e1v reduced to 8 km' 181 189 182 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v = 13 km)183 ij0 = 1 24 ; ij1 = 125; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3190 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v = 13 km) 191 ij0 = 164 - isrow ; ij1 = 165 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3 184 192 IF(lwp) WRITE(numout,*) 185 193 IF(lwp) WRITE(numout,*) ' orca_r1: Ombai : e1v reduced to 13 km' 186 194 187 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v = 20 km)188 ij0 = 1 24 ; ij1 = 125; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3195 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v = 20 km) 196 ij0 = 164 - isrow ; ij1 = 145 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 189 197 IF(lwp) WRITE(numout,*) 190 198 IF(lwp) WRITE(numout,*) ' orca_r1: Timor Passage : e1v reduced to 20 km' 191 199 192 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v = 30 km)193 ij0 = 1 41 ; ij1 = 142; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 30.e3200 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v = 30 km) 201 ij0 = 181 - isrow ; ij1 = 182 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 30.e3 194 202 IF(lwp) WRITE(numout,*) 195 203 IF(lwp) WRITE(numout,*) ' orca_r1: W Halmahera : e1v reduced to 30 km' 196 204 197 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v = 50 km)198 ij0 = 1 41 ; ij1 = 142; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 50.e3205 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v = 50 km) 206 ij0 = 181 - isrow ; ij1 = 182 - isrow ; e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 50.e3 199 207 IF(lwp) WRITE(numout,*) 200 208 IF(lwp) WRITE(numout,*) ' orca_r1: E Halmahera : e1v reduced to 50 km' 201 202 !203 204 !205 !206 209 ! 207 210 ! … … 471 474 re2u_e1u(:,:) = e2u(:,:) / e1u(:,:) 472 475 re1v_e2v(:,:) = e1v(:,:) / e2v(:,:) 476 r1_e1t (:,:) = 1._wp / e1t(:,:) 477 r1_e1u (:,:) = 1._wp / e1u(:,:) 478 r1_e1v (:,:) = 1._wp / e1v(:,:) 479 r1_e1f (:,:) = 1._wp / e1f(:,:) 480 r1_e2t (:,:) = 1._wp / e2t(:,:) 481 r1_e2u (:,:) = 1._wp / e2u(:,:) 482 r1_e2v (:,:) = 1._wp / e2v(:,:) 483 r1_e2f (:,:) = 1._wp / e2f(:,:) 473 484 474 485 ! Control printing : Grid informations (if not restart) … … 616 627 CALL iom_open( 'coordinates', inum ) 617 628 618 CALL iom_get( inum, jpdom_data, 'glamt', glamt )619 CALL iom_get( inum, jpdom_data, 'glamu', glamu )620 CALL iom_get( inum, jpdom_data, 'glamv', glamv )621 CALL iom_get( inum, jpdom_data, 'glamf', glamf )622 623 CALL iom_get( inum, jpdom_data, 'gphit', gphit )624 CALL iom_get( inum, jpdom_data, 'gphiu', gphiu )625 CALL iom_get( inum, jpdom_data, 'gphiv', gphiv )626 CALL iom_get( inum, jpdom_data, 'gphif', gphif )627 628 CALL iom_get( inum, jpdom_data, 'e1t', e1t )629 CALL iom_get( inum, jpdom_data, 'e1u', e1u )630 CALL iom_get( inum, jpdom_data, 'e1v', e1v )631 CALL iom_get( inum, jpdom_data, 'e1f', e1f )632 633 CALL iom_get( inum, jpdom_data, 'e2t', e2t )634 CALL iom_get( inum, jpdom_data, 'e2u', e2u )635 CALL iom_get( inum, jpdom_data, 'e2v', e2v )636 CALL iom_get( inum, jpdom_data, 'e2f', e2f )629 CALL iom_get( inum, jpdom_data, 'glamt', glamt, lrowattr=ln_use_jattr ) 630 CALL iom_get( inum, jpdom_data, 'glamu', glamu, lrowattr=ln_use_jattr ) 631 CALL iom_get( inum, jpdom_data, 'glamv', glamv, lrowattr=ln_use_jattr ) 632 CALL iom_get( inum, jpdom_data, 'glamf', glamf, lrowattr=ln_use_jattr ) 633 634 CALL iom_get( inum, jpdom_data, 'gphit', gphit, lrowattr=ln_use_jattr ) 635 CALL iom_get( inum, jpdom_data, 'gphiu', gphiu, lrowattr=ln_use_jattr ) 636 CALL iom_get( inum, jpdom_data, 'gphiv', gphiv, lrowattr=ln_use_jattr ) 637 CALL iom_get( inum, jpdom_data, 'gphif', gphif, lrowattr=ln_use_jattr ) 638 639 CALL iom_get( inum, jpdom_data, 'e1t', e1t, lrowattr=ln_use_jattr ) 640 CALL iom_get( inum, jpdom_data, 'e1u', e1u, lrowattr=ln_use_jattr ) 641 CALL iom_get( inum, jpdom_data, 'e1v', e1v, lrowattr=ln_use_jattr ) 642 CALL iom_get( inum, jpdom_data, 'e1f', e1f, lrowattr=ln_use_jattr ) 643 644 CALL iom_get( inum, jpdom_data, 'e2t', e2t, lrowattr=ln_use_jattr ) 645 CALL iom_get( inum, jpdom_data, 'e2u', e2u, lrowattr=ln_use_jattr ) 646 CALL iom_get( inum, jpdom_data, 'e2v', e2v, lrowattr=ln_use_jattr ) 647 CALL iom_get( inum, jpdom_data, 'e2f', e2f, lrowattr=ln_use_jattr ) 637 648 638 649 CALL iom_close( inum ) -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90
r5038 r5620 134 134 INTEGER :: ijf, ijl, ij0, ij1 ! - - 135 135 INTEGER :: ios 136 INTEGER :: isrow ! index for ORCA1 starting row 136 137 INTEGER , POINTER, DIMENSION(:,:) :: imsk 137 138 REAL(wp), POINTER, DIMENSION(:,:) :: zwf … … 281 282 CALL lbc_lnk( fmask_i, 'F', 1._wp ) 282 283 284 ! 3. Ocean/land mask at wu-, wv- and w points 285 !---------------------------------------------- 286 wmask (:,:,1) = tmask(:,:,1) ! ???????? 287 wumask(:,:,1) = umask(:,:,1) ! ???????? 288 wvmask(:,:,1) = vmask(:,:,1) ! ???????? 289 DO jk=2,jpk 290 wmask (:,:,jk)=tmask(:,:,jk) * tmask(:,:,jk-1) 291 wumask(:,:,jk)=umask(:,:,jk) * umask(:,:,jk-1) 292 wvmask(:,:,jk)=vmask(:,:,jk) * vmask(:,:,jk-1) 293 END DO 283 294 284 295 ! 4. ocean/land mask for the elliptic equation … … 391 402 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration 392 403 ! ! Increased lateral friction near of some straits 404 ! This dirty section will be suppressed by simplification process: 405 ! all this will come back in input files 406 ! Currently these hard-wired indices relate to configuration with 407 ! extend grid (jpjglo=332) 408 ! 409 isrow = 332 - jpjglo 410 ! 393 411 IF(lwp) WRITE(numout,*) 394 412 IF(lwp) WRITE(numout,*) ' orca_r1: increase friction near the following straits : ' 395 413 IF(lwp) WRITE(numout,*) ' Gibraltar ' 396 ii0 = 28 3 ; ii1 = 284! Gibraltar Strait397 ij0 = 2 00 ; ij1 = 200 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) =2._wp414 ii0 = 282 ; ii1 = 283 ! Gibraltar Strait 415 ij0 = 241 - isrow ; ij1 = 241 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp 398 416 399 417 IF(lwp) WRITE(numout,*) ' Bhosporus ' 400 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait401 ij0 = 2 08 ; ij1 = 208 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) =2._wp418 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait 419 ij0 = 248 - isrow ; ij1 = 248 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp 402 420 403 421 IF(lwp) WRITE(numout,*) ' Makassar (Top) ' 404 ii0 = 48 ; ii1 = 48 ! Makassar Strait (Top)405 ij0 = 1 49 ; ij1 = 150 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) =3._wp422 ii0 = 48 ; ii1 = 48 ! Makassar Strait (Top) 423 ij0 = 189 - isrow ; ij1 = 190 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp 406 424 407 425 IF(lwp) WRITE(numout,*) ' Lombok ' 408 ii0 = 44 ; ii1 = 44 ! Lombok Strait409 ij0 = 1 24 ; ij1 = 125 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) =2._wp426 ii0 = 44 ; ii1 = 44 ! Lombok Strait 427 ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp 410 428 411 429 IF(lwp) WRITE(numout,*) ' Ombai ' 412 ii0 = 53 ; ii1 = 53 ! Ombai Strait413 ij0 = 1 24 ; ij1 = 125 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1),1:jpk ) = 2._wp430 ii0 = 53 ; ii1 = 53 ! Ombai Strait 431 ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp 414 432 415 433 IF(lwp) WRITE(numout,*) ' Timor Passage ' 416 ii0 = 56 ; ii1 = 56 ! Timor Passage417 ij0 = 1 24 ; ij1 = 125 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1),1:jpk ) = 2._wp434 ii0 = 56 ; ii1 = 56 ! Timor Passage 435 ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp 418 436 419 437 IF(lwp) WRITE(numout,*) ' West Halmahera ' 420 ii0 = 58 ; ii1 = 58 ! West Halmahera Strait421 ij0 = 1 41 ; ij1 = 142 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1),1:jpk ) = 3._wp438 ii0 = 58 ; ii1 = 58 ! West Halmahera Strait 439 ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp 422 440 423 441 IF(lwp) WRITE(numout,*) ' East Halmahera ' 424 ii0 = 55 ; ii1 = 55 ! East Halmahera Strait425 ij0 = 1 41 ; ij1 = 142 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1),1:jpk ) = 3._wp442 ii0 = 55 ; ii1 = 55 ! East Halmahera Strait 443 ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp 426 444 ! 427 445 ENDIF -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r5038 r5620 8 8 !! 3.3 ! 2011-10 (M. Leclair) totally rewrote domvvl: 9 9 !! vvl option includes z_star and z_tilde coordinates 10 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 10 11 !!---------------------------------------------------------------------- 11 12 !! 'key_vvl' variable volume … … 125 126 INTEGER :: ji,jj,jk 126 127 INTEGER :: ii0, ii1, ij0, ij1 128 REAL(wp):: zcoef 127 129 !!---------------------------------------------------------------------- 128 130 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_init') … … 164 166 ! t- and w- points depth 165 167 ! ---------------------- 168 ! set the isf depth as it is in the initial step 166 169 fsdept_n(:,:,1) = 0.5_wp * fse3w_n(:,:,1) 167 170 fsdepw_n(:,:,1) = 0.0_wp … … 169 172 fsdept_b(:,:,1) = 0.5_wp * fse3w_b(:,:,1) 170 173 fsdepw_b(:,:,1) = 0.0_wp 171 DO jj = 1,jpj 172 DO ji = 1,jpi 173 DO jk = 2,mikt(ji,jj)-1 174 fsdept_n(ji,jj,jk) = gdept_0(ji,jj,jk) 175 fsdepw_n(ji,jj,jk) = gdepw_0(ji,jj,jk) 176 fsde3w_n(ji,jj,jk) = gdept_0(ji,jj,jk) - sshn(ji,jj) 177 fsdept_b(ji,jj,jk) = gdept_0(ji,jj,jk) 178 fsdepw_b(ji,jj,jk) = gdepw_0(ji,jj,jk) 179 END DO 180 IF (mikt(ji,jj) .GT. 1) THEN 181 jk = mikt(ji,jj) 182 fsdept_n(ji,jj,jk) = gdepw_0(ji,jj,jk) + 0.5_wp * fse3w_n(ji,jj,jk) 183 fsdepw_n(ji,jj,jk) = gdepw_0(ji,jj,jk) 184 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk ) - sshn (ji,jj) 185 fsdept_b(ji,jj,jk) = gdepw_0(ji,jj,jk) + 0.5_wp * fse3w_b(ji,jj,jk) 186 fsdepw_b(ji,jj,jk) = gdepw_0(ji,jj,jk) 187 END IF 188 DO jk = mikt(ji,jj)+1, jpk 189 fsdept_n(ji,jj,jk) = fsdept_n(ji,jj,jk-1) + fse3w_n(ji,jj,jk) 174 175 DO jk = 2, jpk 176 DO jj = 1,jpj 177 DO ji = 1,jpi 178 ! zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 179 ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 180 ! 0.5 where jk = mikt 181 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 190 182 fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1) 191 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk ) - sshn (ji,jj) 192 fsdept_b(ji,jj,jk) = fsdept_b(ji,jj,jk-1) + fse3w_b(ji,jj,jk) 183 fsdept_n(ji,jj,jk) = zcoef * ( fsdepw_n(ji,jj,jk ) + 0.5 * fse3w_n(ji,jj,jk)) & 184 & + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) + fse3w_n(ji,jj,jk)) 185 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - sshn(ji,jj) 193 186 fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1) 187 fsdept_b(ji,jj,jk) = zcoef * ( fsdepw_b(ji,jj,jk ) + 0.5 * fse3w_b(ji,jj,jk)) & 188 & + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) + fse3w_b(ji,jj,jk)) 194 189 END DO 195 190 END DO … … 588 583 INTEGER, INTENT( in ) :: kt ! time step 589 584 !! * Local declarations 590 REAL(wp), POINTER, DIMENSION(:,:,:) :: z_e3t_def591 585 INTEGER :: ji,jj,jk ! dummy loop indices 586 REAL(wp) :: zcoef 592 587 !!---------------------------------------------------------------------- 593 588 594 589 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_sf_swp') 595 !596 CALL wrk_alloc( jpi, jpj, jpk, z_e3t_def )597 590 ! 598 591 IF( kt == nit000 ) THEN … … 638 631 ! t- and w- points depth 639 632 ! ---------------------- 633 ! set the isf depth as it is in the initial step 640 634 fsdept_n(:,:,1) = 0.5_wp * fse3w_n(:,:,1) 641 635 fsdepw_n(:,:,1) = 0.0_wp 642 636 fsde3w_n(:,:,1) = fsdept_n(:,:,1) - sshn(:,:) 643 DO jj = 1,jpj 644 DO ji = 1,jpi 645 DO jk = 2,mikt(ji,jj)-1 646 fsdept_n(ji,jj,jk) = gdept_0(ji,jj,jk) 647 fsdepw_n(ji,jj,jk) = gdepw_0(ji,jj,jk) 648 fsde3w_n(ji,jj,jk) = gdept_0(ji,jj,jk) - sshn(ji,jj) 649 END DO 650 IF (mikt(ji,jj) .GT. 1) THEN 651 jk = mikt(ji,jj) 652 fsdept_n(ji,jj,jk) = gdepw_0(ji,jj,jk) + 0.5_wp * fse3w_n(ji,jj,jk) 653 fsdepw_n(ji,jj,jk) = gdepw_0(ji,jj,jk) 654 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk ) - sshn (ji,jj) 655 END IF 656 DO jk = mikt(ji,jj)+1, jpk 657 fsdept_n(ji,jj,jk) = fsdept_n(ji,jj,jk-1) + fse3w_n(ji,jj,jk) 637 638 DO jk = 2, jpk 639 DO jj = 1,jpj 640 DO ji = 1,jpi 641 ! zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 642 ! 1 for jk = mikt 643 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 658 644 fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1) 659 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk ) - sshn (ji,jj) 645 fsdept_n(ji,jj,jk) = zcoef * ( fsdepw_n(ji,jj,jk ) + 0.5 * fse3w_n(ji,jj,jk)) & 646 & + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) + fse3w_n(ji,jj,jk)) 647 fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - sshn(ji,jj) 660 648 END DO 661 649 END DO 662 650 END DO 651 663 652 ! Local depth and Inverse of the local depth of the water column at u- and v- points 664 653 ! ---------------------------------------------------------------------------------- … … 679 668 ! Write outputs 680 669 ! ============= 681 z_e3t_def(:,:,:) = ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 682 CALL iom_put( "cellthc" , fse3t_n (:,:,:) ) 670 CALL iom_put( "e3t" , fse3t_n (:,:,:) ) 671 CALL iom_put( "e3u" , fse3u_n (:,:,:) ) 672 CALL iom_put( "e3v" , fse3v_n (:,:,:) ) 673 CALL iom_put( "e3w" , fse3w_n (:,:,:) ) 683 674 CALL iom_put( "tpt_dep" , fsde3w_n (:,:,:) ) 684 CALL iom_put( "e3tdef" , z_e3t_def(:,:,:) ) 675 IF( iom_use("e3tdef") ) & 676 CALL iom_put( "e3tdef" , ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 685 677 686 678 ! write restart file 687 679 ! ================== 688 680 IF( lrst_oce ) CALL dom_vvl_rst( kt, 'WRITE' ) 689 !690 CALL wrk_dealloc( jpi, jpj, jpk, z_e3t_def )691 681 ! 692 682 IF( nn_timing == 1 ) CALL timing_stop('dom_vvl_sf_swp') … … 1049 1039 INTEGER :: ji, jj, jk ! dummy loop indices 1050 1040 INTEGER :: ij0, ij1, ii0, ii1 ! dummy loop indices 1041 INTEGER :: isrow ! index for ORCA1 starting row 1051 1042 !! acc 1052 1043 !! Hmm with the time splitting these "fixes" seem to do more harm than good. Temporarily disabled for … … 1132 1123 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration 1133 1124 ! ! ===================== 1134 ! 1135 ii0 = 281 ; ii1 = 282 ! Gibraltar Strait (e2u was modified) 1136 ij0 = 200 ; ij1 = 200 1125 ! This dirty section will be suppressed by simplification process: 1126 ! all this will come back in input files 1127 ! Currently these hard-wired indices relate to configuration with 1128 ! extend grid (jpjglo=332) 1129 ! which had a grid-size of 362x292. 1130 isrow = 332 - jpjglo 1131 ! 1132 ii0 = 282 ; ii1 = 283 ! Gibraltar Strait (e2u was modified) 1133 ij0 = 241 - isrow ; ij1 = 241 - isrow 1137 1134 DO jk = 1, jpkm1 1138 1135 DO jj = mj0(ij0), mj1(ij1) … … 1154 1151 END DO 1155 1152 ! 1156 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u was modified)1157 ij0 = 2 08 ; ij1 = 2081153 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u was modified) 1154 ij0 = 248 - isrow ; ij1 = 248 - isrow 1158 1155 DO jk = 1, jpkm1 1159 1156 DO jj = mj0(ij0), mj1(ij1) … … 1175 1172 END DO 1176 1173 ! 1177 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v was modified)1178 ij0 = 1 24 ; ij1 = 1251174 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v was modified) 1175 ij0 = 164 - isrow ; ij1 = 165 - isrow 1179 1176 DO jk = 1, jpkm1 1180 1177 DO jj = mj0(ij0), mj1(ij1) … … 1191 1188 END DO 1192 1189 ! 1193 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v was modified) [closed from bathy_11 on]1194 ij0 = 1 24 ; ij1 = 1251190 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v was modified) [closed from bathy_11 on] 1191 ij0 = 164 - isrow ; ij1 = 165 - isrow 1195 1192 DO jk = 1, jpkm1 1196 1193 DO jj = mj0(ij0), mj1(ij1) … … 1207 1204 END DO 1208 1205 ! 1209 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v was modified)1210 ij0 = 1 24 ; ij1 = 1251206 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v was modified) 1207 ij0 = 164 - isrow ; ij1 = 165 - isrow 1211 1208 DO jk = 1, jpkm1 1212 1209 DO jj = mj0(ij0), mj1(ij1) … … 1223 1220 END DO 1224 1221 ! 1225 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v was modified)1226 ij0 = 1 24 ; ij1 = 1251222 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v was modified) 1223 ij0 = 164 - isrow ; ij1 = 165 - isrow 1227 1224 DO jk = 1, jpkm1 1228 1225 DO jj = mj0(ij0), mj1(ij1) … … 1239 1236 END DO 1240 1237 ! 1241 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v was modified)1242 ij0 = 1 41 ; ij1 = 1421238 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v was modified) 1239 ij0 = 181 - isrow ; ij1 = 182 - isrow 1243 1240 DO jk = 1, jpkm1 1244 1241 DO jj = mj0(ij0), mj1(ij1) … … 1255 1252 END DO 1256 1253 ! 1257 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v was modified)1258 ij0 = 1 41 ; ij1 = 1421254 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v was modified) 1255 ij0 = 181 - isrow ; ij1 = 182 - isrow 1259 1256 DO jk = 1, jpkm1 1260 1257 DO jj = mj0(ij0), mj1(ij1) -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90
r5038 r5620 215 215 CALL iom_rstput( 0, 0, inum4, 'gdept_1d' , gdept_1d ) ! ! stretched system 216 216 CALL iom_rstput( 0, 0, inum4, 'gdepw_1d' , gdepw_1d ) 217 CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0, ktype = jp_r4 ) 218 CALL iom_rstput( 0, 0, inum4, 'gdepw_0', gdepw_0, ktype = jp_r4 ) 217 219 ENDIF 218 220 -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90
r5038 r5620 17 17 !! 3.4 ! 2012-08 (J. Siddorn) added Siddorn and Furner stretching function 18 18 !! 3.4 ! 2012-12 (R. Bourdalle-Badie and G. Reffray) modify C1D case 19 !! 3.6 ! 2014-11 (P. Mathiot and C. Harris) add ice shelf capabilitye 19 20 !!---------------------------------------------------------------------- 20 21 … … 35 36 USE oce ! ocean variables 36 37 USE dom_oce ! ocean domain 37 USE sbc_oce ! surface variable (isf)38 38 USE closea ! closed seas 39 39 USE c1d ! 1D vertical configuration … … 298 298 ENDIF 299 299 300 IF ( ln_isfcav ) THEN 300 301 ! need to be like this to compute the pressure gradient with ISF. If not, level beneath the ISF are not aligned (sum(e3t) /= depth) 301 302 ! define e3t_0 and e3w_0 as the differences between gdept and gdepw respectively 302 DO jk = 1, jpkm1 303 e3t_1d(jk) = gdepw_1d(jk+1)-gdepw_1d(jk) 304 END DO 305 e3t_1d(jpk) = e3t_1d(jpk-1) ! we don't care because this level is masked in NEMO 306 307 DO jk = 2, jpk 308 e3w_1d(jk) = gdept_1d(jk) - gdept_1d(jk-1) 309 END DO 310 e3w_1d(1 ) = 2._wp * (gdept_1d(1) - gdepw_1d(1)) 303 DO jk = 1, jpkm1 304 e3t_1d(jk) = gdepw_1d(jk+1)-gdepw_1d(jk) 305 END DO 306 e3t_1d(jpk) = e3t_1d(jpk-1) ! we don't care because this level is masked in NEMO 307 308 DO jk = 2, jpk 309 e3w_1d(jk) = gdept_1d(jk) - gdept_1d(jk-1) 310 END DO 311 e3w_1d(1 ) = 2._wp * (gdept_1d(1) - gdepw_1d(1)) 312 END IF 311 313 312 314 !!gm BUG in s-coordinate this does not work! … … 365 367 INTEGER :: ji, jj, jl, jk ! dummy loop indices 366 368 INTEGER :: inum ! temporary logical unit 369 INTEGER :: ierror ! error flag 367 370 INTEGER :: ii_bump, ij_bump, ih ! bump center position 368 371 INTEGER :: ii0, ii1, ij0, ij1, ik ! local indices 369 372 REAL(wp) :: r_bump , h_bump , h_oce ! bump characteristics 370 373 REAL(wp) :: zi, zj, zh, zhmin ! local scalars 371 INTEGER , POINTER, DIMENSION(:,:) :: idta ! global domain integer data372 REAL(wp), POINTER, DIMENSION(:,:) :: zdta ! global domain scalar data374 INTEGER , ALLOCATABLE, DIMENSION(:,:) :: idta ! global domain integer data 375 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zdta ! global domain scalar data 373 376 !!---------------------------------------------------------------------- 374 377 ! 375 378 IF( nn_timing == 1 ) CALL timing_start('zgr_bat') 376 !377 CALL wrk_alloc( jpidta, jpjdta, idta )378 CALL wrk_alloc( jpidta, jpjdta, zdta )379 379 ! 380 380 IF(lwp) WRITE(numout,*) … … 385 385 ! ! ================== ! 386 386 ! ! global domain level and meter bathymetry (idta,zdta) 387 ! 388 ALLOCATE( idta(jpidta,jpjdta), STAT=ierror ) 389 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'zgr_bat: unable to allocate idta array' ) 390 ALLOCATE( zdta(jpidta,jpjdta), STAT=ierror ) 391 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'zgr_bat: unable to allocate zdta array' ) 387 392 ! 388 393 IF( ntopo == 0 ) THEN ! flat basin … … 468 473 misfdep(:,:)=1 469 474 ! 470 ! (ISF) TODO build ice draft netcdf file for isomip and build the corresponding part of code 471 IF( cp_cfg == "isomip" ) THEN 472 ! 473 risfdep(:,:)=200.e0 474 misfdep(:,:)=1 475 ij0 = 1 ; ij1 = 40 476 DO jj = mj0(ij0), mj1(ij1) 477 risfdep(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp 478 END DO 479 WHERE( bathy(:,:) <= 0._wp ) risfdep(:,:) = 0._wp 480 ! 481 ELSEIF ( cp_cfg == "isomip2" ) THEN 482 ! 483 risfdep(:,:)=0.e0 484 misfdep(:,:)=1 485 ij0 = 1 ; ij1 = 40 486 DO jj = mj0(ij0), mj1(ij1) 487 risfdep(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp 488 END DO 489 WHERE( bathy(:,:) <= 0._wp ) risfdep(:,:) = 0._wp 490 END IF 475 DEALLOCATE( idta, zdta ) 491 476 ! 492 477 ! ! ================ ! … … 529 514 IF( ln_zps .OR. ln_sco ) THEN ! zps or sco : read meter bathymetry 530 515 CALL iom_open ( 'bathy_meter.nc', inum ) 531 CALL iom_get ( inum, jpdom_data, 'Bathymetry', bathy ) 516 IF ( ln_isfcav ) THEN 517 CALL iom_get ( inum, jpdom_data, 'Bathymetry_isf', bathy, lrowattr=.false. ) 518 ELSE 519 CALL iom_get ( inum, jpdom_data, 'Bathymetry' , bathy, lrowattr=ln_use_jattr ) 520 END IF 532 521 CALL iom_close( inum ) 533 ! 522 ! 534 523 risfdep(:,:)=0._wp 535 524 misfdep(:,:)=1 … … 579 568 IF ( .not. ln_sco ) THEN !== set a minimum depth ==! 580 569 ! patch to avoid case bathy = ice shelf draft and bathy between 0 and zhmin 581 WHERE (bathy == risfdep) 582 bathy = 0.0_wp ; risfdep = 0.0_wp 583 END WHERE 570 IF ( ln_isfcav ) THEN 571 WHERE (bathy == risfdep) 572 bathy = 0.0_wp ; risfdep = 0.0_wp 573 END WHERE 574 END IF 584 575 ! end patch 585 576 IF( rn_hmin < 0._wp ) THEN ; ik = - INT( rn_hmin ) ! from a nb of level … … 592 583 IF(lwp) write(numout,*) 'Minimum ocean depth: ', zhmin, ' minimum number of ocean levels : ', ik 593 584 ENDIF 594 !595 CALL wrk_dealloc( jpidta, jpjdta, idta )596 CALL wrk_dealloc( jpidta, jpjdta, zdta )597 585 ! 598 586 IF( nn_timing == 1 ) CALL timing_stop('zgr_bat') … … 959 947 !!---------------------------------------------------------------------- 960 948 !! 949 INTEGER :: ji, jj, jk ! dummy loop indices 950 INTEGER :: ik, it, ikb, ikt ! temporary integers 951 LOGICAL :: ll_print ! Allow control print for debugging 952 REAL(wp) :: ze3tp , ze3wp ! Last ocean level thickness at T- and W-points 953 REAL(wp) :: zdepwp, zdepth ! Ajusted ocean depth to avoid too small e3t 954 REAL(wp) :: zmax ! Maximum depth 955 REAL(wp) :: zdiff ! temporary scalar 956 REAL(wp) :: zrefdep ! temporary scalar 957 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprt 958 !!--------------------------------------------------------------------- 959 ! 960 IF( nn_timing == 1 ) CALL timing_start('zgr_zps') 961 ! 962 CALL wrk_alloc( jpi, jpj, jpk, zprt ) 963 ! 964 IF(lwp) WRITE(numout,*) 965 IF(lwp) WRITE(numout,*) ' zgr_zps : z-coordinate with partial steps' 966 IF(lwp) WRITE(numout,*) ' ~~~~~~~ ' 967 IF(lwp) WRITE(numout,*) ' mbathy is recomputed : bathy_level file is NOT used' 968 969 ll_print = .FALSE. ! Local variable for debugging 970 971 IF(lwp .AND. ll_print) THEN ! control print of the ocean depth 972 WRITE(numout,*) 973 WRITE(numout,*) 'dom_zgr_zps: bathy (in hundred of meters)' 974 CALL prihre( bathy, jpi, jpj, 1,jpi, 1, 1, jpj, 1, 1.e-2, numout ) 975 ENDIF 976 977 978 ! bathymetry in level (from bathy_meter) 979 ! =================== 980 zmax = gdepw_1d(jpk) + e3t_1d(jpk) ! maximum depth (i.e. the last ocean level thickness <= 2*e3t_1d(jpkm1) ) 981 bathy(:,:) = MIN( zmax , bathy(:,:) ) ! bounded value of bathy (min already set at the end of zgr_bat) 982 WHERE( bathy(:,:) == 0._wp ) ; mbathy(:,:) = 0 ! land : set mbathy to 0 983 ELSE WHERE ; mbathy(:,:) = jpkm1 ! ocean : initialize mbathy to the max ocean level 984 END WHERE 985 986 ! Compute mbathy for ocean points (i.e. the number of ocean levels) 987 ! find the number of ocean levels such that the last level thickness 988 ! is larger than the minimum of e3zps_min and e3zps_rat * e3t_1d (where 989 ! e3t_1d is the reference level thickness 990 DO jk = jpkm1, 1, -1 991 zdepth = gdepw_1d(jk) + MIN( e3zps_min, e3t_1d(jk)*e3zps_rat ) 992 WHERE( 0._wp < bathy(:,:) .AND. bathy(:,:) <= zdepth ) mbathy(:,:) = jk-1 993 END DO 994 995 IF ( ln_isfcav ) CALL zgr_isf 996 997 ! Scale factors and depth at T- and W-points 998 DO jk = 1, jpk ! intitialization to the reference z-coordinate 999 gdept_0(:,:,jk) = gdept_1d(jk) 1000 gdepw_0(:,:,jk) = gdepw_1d(jk) 1001 e3t_0 (:,:,jk) = e3t_1d (jk) 1002 e3w_0 (:,:,jk) = e3w_1d (jk) 1003 END DO 1004 ! 1005 DO jj = 1, jpj 1006 DO ji = 1, jpi 1007 ik = mbathy(ji,jj) 1008 IF( ik > 0 ) THEN ! ocean point only 1009 ! max ocean level case 1010 IF( ik == jpkm1 ) THEN 1011 zdepwp = bathy(ji,jj) 1012 ze3tp = bathy(ji,jj) - gdepw_1d(ik) 1013 ze3wp = 0.5_wp * e3w_1d(ik) * ( 1._wp + ( ze3tp/e3t_1d(ik) ) ) 1014 e3t_0(ji,jj,ik ) = ze3tp 1015 e3t_0(ji,jj,ik+1) = ze3tp 1016 e3w_0(ji,jj,ik ) = ze3wp 1017 e3w_0(ji,jj,ik+1) = ze3tp 1018 gdepw_0(ji,jj,ik+1) = zdepwp 1019 gdept_0(ji,jj,ik ) = gdept_1d(ik-1) + ze3wp 1020 gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + ze3tp 1021 ! 1022 ELSE ! standard case 1023 IF( bathy(ji,jj) <= gdepw_1d(ik+1) ) THEN ; gdepw_0(ji,jj,ik+1) = bathy(ji,jj) 1024 ELSE ; gdepw_0(ji,jj,ik+1) = gdepw_1d(ik+1) 1025 ENDIF 1026 !gm Bug? check the gdepw_1d 1027 ! ... on ik 1028 gdept_0(ji,jj,ik) = gdepw_1d(ik) + ( gdepw_0(ji,jj,ik+1) - gdepw_1d(ik) ) & 1029 & * ((gdept_1d( ik ) - gdepw_1d(ik) ) & 1030 & / ( gdepw_1d( ik+1) - gdepw_1d(ik) )) 1031 e3t_0 (ji,jj,ik) = e3t_1d (ik) * ( gdepw_0 (ji,jj,ik+1) - gdepw_1d(ik) ) & 1032 & / ( gdepw_1d( ik+1) - gdepw_1d(ik) ) 1033 e3w_0(ji,jj,ik) = 0.5_wp * ( gdepw_0(ji,jj,ik+1) + gdepw_1d(ik+1) - 2._wp * gdepw_1d(ik) ) & 1034 & * ( e3w_1d(ik) / ( gdepw_1d(ik+1) - gdepw_1d(ik) ) ) 1035 ! ... on ik+1 1036 e3w_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik) 1037 e3t_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik) 1038 gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik) 1039 ENDIF 1040 ENDIF 1041 END DO 1042 END DO 1043 ! 1044 it = 0 1045 DO jj = 1, jpj 1046 DO ji = 1, jpi 1047 ik = mbathy(ji,jj) 1048 IF( ik > 0 ) THEN ! ocean point only 1049 e3tp (ji,jj) = e3t_0(ji,jj,ik) 1050 e3wp (ji,jj) = e3w_0(ji,jj,ik) 1051 ! test 1052 zdiff= gdepw_0(ji,jj,ik+1) - gdept_0(ji,jj,ik ) 1053 IF( zdiff <= 0._wp .AND. lwp ) THEN 1054 it = it + 1 1055 WRITE(numout,*) ' it = ', it, ' ik = ', ik, ' (i,j) = ', ji, jj 1056 WRITE(numout,*) ' bathy = ', bathy(ji,jj) 1057 WRITE(numout,*) ' gdept_0 = ', gdept_0(ji,jj,ik), ' gdepw_0 = ', gdepw_0(ji,jj,ik+1), ' zdiff = ', zdiff 1058 WRITE(numout,*) ' e3tp = ', e3t_0 (ji,jj,ik), ' e3wp = ', e3w_0 (ji,jj,ik ) 1059 ENDIF 1060 ENDIF 1061 END DO 1062 END DO 1063 ! 1064 IF ( ln_isfcav ) THEN 1065 ! (ISF) Definition of e3t, u, v, w for ISF case 1066 DO jj = 1, jpj 1067 DO ji = 1, jpi 1068 ik = misfdep(ji,jj) 1069 IF( ik > 1 ) THEN ! ice shelf point only 1070 IF( risfdep(ji,jj) < gdepw_1d(ik) ) risfdep(ji,jj)= gdepw_1d(ik) 1071 gdepw_0(ji,jj,ik) = risfdep(ji,jj) 1072 !gm Bug? check the gdepw_0 1073 ! ... on ik 1074 gdept_0(ji,jj,ik) = gdepw_1d(ik+1) - ( gdepw_1d(ik+1) - gdepw_0(ji,jj,ik) ) & 1075 & * ( gdepw_1d(ik+1) - gdept_1d(ik) ) & 1076 & / ( gdepw_1d(ik+1) - gdepw_1d(ik) ) 1077 e3t_0 (ji,jj,ik ) = gdepw_1d(ik+1) - gdepw_0(ji,jj,ik) 1078 e3w_0 (ji,jj,ik+1) = gdept_1d(ik+1) - gdept_0(ji,jj,ik) 1079 1080 IF( ik + 1 == mbathy(ji,jj) ) THEN ! ice shelf point only (2 cell water column) 1081 e3w_0 (ji,jj,ik+1) = gdept_0(ji,jj,ik+1) - gdept_0(ji,jj,ik) 1082 ENDIF 1083 ! ... on ik / ik-1 1084 e3w_0 (ji,jj,ik ) = 2._wp * (gdept_0(ji,jj,ik) - gdepw_0(ji,jj,ik)) 1085 e3t_0 (ji,jj,ik-1) = gdepw_0(ji,jj,ik) - gdepw_1d(ik-1) 1086 ! The next line isn't required and doesn't affect results - included for consistency with bathymetry code 1087 gdept_0(ji,jj,ik-1) = gdept_1d(ik-1) 1088 ENDIF 1089 END DO 1090 END DO 1091 ! 1092 it = 0 1093 DO jj = 1, jpj 1094 DO ji = 1, jpi 1095 ik = misfdep(ji,jj) 1096 IF( ik > 1 ) THEN ! ice shelf point only 1097 e3tp (ji,jj) = e3t_0(ji,jj,ik ) 1098 e3wp (ji,jj) = e3w_0(ji,jj,ik+1 ) 1099 ! test 1100 zdiff= gdept_0(ji,jj,ik) - gdepw_0(ji,jj,ik ) 1101 IF( zdiff <= 0. .AND. lwp ) THEN 1102 it = it + 1 1103 WRITE(numout,*) ' it = ', it, ' ik = ', ik, ' (i,j) = ', ji, jj 1104 WRITE(numout,*) ' risfdep = ', risfdep(ji,jj) 1105 WRITE(numout,*) ' gdept = ', gdept_0(ji,jj,ik), ' gdepw = ', gdepw_0(ji,jj,ik+1), ' zdiff = ', zdiff 1106 WRITE(numout,*) ' e3tp = ', e3tp(ji,jj), ' e3wp = ', e3wp(ji,jj) 1107 ENDIF 1108 ENDIF 1109 END DO 1110 END DO 1111 END IF 1112 ! END (ISF) 1113 1114 ! Scale factors and depth at U-, V-, UW and VW-points 1115 DO jk = 1, jpk ! initialisation to z-scale factors 1116 e3u_0 (:,:,jk) = e3t_1d(jk) 1117 e3v_0 (:,:,jk) = e3t_1d(jk) 1118 e3uw_0(:,:,jk) = e3w_1d(jk) 1119 e3vw_0(:,:,jk) = e3w_1d(jk) 1120 END DO 1121 DO jk = 1,jpk ! Computed as the minimum of neighbooring scale factors 1122 DO jj = 1, jpjm1 1123 DO ji = 1, fs_jpim1 ! vector opt. 1124 e3u_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) ) 1125 e3v_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) ) 1126 e3uw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji+1,jj,jk) ) 1127 e3vw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji,jj+1,jk) ) 1128 END DO 1129 END DO 1130 END DO 1131 IF ( ln_isfcav ) THEN 1132 ! (ISF) define e3uw (adapted for 2 cells in the water column) 1133 DO jj = 2, jpjm1 1134 DO ji = 2, fs_jpim1 ! vector opt. 1135 ikb = MAX(mbathy (ji,jj),mbathy (ji+1,jj)) 1136 ikt = MAX(misfdep(ji,jj),misfdep(ji+1,jj)) 1137 IF (ikb == ikt+1) e3uw_0(ji,jj,ikb) = MIN( gdept_0(ji,jj,ikb ), gdept_0(ji+1,jj ,ikb ) ) & 1138 & - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji+1,jj ,ikb-1) ) 1139 ikb = MAX(mbathy (ji,jj),mbathy (ji,jj+1)) 1140 ikt = MAX(misfdep(ji,jj),misfdep(ji,jj+1)) 1141 IF (ikb == ikt+1) e3vw_0(ji,jj,ikb) = MIN( gdept_0(ji,jj,ikb ), gdept_0(ji ,jj+1,ikb ) ) & 1142 & - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji ,jj+1,ikb-1) ) 1143 END DO 1144 END DO 1145 END IF 1146 1147 CALL lbc_lnk( e3u_0 , 'U', 1._wp ) ; CALL lbc_lnk( e3uw_0, 'U', 1._wp ) ! lateral boundary conditions 1148 CALL lbc_lnk( e3v_0 , 'V', 1._wp ) ; CALL lbc_lnk( e3vw_0, 'V', 1._wp ) 1149 ! 1150 DO jk = 1, jpk ! set to z-scale factor if zero (i.e. along closed boundaries) 1151 WHERE( e3u_0 (:,:,jk) == 0._wp ) e3u_0 (:,:,jk) = e3t_1d(jk) 1152 WHERE( e3v_0 (:,:,jk) == 0._wp ) e3v_0 (:,:,jk) = e3t_1d(jk) 1153 WHERE( e3uw_0(:,:,jk) == 0._wp ) e3uw_0(:,:,jk) = e3w_1d(jk) 1154 WHERE( e3vw_0(:,:,jk) == 0._wp ) e3vw_0(:,:,jk) = e3w_1d(jk) 1155 END DO 1156 1157 ! Scale factor at F-point 1158 DO jk = 1, jpk ! initialisation to z-scale factors 1159 e3f_0(:,:,jk) = e3t_1d(jk) 1160 END DO 1161 DO jk = 1, jpk ! Computed as the minimum of neighbooring V-scale factors 1162 DO jj = 1, jpjm1 1163 DO ji = 1, fs_jpim1 ! vector opt. 1164 e3f_0(ji,jj,jk) = MIN( e3v_0(ji,jj,jk), e3v_0(ji+1,jj,jk) ) 1165 END DO 1166 END DO 1167 END DO 1168 CALL lbc_lnk( e3f_0, 'F', 1._wp ) ! Lateral boundary conditions 1169 ! 1170 DO jk = 1, jpk ! set to z-scale factor if zero (i.e. along closed boundaries) 1171 WHERE( e3f_0(:,:,jk) == 0._wp ) e3f_0(:,:,jk) = e3t_1d(jk) 1172 END DO 1173 !!gm bug ? : must be a do loop with mj0,mj1 1174 ! 1175 e3t_0(:,mj0(1),:) = e3t_0(:,mj0(2),:) ! we duplicate factor scales for jj = 1 and jj = 2 1176 e3w_0(:,mj0(1),:) = e3w_0(:,mj0(2),:) 1177 e3u_0(:,mj0(1),:) = e3u_0(:,mj0(2),:) 1178 e3v_0(:,mj0(1),:) = e3v_0(:,mj0(2),:) 1179 e3f_0(:,mj0(1),:) = e3f_0(:,mj0(2),:) 1180 1181 ! Control of the sign 1182 IF( MINVAL( e3t_0 (:,:,:) ) <= 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r e3t_0 <= 0' ) 1183 IF( MINVAL( e3w_0 (:,:,:) ) <= 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r e3w_0 <= 0' ) 1184 IF( MINVAL( gdept_0(:,:,:) ) < 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r gdept_0 < 0' ) 1185 IF( MINVAL( gdepw_0(:,:,:) ) < 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r gdepw_0 < 0' ) 1186 1187 ! Compute gdep3w_0 (vertical sum of e3w) 1188 IF ( ln_isfcav ) THEN ! if cavity 1189 WHERE (misfdep == 0) misfdep = 1 1190 DO jj = 1,jpj 1191 DO ji = 1,jpi 1192 gdep3w_0(ji,jj,1) = 0.5_wp * e3w_0(ji,jj,1) 1193 DO jk = 2, misfdep(ji,jj) 1194 gdep3w_0(ji,jj,jk) = gdep3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk) 1195 END DO 1196 IF (misfdep(ji,jj) .GE. 2) gdep3w_0(ji,jj,misfdep(ji,jj)) = risfdep(ji,jj) + 0.5_wp * e3w_0(ji,jj,misfdep(ji,jj)) 1197 DO jk = misfdep(ji,jj) + 1, jpk 1198 gdep3w_0(ji,jj,jk) = gdep3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk) 1199 END DO 1200 END DO 1201 END DO 1202 ELSE ! no cavity 1203 gdep3w_0(:,:,1) = 0.5_wp * e3w_0(:,:,1) 1204 DO jk = 2, jpk 1205 gdep3w_0(:,:,jk) = gdep3w_0(:,:,jk-1) + e3w_0(:,:,jk) 1206 END DO 1207 END IF 1208 ! ! ================= ! 1209 IF(lwp .AND. ll_print) THEN ! Control print ! 1210 ! ! ================= ! 1211 DO jj = 1,jpj 1212 DO ji = 1, jpi 1213 ik = MAX( mbathy(ji,jj), 1 ) 1214 zprt(ji,jj,1) = e3t_0 (ji,jj,ik) 1215 zprt(ji,jj,2) = e3w_0 (ji,jj,ik) 1216 zprt(ji,jj,3) = e3u_0 (ji,jj,ik) 1217 zprt(ji,jj,4) = e3v_0 (ji,jj,ik) 1218 zprt(ji,jj,5) = e3f_0 (ji,jj,ik) 1219 zprt(ji,jj,6) = gdep3w_0(ji,jj,ik) 1220 END DO 1221 END DO 1222 WRITE(numout,*) 1223 WRITE(numout,*) 'domzgr e3t(mbathy)' ; CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1224 WRITE(numout,*) 1225 WRITE(numout,*) 'domzgr e3w(mbathy)' ; CALL prihre(zprt(:,:,2),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1226 WRITE(numout,*) 1227 WRITE(numout,*) 'domzgr e3u(mbathy)' ; CALL prihre(zprt(:,:,3),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1228 WRITE(numout,*) 1229 WRITE(numout,*) 'domzgr e3v(mbathy)' ; CALL prihre(zprt(:,:,4),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1230 WRITE(numout,*) 1231 WRITE(numout,*) 'domzgr e3f(mbathy)' ; CALL prihre(zprt(:,:,5),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1232 WRITE(numout,*) 1233 WRITE(numout,*) 'domzgr gdep3w(mbathy)' ; CALL prihre(zprt(:,:,6),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1234 ENDIF 1235 ! 1236 CALL wrk_dealloc( jpi, jpj, jpk, zprt ) 1237 ! 1238 IF( nn_timing == 1 ) CALL timing_stop('zgr_zps') 1239 ! 1240 END SUBROUTINE zgr_zps 1241 1242 SUBROUTINE zgr_isf 1243 !!---------------------------------------------------------------------- 1244 !! *** ROUTINE zgr_isf *** 1245 !! 1246 !! ** Purpose : check the bathymetry in levels 1247 !! 1248 !! ** Method : THe water column have to contained at least 2 cells 1249 !! Bathymetry and isfdraft are modified (dig/close) to respect 1250 !! this criterion. 1251 !! 1252 !! 1253 !! ** Action : - test compatibility between isfdraft and bathy 1254 !! - bathy and isfdraft are modified 1255 !!---------------------------------------------------------------------- 1256 !! 961 1257 INTEGER :: ji, jj, jk, jl ! dummy loop indices 962 1258 INTEGER :: ik, it ! temporary integers … … 969 1265 REAL(wp) :: zdiff ! temporary scalar 970 1266 REAL(wp) :: zrefdep ! temporary scalar 971 REAL(wp) :: zbathydiff, zrisfdepdiff 972 REAL(wp), POINTER, DIMENSION(:,:) :: zrisfdep, zbathy, zmask ! 3D workspace (ISH) 973 INTEGER , POINTER, DIMENSION(:,:) :: zmbathy, zmisfdep ! 3D workspace (ISH) 974 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprt 1267 REAL(wp) :: zbathydiff, zrisfdepdiff ! isf temporary scalar 1268 REAL(wp), POINTER, DIMENSION(:,:) :: zrisfdep, zbathy, zmask ! 2D workspace (ISH) 1269 INTEGER , POINTER, DIMENSION(:,:) :: zmbathy, zmisfdep ! 2D workspace (ISH) 975 1270 !!--------------------------------------------------------------------- 976 1271 ! 977 IF( nn_timing == 1 ) CALL timing_start('zgr_zps') 978 ! 979 CALL wrk_alloc( jpi, jpj, jpk, zprt ) 1272 IF( nn_timing == 1 ) CALL timing_start('zgr_isf') 1273 ! 980 1274 CALL wrk_alloc( jpi, jpj, zbathy, zmask, zrisfdep) 981 CALL wrk_alloc( jpi, jpj, zmbathy, zmisfdep) 982 ! 983 IF(lwp) WRITE(numout,*) 984 IF(lwp) WRITE(numout,*) ' zgr_zps : z-coordinate with partial steps' 985 IF(lwp) WRITE(numout,*) ' ~~~~~~~ ' 986 IF(lwp) WRITE(numout,*) ' mbathy is recomputed : bathy_level file is NOT used' 987 988 ll_print = .FALSE. ! Local variable for debugging 989 990 IF(lwp .AND. ll_print) THEN ! control print of the ocean depth 991 WRITE(numout,*) 992 WRITE(numout,*) 'dom_zgr_zps: bathy (in hundred of meters)' 993 CALL prihre( bathy, jpi, jpj, 1,jpi, 1, 1, jpj, 1, 1.e-2, numout ) 994 ENDIF 995 996 ! bathymetry in level (from bathy_meter) 997 ! =================== 998 zmax = gdepw_1d(jpk) + e3t_1d(jpk) ! maximum depth (i.e. the last ocean level thickness <= 2*e3t_1d(jpkm1) ) 999 bathy(:,:) = MIN( zmax , bathy(:,:) ) ! bounded value of bathy (min already set at the end of zgr_bat) 1000 WHERE( bathy(:,:) == 0._wp ) ; mbathy(:,:) = 0 ! land : set mbathy to 0 1001 ELSE WHERE ; mbathy(:,:) = jpkm1 ! ocean : initialize mbathy to the max ocean level 1002 END WHERE 1003 1004 ! Compute mbathy for ocean points (i.e. the number of ocean levels) 1005 ! find the number of ocean levels such that the last level thickness 1006 ! is larger than the minimum of e3zps_min and e3zps_rat * e3t_1d (where 1007 ! e3t_1d is the reference level thickness 1008 DO jk = jpkm1, 1, -1 1009 zdepth = gdepw_1d(jk) + MIN( e3zps_min, e3t_1d(jk)*e3zps_rat ) 1010 WHERE( 0._wp < bathy(:,:) .AND. bathy(:,:) <= zdepth ) mbathy(:,:) = jk-1 1011 END DO 1275 CALL wrk_alloc( jpi, jpj, zmisfdep, zmbathy ) 1276 1277 1012 1278 ! (ISF) compute misfdep 1013 1279 WHERE( risfdep(:,:) == 0._wp .AND. bathy(:,:) .NE. 0) ; misfdep(:,:) = 1 ! open water : set misfdep to 1 … … 1053 1319 misfdep(jpi,:) = misfdep( 2 ,:) 1054 1320 ENDIF 1055 1321 1056 1322 IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN 1057 1323 mbathy( 1 ,:) = mbathy(jpim1,:) ! local domain is cyclic east-west 1058 1324 mbathy(jpi,:) = mbathy( 2 ,:) 1059 1325 ENDIF 1060 1326 1061 1327 ! split last cell if possible (only where water column is 2 cell or less) 1062 1328 DO jk = jpkm1, 1, -1 … … 1076 1342 END WHERE 1077 1343 END DO 1078 1344 1079 1345 1080 1346 ! Case where bathy and risfdep compatible but not the level variable mbathy/misfdep because of partial cell condition … … 1252 1518 1253 1519 ! remove single point "bay" on isf coast line in the ice shelf draft' 1254 DO jk = 1, jpk1520 DO jk = 2, jpk 1255 1521 WHERE (misfdep==0) misfdep=jpk 1256 1522 zmask=0 … … 1357 1623 IF( zmbathy(ji,jj) .LT. misfdep(ji ,jj+1) ) ibtestjp1 = 0 1358 1624 ibtest=MAX(ibtestim1, ibtestip1, ibtestjm1, ibtestjp1) 1359 IF( ibtest == 0 ) THEN1625 IF( ibtest == 0 .AND. misfdep(ji,jj) .GE. 2) THEN 1360 1626 mbathy(ji,jj) = 0 ; bathy(ji,jj) = 0.0_wp ; misfdep(ji,jj) = 0 ; risfdep(ji,jj) = 0.0_wp ; 1361 1627 END IF … … 1473 1739 ENDIF 1474 1740 1475 ! Scale factors and depth at T- and W-points1476 DO jk = 1, jpk ! intitialization to the reference z-coordinate1477 gdept_0(:,:,jk) = gdept_1d(jk)1478 gdepw_0(:,:,jk) = gdepw_1d(jk)1479 e3t_0 (:,:,jk) = e3t_1d (jk)1480 e3w_0 (:,:,jk) = e3w_1d (jk)1481 END DO1482 !1483 DO jj = 1, jpj1484 DO ji = 1, jpi1485 ik = mbathy(ji,jj)1486 IF( ik > 0 ) THEN ! ocean point only1487 ! max ocean level case1488 IF( ik == jpkm1 ) THEN1489 zdepwp = bathy(ji,jj)1490 ze3tp = bathy(ji,jj) - gdepw_1d(ik)1491 ze3wp = 0.5_wp * e3w_1d(ik) * ( 1._wp + ( ze3tp/e3t_1d(ik) ) )1492 e3t_0(ji,jj,ik ) = ze3tp1493 e3t_0(ji,jj,ik+1) = ze3tp1494 e3w_0(ji,jj,ik ) = ze3wp1495 e3w_0(ji,jj,ik+1) = ze3tp1496 gdepw_0(ji,jj,ik+1) = zdepwp1497 gdept_0(ji,jj,ik ) = gdept_1d(ik-1) + ze3wp1498 gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + ze3tp1499 !1500 ELSE ! standard case1501 IF( bathy(ji,jj) <= gdepw_1d(ik+1) ) THEN ; gdepw_0(ji,jj,ik+1) = bathy(ji,jj)1502 ELSE ; gdepw_0(ji,jj,ik+1) = gdepw_1d(ik+1)1503 ENDIF1504 !gm Bug? check the gdepw_1d1505 ! ... on ik1506 gdept_0(ji,jj,ik) = gdepw_1d(ik) + ( gdepw_0(ji,jj,ik+1) - gdepw_1d(ik) ) &1507 & * ((gdept_1d( ik ) - gdepw_1d(ik) ) &1508 & / ( gdepw_1d( ik+1) - gdepw_1d(ik) ))1509 e3t_0(ji,jj,ik) = e3t_1d (ik) * ( gdepw_0 (ji,jj,ik+1) - gdepw_1d(ik) ) &1510 & / ( gdepw_1d( ik+1) - gdepw_1d(ik) )1511 e3w_0(ji,jj,ik) = 0.5_wp * ( gdepw_0(ji,jj,ik+1) + gdepw_1d(ik+1) - 2._wp * gdepw_1d(ik) ) &1512 & * ( e3w_1d(ik) / ( gdepw_1d(ik+1) - gdepw_1d(ik) ) )1513 ! ... on ik+11514 e3w_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik)1515 e3t_0 (ji,jj,ik+1) = e3t_0 (ji,jj,ik)1516 gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik)1517 ENDIF1518 ENDIF1519 END DO1520 END DO1521 !1522 it = 01523 DO jj = 1, jpj1524 DO ji = 1, jpi1525 ik = mbathy(ji,jj)1526 IF( ik > 0 ) THEN ! ocean point only1527 e3tp (ji,jj) = e3t_0(ji,jj,ik)1528 e3wp (ji,jj) = e3w_0(ji,jj,ik)1529 ! test1530 zdiff= gdepw_0(ji,jj,ik+1) - gdept_0(ji,jj,ik )1531 IF( zdiff <= 0._wp .AND. lwp ) THEN1532 it = it + 11533 WRITE(numout,*) ' it = ', it, ' ik = ', ik, ' (i,j) = ', ji, jj1534 WRITE(numout,*) ' bathy = ', bathy(ji,jj)1535 WRITE(numout,*) ' gdept_0 = ', gdept_0(ji,jj,ik), ' gdepw_0 = ', gdepw_0(ji,jj,ik+1), ' zdiff = ', zdiff1536 WRITE(numout,*) ' e3tp = ', e3t_0 (ji,jj,ik), ' e3wp = ', e3w_0 (ji,jj,ik )1537 ENDIF1538 ENDIF1539 END DO1540 END DO1541 !1542 ! (ISF) Definition of e3t, u, v, w for ISF case1543 DO jj = 1, jpj1544 DO ji = 1, jpi1545 ik = misfdep(ji,jj)1546 IF( ik > 1 ) THEN ! ice shelf point only1547 IF( risfdep(ji,jj) < gdepw_1d(ik) ) risfdep(ji,jj)= gdepw_1d(ik)1548 gdepw_0(ji,jj,ik) = risfdep(ji,jj)1549 !gm Bug? check the gdepw_01550 ! ... on ik1551 gdept_0(ji,jj,ik) = gdepw_1d(ik+1) - ( gdepw_1d(ik+1) - gdepw_0(ji,jj,ik) ) &1552 & * ( gdepw_1d(ik+1) - gdept_1d(ik) ) &1553 & / ( gdepw_1d(ik+1) - gdepw_1d(ik) )1554 e3t_0 (ji,jj,ik ) = gdepw_1d(ik+1) - gdepw_0(ji,jj,ik)1555 e3w_0 (ji,jj,ik+1) = gdept_1d(ik+1) - gdept_0(ji,jj,ik)1556 1557 IF( ik + 1 == mbathy(ji,jj) ) THEN ! ice shelf point only (2 cell water column)1558 e3w_0 (ji,jj,ik+1) = gdept_0(ji,jj,ik+1) - gdept_0(ji,jj,ik)1559 ENDIF1560 ! ... on ik / ik-11561 e3w_0 (ji,jj,ik ) = 2._wp * (gdept_0(ji,jj,ik) - gdepw_0(ji,jj,ik))1562 e3t_0 (ji,jj,ik-1) = gdepw_0(ji,jj,ik) - gdepw_1d(ik-1)1563 ! The next line isn't required and doesn't affect results - included for consistency with bathymetry code1564 gdept_0(ji,jj,ik-1) = gdept_1d(ik-1)1565 ENDIF1566 END DO1567 END DO1568 !1569 it = 01570 DO jj = 1, jpj1571 DO ji = 1, jpi1572 ik = misfdep(ji,jj)1573 IF( ik > 1 ) THEN ! ice shelf point only1574 e3tp (ji,jj) = e3t_0(ji,jj,ik )1575 e3wp (ji,jj) = e3w_0(ji,jj,ik+1 )1576 ! test1577 zdiff= gdept_0(ji,jj,ik) - gdepw_0(ji,jj,ik )1578 IF( zdiff <= 0. .AND. lwp ) THEN1579 it = it + 11580 WRITE(numout,*) ' it = ', it, ' ik = ', ik, ' (i,j) = ', ji, jj1581 WRITE(numout,*) ' risfdep = ', risfdep(ji,jj)1582 WRITE(numout,*) ' gdept = ', gdept_0(ji,jj,ik), ' gdepw = ', gdepw_0(ji,jj,ik+1), ' zdiff = ', zdiff1583 WRITE(numout,*) ' e3tp = ', e3tp(ji,jj), ' e3wp = ', e3wp(ji,jj)1584 ENDIF1585 ENDIF1586 END DO1587 END DO1588 ! END (ISF)1589 1590 ! Scale factors and depth at U-, V-, UW and VW-points1591 DO jk = 1, jpk ! initialisation to z-scale factors1592 e3u_0 (:,:,jk) = e3t_1d(jk)1593 e3v_0 (:,:,jk) = e3t_1d(jk)1594 e3uw_0(:,:,jk) = e3w_1d(jk)1595 e3vw_0(:,:,jk) = e3w_1d(jk)1596 END DO1597 DO jk = 1,jpk ! Computed as the minimum of neighbooring scale factors1598 DO jj = 1, jpjm11599 DO ji = 1, fs_jpim1 ! vector opt.1600 e3u_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) )1601 e3v_0 (ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) )1602 e3uw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji+1,jj,jk) )1603 e3vw_0(ji,jj,jk) = MIN( e3w_0(ji,jj,jk), e3w_0(ji,jj+1,jk) )1604 END DO1605 END DO1606 END DO1607 ! (ISF) define e3uw1608 DO jk = 2,jpk1609 DO jj = 1, jpjm11610 DO ji = 1, fs_jpim1 ! vector opt.1611 e3uw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji+1,jj ,jk) ) &1612 & - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji+1,jj ,jk-1) )1613 e3vw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji ,jj+1,jk) ) &1614 & - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji ,jj+1,jk-1) )1615 END DO1616 END DO1617 END DO1618 !End (ISF)1619 1620 CALL lbc_lnk( e3u_0 , 'U', 1._wp ) ; CALL lbc_lnk( e3uw_0, 'U', 1._wp ) ! lateral boundary conditions1621 CALL lbc_lnk( e3v_0 , 'V', 1._wp ) ; CALL lbc_lnk( e3vw_0, 'V', 1._wp )1622 !1623 DO jk = 1, jpk ! set to z-scale factor if zero (i.e. along closed boundaries)1624 WHERE( e3u_0 (:,:,jk) == 0._wp ) e3u_0 (:,:,jk) = e3t_1d(jk)1625 WHERE( e3v_0 (:,:,jk) == 0._wp ) e3v_0 (:,:,jk) = e3t_1d(jk)1626 WHERE( e3uw_0(:,:,jk) == 0._wp ) e3uw_0(:,:,jk) = e3w_1d(jk)1627 WHERE( e3vw_0(:,:,jk) == 0._wp ) e3vw_0(:,:,jk) = e3w_1d(jk)1628 END DO1629 1630 ! Scale factor at F-point1631 DO jk = 1, jpk ! initialisation to z-scale factors1632 e3f_0(:,:,jk) = e3t_1d(jk)1633 END DO1634 DO jk = 1, jpk ! Computed as the minimum of neighbooring V-scale factors1635 DO jj = 1, jpjm11636 DO ji = 1, fs_jpim1 ! vector opt.1637 e3f_0(ji,jj,jk) = MIN( e3v_0(ji,jj,jk), e3v_0(ji+1,jj,jk) )1638 END DO1639 END DO1640 END DO1641 CALL lbc_lnk( e3f_0, 'F', 1._wp ) ! Lateral boundary conditions1642 !1643 DO jk = 1, jpk ! set to z-scale factor if zero (i.e. along closed boundaries)1644 WHERE( e3f_0(:,:,jk) == 0._wp ) e3f_0(:,:,jk) = e3t_1d(jk)1645 END DO1646 !!gm bug ? : must be a do loop with mj0,mj11647 !1648 e3t_0(:,mj0(1),:) = e3t_0(:,mj0(2),:) ! we duplicate factor scales for jj = 1 and jj = 21649 e3w_0(:,mj0(1),:) = e3w_0(:,mj0(2),:)1650 e3u_0(:,mj0(1),:) = e3u_0(:,mj0(2),:)1651 e3v_0(:,mj0(1),:) = e3v_0(:,mj0(2),:)1652 e3f_0(:,mj0(1),:) = e3f_0(:,mj0(2),:)1653 1654 ! Control of the sign1655 IF( MINVAL( e3t_0 (:,:,:) ) <= 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r e3t_0 <= 0' )1656 IF( MINVAL( e3w_0 (:,:,:) ) <= 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r e3w_0 <= 0' )1657 IF( MINVAL( gdept_0(:,:,:) ) < 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r gdept_0 < 0' )1658 IF( MINVAL( gdepw_0(:,:,:) ) < 0._wp ) CALL ctl_stop( ' zgr_zps : e r r o r gdepw_0 < 0' )1659 1660 ! Compute gdep3w_0 (vertical sum of e3w)1661 WHERE (misfdep == 0) misfdep = 11662 DO jj = 1,jpj1663 DO ji = 1,jpi1664 gdep3w_0(ji,jj,1) = 0.5_wp * e3w_0(ji,jj,1)1665 DO jk = 2, misfdep(ji,jj)1666 gdep3w_0(ji,jj,jk) = gdep3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk)1667 END DO1668 IF (misfdep(ji,jj) .GE. 2) gdep3w_0(ji,jj,misfdep(ji,jj)) = risfdep(ji,jj) + 0.5_wp * e3w_0(ji,jj,misfdep(ji,jj))1669 DO jk = misfdep(ji,jj) + 1, jpk1670 gdep3w_0(ji,jj,jk) = gdep3w_0(ji,jj,jk-1) + e3w_0(ji,jj,jk)1671 END DO1672 END DO1673 END DO1674 ! ! ================= !1675 IF(lwp .AND. ll_print) THEN ! Control print !1676 ! ! ================= !1677 DO jj = 1,jpj1678 DO ji = 1, jpi1679 ik = MAX( mbathy(ji,jj), 1 )1680 zprt(ji,jj,1) = e3t_0 (ji,jj,ik)1681 zprt(ji,jj,2) = e3w_0 (ji,jj,ik)1682 zprt(ji,jj,3) = e3u_0 (ji,jj,ik)1683 zprt(ji,jj,4) = e3v_0 (ji,jj,ik)1684 zprt(ji,jj,5) = e3f_0 (ji,jj,ik)1685 zprt(ji,jj,6) = gdep3w_0(ji,jj,ik)1686 END DO1687 END DO1688 WRITE(numout,*)1689 WRITE(numout,*) 'domzgr e3t(mbathy)' ; CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1690 WRITE(numout,*)1691 WRITE(numout,*) 'domzgr e3w(mbathy)' ; CALL prihre(zprt(:,:,2),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1692 WRITE(numout,*)1693 WRITE(numout,*) 'domzgr e3u(mbathy)' ; CALL prihre(zprt(:,:,3),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1694 WRITE(numout,*)1695 WRITE(numout,*) 'domzgr e3v(mbathy)' ; CALL prihre(zprt(:,:,4),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1696 WRITE(numout,*)1697 WRITE(numout,*) 'domzgr e3f(mbathy)' ; CALL prihre(zprt(:,:,5),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1698 WRITE(numout,*)1699 WRITE(numout,*) 'domzgr gdep3w(mbathy)' ; CALL prihre(zprt(:,:,6),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1700 ENDIF1701 !1702 CALL wrk_dealloc( jpi, jpj, jpk, zprt )1703 1741 CALL wrk_dealloc( jpi, jpj, zmask, zbathy, zrisfdep ) 1704 1742 CALL wrk_dealloc( jpi, jpj, zmisfdep, zmbathy ) 1705 ! 1706 IF( nn_timing == 1 ) CALL timing_stop('zgr_ zps')1707 !1708 END SUBROUTINE zgr_zps1743 1744 IF( nn_timing == 1 ) CALL timing_stop('zgr_isf') 1745 1746 END SUBROUTINE 1709 1747 1710 1748 SUBROUTINE zgr_sco -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/dtatsd.F90
- Property svn:keywords set to Id
r5038 r5620 39 39 !!---------------------------------------------------------------------- 40 40 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 41 !! $Id : dtatem.F90 2392 2010-11-15 21:20:05Z gm$41 !! $Id$ 42 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 43 43 !!---------------------------------------------------------------------- -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90
r5038 r5620 69 69 !! ** Purpose : Initialization of the dynamics and tracer fields. 70 70 !!---------------------------------------------------------------------- 71 ! - ML - needed for initialization of e3t_b 72 INTEGER :: ji,jj,jk ! dummy loop indices 73 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zuvd ! U & V data workspace 71 INTEGER :: ji, jj, jk ! dummy loop indices 72 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zuvd ! U & V data workspace 74 73 !!---------------------------------------------------------------------- 75 74 ! … … 84 83 IF( lk_c1d ) CALL dta_uvd_init ! Initialization of U & V input data 85 84 86 rhd (:,:,: ) = 0._wp 87 rhop (:,:,: ) = 0._wp 88 rn2 (:,:,: ) = 0._wp 89 tsa (:,:,:,:) = 0._wp 90 rab_b(:,:,:,:) = 0._wp 91 rab_n(:,:,:,:) = 0._wp 85 rhd (:,:,: ) = 0._wp ; rhop (:,:,: ) = 0._wp ! set one for all to 0 at level jpk 86 rn2b (:,:,: ) = 0._wp ; rn2 (:,:,: ) = 0._wp ! set one for all to 0 at levels 1 and jpk 87 tsa (:,:,:,:) = 0._wp ! set one for all to 0 at level jpk 88 rab_b(:,:,:,:) = 0._wp ; rab_n(:,:,:,:) = 0._wp ! set one for all to 0 at level jpk 92 89 93 90 IF( ln_rstart ) THEN ! Restart from a file … … 113 110 ELSEIF( cp_cfg == 'gyre' ) THEN 114 111 CALL istate_gyre ! GYRE configuration : start from pre-defined T-S fields 115 ELSEIF( cp_cfg == 'isomip' .OR. cp_cfg == 'isomip2') THEN116 IF(lwp) WRITE(numout,*) 'Initialization of T+S for ISOMIP domain'117 tsn(:,:,:,jp_tem)=-1.9*tmask(:,:,:) ! ISOMIP configuration : start from constant T+S fields118 tsn(:,:,:,jp_sal)=34.4*tmask(:,:,:)119 tsb(:,:,:,:)=tsn(:,:,:,:)120 112 ELSE ! Initial T-S, U-V fields read in files 121 113 IF ( ln_tsd_init ) THEN ! read 3D T and S data at nit000 … … 137 129 CALL eos( tsb, rhd, rhop, gdept_0(:,:,:) ) ! before potential and in situ densities 138 130 #if ! defined key_c1d 139 IF( ln_zps ) CALL zps_hde( nit000, jpts, tsb, gtsu, gtsv, & ! Partial steps: before horizontal gradient 140 & rhd, gru , grv, aru, arv, gzu, gzv, ge3ru, ge3rv, & ! 141 & gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi ) ! of t, s, rd at the last ocean level 131 IF( ln_zps .AND. .NOT. ln_isfcav) & 132 & CALL zps_hde ( nit000, jpts, tsb, gtsu, gtsv, & ! Partial steps: before horizontal gradient 133 & rhd, gru , grv ) ! of t, s, rd at the last ocean level 134 IF( ln_zps .AND. ln_isfcav) & 135 & CALL zps_hde_isf( nit000, jpts, tsb, gtsu, gtsv, & ! Partial steps for top cell (ISF) 136 & rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv , & 137 & gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi ) ! of t, s, rd at the last ocean level 142 138 #endif 143 139 ! -
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90
r5038 r5620 41 41 REAL(wp), PUBLIC :: rt0 = 273.15_wp !: freezing point of fresh water [Kelvin] 42 42 #if defined key_lim3 43 REAL(wp), PUBLIC :: rt0_snow = 273.1 6_wp !: melting point of snow [Kelvin]44 REAL(wp), PUBLIC :: rt0_ice = 273.1 6_wp !: melting point of ice [Kelvin]43 REAL(wp), PUBLIC :: rt0_snow = 273.15_wp !: melting point of snow [Kelvin] 44 REAL(wp), PUBLIC :: rt0_ice = 273.15_wp !: melting point of ice [Kelvin] 45 45 #else 46 46 REAL(wp), PUBLIC :: rt0_snow = 273.15_wp !: melting point of snow [Kelvin] … … 51 51 REAL(wp), PUBLIC :: rcp !: ocean specific heat [J/Kelvin] 52 52 REAL(wp), PUBLIC :: r1_rcp !: = 1. / rcp [Kelvin/J] 53 REAL(wp), PUBLIC :: rau0_rcp !: = rau0 * rcp 53 54 REAL(wp), PUBLIC :: r1_rau0_rcp !: = 1. / ( rau0 * rcp ) 54 55 … … 82 83 REAL(wp), PUBLIC :: xlic = 300.33e+6_wp !: volumetric latent heat fusion of ice [J/m3] 83 84 REAL(wp), PUBLIC :: xsn = 2.8e+6_wp !: volumetric latent heat of sublimation of snow [J/m3] 85 #endif 86 #if defined key_lim3 87 REAL(wp), PUBLIC :: r1_rhoic !: 1 / rhoic 88 REAL(wp), PUBLIC :: r1_rhosn !: 1 / rhosn 84 89 #endif 85 90 !!---------------------------------------------------------------------- … … 166 171 lfus = xlsn / rhosn ! latent heat of fusion of fresh ice 167 172 #endif 168 173 #if defined key_lim3 174 r1_rhoic = 1._wp / rhoic 175 r1_rhosn = 1._wp / rhosn 176 #endif 169 177 IF(lwp) THEN 170 178 WRITE(numout,*)
Note: See TracChangeset
for help on using the changeset viewer.