Changeset 6491 for branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO
- Timestamp:
- 2016-04-21T18:15:17+02:00 (8 years ago)
- Location:
- branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 17 edited
- 1 copied
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r6487 r6491 39 39 USE zdfmxl ! mixed layer 40 40 USE dianam ! build name of file (routine) 41 USE zdftke ! vertical physics: one-equation scheme 41 42 USE zdfddm ! vertical physics: double diffusion 42 43 USE diahth ! thermocline diagnostics … … 242 243 CALL iom_put( "avt" , avt ) ! T vert. eddy diff. coef. 243 244 CALL iom_put( "avm" , avmu ) ! T vert. eddy visc. coef. 245 IF( lk_zdftke ) THEN 246 CALL iom_put( "tke" , en ) ! TKE budget: Turbulent Kinetic Energy 247 CALL iom_put( "tke_niw" , e_niw ) ! TKE budget: Near-inertial waves 248 ENDIF 244 249 CALL iom_put( "avs" , fsavs(:,:,:) ) ! S vert. eddy diff. coef. (useful only with key_zdfddm) 245 250 -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90
r6486 r6491 36 36 PUBLIC clo_bat ! routine called in domzgr module 37 37 38 INTEGER, PUBLIC, PARAMETER :: jpncs = 4!: number of closed sea38 INTEGER, PUBLIC, PARAMETER :: jpncs = 10 !: number of closed sea 39 39 INTEGER, PUBLIC, DIMENSION(jpncs) :: ncstt !: Type of closed sea 40 40 INTEGER, PUBLIC, DIMENSION(jpncs) :: ncsi1, ncsj1 !: south-west closed sea limits (i,j) … … 155 155 ncsi2(4) = 76 ; ncsj2(4) = 61 156 156 ncsir(4,1) = 84 ; ncsjr(4,1) = 59 157 ! ! ======================= 158 CASE ( 025 ) ! ORCA_R025 configuration159 ! ! ======================= 160 ncsnr(1) = 1 ; ncstt(1) = 0 ! Caspian + Aralsea161 ncsi1(1) = 1330 ; ncsj1(1) = 645162 ncsi2(1) = 1 400 ; ncsj2(1) = 795157 ! ! ================================ 158 CASE ( 025 ) ! ORCA_R025 extended configuration 159 ! ! ================================ 160 ncsnr(1) = 1 ; ncstt(1) = 0 ! Caspian sea 161 ncsi1(1) = 1330 ; ncsj1(1) = 831 162 ncsi2(1) = 1375 ; ncsj2(1) = 981 163 163 ncsir(1,1) = 1 ; ncsjr(1,1) = 1 164 164 ! 165 ncsnr(2) = 1 ; ncstt(2) = 0 ! A zov Sea166 ncsi1(2) = 1 284 ; ncsj1(2) = 722167 ncsi2(2) = 1 304 ; ncsj2(2) = 747165 ncsnr(2) = 1 ; ncstt(2) = 0 ! Aral sea 166 ncsi1(2) = 1376 ; ncsj1(2) = 900 167 ncsi2(2) = 1400 ; ncsj2(2) = 981 168 168 ncsir(2,1) = 1 ; ncsjr(2,1) = 1 169 ! 170 ncsnr(3) = 1 ; ncstt(3) = 0 ! Azov Sea 171 ncsi1(3) = 1284 ; ncsj1(3) = 908 172 ncsi2(3) = 1304 ; ncsj2(3) = 933 173 ncsir(3,1) = 1 ; ncsjr(3,1) = 1 174 ! 175 ncsnr(4) = 1 ; ncstt(4) = 0 ! Lake Superior 176 ncsi1(4) = 781 ; ncsj1(4) = 904 177 ncsi2(4) = 815 ; ncsj2(4) = 926 178 ncsir(4,1) = 1 ; ncsjr(4,1) = 1 179 ! 180 ncsnr(5) = 1 ; ncstt(5) = 0 ! Lake Michigan 181 ncsi1(5) = 795 ; ncsj1(5) = 871 182 ncsi2(5) = 813 ; ncsj2(5) = 905 183 ncsir(5,1) = 1 ; ncsjr(5,1) = 1 184 ! 185 ncsnr(6) = 1 ; ncstt(6) = 0 ! Lake Huron part 1 186 ncsi1(6) = 814 ; ncsj1(6) = 882 187 ncsi2(6) = 825 ; ncsj2(6) = 905 188 ncsir(6,1) = 1 ; ncsjr(6,1) = 1 189 ! 190 ncsnr(7) = 1 ; ncstt(7) = 0 ! Lake Huron part 2 191 ncsi1(7) = 826 ; ncsj1(7) = 889 192 ncsi2(7) = 833 ; ncsj2(7) = 905 193 ncsir(7,1) = 1 ; ncsjr(7,1) = 1 194 ! 195 ncsnr(8) = 1 ; ncstt(8) = 0 ! Lake Erie 196 ncsi1(8) = 816 ; ncsj1(8) = 871 197 ncsi2(8) = 837 ; ncsj2(8) = 881 198 ncsir(8,1) = 1 ; ncsjr(8,1) = 1 199 ! 200 ncsnr(9) = 1 ; ncstt(9) = 0 ! Lake Ontario 201 ncsi1(9) = 831 ; ncsj1(9) = 882 202 ncsi2(9) = 847 ; ncsj2(9) = 889 203 ncsir(9,1) = 1 ; ncsjr(9,1) = 1 204 ! 205 ncsnr(10) = 1 ; ncstt(10) = 0 ! Lake Victoria 206 ncsi1(10) = 1274 ; ncsj1(10) = 672 207 ncsi2(10) = 1289 ; ncsj2(10) = 687 208 ncsir(10,1) = 1 ; ncsjr(10,1) = 1 169 209 ! 170 210 END SELECT -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r6486 r6491 136 136 USE ioipsl 137 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 & nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , ln_rstdate, nn_rstctl, & 139 139 & nn_it000, nn_itend , nn_date0 , nn_leapy , nn_istate , nn_stock , & 140 140 & nn_write, ln_dimgnnn, ln_mskland , ln_cfmeta , ln_clobber, nn_chunksz, nn_euler … … 173 173 WRITE(numout,*) ' file prefix restart output cn_ocerst_out= ', cn_ocerst_out 174 174 WRITE(numout,*) ' restart output directory cn_ocerst_outdir= ', cn_ocerst_outdir 175 WRITE(numout,*) ' restart logical ln_rstart = ', ln_rstart 175 WRITE(numout,*) ' restart logical ln_rstart = ' , ln_rstart 176 WRITE(numout,*) ' datestamping of restarts ln_rstdate = ', ln_rstdate 176 177 WRITE(numout,*) ' start with forward time step nn_euler = ', nn_euler 177 178 WRITE(numout,*) ' control of time step nn_rstctl = ', nn_rstctl -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90
r6487 r6491 136 136 INTEGER :: isrow ! index for ORCA1 starting row 137 137 INTEGER , POINTER, DIMENSION(:,:) :: imsk 138 REAL(wp) :: zphi_drake_passage, zshlat_antarc 138 139 REAL(wp), POINTER, DIMENSION(:,:) :: zwf 139 140 !! … … 445 446 ENDIF 446 447 ! 448 IF( cp_cfg == "orca" .AND. jp_cfg == 025 .AND. rn_shlat == 0.0 ) THEN 449 ! ! ORCA_R025 configuration 450 ! ! Increased lateral friction on parts of Antarctic coastline 451 ! ! for increased stability 452 ! ! NB. This only works to do this here if we have free slip 453 ! ! generally, so fmask is zero at coast points. 454 IF(lwp) WRITE(numout,*) 455 IF(lwp) WRITE(numout,*) ' orca_r025: increase friction in following regions : ' 456 IF(lwp) WRITE(numout,*) ' whole Antarctic coastline: partial slip shlat=1 ' 457 458 zphi_drake_passage = -58.0_wp 459 zshlat_antarc = 1.0_wp 460 zwf(:,:) = fmask(:,:,1) 461 DO jj = 2, jpjm1 462 DO ji = fs_2, fs_jpim1 ! vector opt. 463 IF( gphif(ji,jj) .lt. zphi_drake_passage .and. fmask(ji,jj,1) == 0._wp ) THEN 464 fmask(ji,jj,:) = zshlat_antarc * MIN( 1._wp , MAX( zwf(ji+1,jj), zwf(ji,jj+1), & 465 & zwf(ji-1,jj), zwf(ji,jj-1) ) ) 466 ENDIF 467 END DO 468 END DO 469 END IF 470 ! 447 471 CALL lbc_lnk( fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask 448 472 -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r6486 r6491 594 594 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~ - interpolate scale factors and compute depths for next time step' 595 595 ENDIF 596 597 ! Write outputs 598 ! ============= 599 CALL iom_put( "e3t" , fse3t_n (:,:,:) ) 600 CALL iom_put( "e3u" , fse3u_n (:,:,:) ) 601 CALL iom_put( "e3v" , fse3v_n (:,:,:) ) 602 CALL iom_put( "e3w" , fse3w_n (:,:,:) ) 603 CALL iom_put( "tpt_dep" , fsde3w_n (:,:,:) ) 604 IF( iom_use("e3tdef") ) & 605 CALL iom_put( "e3tdef" , ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 606 596 607 ! 597 608 ! Time filter and swap of scale factors … … 666 677 END DO 667 678 668 ! Write outputs669 ! =============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 (:,:,:) )674 CALL iom_put( "tpt_dep" , fsde3w_n (:,:,:) )675 IF( iom_use("e3tdef") ) &676 CALL iom_put( "e3tdef" , ( ( fse3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )677 679 678 680 ! write restart file -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ICB/icbdia.F90
r6486 r6491 371 371 IF( .NOT. ln_bergdia ) RETURN !!gm useless iom will control whether it is output or not 372 372 ! 373 CALL iom_put( "berg_total_melt" , berg_grid%floating_melt(:,:) ) ! Total melt flux to ocean [kg/m2/s] 374 CALL iom_put( "berg_total_heat_flux" , berg_grid%calving_hflx(:,:) ) ! Total iceberg-ocean heat flux [W/m2] 373 375 CALL iom_put( "berg_melt" , berg_melt (:,:) ) ! Melt rate of icebergs [kg/m2/s] 374 376 CALL iom_put( "berg_buoy_melt" , buoy_melt (:,:) ) ! Buoyancy component of iceberg melt rate [kg/m2/s] -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ICB/icbrst.F90
r6486 r6491 12 12 !! - ! Currently needs a fixed processor 13 13 !! - ! layout between restarts 14 !! - ! 2015-11 Dave Storkey Convert icb_rst_read to use IOM so can 15 !! read single restart files 14 16 !!---------------------------------------------------------------------- 15 17 !!---------------------------------------------------------------------- … … 18 20 !!---------------------------------------------------------------------- 19 21 USE par_oce ! NEMO parameters 22 USE phycst ! for rday 20 23 USE dom_oce ! NEMO domain 21 24 USE in_out_manager ! NEMO IO routines 25 USE ioipsl, ONLY : ju2ymds ! for calendar 22 26 USE lib_mpp ! NEMO MPI library, lk_mpp in particular 23 27 USE netcdf ! netcdf routines for IO 28 USE iom 24 29 USE icb_oce ! define iceberg arrays 25 30 USE icbutl ! iceberg utility routines … … 57 62 INTEGER :: idim, ivar, iatt 58 63 INTEGER :: jn, iunlim_dim, ibergs_in_file 59 INTEGER :: iclass 60 INTEGER, DIMENSION(1) :: istrt, ilngth, idata 61 INTEGER, DIMENSION(2) :: istrt2, ilngth2 62 INTEGER, DIMENSION(nkounts) :: idata2 63 REAL(wp), DIMENSION(1) :: zdata ! need 1d array to read in with 64 ! start and count arrays 64 INTEGER :: ii,ij,iclass 65 REAL(wp), DIMENSION(nkounts) :: zdata 65 66 LOGICAL :: ll_found_restart 66 67 CHARACTER(len=256) :: cl_path … … 71 72 !!---------------------------------------------------------------------- 72 73 73 ! Find a restart file. Assume iceberg restarts in same directory as ocean restarts. 74 ! Find a restart file. Assume iceberg restarts in same directory as ocean restarts 75 ! and are called TRIM(cn_ocerst)//'_icebergs' 74 76 cl_path = TRIM(cn_ocerst_indir) 75 77 IF( cl_path(LEN_TRIM(cl_path):) /= '/' ) cl_path = TRIM(cl_path) // '/' 76 cl_filename = ' ' 77 IF ( lk_mpp ) THEN 78 cl_filename = ' ' 79 WRITE( cl_filename, '("restart_icebergs_",I4.4,".nc")' ) narea-1 80 INQUIRE( file=TRIM(cl_path)//TRIM(cl_filename), exist=ll_found_restart ) 81 ELSE 82 cl_filename = 'restart_icebergs.nc' 83 INQUIRE( file=TRIM(cl_path)//TRIM(cl_filename), exist=ll_found_restart ) 84 ENDIF 85 86 IF ( .NOT. ll_found_restart) THEN ! only do the following if a file was found 87 CALL ctl_stop('icebergs: no restart file found') 88 ENDIF 89 90 IF (nn_verbose_level >= 0 .AND. lwp) & 91 WRITE(numout,'(2a)') 'icebergs, read_restart_bergs: found restart file = ',TRIM(cl_path)//TRIM(cl_filename) 92 93 nret = NF90_OPEN(TRIM(cl_path)//TRIM(cl_filename), NF90_NOWRITE, ncid) 94 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart_bergs: nf_open failed') 95 96 nret = nf90_inquire(ncid, idim, ivar, iatt, iunlim_dim) 97 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart_bergs: nf_inquire failed') 98 99 IF( iunlim_dim .NE. -1) THEN 100 101 nret = nf90_inquire_dimension(ncid, iunlim_dim, cl_dname, ibergs_in_file) 102 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart_bergs: nf_inq_dimlen failed') 103 104 nret = NF90_INQ_VARID(ncid, 'number', numberid) 105 nret = NF90_INQ_VARID(ncid, 'mass_scaling', nscaling_id) 106 nret = NF90_INQ_VARID(ncid, 'xi', nxid) 107 nret = NF90_INQ_VARID(ncid, 'yj', nyid) 108 nret = NF90_INQ_VARID(ncid, 'lon', nlonid) 109 nret = NF90_INQ_VARID(ncid, 'lat', nlatid) 110 nret = NF90_INQ_VARID(ncid, 'uvel', nuvelid) 111 nret = NF90_INQ_VARID(ncid, 'vvel', nvvelid) 112 nret = NF90_INQ_VARID(ncid, 'mass', nmassid) 113 nret = NF90_INQ_VARID(ncid, 'thickness', nthicknessid) 114 nret = NF90_INQ_VARID(ncid, 'width', nwidthid) 115 nret = NF90_INQ_VARID(ncid, 'length', nlengthid) 116 nret = NF90_INQ_VARID(ncid, 'year', nyearid) 117 nret = NF90_INQ_VARID(ncid, 'day', ndayid) 118 nret = NF90_INQ_VARID(ncid, 'mass_of_bits', nmass_of_bits_id) 119 nret = NF90_INQ_VARID(ncid, 'heat_density', nheat_density_id) 120 121 ilngth(1) = 1 122 istrt2(1) = 1 123 ilngth2(1) = nkounts 124 ilngth2(2) = 1 125 DO jn=1, ibergs_in_file 126 127 istrt(1) = jn 128 istrt2(2) = jn 129 130 nret = NF90_GET_VAR(ncid, numberid, idata2, istrt2, ilngth2 ) 131 localberg%number(:) = idata2(:) 132 133 nret = NF90_GET_VAR(ncid, nscaling_id, zdata, istrt, ilngth ) 134 localberg%mass_scaling = zdata(1) 135 136 nret = NF90_GET_VAR(ncid, nlonid, zdata, istrt, ilngth) 137 localpt%lon = zdata(1) 138 nret = NF90_GET_VAR(ncid, nlatid, zdata, istrt, ilngth) 139 localpt%lat = zdata(1) 140 IF (nn_verbose_level >= 2 .AND. lwp) THEN 141 WRITE(numout,'(a,i5,a,2f10.4,a,i5)') 'icebergs, read_restart_bergs: berg ',jn,' is at ', & 142 localpt%lon,localpt%lat,' on PE ',narea-1 78 cl_filename = TRIM(cn_ocerst_in)//'_icebergs' 79 CALL iom_open( TRIM(cl_path)//cl_filename, ncid ) 80 81 IF( iom_file(ncid)%iduld .GE. 0) THEN 82 83 ibergs_in_file = iom_file(ncid)%lenuld 84 DO jn = 1,ibergs_in_file 85 86 ! iom_get treats the unlimited dimension as time. Here the unlimited dimension 87 ! is the iceberg index, but we can still use the ktime keyword to get the iceberg we want. 88 89 CALL iom_get( ncid, 'xi' ,localpt%xi , ktime=jn ) 90 CALL iom_get( ncid, 'yj' ,localpt%yj , ktime=jn ) 91 92 ii = INT( localpt%xi + 0.5 ) 93 ij = INT( localpt%yj + 0.5 ) 94 ! Only proceed if this iceberg is on the local processor (excluding halos). 95 IF ( ii .GE. nldi+nimpp-1 .AND. ii .LE. nlei+nimpp-1 .AND. & 96 & ij .GE. nldj+njmpp-1 .AND. ij .LE. nlej+njmpp-1 ) THEN 97 98 CALL iom_get( ncid, jpdom_unknown, 'number' , zdata(:) , ktime=jn, kstart=(/1/), kcount=(/nkounts/) ) 99 localberg%number(:) = INT(zdata(:)) 100 CALL iom_get( ncid, 'mass_scaling' , localberg%mass_scaling, ktime=jn ) 101 CALL iom_get( ncid, 'lon' , localpt%lon , ktime=jn ) 102 CALL iom_get( ncid, 'lat' , localpt%lat , ktime=jn ) 103 CALL iom_get( ncid, 'uvel' , localpt%uvel , ktime=jn ) 104 CALL iom_get( ncid, 'vvel' , localpt%vvel , ktime=jn ) 105 CALL iom_get( ncid, 'mass' , localpt%mass , ktime=jn ) 106 CALL iom_get( ncid, 'thickness' , localpt%thickness , ktime=jn ) 107 CALL iom_get( ncid, 'width' , localpt%width , ktime=jn ) 108 CALL iom_get( ncid, 'length' , localpt%length , ktime=jn ) 109 CALL iom_get( ncid, 'year' , zdata(1) , ktime=jn ) 110 localpt%year = INT(zdata(1)) 111 CALL iom_get( ncid, 'day' , localpt%day , ktime=jn ) 112 CALL iom_get( ncid, 'mass_of_bits' , localpt%mass_of_bits , ktime=jn ) 113 CALL iom_get( ncid, 'heat_density' , localpt%heat_density , ktime=jn ) 114 115 ! 116 CALL icb_utl_add( localberg, localpt ) 117 143 118 ENDIF 144 nret = NF90_GET_VAR(ncid, nxid, zdata, istrt, ilngth) 145 localpt%xi = zdata(1) 146 nret = NF90_GET_VAR(ncid, nyid, zdata, istrt, ilngth) 147 localpt%yj = zdata(1) 148 nret = NF90_GET_VAR(ncid, nuvelid, zdata, istrt, ilngth ) 149 localpt%uvel = zdata(1) 150 nret = NF90_GET_VAR(ncid, nvvelid, zdata, istrt, ilngth ) 151 localpt%vvel = zdata(1) 152 nret = NF90_GET_VAR(ncid, nmassid, zdata, istrt, ilngth ) 153 localpt%mass = zdata(1) 154 nret = NF90_GET_VAR(ncid, nthicknessid, zdata, istrt, ilngth ) 155 localpt%thickness = zdata(1) 156 nret = NF90_GET_VAR(ncid, nwidthid, zdata, istrt, ilngth ) 157 localpt%width = zdata(1) 158 nret = NF90_GET_VAR(ncid, nlengthid, zdata, istrt, ilngth ) 159 localpt%length = zdata(1) 160 nret = NF90_GET_VAR(ncid, nyearid, idata, istrt, ilngth ) 161 localpt%year = idata(1) 162 nret = NF90_GET_VAR(ncid, ndayid, zdata, istrt, ilngth ) 163 localpt%day = zdata(1) 164 nret = NF90_GET_VAR(ncid, nmass_of_bits_id, zdata, istrt, ilngth ) 165 localpt%mass_of_bits = zdata(1) 166 nret = NF90_GET_VAR(ncid, nheat_density_id, zdata, istrt, ilngth ) 167 localpt%heat_density = zdata(1) 168 ! 169 CALL icb_utl_add( localberg, localpt ) 119 170 120 END DO 171 ! 172 ENDIF 173 174 nret = NF90_INQ_DIMID( ncid, 'c', nc_dim ) 175 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart: nf_inq_dimid c failed') 176 177 nret = NF90_INQUIRE_DIMENSION( ncid, nc_dim, cl_dname, iclass ) 178 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart: nf_inquire_dimension failed') 179 180 nret = NF90_INQ_VARID(ncid, 'kount' , nkountid) 181 nret = NF90_INQ_VARID(ncid, 'calving' , ncalvid) 182 nret = NF90_INQ_VARID(ncid, 'calving_hflx', ncalvhid) 183 nret = NF90_INQ_VARID(ncid, 'stored_ice' , nsiceid) 184 nret = NF90_INQ_VARID(ncid, 'stored_heat' , nsheatid) 185 186 nstrt3(1) = 1 187 nstrt3(2) = 1 188 nlngth3(1) = jpi 189 nlngth3(2) = jpj 190 nlngth3(3) = 1 191 192 DO jn = 1, iclass 193 nstrt3(3) = jn 194 nret = NF90_GET_VAR( ncid, nsiceid , griddata, nstrt3, nlngth3 ) 195 berg_grid%stored_ice(:,:,jn) = griddata(:,:,1) 196 END DO 197 198 nret = NF90_GET_VAR( ncid, ncalvid , src_calving (:,:) ) 199 nret = NF90_GET_VAR( ncid, ncalvhid, src_calving_hflx (:,:) ) 200 nret = NF90_GET_VAR( ncid, nsheatid, berg_grid%stored_heat(:,:) ) 201 nret = NF90_GET_VAR( ncid, nkountid, idata2(:) ) 202 num_bergs(:) = idata2(:) 203 204 ! Finish up 205 nret = NF90_CLOSE(ncid) 206 IF (nret .ne. NF90_NOERR) CALL ctl_stop('icebergs, read_restart: nf_close failed') 121 122 ENDIF 123 124 ! Gridded variables 125 CALL iom_get( ncid, jpdom_autoglo, 'calving' , src_calving ) 126 CALL iom_get( ncid, jpdom_autoglo, 'calving_hflx', src_calving_hflx ) 127 CALL iom_get( ncid, jpdom_autoglo, 'stored_heat' , berg_grid%stored_heat ) 128 CALL iom_get( ncid, jpdom_autoglo_xy, 'stored_ice' , berg_grid%stored_ice, kstart=(/1,1,1/), kcount=(/1,1,nclasses/) ) 129 130 CALL iom_get( ncid, jpdom_unknown, 'kount' , zdata(:) ) 131 num_bergs(:) = INT(zdata(:)) 207 132 208 133 ! Sanity check … … 211 136 WRITE(numout,'(2(a,i5))') 'icebergs, read_restart_bergs: # bergs =',jn,' on PE',narea-1 212 137 IF( lk_mpp ) THEN 213 CALL mpp_sum(ibergs_in_file) 138 ! Only mpp_sum ibergs_in_file if we are reading from multiple restart files. 139 IF( INDEX(iom_file(ncid)%name,'icebergs.nc' ) .EQ. 0 ) CALL mpp_sum(ibergs_in_file) 214 140 CALL mpp_sum(jn) 215 141 ENDIF … … 217 143 & ' bergs in the restart file and', jn,' bergs have been read' 218 144 ! 145 ! Finish up 146 CALL iom_close( ncid ) 147 ! 219 148 IF( lwp .and. nn_verbose_level >= 0) WRITE(numout,'(a)') 'icebergs, read_restart_bergs: completed' 220 149 ! … … 231 160 INTEGER :: jn ! dummy loop index 232 161 INTEGER :: ix_dim, iy_dim, ik_dim, in_dim 233 CHARACTER(len=256) :: cl_path 234 CHARACTER(len=256) :: cl_filename 162 INTEGER :: iyear, imonth, iday 163 REAL (wp) :: zsec 164 CHARACTER(len=256) :: cl_path 165 CHARACTER(len=256) :: cl_filename 166 CHARACTER(LEN=20) :: clkt ! ocean time-step deine as a character 235 167 TYPE(iceberg), POINTER :: this 236 168 TYPE(point) , POINTER :: pt … … 240 172 cl_path = TRIM(cn_ocerst_outdir) 241 173 IF( cl_path(LEN_TRIM(cl_path):) /= '/' ) cl_path = TRIM(cl_path) // '/' 174 IF ( ln_rstdate ) THEN 175 CALL ju2ymds( fjulday + rdttra(1) / rday, iyear, imonth, iday, zsec ) 176 WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday 177 ELSE 178 IF( kt > 999999999 ) THEN ; WRITE(clkt, * ) kt 179 ELSE ; WRITE(clkt, '(i8.8)') kt 180 ENDIF 181 ENDIF 242 182 IF( lk_mpp ) THEN 243 WRITE(cl_filename,'(A,"_icebergs_", I8.8,"_restart_",I4.4,".nc")') TRIM(cexper), kt, narea-1183 WRITE(cl_filename,'(A,"_icebergs_",A,"_restart_",I4.4,".nc")') TRIM(cexper), TRIM(ADJUSTL(clkt)), narea-1 244 184 ELSE 245 WRITE(cl_filename,'(A,"_icebergs_", I8.8,"_restart.nc")') TRIM(cexper), kt185 WRITE(cl_filename,'(A,"_icebergs_",A,"_restart.nc")') TRIM(cexper), TRIM(ADJUSTL(clkt)) 246 186 ENDIF 247 187 IF (nn_verbose_level >= 0) WRITE(numout,'(2a)') 'icebergs, write_restart: creating ',TRIM(cl_path)//TRIM(cl_filename) -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ICB/icbthm.F90
r6486 r6491 18 18 USE dom_oce ! NEMO domain 19 19 USE in_out_manager ! NEMO IO routines, numout in particular 20 USE iom 20 21 USE lib_mpp ! NEMO MPI routines, ctl_stop in particular 21 22 USE phycst ! NEMO physical constants … … 160 161 zmelt = ( zdM - ( zdMbitsE - zdMbitsM ) ) * z1_dt ! kg/s 161 162 berg_grid%floating_melt(ii,ij) = berg_grid%floating_melt(ii,ij) + zmelt * z1_e1e2 ! kg/m2/s 162 zheat = zmelt * pt%heat_density ! kg/s x J/kg = J/s 163 ! zheat = zmelt * pt%heat_density ! kg/s x J/kg = J/s 164 zheat = zmelt * lfus !rma kg/s x J/kg (latent heat of fusion) = J/s 163 165 berg_grid%calving_hflx (ii,ij) = berg_grid%calving_hflx (ii,ij) + zheat * z1_e1e2 ! W/m2 164 166 CALL icb_dia_melt( ii, ij, zMnew, zheat, this%mass_scaling, & … … 208 210 IF(.NOT. ln_passive_mode ) THEN 209 211 emp (:,:) = emp (:,:) - berg_grid%floating_melt(:,:) 210 !! qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:) !!gm heat flux not yet properly coded ==>> need it, SOLVE that! 212 qns (:,:) = qns (:,:) - berg_grid%calving_hflx (:,:) 211 213 ENDIF 212 214 ! -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/in_out_manager.F90
r6486 r6491 30 30 CHARACTER(lc) :: cn_ocerst_outdir !: restart output directory 31 31 LOGICAL :: ln_rstart !: start from (F) rest or (T) a restart file 32 LOGICAL :: ln_rstdate !: datestamping of restarts 32 33 LOGICAL :: ln_rst_list !: output restarts at list of times (T) or by frequency (F) 33 34 INTEGER :: nn_no !: job number -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90
r6487 r6491 710 710 CHARACTER(LEN=256) :: clname ! file name 711 711 CHARACTER(LEN=1) :: clrankpv, cldmspc ! 712 LOGICAL :: ll_depth_spec ! T => if kstart, kcount present then *only* use values for 3rd spatial dimension. 712 713 !--------------------------------------------------------------------- 713 714 ! … … 722 723 IF( PRESENT(kcount) .AND. (.NOT. PRESENT(kstart)) ) CALL ctl_stop(trim(clinfo), 'kcount present needs kstart present') 723 724 IF( PRESENT(kstart) .AND. (.NOT. PRESENT(kcount)) ) CALL ctl_stop(trim(clinfo), 'kstart present needs kcount present') 724 IF( PRESENT(kstart) .AND. idom /= jpdom_unknown ) CALL ctl_stop(trim(clinfo), 'kstart present needs kdom = jpdom_unknown') 725 IF( PRESENT(kstart) .AND. idom /= jpdom_unknown .AND. idom /= jpdom_autoglo_xy ) & 726 & CALL ctl_stop(trim(clinfo), 'kstart present needs kdom = jpdom_unknown or kdom = jpdom_autoglo_xy') 725 727 726 728 luse_jattr = .false. … … 755 757 ! update idom definition... 756 758 ! Identify the domain in case of jpdom_auto(glo/dta) definition 759 IF( idom == jpdom_autoglo_xy ) THEN 760 ll_depth_spec = .TRUE. 761 idom = jpdom_autoglo 762 ELSE 763 ll_depth_spec = .FALSE. 764 ENDIF 757 765 IF( idom == jpdom_autoglo .OR. idom == jpdom_autodta ) THEN 758 766 IF( idom == jpdom_autoglo ) THEN ; idom = jpdom_global … … 808 816 istart(idmspc+1) = itime 809 817 810 IF( PRESENT(kstart)) THEN ; istart(1:idmspc) = kstart(1:idmspc) ; icnt(1:idmspc) = kcount(1:idmspc)818 IF( PRESENT(kstart) .AND. .NOT. ll_depth_spec ) THEN ; istart(1:idmspc) = kstart(1:idmspc) ; icnt(1:idmspc) = kcount(1:idmspc) 811 819 ELSE 812 IF( idom == jpdom_unknown ) THEN ; icnt(1:idmspc) = idimsz(1:idmspc)820 IF( idom == jpdom_unknown ) THEN ; icnt(1:idmspc) = idimsz(1:idmspc) 813 821 ELSE 814 822 IF( .NOT. PRESENT(pv_r1d) ) THEN ! not a 1D array … … 833 841 ENDIF 834 842 IF( PRESENT(pv_r3d) ) THEN 835 IF( idom == jpdom_data ) THEN ; icnt(3) = jpkdta 836 ELSE ; icnt(3) = jpk 843 IF( idom == jpdom_data ) THEN ; icnt(3) = jpkdta 844 ELSE IF( ll_depth_spec .AND. PRESENT(kstart) ) THEN ; istart(3) = kstart(3); icnt(3) = kcount(3) 845 ELSE ; icnt(3) = jpk 837 846 ENDIF 838 847 ENDIF -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/iom_def.F90
r6486 r6491 26 26 INTEGER, PARAMETER, PUBLIC :: jpdom_unknown = 7 !: No dimension checking 27 27 INTEGER, PARAMETER, PUBLIC :: jpdom_autoglo = 8 !: 28 INTEGER, PARAMETER, PUBLIC :: jpdom_autodta = 9 !: 28 INTEGER, PARAMETER, PUBLIC :: jpdom_autoglo_xy = 9 !: Automatically set horizontal dimensions only 29 INTEGER, PARAMETER, PUBLIC :: jpdom_autodta = 10 !: 29 30 30 31 INTEGER, PARAMETER, PUBLIC :: jpioipsl = 100 !: Use ioipsl (fliocom only) library … … 57 58 INTEGER :: nvars !: number of identified varibles in the file 58 59 INTEGER :: iduld !: id of the unlimited dimension 60 INTEGER :: lenuld !: length of the unlimited dimension (number of records in file) 59 61 INTEGER :: irec !: writing record position 60 62 CHARACTER(LEN=32) :: uldname !: name of the unlimited dimension -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/iom_nf90.F90
r6486 r6491 154 154 CALL iom_nf90_check(NF90_Inquire(if90id, unlimitedDimId = iom_file(kiomid)%iduld), clinfo) 155 155 IF ( iom_file(kiomid)%iduld .GE. 0 ) THEN 156 CALL iom_nf90_check(NF90_Inquire_Dimension(if90id, iom_file(kiomid)%iduld, & 157 & name = iom_file(kiomid)%uldname), clinfo) 156 CALL iom_nf90_check(NF90_Inquire_Dimension(if90id, iom_file(kiomid)%iduld, & 157 & name = iom_file(kiomid)%uldname, & 158 & len = iom_file(kiomid)%lenuld ), clinfo ) 158 159 ENDIF 159 160 IF(lwp) WRITE(numout,*) ' ---> '//TRIM(cdname)//' OK' -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/restart.F90
r6488 r6491 21 21 USE in_out_manager ! I/O manager 22 22 USE iom ! I/O module 23 USE ioipsl, ONLY : ju2ymds ! for calendar 23 24 USE eosbn2 ! equation of state (eos bn2 routine) 24 25 USE trdmxl_oce ! ocean active mixed layer tracers trends variables … … 55 56 !!---------------------------------------------------------------------- 56 57 INTEGER, INTENT(in) :: kt ! ocean time-step 58 INTEGER :: iyear, imonth, iday 59 REAL (wp) :: zsec 57 60 !! 58 61 CHARACTER(LEN=20) :: clkt ! ocean time-step deine as a character 59 62 CHARACTER(LEN=50) :: clname ! ocean output restart file name 60 CHARACTER( lc):: clpath ! full path to ocean output restart file63 CHARACTER(LEN=150) :: clpath ! full path to ocean output restart file 61 64 !!---------------------------------------------------------------------- 62 65 ! … … 82 85 IF( kt == nitrst - 1 .OR. nstock == 1 .OR. ( kt == nitend .AND. .NOT. lrst_oce ) ) THEN 83 86 IF( nitrst <= nitend .AND. nitrst > 0 ) THEN 84 ! beware of the format used to write kt (default is i8.8, that should be large enough...) 85 IF( nitrst > 999999999 ) THEN ; WRITE(clkt, * ) nitrst 86 ELSE ; WRITE(clkt, '(i8.8)') nitrst 87 IF ( ln_rstdate ) THEN 88 CALL ju2ymds( fjulday + rdttra(1) / rday, iyear, imonth, iday, zsec ) 89 WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday 90 ELSE 91 ! beware of the format used to write kt (default is i8.8, that should be large enough...) 92 IF( nitrst > 999999999 ) THEN ; WRITE(clkt, * ) nitrst 93 ELSE ; WRITE(clkt, '(i8.8)') nitrst 94 ENDIF 87 95 ENDIF 88 96 ! create the file -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90
r6487 r6491 18 18 USE phycst ! physical constants 19 19 USE iom ! I/O library 20 USE eosbn2 ! for zdf_mxl_zint 20 21 USE lib_mpp ! MPP library 21 22 USE wrk_nemo ! work arrays … … 33 34 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmlp !: mixed layer depth (rho=rho0+zdcrit) [m] 34 35 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmlpt !: mixed layer depth at t-points [m] 36 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmld_zint !: vertically-interpolated mixed layer depth [m] 37 LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ll_found ! Is T_b to be found by interpolation ? 38 LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ll_belowml ! Flag points below mixed layer when ll_found=F 35 39 36 40 REAL(wp), PUBLIC :: rho_c = 0.01_wp !: density criterion for mixed layer depth 37 41 REAL(wp) :: avt_c = 5.e-4_wp ! Kz criterion for the turbocline depth 42 43 ! Namelist variables for namzdf_mldzint 44 INTEGER :: nn_mld_type ! mixed layer type 45 REAL(wp) :: rn_zref ! depth of initial T_ref 46 REAL(wp) :: rn_dT_crit ! Critical temp diff 47 REAL(wp) :: rn_iso_frac ! Fraction of rn_dT_crit used 38 48 39 49 !! * Substitutions … … 52 62 zdf_mxl_alloc = 0 ! set to zero if no array to be allocated 53 63 IF( .NOT. ALLOCATED( nmln ) ) THEN 54 ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), STAT= zdf_mxl_alloc ) 64 ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), hmld_zint(jpi,jpj), & 65 & ll_found(jpi,jpj), ll_belowml(jpi,jpj,jpk), STAT= zdf_mxl_alloc ) 55 66 ! 56 67 IF( lk_mpp ) CALL mpp_sum ( zdf_mxl_alloc ) … … 90 101 CALL wrk_alloc( jpi,jpj, imld ) 91 102 92 IF( kt == nit000 ) THEN103 IF( kt <= nit000 ) THEN 93 104 IF(lwp) WRITE(numout,*) 94 105 IF(lwp) WRITE(numout,*) 'zdf_mxl : mixed layer depth' … … 134 145 END DO 135 146 IF( .NOT.lk_offline ) THEN ! no need to output in offline mode 147 IF( kt >= nit000 ) THEN ! workaround for calls before SOMETHING reads the XIOS namelist 136 148 CALL iom_put( "mldr10_1", hmlp ) ! mixed layer depth 137 149 CALL iom_put( "mldkz5" , hmld ) ! turbocline depth 138 150 ENDIF 151 ENDIF 139 152 153 ! Vertically-interpolated mixed-layer depth diagnostic 154 IF( iom_use( "mldzint" ) ) THEN 155 CALL zdf_mxl_zint( kt ) 156 CALL iom_put( "mldzint" , hmld_zint ) 157 ENDIF 158 140 159 IF(ln_ctl) CALL prt_ctl( tab2d_1=REAL(nmln,wp), clinfo1=' nmln : ', tab2d_2=hmlp, clinfo2=' hmlp : ', ovlap=1 ) 141 160 ! … … 145 164 ! 146 165 END SUBROUTINE zdf_mxl 166 167 SUBROUTINE zdf_mxl_zint( kt ) 168 !!---------------------------------------------------------------------------------- 169 !! *** ROUTINE zdf_mxl_zint *** 170 ! 171 ! Calculate vertically-interpolated mixed layer depth diagnostic. 172 ! 173 ! This routine can calculate the mixed layer depth diagnostic suggested by 174 ! Kara et al, 2000, JGR, 105, 16803, but is more general and can calculate 175 ! vertically-interpolated mixed-layer depth diagnostics with other parameter 176 ! settings set in the namzdf_mldzint namelist. 177 ! 178 ! If mld_type=1 the mixed layer depth is calculated as the depth at which the 179 ! density has increased by an amount equivalent to a temperature difference of 180 ! 0.8C at the surface. 181 ! 182 ! For other values of mld_type the mixed layer is calculated as the depth at 183 ! which the temperature differs by 0.8C from the surface temperature. 184 ! 185 ! David Acreman, Daley Calvert 186 ! 187 !!----------------------------------------------------------------------------------- 188 189 INTEGER, INTENT(in) :: kt ! ocean time-step index 190 ! 191 ! Local variables 192 INTEGER, POINTER, DIMENSION(:,:) :: ikmt ! number of active tracer levels 193 INTEGER, POINTER, DIMENSION(:,:) :: ik_ref ! index of reference level 194 INTEGER, POINTER, DIMENSION(:,:) :: ik_iso ! index of last uniform temp level 195 REAL, POINTER, DIMENSION(:,:,:) :: zT ! Temperature or density 196 REAL, POINTER, DIMENSION(:,:) :: ppzdep ! depth for use in calculating d(rho) 197 REAL, POINTER, DIMENSION(:,:) :: zT_ref ! reference temperature 198 REAL :: zT_b ! base temperature 199 REAL, POINTER, DIMENSION(:,:,:) :: zdTdz ! gradient of zT 200 REAL, POINTER, DIMENSION(:,:,:) :: zmoddT ! Absolute temperature difference 201 REAL :: zdz ! depth difference 202 REAL :: zdT ! temperature difference 203 REAL, POINTER, DIMENSION(:,:) :: zdelta_T ! difference critereon 204 REAL, POINTER, DIMENSION(:,:) :: zRHO1, zRHO2 ! Densities 205 INTEGER :: ji, jj, jk ! loop counter 206 INTEGER :: ios 207 208 NAMELIST/namzdf_mldzint/ nn_mld_type, rn_zref, rn_dT_crit, rn_iso_frac 209 210 !!------------------------------------------------------------------------------------- 211 ! 212 CALL wrk_alloc( jpi, jpj, ikmt, ik_ref, ik_iso) 213 CALL wrk_alloc( jpi, jpj, ppzdep, zT_ref, zdelta_T, zRHO1, zRHO2 ) 214 CALL wrk_alloc( jpi, jpj, jpk, zT, zdTdz, zmoddT ) 215 216 IF( kt == nit000 ) THEN 217 REWIND( numnam_ref ) ! Namelist namzdf_mldzint in reference namelist 218 READ ( numnam_ref, namzdf_mldzint, IOSTAT = ios, ERR = 901) 219 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in reference namelist', lwp ) 220 221 REWIND( numnam_cfg ) ! Namelist namzdf_mldzint in configuration namelist 222 READ ( numnam_cfg, namzdf_mldzint, IOSTAT = ios, ERR = 902 ) 223 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in configuration namelist', lwp ) 224 IF(lwm) WRITE ( numond, namzdf_mldzint ) 225 226 WRITE(numout,*) '===== Vertically-interpolated mixed layer =====' 227 WRITE(numout,*) 'nn_mld_type = ',nn_mld_type 228 WRITE(numout,*) 'rn_zref = ',rn_zref 229 WRITE(numout,*) 'rn_dT_crit = ',rn_dT_crit 230 WRITE(numout,*) 'rn_iso_frac = ',rn_iso_frac 231 WRITE(numout,*) '===============================================' 232 ENDIF 233 234 ! Set the mixed layer depth criterion at each grid point 235 IF (nn_mld_type == 1) THEN 236 ppzdep(:,:)=0.0 237 call eos ( tsn(:,:,1,:), ppzdep(:,:), zRHO1(:,:) ) 238 ! Use zT temporarily as a copy of tsn with rn_dT_crit added to SST 239 ! [assumes number of tracers less than number of vertical levels] 240 zT(:,:,1:jpts)=tsn(:,:,1,1:jpts) 241 zT(:,:,jp_tem)=zT(:,:,1)+rn_dT_crit 242 CALL eos( zT(:,:,1:jpts), ppzdep(:,:), zRHO2(:,:) ) 243 zdelta_T(:,:) = abs( zRHO1(:,:) - zRHO2(:,:) ) * rau0 244 ! RHO from eos (2d version) doesn't calculate north or east halo: 245 CALL lbc_lnk( zdelta_T, 'T', 1. ) 246 zT(:,:,:) = rhop(:,:,:) 247 ELSE 248 zdelta_T(:,:) = rn_dT_crit 249 zT(:,:,:) = tsn(:,:,:,jp_tem) 250 END IF 251 252 ! Calculate the gradient of zT and absolute difference for use later 253 DO jk = 1 ,jpk-2 254 zdTdz(:,:,jk) = ( zT(:,:,jk+1) - zT(:,:,jk) ) / fse3w(:,:,jk+1) 255 zmoddT(:,:,jk) = abs( zT(:,:,jk+1) - zT(:,:,jk) ) 256 END DO 257 258 ! Find density/temperature at the reference level (Kara et al use 10m). 259 ! ik_ref is the index of the box centre immediately above or at the reference level 260 ! Find rn_zref in the array of model level depths and find the ref 261 ! density/temperature by linear interpolation. 262 DO jk = jpkm1, 2, -1 263 WHERE ( fsdept(:,:,jk) > rn_zref ) 264 ik_ref(:,:) = jk - 1 265 zT_ref(:,:) = zT(:,:,jk-1) + zdTdz(:,:,jk-1) * ( rn_zref - fsdept(:,:,jk-1) ) 266 END WHERE 267 END DO 268 269 ! If the first grid box centre is below the reference level then use the 270 ! top model level to get zT_ref 271 WHERE ( fsdept(:,:,1) > rn_zref ) 272 zT_ref = zT(:,:,1) 273 ik_ref = 1 274 END WHERE 275 276 ! The number of active tracer levels is 1 less than the number of active w levels 277 ikmt(:,:) = mbathy(:,:) - 1 278 279 ! Search for a uniform density/temperature region where adjacent levels 280 ! differ by less than rn_iso_frac * deltaT. 281 ! ik_iso is the index of the last level in the uniform layer 282 ! ll_found indicates whether the mixed layer depth can be found by interpolation 283 ik_iso(:,:) = ik_ref(:,:) 284 ll_found(:,:) = .false. 285 DO jj = 1, nlcj 286 DO ji = 1, nlci 287 !CDIR NOVECTOR 288 DO jk = ik_ref(ji,jj), ikmt(ji,jj)-1 289 IF ( zmoddT(ji,jj,jk) > ( rn_iso_frac * zdelta_T(ji,jj) ) ) THEN 290 ik_iso(ji,jj) = jk 291 ll_found(ji,jj) = ( zmoddT(ji,jj,jk) > zdelta_T(ji,jj) ) 292 EXIT 293 END IF 294 END DO 295 END DO 296 END DO 297 298 ! Use linear interpolation to find depth of mixed layer base where possible 299 hmld_zint(:,:) = rn_zref 300 DO jj = 1, jpj 301 DO ji = 1, jpi 302 IF (ll_found(ji,jj) .and. tmask(ji,jj,1) == 1.0) THEN 303 zdz = abs( zdelta_T(ji,jj) / zdTdz(ji,jj,ik_iso(ji,jj)) ) 304 hmld_zint(ji,jj) = fsdept(ji,jj,ik_iso(ji,jj)) + zdz 305 END IF 306 END DO 307 END DO 308 309 ! If ll_found = .false. then calculate MLD using difference of zdelta_T 310 ! from the reference density/temperature 311 312 ! Prevent this section from working on land points 313 WHERE ( tmask(:,:,1) /= 1.0 ) 314 ll_found = .true. 315 END WHERE 316 317 DO jk=1, jpk 318 ll_belowml(:,:,jk) = abs( zT(:,:,jk) - zT_ref(:,:) ) >= zdelta_T(:,:) 319 END DO 320 321 ! Set default value where interpolation cannot be used (ll_found=false) 322 DO jj = 1, jpj 323 DO ji = 1, jpi 324 IF ( .not. ll_found(ji,jj) ) hmld_zint(ji,jj) = fsdept(ji,jj,ikmt(ji,jj)) 325 END DO 326 END DO 327 328 DO jj = 1, jpj 329 DO ji = 1, jpi 330 !CDIR NOVECTOR 331 DO jk = ik_ref(ji,jj)+1, ikmt(ji,jj) 332 IF ( ll_found(ji,jj) ) EXIT 333 IF ( ll_belowml(ji,jj,jk) ) THEN 334 zT_b = zT_ref(ji,jj) + zdelta_T(ji,jj) * SIGN(1.0, zdTdz(ji,jj,jk-1) ) 335 zdT = zT_b - zT(ji,jj,jk-1) 336 zdz = zdT / zdTdz(ji,jj,jk-1) 337 hmld_zint(ji,jj) = fsdept(ji,jj,jk-1) + zdz 338 EXIT 339 END IF 340 END DO 341 END DO 342 END DO 343 344 hmld_zint(:,:) = hmld_zint(:,:)*tmask(:,:,1) 345 ! 346 CALL wrk_dealloc( jpi, jpj, ikmt, ik_ref, ik_iso) 347 CALL wrk_dealloc( jpi, jpj, ppzdep, zT_ref, zdelta_T, zRHO1, zRHO2 ) 348 CALL wrk_dealloc( jpi,jpj, jpk, zT, zdTdz, zmoddT ) 349 ! 350 END SUBROUTINE zdf_mxl_zint 147 351 148 352 !!====================================================================== -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
r6487 r6491 83 83 INTEGER :: nn_htau ! type of tke profile of penetration (=0/1) 84 84 REAL(wp) :: rn_efr ! fraction of TKE surface value which penetrates in the ocean 85 REAL(wp) :: rn_c ! fraction of TKE added within the mixed layer by nn_etau 85 86 LOGICAL :: ln_lc ! Langmuir cells (LC) as a source term of TKE or not 86 87 REAL(wp) :: rn_lc ! coef to compute vertical velocity of Langmuir cells … … 88 89 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) 89 90 REAL(wp) :: rmxl_min ! minimum mixing length value (deduced from rn_ediff and rn_emin values) [m] 91 REAL(wp) :: rhtau ! coefficient to relate MLD to htau when nn_htau == 2 90 92 REAL(wp) :: rhftau_add = 1.e-3_wp ! add offset applied to HF part of taum (nn_etau=3) 91 93 REAL(wp) :: rhftau_scl = 1.0_wp ! scale factor applied to HF part of taum (nn_etau=3) 92 94 93 95 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: htau ! depth of tke penetration (nn_htau) 96 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e_niw !: TKE budget- near-inertial waves term 97 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: efr ! surface boundary condition for nn_etau = 4 94 98 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dissl ! now mixing lenght of dissipation 95 99 #if defined key_c1d … … 114 118 !!---------------------------------------------------------------------- 115 119 ALLOCATE( & 120 & efr (jpi,jpj) , e_niw(jpi,jpj,jpk) , & 116 121 #if defined key_c1d 117 122 & e_dis(jpi,jpj,jpk) , e_mix(jpi,jpj,jpk) , & … … 420 425 END DO 421 426 427 ! ! Save TKE prior to nn_etau addition 428 e_niw(:,:,:) = en(:,:,:) 429 ! 422 430 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 423 431 ! ! TKE due to surface and internal wave breaking 424 432 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 433 IF( nn_htau == 2 ) THEN !* mixed-layer depth dependant length scale 434 DO jj = 2, jpjm1 435 DO ji = fs_2, fs_jpim1 ! vector opt. 436 htau(ji,jj) = rhtau * hmlp(ji,jj) 437 END DO 438 END DO 439 ENDIF 440 #if defined key_iomput 441 ! 442 CALL iom_put( "htau", htau(:,:) ) ! Check htau (even if constant in time) 443 #endif 444 ! 425 445 IF( nn_etau == 1 ) THEN !* penetration below the mixed layer (rn_efr fraction) 426 446 DO jk = 2, jpkm1 … … 457 477 END DO 458 478 END DO 479 ELSEIF( nn_etau == 4 ) THEN !* column integral independant of htau (rn_efr must be scaled up) 480 IF( nn_htau == 2 ) THEN ! efr dependant on time-varying htau 481 DO jj = 2, jpjm1 482 DO ji = fs_2, fs_jpim1 ! vector opt. 483 efr(ji,jj) = rn_efr / ( htau(ji,jj) * ( 1._wp - EXP( -bathy(ji,jj) / htau(ji,jj) ) ) ) 484 END DO 485 END DO 486 ENDIF 487 DO jk = 2, jpkm1 488 DO jj = 2, jpjm1 489 DO ji = fs_2, fs_jpim1 ! vector opt. 490 en(ji,jj,jk) = en(ji,jj,jk) + efr(ji,jj) * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 491 & * ( 1._wp - fr_i(ji,jj) ) * tmask(ji,jj,jk) 492 END DO 493 END DO 494 END DO 459 495 ENDIF 460 496 CALL lbc_lnk( en, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) 497 ! 498 DO jk = 2, jpkm1 ! TKE budget: near-inertial waves term 499 DO jj = 2, jpjm1 500 DO ji = fs_2, fs_jpim1 ! vector opt. 501 e_niw(ji,jj,jk) = en(ji,jj,jk) - e_niw(ji,jj,jk) 502 END DO 503 END DO 504 END DO 505 ! 506 CALL lbc_lnk( e_niw, 'W', 1. ) 461 507 ! 462 508 CALL wrk_dealloc( jpi,jpj, imlc ) ! integer … … 722 768 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 723 769 & rn_mxl0 , nn_pdl , ln_lc , rn_lc , & 724 & nn_etau , nn_htau , rn_efr 725 !!---------------------------------------------------------------------- 726 ! 770 & nn_etau , nn_htau , rn_efr , rn_c 771 !!---------------------------------------------------------------------- 772 727 773 REWIND( numnam_ref ) ! Namelist namzdf_tke in reference namelist : Turbulent Kinetic Energy 728 774 READ ( numnam_ref, namzdf_tke, IOSTAT = ios, ERR = 901) … … 757 803 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 758 804 WRITE(numout,*) ' fraction of en which pene. the thermocline rn_efr = ', rn_efr 805 WRITE(numout,*) ' fraction of TKE added within the mixed layer by nn_etau rn_c = ', rn_c 759 806 WRITE(numout,*) 760 807 WRITE(numout,*) ' critical Richardson nb with your parameters ri_cri = ', ri_cri … … 767 814 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 768 815 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 769 IF( nn_htau < 0 .OR. nn_htau > 1 ) CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2' )816 IF( nn_htau < 0 .OR. nn_htau > 5 ) CALL ctl_stop( 'bad flag: nn_htau is 0 to 5 ' ) 770 817 IF( nn_etau == 3 .AND. .NOT. ln_cpl ) CALL ctl_stop( 'nn_etau == 3 : HF taum only known in coupled mode' ) 771 818 … … 780 827 ENDIF 781 828 829 IF( nn_etau /= 0 .and. nn_htau == 2 ) THEN 830 ierr = zdf_mxl_alloc() 831 nmln(:,:) = nlb10 ! Initialization of nmln 832 ENDIF 833 782 834 ! !* depth of penetration of surface tke 783 835 IF( nn_etau /= 0 ) THEN 836 htau(:,:) = 0._wp 784 837 SELECT CASE( nn_htau ) ! Choice of the depth of penetration 785 838 CASE( 0 ) ! constant depth penetration (here 10 meters) … … 787 840 CASE( 1 ) ! F(latitude) : 0.5m to 30m poleward of 40 degrees 788 841 htau(:,:) = MAX( 0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) ) ) 842 CASE( 2 ) ! fraction of depth-integrated TKE within mixed-layer 843 rhtau = -1._wp / LOG( 1._wp - rn_c ) 844 CASE( 3 ) ! F(latitude) : 0.5m to 15m poleward of 20 degrees 845 htau(:,:) = MAX( 0.5_wp, MIN( 15._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) ) ) 846 CASE( 4 ) ! F(latitude) : 0.5m to 10m/30m poleward of 13/40 degrees north/south 847 DO jj = 2, jpjm1 848 DO ji = fs_2, fs_jpim1 ! vector opt. 849 IF( gphit(ji,jj) <= 0._wp ) THEN 850 htau(ji,jj) = MAX( 0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) ) ) 851 ELSE 852 htau(ji,jj) = MAX( 0.5_wp, MIN( 10._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) ) ) 853 ENDIF 854 END DO 855 END DO 856 CASE ( 5 ) ! F(latitude) : 0.5m to 10m poleward of 13 degrees north/south, 857 DO jj = 2, jpjm1 ! 10m to 30m between 30/45 degrees south 858 DO ji = fs_2, fs_jpim1 ! vector opt. 859 IF( gphit(ji,jj) <= -30._wp ) THEN 860 htau(ji,jj) = MAX( 10._wp, MIN( 30._wp, 55._wp* ABS( SIN( rpi/120._wp * ( gphit(ji,jj) + 23._wp ) ) ) ) ) 861 ELSE 862 htau(ji,jj) = MAX( 0.5_wp, MIN( 10._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) ) ) 863 ENDIF 864 END DO 865 END DO 789 866 END SELECT 867 ! 868 IF( nn_etau == 4 .AND. nn_htau /= 2 ) THEN ! efr dependant on constant htau 869 DO jj = 2, jpjm1 870 DO ji = fs_2, fs_jpim1 ! vector opt. 871 efr(ji,jj) = rn_efr / ( htau(ji,jj) * ( 1._wp - EXP( -bathy(ji,jj) / htau(ji,jj) ) ) ) 872 END DO 873 END DO 874 ENDIF 790 875 ENDIF 791 876 ! !* set vertical eddy coef. to the background value -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/step.F90
r6487 r6491 234 234 IF( lk_diaar5 ) CALL dia_ar5( kstp ) ! ar5 diag 235 235 IF( lk_diaharm ) CALL dia_harm( kstp ) ! Tidal harmonic analysis 236 CALL dia_prod( kstp ) ! ocean model: product diagnostics 236 237 CALL dia_wri( kstp ) ! ocean model: outputs 237 238 ! -
branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/step_oce.F90
r6487 r6491 88 88 89 89 USE diawri ! Standard run outputs (dia_wri routine) 90 USE diaprod ! Product diagnostics (dia_prod routine) 90 91 USE diaptr ! poleward transports (dia_ptr routine) 91 92 USE diadct ! sections transports (dia_dct routine)
Note: See TracChangeset
for help on using the changeset viewer.