Changeset 14995
- Timestamp:
- 2021-06-15T19:15:26+02:00 (3 years ago)
- Location:
- NEMO/trunk
- Files:
-
- 16 edited
-
cfgs/SHARED/field_def_nemo-oce.xml (modified) (2 diffs)
-
cfgs/SHARED/namelist_ref (modified) (3 diffs)
-
src/OCE/ISF/isf_oce.F90 (modified) (1 diff)
-
src/OCE/ISF/isfcav.F90 (modified) (3 diffs)
-
src/OCE/ISF/isfcavmlt.F90 (modified) (3 diffs)
-
src/OCE/ISF/isfpar.F90 (modified) (6 diffs)
-
src/OCE/ISF/isfparmlt.F90 (modified) (6 diffs)
-
src/OCE/ISF/isfstp.F90 (modified) (2 diffs)
-
src/OCE/SBC/sbccpl.F90 (modified) (2 diffs)
-
tests/ISOMIP+/MY_SRC/dtatsd.F90 (modified) (10 diffs)
-
tests/ISOMIP+/MY_SRC/eosbn2.F90 (modified) (29 diffs)
-
tests/ISOMIP+/MY_SRC/isf_oce.F90 (modified) (10 diffs)
-
tests/ISOMIP+/MY_SRC/isfstp.F90 (modified) (10 diffs)
-
tests/ISOMIP+/MY_SRC/istate.F90 (modified) (5 diffs)
-
tests/ISOMIP+/MY_SRC/sbcfwb.F90 (modified) (8 diffs)
-
tests/ISOMIP+/MY_SRC/tradmp.F90 (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/cfgs/SHARED/field_def_nemo-oce.xml
r14983 r14995 245 245 <field id="hflx_dmp_cea" long_name="heat flux due to damping" standard_name="heat_flux_due_to_damping" unit="W/m2" /> 246 246 247 <!-- * variable related to ice shelf forcing * --> 248 249 <!-- * fwf * --> 250 <field id="fwfisf_cav" long_name="Ice shelf fresh water flux ( from oce to isf )" unit="kg/m2/s" /> 251 <field id="fwfisf_par" long_name="Ice shelf fresh water flux ( from oce to isf )" unit="kg/m2/s" /> 252 <field id="fwfisf3d_cav" long_name="3d Ice shelf fresh water flux ( from oce to isf )" unit="kg/m2/s" grid_ref="grid_T_3D" /> 253 <field id="fwfisf3d_par" long_name="3d Ice shelf fresh water flux ( from oce to isf )" unit="kg/m2/s" grid_ref="grid_T_3D" /> 254 255 <!-- * heat fluxes * --> 256 <field id="qoceisf_cav" long_name="Ice shelf ocean heat flux ( from oce to isf )" unit="W/m2" /> 257 <field id="qoceisf_par" long_name="Ice shelf ocean heat flux ( from oce to isf )" unit="W/m2" /> 258 <field id="qlatisf_cav" long_name="Ice shelf latent heat flux ( from oce to isf )" unit="W/m2" /> 259 <field id="qlatisf_par" long_name="Ice shelf latent heat flux ( from oce to isf )" unit="W/m2" /> 260 <field id="qhcisf_cav" long_name="Ice shelf heat content flux of injected water ( from oce to isf )" unit="W/m2" /> 261 <field id="qhcisf_par" long_name="Ice shelf heat content flux of injected water ( from oce to isf )" unit="W/m2" /> 262 <field id="qoceisf3d_cav" long_name="Ice shelf ocean heat flux ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 263 <field id="qoceisf3d_par" long_name="Ice shelf ocean heat flux ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 264 <field id="qlatisf3d_cav" long_name="Ice shelf latent heat flux ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 265 <field id="qlatisf3d_par" long_name="Ice shelf latent heat flux ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 266 <field id="qhcisf3d_cav" long_name="Ice shelf heat content flux of injected water ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 267 <field id="qhcisf3d_par" long_name="Ice shelf heat content flux of injected water ( from oce to isf )" unit="W/m2" grid_ref="grid_T_3D" /> 268 <field id="qconisf" long_name="Conductive heat flux through the ice shelf ( from oce to isf )" unit="W/m2" /> 269 270 <!-- top boundary layer properties --> 271 <field id="isftfrz_cav" long_name="freezing point temperature at ocean/isf interface" unit="degC" /> 272 <field id="isftfrz_par" long_name="freezing point temperature in the parametrization boundary layer" unit="degC" /> 273 <field id="isfthermald_cav" long_name="thermal driving of ice shelf melting" unit="degC" /> 274 <field id="isfthermald_par" long_name="thermal driving of ice shelf melting" unit="degC" /> 275 <field id="isfgammat" long_name="Ice shelf heat-transfert velocity" unit="m/s" /> 276 <field id="isfgammas" long_name="Ice shelf salt-transfert velocity" unit="m/s" /> 277 <field id="ttbl_cav" long_name="temperature in Losch tbl" unit="degC" /> 278 <field id="ttbl_par" long_name="temperature in the parametrisation boundary layer" unit="degC" /> 279 <field id="stbl" long_name="salinity in the Losh tbl" unit="1e-3" /> 280 <field id="utbl" long_name="zonal current in the Losh tbl at T point" unit="m/s" /> 281 <field id="vtbl" long_name="merid current in the Losh tbl at T point" unit="m/s" /> 282 <field id="isfustar" long_name="ustar at T point used in ice shelf melting" unit="m/s" /> 283 247 284 </field_group> <!-- grid_T --> 248 285 … … 403 440 <!-- * variable relative to atmospheric pressure forcing : available with ln_apr_dyn --> 404 441 <field id="ssh_ib" long_name="Inverse barometer sea surface height" standard_name="sea_surface_height_correction_due_to_air_pressure_at_low_frequency" unit="m" /> 405 406 <!-- * variable related to ice shelf forcing * -->407 <field id="isftfrz_cav" long_name="freezing point temperature at ocean/isf interface" unit="degC" />408 <field id="isftfrz_par" long_name="freezing point temperature in the parametrization boundary layer" unit="degC" />409 <field id="fwfisf_cav" long_name="Ice shelf melt rate" unit="kg/m2/s" />410 <field id="fwfisf_par" long_name="Ice shelf melt rate" unit="kg/m2/s" />411 <field id="qoceisf_cav" long_name="Ice shelf ocean heat flux" unit="W/m2" />412 <field id="qoceisf_par" long_name="Ice shelf ocean heat flux" unit="W/m2" />413 <field id="qlatisf_cav" long_name="Ice shelf latent heat flux" unit="W/m2" />414 <field id="qlatisf_par" long_name="Ice shelf latent heat flux" unit="W/m2" />415 <field id="qhcisf_cav" long_name="Ice shelf heat content flux of injected water" unit="W/m2" />416 <field id="qhcisf_par" long_name="Ice shelf heat content flux of injected water" unit="W/m2" />417 <field id="fwfisf3d_cav" long_name="Ice shelf melt rate" unit="kg/m2/s" grid_ref="grid_T_3D" />418 <field id="fwfisf3d_par" long_name="Ice shelf melt rate" unit="kg/m2/s" grid_ref="grid_T_3D" />419 <field id="qoceisf3d_cav" long_name="Ice shelf ocean heat flux" unit="W/m2" grid_ref="grid_T_3D" />420 <field id="qoceisf3d_par" long_name="Ice shelf ocean heat flux" unit="W/m2" grid_ref="grid_T_3D" />421 <field id="qlatisf3d_cav" long_name="Ice shelf latent heat flux" unit="W/m2" grid_ref="grid_T_3D" />422 <field id="qlatisf3d_par" long_name="Ice shelf latent heat flux" unit="W/m2" grid_ref="grid_T_3D" />423 <field id="qhcisf3d_cav" long_name="Ice shelf heat content flux of injected water" unit="W/m2" grid_ref="grid_T_3D" />424 <field id="qhcisf3d_par" long_name="Ice shelf heat content flux of injected water" unit="W/m2" grid_ref="grid_T_3D" />425 <field id="ttbl_cav" long_name="temperature in Losch tbl" unit="degC" />426 <field id="ttbl_par" long_name="temperature in the parametrisation boundary layer" unit="degC" />427 <field id="isfthermald_cav" long_name="thermal driving of ice shelf melting" unit="degC" />428 <field id="isfthermald_par" long_name="thermal driving of ice shelf melting" unit="degC" />429 <field id="isfgammat" long_name="Ice shelf heat-transfert velocity" unit="m/s" />430 <field id="isfgammas" long_name="Ice shelf salt-transfert velocity" unit="m/s" />431 <field id="stbl" long_name="salinity in the Losh tbl" unit="1e-3" />432 <field id="utbl" long_name="zonal current in the Losh tbl at T point" unit="m/s" />433 <field id="vtbl" long_name="merid current in the Losh tbl at T point" unit="m/s" />434 <field id="isfustar" long_name="ustar at T point used in ice shelf melting" unit="m/s" />435 <field id="qconisf" long_name="Conductive heat flux through the ice shelf" unit="W/m2" />436 442 437 443 <!-- *_oce variables available with ln_blk_clio or ln_blk_core --> -
NEMO/trunk/cfgs/SHARED/namelist_ref
r14966 r14995 531 531 ln_isfcav_mlt = .false. ! ice shelf melting into the cavity (need ln_isfcav = .true. in domain_cfg.nc) 532 532 cn_isfcav_mlt = '3eq' ! ice shelf melting formulation (spe/2eq/3eq/oasis) 533 ! ! spe = fwfisf is read from a forcing field 533 ! ! spe = fwfisf is read from a forcing field ( melt > 0; freezing < 0 ) 534 534 ! ! 2eq = ISOMIP like: 2 equations formulation (Hunter et al., 2006 for a short description) 535 535 ! ! 3eq = ISOMIP+ like: 3 equations formulation (Asay-Davis et al., 2016 for a short description) … … 556 556 ln_isfpar_mlt = .false. ! ice shelf melting parametrised 557 557 cn_isfpar_mlt = 'spe' ! ice shelf melting parametrisation (spe/bg03/oasis) 558 ! ! spe = fwfisf is read from a forcing field 558 ! ! spe = fwfisf is read from a forcing field ( melt > 0; freezing < 0 ) 559 559 ! ! bg03 = melt computed using Beckmann and Goosse parametrisation 560 560 ! ! oasis = fwfisf is given by oasis and pattern by file sn_isfpar_fwf 561 ! 562 !* bg03 case 563 rn_isfpar_bg03_gt0 = 1.0e-4 ! gamma coeficient used in bg03 paper [m/s] 564 ! 565 !*** File definition *** 561 566 ! 562 567 !* all cases … … 566 571 sn_isfpar_zmax = 'isfmlt_par', 0. ,'sozisfmax', .false. , .true. , 'yearly' , '' , '' , '' 567 572 sn_isfpar_zmin = 'isfmlt_par', 0. ,'sozisfmin', .false. , .true. , 'yearly' , '' , '' , '' 573 ! 568 574 !* 'spe' and 'oasis' case 569 575 sn_isfpar_fwf = 'isfmlt_par' , -12. ,'sofwfisf' , .false. , .true. , 'yearly' , '' , '' , '' 576 ! 570 577 !* 'bg03' case 578 !* Leff is in [km] 571 579 sn_isfpar_Leff = 'isfmlt_par', 0. ,'Leff' , .false. , .true. , 'yearly' , '' , '' , '' 572 580 ! -
NEMO/trunk/src/OCE/ISF/isf_oce.F90
r14064 r14995 47 47 ! 48 48 ! 0.3 -------- ice shelf cavity parametrised namelist parameter ------------- 49 LOGICAL , PUBLIC :: ln_isfpar_mlt !: logical for the computation of melt inside the cavity 50 CHARACTER(LEN=256), PUBLIC :: cn_isfpar_mlt !: melt formulation (cavity/param) 51 TYPE(FLD_N) , PUBLIC :: sn_isfpar_fwf !: information about the isf melting file to be read 52 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmax !: information about the grounding line depth file to be read 53 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmin !: information about the calving line depth file to be read 54 TYPE(FLD_N) , PUBLIC :: sn_isfpar_Leff !: information about the effective length file to be read 49 LOGICAL , PUBLIC :: ln_isfpar_mlt !: logical for the computation of melt inside the cavity 50 REAL(wp) , PUBLIC :: rn_isfpar_bg03_gt0 !: temperature exchange coeficient [m/s] 51 CHARACTER(LEN=256), PUBLIC :: cn_isfpar_mlt !: melt formulation (cavity/param) 52 TYPE(FLD_N) , PUBLIC :: sn_isfpar_fwf !: information about the isf melting file to be read 53 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmax !: information about the grounding line depth file to be read 54 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmin !: information about the calving line depth file to be read 55 TYPE(FLD_N) , PUBLIC :: sn_isfpar_Leff !: information about the effective length file to be read 55 56 ! 56 57 ! 0.4 -------- coupling namelist parameter ------------- -
NEMO/trunk/src/OCE/ISF/isfcav.F90
r14433 r14995 59 59 !! - compute heat and fwf fluxes 60 60 !! - output 61 !! 62 !! ** Convention : all fluxes are from oce to isf ( > 0 out of the ocean ) 63 !! 61 64 !!--------------------------------------------------------------------- 62 65 !!-------------------------- OUT -------------------------------------- … … 132 135 zqlat(:,:) = - pqfwf(:,:) * rLfusisf ! 2d latent heat flux (W/m2) 133 136 ! 134 ! total heat flux ( > 0 out )137 ! total heat flux ( > 0 out ) 135 138 zqh(:,:) = ( zqhc (:,:) + zqoce(:,:) ) 136 139 ! … … 148 151 ! 149 152 IF ( ln_isfdebug ) THEN 153 IF(lwp) WRITE(numout,*) '' 150 154 CALL debug('isf_cav: ptsc T',ptsc(:,:,1)) 151 155 CALL debug('isf_cav: ptsc S',ptsc(:,:,2)) 152 156 CALL debug('isf_cav: pqfwf fwf',pqfwf(:,:)) 157 IF(lwp) WRITE(numout,*) '' 153 158 END IF 154 159 ! -
NEMO/trunk/src/OCE/ISF/isfcavmlt.F90
r13472 r14995 79 79 ! 80 80 IF (ln_isfdebug) THEN 81 IF(lwp) WRITE(numout,*) '' 81 82 CALL debug( 'isfcav_mlt qhc :', pqhc (:,:) ) 82 83 CALL debug( 'isfcav_mlt qoce :', pqoce(:,:) ) 83 84 CALL debug( 'isfcav_mlt qfwf :', pqfwf(:,:) ) 85 IF(lwp) WRITE(numout,*) '' 84 86 END IF 85 87 ! … … 122 124 ! 123 125 ! output freezing point at the interface 124 CALL iom_put('isftfrz_cav', ztfrz )126 CALL iom_put('isftfrz_cav', ztfrz(:,:) * mskisf_cav(:,:) ) 125 127 ! 126 128 END SUBROUTINE isfcav_mlt_spe … … 167 169 ! output thermal driving and freezinpoint at the ice shelf interface 168 170 CALL iom_put('isfthermald_cav', zthd ) 169 CALL iom_put('isftfrz_cav' , ztfrz )171 CALL iom_put('isftfrz_cav' , ztfrz(:,:) * mskisf_cav(:,:) ) 170 172 ! 171 173 END SUBROUTINE isfcav_mlt_2eq -
NEMO/trunk/src/OCE/ISF/isfpar.F90
r14433 r14995 9 9 !! 3.4 ! 2013-03 (P. Mathiot) Merging + parametrization 10 10 !! 4.1 ! 2019-09 (P. Mathiot) Restructuration 11 !! 4.2 ! 2021-05 (C. Ethe ) Test and fix oasis case 11 12 !!---------------------------------------------------------------------- 12 13 … … 56 57 !! the name tbl was kept. 57 58 !! 59 !! ** Convention : all fluxes are from oce to isf ( > 0 out of the ocean ) 60 !! 58 61 !!--------------------------------------------------------------------- 59 62 !!-------------------------- OUT -------------------------------------- … … 75 78 zqhc (:,:) = zqhc(:,:) * mskisf_par(:,:) 76 79 ! 77 ! compute heat content flux ( > 0 out )78 zqlat(:,:) = pqfwf(:,:) * rLfusisf ! 2d latent heat flux (W/m2)80 ! compute latent heat flux ( > 0 out ) 81 zqlat(:,:) = - pqfwf(:,:) * rLfusisf ! 2d latent heat flux (W/m2) 79 82 ! 80 83 ! total heat flux ( > 0 out ) … … 88 91 ! 89 92 ! set temperature content 90 ptsc(:,:,jp_tem) = zqh(:,:) * r1_rho0_rcp93 ptsc(:,:,jp_tem) = - zqh(:,:) * r1_rho0_rcp 91 94 ! 92 95 ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) … … 94 97 ! 95 98 IF ( ln_isfdebug ) THEN 99 IF(lwp) WRITE(numout,*) 96 100 CALL debug('isf_par: ptsc T',ptsc(:,:,1)) 97 101 CALL debug('isf_par: ptsc S',ptsc(:,:,2)) 98 102 CALL debug('isf_par: pqfwf fwf',pqfwf(:,:)) 103 IF(lwp) WRITE(numout,*) 99 104 END IF 100 105 ! … … 175 180 CASE ( 'oasis' ) 176 181 ! 182 ALLOCATE( sf_isfpar_fwf(1), STAT=ierr ) 183 ALLOCATE( sf_isfpar_fwf(1)%fnow(jpi,jpj,1), sf_isfpar_fwf(1)%fdta(jpi,jpj,1,2) ) 184 CALL fld_fill( sf_isfpar_fwf, (/ sn_isfpar_fwf /), cn_isfdir, 'isf_par_init', 'read fresh water flux isf data', 'namisf' ) 185 ! 177 186 IF(lwp) WRITE(numout,*) 178 187 IF(lwp) WRITE(numout,*) ' ==>>> isf melt provided by OASIS (cn_isfmlt_par = oasis)' -
NEMO/trunk/src/OCE/ISF/isfparmlt.F90
r12489 r14995 10 10 USE isf_oce ! ice shelf 11 11 USE isftbl , ONLY: isf_tbl ! ice shelf depth average 12 USE isfutils,ONLY: debug ! debug subroutine 12 13 13 14 USE dom_oce ! ocean space and time domain … … 39 40 ! ------------------------------------------------------------------------------------------------------- 40 41 41 SUBROUTINE isfpar_mlt( kt, Kmm, pqhc, pqoce, pqfwf )42 SUBROUTINE isfpar_mlt( kt, Kmm, pqhc, pqoce, pqfwf ) 42 43 !!--------------------------------------------------------------------- 43 44 !! *** ROUTINE isfpar_mlt *** … … 69 70 END SELECT 70 71 ! 72 IF (ln_isfdebug) THEN 73 IF(lwp) WRITE(numout,*) '' 74 CALL debug( 'isfpar_mlt qhc :', pqhc (:,:) ) 75 CALL debug( 'isfpar_mlt qoce :', pqoce(:,:) ) 76 CALL debug( 'isfpar_mlt qfwf :', pqfwf(:,:) ) 77 IF(lwp) WRITE(numout,*) '' 78 END IF 79 ! 71 80 END SUBROUTINE isfpar_mlt 72 81 … … 103 112 CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par ) 104 113 ! 105 pqfwf(:,:) = - sf_isfpar_fwf(1)%fnow(:,:,1) ! fresh water flux from the isf (fwfisf <0 mean melting) 106 pqoce(:,:) = pqfwf(:,:) * rLfusisf ! ocean/ice shelf flux assume to be equal to latent heat flux107 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux 108 ! 109 CALL iom_put('isftfrz_par', ztfrz )114 pqfwf(:,:) = - sf_isfpar_fwf(1)%fnow(:,:,1) ! fresh water flux from the isf (fwfisf <0 mean melting) ( > 0 out ) 115 pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocean/ice shelf flux assume to be equal to latent heat flux ( > 0 out ) 116 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( > 0 out ) 117 ! 118 CALL iom_put('isftfrz_par', ztfrz(:,:) * mskisf_par(:,:) ) 110 119 ! 111 120 END SUBROUTINE isfpar_mlt_spe … … 145 154 ! 146 155 ! 1. ------------Mean temperature 147 CALL isf_tbl(Kmm, ts(:,:, jk,jp_tem,Kmm), ztavg, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )156 CALL isf_tbl(Kmm, ts(:,:,:,jp_tem,Kmm), ztavg, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par ) 148 157 ! 149 158 ! 2. ------------Net heat flux and fresh water flux due to the ice shelf 150 pq oce(:,:) = rho0 * rcp * rn_gammat0 * risfLeff(:,:) * e1t(:,:) * ( ztavg(:,:) - ztfrz(:,:) ) * r1_e1e2t(:,:)151 pq fwf(:,:) = - pqoce(:,:) / rLfusisf ! derived from the latent heat flux152 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux 159 pqfwf(:,:) = - rho0 * rcp * rn_isfpar_bg03_gt0 * risfLeff(:,:) * e1t(:,:) * (ztavg(:,:) - ztfrz(:,:) ) * r1_e1e2t(:,:) / rLfusisf ! ( > 0 out ) 160 pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocean/ice shelf flux assume to be equal to latent heat flux ( > 0 out ) 161 pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( > 0 out ) 153 162 ! 154 163 ! 3. ------------BG03 output … … 157 166 ! 158 167 ! output thermal driving 159 CALL iom_put('isfthermald_par',( zt frz(:,:) - ztavg(:,:) ) * mskisf_par(:,:))168 CALL iom_put('isfthermald_par',( ztavg(:,:) - ztfrz(:,:) ) * mskisf_par(:,:)) 160 169 ! 161 170 ! output freezing point used to define the thermal driving and heat content fluxes 162 CALL iom_put('isftfrz_par', ztfrz )171 CALL iom_put('isftfrz_par', ztfrz(:,:) * mskisf_par(:,:) ) 163 172 ! 164 173 END SUBROUTINE isfpar_mlt_bg03 -
NEMO/trunk/src/OCE/ISF/isfstp.F90
r14143 r14995 262 262 IF ( l_isfoasis .AND. ln_isf ) THEN 263 263 ! 264 CALL ctl_stop( 'namelist combination ln_cpl and ln_isf not tested' )265 !266 264 ! NEMO coupled to ATMO model with isf cavity need oasis method for melt computation 267 265 IF ( ln_isfcav_mlt .AND. TRIM(cn_isfcav_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis is the only option availble if fwf send by oasis' ) 268 266 IF ( ln_isfpar_mlt .AND. TRIM(cn_isfpar_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis is the only option availble if fwf send by oasis' ) 269 !270 ! oasis melt computation not tested (coded but not tested)271 IF ( ln_isfcav_mlt .OR. ln_isfpar_mlt ) THEN272 IF ( TRIM(cn_isfcav_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis not tested' )273 IF ( TRIM(cn_isfpar_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis not tested' )274 END IF275 267 ! 276 268 ! oasis melt computation with cavity open and cavity parametrised (not coded) … … 301 293 & sn_isfpar_zmin, sn_isfpar_zmax, sn_isfpar_Leff, & 302 294 & ln_isfcpl , nn_drown , ln_isfcpl_cons, ln_isfdebug, & 303 & cn_isfload , rn_isfload_T , rn_isfload_S , cn_isfdir 295 & cn_isfload , rn_isfload_T , rn_isfload_S , cn_isfdir , & 296 & rn_isfpar_bg03_gt0 304 297 !!---------------------------------------------------------------------- 305 298 ! -
NEMO/trunk/src/OCE/SBC/sbccpl.F90
r14834 r14995 285 285 & sn_rcv_charn , sn_rcv_taw , sn_rcv_bhd , sn_rcv_tusd , sn_rcv_tvsd, & 286 286 & sn_rcv_wdrag , sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , & 287 & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_icb , sn_rcv_isf , sn_rcv_ts_ice !!, sn_rcv_qtrice 287 & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_icb , sn_rcv_isf , sn_rcv_ts_ice, & !!, sn_rcv_qtrice 288 & sn_rcv_mslp 288 289 289 290 !!--------------------------------------------------------------------- … … 528 529 IF(lwp) WRITE(numout,*) 529 530 IF(lwp) WRITE(numout,*) ' iceshelf received from oasis ' 530 CALL ctl_stop('STOP','not coded')531 ENDIF531 ENDIF 532 ! 532 533 ! 533 534 ! ! ------------------------- ! -
NEMO/trunk/tests/ISOMIP+/MY_SRC/dtatsd.F90
r14857 r14995 6 6 !! History : OPA ! 1991-03 () Original code 7 7 !! - ! 1992-07 (M. Imbard) 8 !! 8.0 ! 1999-10 (M.A. Foujols, M. Imbard) NetCDF FORMAT 9 !! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module 8 !! 8.0 ! 1999-10 (M.A. Foujols, M. Imbard) NetCDF FORMAT 9 !! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module 10 10 !! 3.3 ! 2010-10 (C. Bricaud, S. Masson) use of fldread 11 11 !! 3.4 ! 2010-11 (G. Madec, C. Ethe) Merge of dtatem and dtasal + remove CPP keys … … 49 49 !!---------------------------------------------------------------------- 50 50 !! *** ROUTINE dta_tsd_init *** 51 !! 52 !! ** Purpose : initialisation of T & S input data 53 !! 51 !! 52 !! ** Purpose : initialisation of T & S input data 53 !! 54 54 !! ** Method : - Read namtsd namelist 55 !! - allocates T & S data structure 55 !! - allocates T & S data structure 56 56 !!---------------------------------------------------------------------- 57 57 LOGICAL, INTENT(in), OPTIONAL :: ld_tradmp ! force the initialization when tradp is used … … 77 77 78 78 IF( PRESENT( ld_tradmp ) ) ln_tsd_dmp = .TRUE. ! forces the initialization when tradmp is used 79 79 80 80 IF(lwp) THEN ! control print 81 81 WRITE(numout,*) … … 114 114 CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' ) ; RETURN 115 115 ENDIF 116 !117 116 ! ! fill sf_tsd with sn_tem & sn_sal and control print 118 117 slf_i(jp_tem) = sn_tem ; slf_i(jp_sal) = sn_sal … … 150 149 !!---------------------------------------------------------------------- 151 150 !! *** ROUTINE dta_tsd *** 152 !! 151 !! 153 152 !! ** Purpose : provides T and S data at kt 154 !! 153 !! 155 154 !! ** Method : - call fldread routine 156 !! - ORCA_R2: add some hand made alteration to read data 155 !! - ORCA_R2: add some hand made alteration to read data 157 156 !! - 'key_orca_lev10' interpolates on 10 times more levels 158 157 !! - s- or mixed z-s coordinate: vertical interpolation on model mesh … … 162 161 !! ** Action : ptsd T-S data on medl mesh and interpolated at time-step kt 163 162 !!---------------------------------------------------------------------- 164 INTEGER 165 CHARACTER(LEN=3) 163 INTEGER , INTENT(in ) :: kt ! ocean time-step 164 CHARACTER(LEN=3) , INTENT(in ) :: cddta ! dmp or ini 166 165 REAL(wp), DIMENSION(A2D(nn_hls),jpk,jpts), INTENT( out) :: ptsd ! T & S data 167 166 ! … … 224 223 IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 225 224 zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 226 ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi 225 ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi 227 226 zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi 228 227 ENDIF … … 237 236 ptsd(ji,jj,jpk,jp_sal) = 0._wp 238 237 END_2D 239 ! 238 ! 240 239 ELSE !== z- or zps- coordinate ==! 241 ! 240 ! 242 241 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk ) 243 242 ptsd(ji,jj,jk,jp_tem) = ptsd(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) ! Mask … … 247 246 IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level 248 247 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 249 ik = mbkt(ji,jj) 248 ik = mbkt(ji,jj) 250 249 IF( ik > 1 ) THEN 251 250 zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) … … 255 254 ik = mikt(ji,jj) 256 255 IF( ik > 1 ) THEN 257 zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 256 zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 258 257 ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik+1,jp_tem) 259 258 ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik+1,jp_sal) -
NEMO/trunk/tests/ISOMIP+/MY_SRC/eosbn2.F90
r14857 r14995 31 31 !! bn2 : compute the Brunt-Vaisala frequency 32 32 !! eos_pt_from_ct: compute the potential temperature from the Conservative Temperature 33 !! eos_rab : generic interface of in situ thermal/haline expansion ratio 33 !! eos_rab : generic interface of in situ thermal/haline expansion ratio 34 34 !! eos_rab_3d : compute in situ thermal/haline expansion ratio 35 35 !! eos_rab_2d : compute in situ thermal/haline expansion ratio for 2d fields … … 46 46 USE in_out_manager ! I/O manager 47 47 USE lib_mpp ! MPP library 48 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 48 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 49 49 USE prtctl ! Print control 50 50 USE lbclnk ! ocean lateral boundary conditions … … 63 63 END INTERFACE 64 64 ! 65 INTERFACE eos_fzp 65 INTERFACE eos_fzp 66 66 MODULE PROCEDURE eos_fzp_2d, eos_fzp_0d 67 67 END INTERFACE … … 91 91 92 92 ! !!! simplified eos coefficients (default value: Vallis 2006) 93 REAL(wp), PUBLIC :: rn_a0 = 1.6550e-1_wp ! thermal expansion coeff. 94 REAL(wp), PUBLIC :: rn_b0 = 7.6554e-1_wp ! saline expansion coeff. 95 REAL(wp) :: rn_lambda1 = 5.9520e-2_wp ! cabbeling coeff. in T^2 96 REAL(wp) :: rn_lambda2 = 5.4914e-4_wp ! cabbeling coeff. in S^2 97 REAL(wp) :: rn_mu1 = 1.4970e-4_wp ! thermobaric coeff. in T 98 REAL(wp) :: rn_mu2 = 1.1090e-5_wp ! thermobaric coeff. in S 99 REAL(wp) :: rn_nu = 2.4341e-3_wp ! cabbeling coeff. in theta*salt 100 93 REAL(wp), PUBLIC :: rn_a0 = 1.6550e-1_wp ! thermal expansion coeff. 94 REAL(wp), PUBLIC :: rn_b0 = 7.6554e-1_wp ! saline expansion coeff. 95 REAL(wp) :: rn_lambda1 = 5.9520e-2_wp ! cabbeling coeff. in T^2 96 REAL(wp) :: rn_lambda2 = 5.4914e-4_wp ! cabbeling coeff. in S^2 97 REAL(wp) :: rn_mu1 = 1.4970e-4_wp ! thermobaric coeff. in T 98 REAL(wp) :: rn_mu2 = 1.1090e-5_wp ! thermobaric coeff. in S 99 REAL(wp) :: rn_nu = 2.4341e-3_wp ! cabbeling coeff. in theta*salt 100 101 101 ! TEOS10/EOS80 parameters 102 102 REAL(wp) :: r1_S0, r1_T0, r1_Z0, rdeltaS 103 103 104 104 ! EOS parameters 105 105 REAL(wp) :: EOS000 , EOS100 , EOS200 , EOS300 , EOS400 , EOS500 , EOS600 … … 119 119 REAL(wp) :: EOS022 120 120 REAL(wp) :: EOS003 , EOS103 121 REAL(wp) :: EOS013 122 121 REAL(wp) :: EOS013 122 123 123 ! ALPHA parameters 124 124 REAL(wp) :: ALP000 , ALP100 , ALP200 , ALP300 , ALP400 , ALP500 … … 135 135 REAL(wp) :: ALP012 136 136 REAL(wp) :: ALP003 137 137 138 138 ! BETA parameters 139 139 REAL(wp) :: BET000 , BET100 , BET200 , BET300 , BET400 , BET500 … … 162 162 REAL(wp) :: PEN002 , PEN102 163 163 REAL(wp) :: PEN012 164 164 165 165 ! ALPHA_PEN parameters 166 166 REAL(wp) :: APE000 , APE100 , APE200 , APE300 … … 241 241 INTEGER , INTENT(in ) :: ktts, ktrd, ktdep 242 242 REAL(wp), DIMENSION(A2D_T(ktts) ,JPK,JPTS), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 243 ! ! 2 : salinity [psu]243 ! ! 2 : salinity [psu] 244 244 REAL(wp), DIMENSION(A2D_T(ktrd) ,JPK ), INTENT( out) :: prd ! in situ density [-] 245 245 REAL(wp), DIMENSION(A2D_T(ktdep),JPK ), INTENT(in ) :: pdep ! depth [m] … … 301 301 & + rn_b0 * ( 1._wp - 0.5_wp*rn_lambda2*zs - rn_mu2*zh ) * zs & 302 302 & - rn_nu * zt * zs 303 ! 303 ! 304 304 prd(ji,jj,jk) = zn * r1_rho0 * ztm ! density anomaly (masked) 305 305 END_3D … … 354 354 INTEGER , INTENT(in ) :: ktts, ktrd, ktrhop, ktdep 355 355 REAL(wp), DIMENSION(A2D_T(ktts) ,JPK,JPTS), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 356 ! ! 2 : salinity [psu]356 ! ! 2 : salinity [psu] 357 357 REAL(wp), DIMENSION(A2D_T(ktrd) ,JPK ), INTENT( out) :: prd ! in situ density [-] 358 358 REAL(wp), DIMENSION(A2D_T(ktrhop),JPK ), INTENT( out) :: prhop ! potential density (surface referenced) … … 467 467 END_3D 468 468 ENDIF 469 469 470 470 CASE( np_seos ) !== simplified EOS ==! 471 471 ! … … 503 503 END SELECT 504 504 ! 505 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-pot: ', tab3d_2=prhop, clinfo2=' pot : ', kdim=jpk ) 505 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-pot: ', & 506 & tab3d_2=prhop, clinfo2=' pot : ', kdim=jpk ) 506 507 ! 507 508 IF( ln_timing ) CALL timing_stop('eos-pot') … … 534 535 INTEGER , INTENT(in ) :: ktts, ktdep, ktrd 535 536 REAL(wp), DIMENSION(A2D_T(ktts),JPTS), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 536 ! ! 2 : salinity [psu]537 ! ! 2 : salinity [psu] 537 538 REAL(wp), DIMENSION(A2D_T(ktdep) ), INTENT(in ) :: pdep ! depth [m] 538 539 REAL(wp), DIMENSION(A2D_T(ktrd) ), INTENT( out) :: prd ! in situ density … … 666 667 ! 667 668 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 669 ! 670 zt = pts (ji,jj,jp_tem) * r1_T0 ! temperature 671 zs = SQRT( ABS( pts(ji,jj,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 672 ztm = tmask(ji,jj,1) ! tmask 673 ! 674 zn0 = (((((EOS060*zt & 675 & + EOS150*zs+EOS050)*zt & 676 & + (EOS240*zs+EOS140)*zs+EOS040)*zt & 677 & + ((EOS330*zs+EOS230)*zs+EOS130)*zs+EOS030)*zt & 678 & + (((EOS420*zs+EOS320)*zs+EOS220)*zs+EOS120)*zs+EOS020)*zt & 679 & + ((((EOS510*zs+EOS410)*zs+EOS310)*zs+EOS210)*zs+EOS110)*zs+EOS010)*zt & 680 & + (((((EOS600*zs+EOS500)*zs+EOS400)*zs+EOS300)*zs+EOS200)*zs+EOS100)*zs+EOS000 681 ! 682 ! 683 prhop(ji,jj) = zn0 * ztm ! potential density referenced at the surface 684 ! 685 END_2D 685 686 686 687 CASE( np_seos ) !== simplified EOS ==! … … 713 714 ! 714 715 END SELECT 716 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=prhop, clinfo1=' pot: ', kdim=1 ) 715 717 ! 716 718 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=prhop, clinfo1=' eos-pot: ' ) … … 741 743 !! ** Action : - pab : thermal/haline expansion ratio at T-points 742 744 !!---------------------------------------------------------------------- 743 INTEGER , INTENT(in ) :: Kmm ! time level index745 INTEGER , INTENT(in ) :: Kmm ! time level index 744 746 INTEGER , INTENT(in ) :: ktts, ktab 745 747 REAL(wp), DIMENSION(A2D_T(ktts),JPK,JPTS), INTENT(in ) :: pts ! pot. temperature & salinity … … 873 875 !! ** Action : - pab : thermal/haline expansion ratio at T-points 874 876 !!---------------------------------------------------------------------- 875 INTEGER 877 INTEGER , INTENT(in ) :: Kmm ! time level index 876 878 INTEGER , INTENT(in ) :: ktts, ktdep, ktab 877 879 REAL(wp), DIMENSION(A2D_T(ktts),JPTS), INTENT(in ) :: pts ! pot. temperature & salinity … … 1113 1115 !! *** ROUTINE bn2 *** 1114 1116 !! 1115 !! ** Purpose : Compute the local Brunt-Vaisala frequency at the 1117 !! ** Purpose : Compute the local Brunt-Vaisala frequency at the 1116 1118 !! time-step of the input arguments 1117 1119 !! … … 1120 1122 !! N.B. N^2 is set one for all to zero at jk=1 in istate module. 1121 1123 !! 1122 !! ** Action : pn2 : square of the brunt-vaisala frequency at w-point 1123 !! 1124 !!---------------------------------------------------------------------- 1125 INTEGER , INTENT(in ) ::Kmm ! time level index1124 !! ** Action : pn2 : square of the brunt-vaisala frequency at w-point 1125 !! 1126 !!---------------------------------------------------------------------- 1127 INTEGER , INTENT(in ) :: Kmm ! time level index 1126 1128 INTEGER , INTENT(in ) :: ktab, ktn2 1127 REAL(wp), DIMENSION(jpi,jpj, jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu]1129 REAL(wp), DIMENSION(jpi,jpj, jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 1128 1130 REAL(wp), DIMENSION(A2D_T(ktab),JPK,JPTS), INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 1129 1131 REAL(wp), DIMENSION(A2D_T(ktn2),JPK ), INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] … … 1137 1139 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 2, jpkm1 ) ! interior points only (2=< jk =< jpkm1 ); surface and bottom value set to zero one for all in istate.F90 1138 1140 zrw = ( gdepw(ji,jj,jk ,Kmm) - gdept(ji,jj,jk,Kmm) ) & 1139 & / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) ) 1140 ! 1141 zaw = pab(ji,jj,jk,jp_tem) * (1. - zrw) + pab(ji,jj,jk-1,jp_tem) * zrw 1141 & / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) ) 1142 ! 1143 zaw = pab(ji,jj,jk,jp_tem) * (1. - zrw) + pab(ji,jj,jk-1,jp_tem) * zrw 1142 1144 zbw = pab(ji,jj,jk,jp_sal) * (1. - zrw) + pab(ji,jj,jk-1,jp_sal) * zrw 1143 1145 ! … … 1211 1213 1212 1214 1213 SUBROUTINE 1215 SUBROUTINE eos_fzp_2d( psal, ptf, pdep ) 1214 1216 !! 1215 1217 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: psal ! salinity [psu] … … 1234 1236 !!---------------------------------------------------------------------- 1235 1237 INTEGER , INTENT(in ) :: kttf 1236 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: psal ! salinity [psu]1237 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ), OPTIONAL :: pdep ! depth [m]1238 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: psal ! salinity [psu] 1239 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ), OPTIONAL :: pdep ! depth [m] 1238 1240 REAL(wp), DIMENSION(A2D_T(kttf)), INTENT(out ) :: ptf ! freezing temperature [Celsius] 1239 1241 ! … … 1267 1269 CALL ctl_stop( 'eos_fzp_2d:', ctmp1 ) 1268 1270 ! 1269 END SELECT 1271 END SELECT 1270 1272 ! 1271 1273 END SUBROUTINE eos_fzp_2d_t … … 1324 1326 !! ** Purpose : Calculates nonlinear anomalies of alpha_PE, beta_PE and PE at T-points 1325 1327 !! 1326 !! ** Method : PE is defined analytically as the vertical 1328 !! ** Method : PE is defined analytically as the vertical 1327 1329 !! primitive of EOS times -g integrated between 0 and z>0. 1328 1330 !! pen is the nonlinear bsq-PE anomaly: pen = ( PE - rho0 gz ) / rho0 gz - rd 1329 !! = 1/z * /int_0^z rd dz - rd 1331 !! = 1/z * /int_0^z rd dz - rd 1330 1332 !! where rd is the density anomaly (see eos_rhd function) 1331 1333 !! ab_pe are partial derivatives of PE anomaly with respect to T and S: … … 1391 1393 ! 1392 1394 zn = ( zn2 * zh + zn1 ) * zh + zn0 1393 ! 1395 ! 1394 1396 pab_pe(ji,jj,jk,jp_tem) = zn * zh * r1_rho0 * ztm 1395 1397 ! … … 1406 1408 ! 1407 1409 zn = ( zn2 * zh + zn1 ) * zh + zn0 1408 ! 1410 ! 1409 1411 pab_pe(ji,jj,jk,jp_sal) = zn / zs * zh * r1_rho0 * ztm 1410 1412 ! … … 1504 1506 IF(lwp) WRITE(numout,*) ' ==>>> use of TEOS-10 equation of state (cons. temp. and abs. salinity)' 1505 1507 ! 1506 l_useCT = .TRUE. ! model temperature is Conservative temperature 1508 l_useCT = .TRUE. ! model temperature is Conservative temperature 1507 1509 ! 1508 1510 rdeltaS = 32._wp … … 1885 1887 1886 1888 r1_S0 = 0.875_wp/35.16504_wp ! Used to convert CT in potential temperature when using bulk formulae (eos_pt_from_ct) 1887 1889 1888 1890 IF(lwp) THEN 1889 1891 WRITE(numout,*) … … 1923 1925 END SELECT 1924 1926 ! 1925 rho0_rcp = rho0 * rcp 1927 rho0_rcp = rho0 * rcp 1926 1928 r1_rho0 = 1._wp / rho0 1927 1929 r1_rcp = 1._wp / rcp 1928 r1_rho0_rcp = 1._wp / rho0_rcp 1930 r1_rho0_rcp = 1._wp / rho0_rcp 1929 1931 ! 1930 1932 IF(lwp) THEN -
NEMO/trunk/tests/ISOMIP+/MY_SRC/isf_oce.F90
r13583 r14995 1 1 MODULE isf_oce 2 2 !!====================================================================== 3 !! *** MODULE sbcisf***4 !! Surface module : compute iceshelf melt and heat flux3 !! *** MODULE isf_oce *** 4 !! Ice shelves : ice shelves variables defined in memory 5 5 !!====================================================================== 6 6 !! History : 3.2 ! 2011-02 (C.Harris ) Original code isf cav … … 48 48 ! 49 49 ! 0.3 -------- ice shelf cavity parametrised namelist parameter ------------- 50 LOGICAL , PUBLIC :: ln_isfpar_mlt !: logical for the computation of melt inside the cavity 51 CHARACTER(LEN=256), PUBLIC :: cn_isfpar_mlt !: melt formulation (cavity/param) 52 TYPE(FLD_N) , PUBLIC :: sn_isfpar_fwf !: information about the isf melting file to be read 53 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmax !: information about the grounding line depth file to be read 54 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmin !: information about the calving line depth file to be read 55 TYPE(FLD_N) , PUBLIC :: sn_isfpar_Leff !: information about the effective length file to be read 50 LOGICAL , PUBLIC :: ln_isfpar_mlt !: logical for the computation of melt inside the cavity 51 REAL(wp) , PUBLIC :: rn_isfpar_bg03_gt0 !: temperature exchange coeficient [m/s] 52 CHARACTER(LEN=256), PUBLIC :: cn_isfpar_mlt !: melt formulation (cavity/param) 53 TYPE(FLD_N) , PUBLIC :: sn_isfpar_fwf !: information about the isf melting file to be read 54 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmax !: information about the grounding line depth file to be read 55 TYPE(FLD_N) , PUBLIC :: sn_isfpar_zmin !: information about the calving line depth file to be read 56 TYPE(FLD_N) , PUBLIC :: sn_isfpar_Leff !: information about the effective length file to be read 56 57 ! 57 58 ! 0.4 -------- coupling namelist parameter ------------- … … 147 148 END SUBROUTINE isf_alloc_par 148 149 150 149 151 SUBROUTINE isf_alloc_cav() 150 152 !!--------------------------------------------------------------------- … … 174 176 END SUBROUTINE isf_alloc_cav 175 177 178 176 179 SUBROUTINE isf_alloc_cpl() 177 180 !!--------------------------------------------------------------------- … … 185 188 ierr = 0 186 189 ! 187 ALLOCATE( risfcpl_ssh(jpi,jpj) , risfcpl_tsc(jpi,jpj,jpk,jpts), risfcpl_vol(jpi,jpj,jpk), STAT=ialloc )188 ierr = ierr + ialloc 189 ! 190 risfcpl_tsc(:,:,:,:) = 0. 0 ; risfcpl_vol(:,:,:) = 0.0 ; risfcpl_ssh(:,:) = 0.0191 192 IF ( ln_isfcpl_cons ) THEN193 ALLOCATE( risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) , risfcpl_cons_ssh(jpi,jpj), STAT=ialloc )190 ALLOCATE( risfcpl_ssh(jpi,jpj) , risfcpl_tsc(jpi,jpj,jpk,jpts) , risfcpl_vol(jpi,jpj,jpk) , STAT=ialloc ) 191 ierr = ierr + ialloc 192 ! 193 risfcpl_tsc(:,:,:,:) = 0._wp ; risfcpl_vol(:,:,:) = 0._wp ; risfcpl_ssh(:,:) = 0._wp 194 195 IF ( ln_isfcpl_cons ) THEN 196 ALLOCATE( risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) , risfcpl_cons_ssh(jpi,jpj) , STAT=ialloc ) 194 197 ierr = ierr + ialloc 195 198 ! 196 risfcpl_cons_tsc(:,:,:,:) = 0. 0 ; risfcpl_cons_vol(:,:,:) = 0.0 ; risfcpl_cons_ssh(:,:) = 0.0199 risfcpl_cons_tsc(:,:,:,:) = 0._wp ; risfcpl_cons_vol(:,:,:) = 0._wp ; risfcpl_cons_ssh(:,:) = 0._wp 197 200 ! 198 201 END IF … … 203 206 END SUBROUTINE isf_alloc_cpl 204 207 208 205 209 SUBROUTINE isf_dealloc_cpl() 206 210 !!--------------------------------------------------------------------- … … 214 218 ierr = 0 215 219 ! 216 DEALLOCATE( risfcpl_ssh , risfcpl_tsc, risfcpl_vol, STAT=ialloc )220 DEALLOCATE( risfcpl_ssh , risfcpl_tsc , risfcpl_vol , STAT=ialloc ) 217 221 ierr = ierr + ialloc 218 222 ! … … 222 226 END SUBROUTINE isf_dealloc_cpl 223 227 228 224 229 SUBROUTINE isf_alloc() 225 230 !!--------------------------------------------------------------------- … … 234 239 ierr = 0 ! set to zero if no array to be allocated 235 240 ! 236 ALLOCATE( fwfisf_par(jpi,jpj) , fwfisf_par_b(jpi,jpj),&237 & fwfisf_cav(jpi,jpj) , fwfisf_cav_b(jpi,jpj),&238 & fwfisf_oasis(jpi,jpj),STAT=ialloc )239 ierr = ierr + ialloc 240 ! 241 ALLOCATE( risf_par_tsc(jpi,jpj,jpts), risf_par_tsc_b(jpi,jpj,jpts), STAT=ialloc )242 ierr = ierr + ialloc 243 ! 244 ALLOCATE( risf_cav_tsc(jpi,jpj,jpts), risf_cav_tsc_b(jpi,jpj,jpts), STAT=ialloc )245 ierr = ierr + ialloc 246 ! 247 ALLOCATE( risfload(jpi,jpj), STAT=ialloc)248 ierr = ierr + ialloc 249 ! 250 ALLOCATE( mskisf_cav(jpi,jpj) , STAT=ialloc)241 ALLOCATE( fwfisf_par (jpi,jpj) , fwfisf_par_b(jpi,jpj) , & 242 & fwfisf_cav (jpi,jpj) , fwfisf_cav_b(jpi,jpj) , & 243 & fwfisf_oasis(jpi,jpj) , STAT=ialloc ) 244 ierr = ierr + ialloc 245 ! 246 ALLOCATE( risf_par_tsc(jpi,jpj,jpts) , risf_par_tsc_b(jpi,jpj,jpts) , STAT=ialloc ) 247 ierr = ierr + ialloc 248 ! 249 ALLOCATE( risf_cav_tsc(jpi,jpj,jpts) , risf_cav_tsc_b(jpi,jpj,jpts) , STAT=ialloc ) 250 ierr = ierr + ialloc 251 ! 252 ALLOCATE( risfload(jpi,jpj) , STAT=ialloc ) 253 ierr = ierr + ialloc 254 ! 255 ALLOCATE( mskisf_cav(jpi,jpj) , STAT=ialloc ) 251 256 ierr = ierr + ialloc 252 257 ! … … 255 260 ! 256 261 ! initalisation of fwf and tsc array to 0 257 risfload(:,:) = 0.0_wp 258 fwfisf_oasis(:,:) = 0.0_wp 259 fwfisf_par(:,:) = 0.0_wp ; fwfisf_par_b(:,:) = 0.0_wp 260 fwfisf_cav(:,:) = 0.0_wp ; fwfisf_cav_b(:,:) = 0.0_wp 261 risf_cav_tsc(:,:,:) = 0.0_wp ; risf_cav_tsc_b(:,:,:) = 0.0_wp 262 risf_par_tsc(:,:,:) = 0.0_wp ; risf_par_tsc_b(:,:,:) = 0.0_wp 263 ! 264 262 risfload (:,:) = 0._wp 263 fwfisf_oasis(:,:) = 0._wp 264 fwfisf_par (:,:) = 0._wp ; fwfisf_par_b (:,:) = 0._wp 265 fwfisf_cav (:,:) = 0._wp ; fwfisf_cav_b (:,:) = 0._wp 266 risf_cav_tsc(:,:,:) = 0._wp ; risf_cav_tsc_b(:,:,:) = 0._wp 267 risf_par_tsc(:,:,:) = 0._wp ; risf_par_tsc_b(:,:,:) = 0._wp 268 ! 265 269 END SUBROUTINE isf_alloc 266 270 271 !!====================================================================== 267 272 END MODULE isf_oce -
NEMO/trunk/tests/ISOMIP+/MY_SRC/isfstp.F90
r13583 r14995 2 2 !!====================================================================== 3 3 !! *** MODULE isfstp *** 4 !! Surface module: compute iceshelf load, melt and heat flux4 !! Ice Shelves : compute iceshelf load, melt and heat flux 5 5 !!====================================================================== 6 6 !! History : 3.2 ! 2011-02 (C.Harris ) Original code isf cav … … 42 42 !! Software governed by the CeCILL license (see ./LICENSE) 43 43 !!---------------------------------------------------------------------- 44 45 44 CONTAINS 46 45 47 SUBROUTINE isf_stp( kt, Kmm )46 SUBROUTINE isf_stp( kt, Kmm ) 48 47 !!--------------------------------------------------------------------- 49 48 !! *** ROUTINE isf_stp *** … … 58 57 !! - compute fluxes 59 58 !! - write restart variables 60 !! 61 !!---------------------------------------------------------------------- 62 INTEGER, INTENT(in) :: kt ! ocean time step 63 INTEGER, INTENT(in) :: Kmm ! ocean time level index 64 !!---------------------------------------------------------------------- 65 INTEGER :: jk ! loop index 66 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t ! e3t 59 !!---------------------------------------------------------------------- 60 INTEGER, INTENT(in) :: kt ! ocean time step 61 INTEGER, INTENT(in) :: Kmm ! ocean time level index 62 ! 63 INTEGER :: jk ! loop index 64 #if defined key_qco 65 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t ! 3D workspace 66 #endif 67 67 !!--------------------------------------------------------------------- 68 68 ! … … 83 83 ! 1.2: compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 84 84 rhisf_tbl_cav(:,:) = rn_htbl * mskisf_cav(:,:) 85 #if defined key_qco 85 86 DO jk = 1, jpk 86 87 ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 87 88 END DO 88 CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav) 89 CALL isf_tbl_lvl( ht(:,:), ze3t , misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav ) 90 #else 91 CALL isf_tbl_lvl( ht(:,:), e3t(:,:,:,Kmm), misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav ) 92 #endif 89 93 ! 90 94 ! 1.3: compute ice shelf melt 91 CALL isf_cav( kt, Kmm, risf_cav_tsc, fwfisf_cav )95 CALL isf_cav( kt, Kmm, risf_cav_tsc, fwfisf_cav ) 92 96 ! 93 97 END IF … … 108 112 ! by simplicity, we assume the top level where param applied do not change with time (done in init part) 109 113 rhisf_tbl_par(:,:) = rhisf0_tbl_par(:,:) 114 #if defined key_qco 110 115 DO jk = 1, jpk 111 116 ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 112 117 END DO 113 CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par) 118 CALL isf_tbl_lvl( ht(:,:), ze3t , misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par ) 119 #else 120 CALL isf_tbl_lvl( ht(:,:), e3t(:,:,:,Kmm), misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par ) 121 #endif 114 122 ! 115 123 ! 2.3: compute ice shelf melt 116 CALL isf_par( kt, Kmm, risf_par_tsc, fwfisf_par )124 CALL isf_par( kt, Kmm, risf_par_tsc, fwfisf_par ) 117 125 ! 118 126 END IF … … 128 136 END SUBROUTINE isf_stp 129 137 130 SUBROUTINE isf_init(Kbb, Kmm, Kaa) 138 139 SUBROUTINE isf_init( Kbb, Kmm, Kaa ) 131 140 !!--------------------------------------------------------------------- 132 141 !! *** ROUTINE isfstp_init *** … … 142 151 !! - call cav/param/isfcpl init routine 143 152 !!---------------------------------------------------------------------- 144 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 153 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 154 !!---------------------------------------------------------------------- 145 155 ! 146 156 ! constrain: l_isfoasis need to be known 147 157 ! 148 ! Read namelist 149 CALL isf_nam() 150 ! 151 ! Allocate public array 152 CALL isf_alloc() 153 ! 154 ! check option compatibility 155 CALL isf_ctl() 156 ! 157 ! compute ice shelf load 158 IF ( ln_isfcav ) CALL isf_load( Kmm, risfload ) 158 CALL isf_nam() ! Read namelist 159 ! 160 CALL isf_alloc() ! Allocate public array 161 ! 162 CALL isf_ctl() ! check option compatibility 163 ! 164 IF( ln_isfcav ) CALL isf_load( Kmm, risfload ) ! compute ice shelf load 159 165 ! 160 166 ! terminate routine now if no ice shelf melt formulation specify 161 IF ( ln_isf ) THEN 162 ! 163 !--------------------------------------------------------------------------------------------------------------------- 164 ! initialisation melt in the cavity 165 IF ( ln_isfcav_mlt ) CALL isf_cav_init() 166 ! 167 !--------------------------------------------------------------------------------------------------------------------- 168 ! initialisation parametrised melt 169 IF ( ln_isfpar_mlt ) CALL isf_par_init() 170 ! 171 !--------------------------------------------------------------------------------------------------------------------- 172 ! initialisation ice sheet coupling 173 IF( ln_isfcpl ) CALL isfcpl_init(Kbb, Kmm, Kaa) 167 IF( ln_isf ) THEN 168 ! 169 IF( ln_isfcav_mlt ) CALL isf_cav_init() ! initialisation melt in the cavity 170 ! 171 IF( ln_isfpar_mlt ) CALL isf_par_init() ! initialisation parametrised melt 172 ! 173 IF( ln_isfcpl ) CALL isfcpl_init( Kbb, Kmm, Kaa ) ! initialisation ice sheet coupling 174 174 ! 175 175 END IF … … 177 177 END SUBROUTINE isf_init 178 178 179 179 180 SUBROUTINE isf_ctl() 180 181 !!--------------------------------------------------------------------- … … 283 284 END IF 284 285 END SUBROUTINE isf_ctl 285 ! 286 287 286 288 SUBROUTINE isf_nam 287 289 !!--------------------------------------------------------------------- … … 299 301 & sn_isfpar_zmin, sn_isfpar_zmax, sn_isfpar_Leff, & 300 302 & ln_isfcpl , nn_drown , ln_isfcpl_cons, ln_isfdebug, rn_vtide, & 301 & cn_isfload , rn_isfload_T , rn_isfload_S , cn_isfdir 303 & cn_isfload , rn_isfload_T , rn_isfload_S , cn_isfdir , & 304 & rn_isfpar_bg03_gt0 302 305 !!---------------------------------------------------------------------- 303 306 ! -
NEMO/trunk/tests/ISOMIP+/MY_SRC/istate.F90
r14857 r14995 34 34 USE lib_mpp ! MPP library 35 35 USE restart ! restart 36 36 37 #if defined key_agrif 38 USE agrif_oce ! initial state interpolation 37 39 USE agrif_oce_interp 38 USE agrif_oce39 40 #endif 40 41 … … 42 43 PRIVATE 43 44 44 PUBLIC istate_init ! routine called by step.F9045 PUBLIC istate_init ! routine called by nemogcm.F90 45 46 46 47 !! * Substitutions … … 63 64 ! 64 65 INTEGER :: ji, jj, jk ! dummy loop indices 65 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgdept ! 3D table !!st patch to use gdept subtitute66 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgdept ! 3D table for qco substitute 66 67 !!gm see comment further down 67 68 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zuvd ! U & V data workspace … … 73 74 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 74 75 75 !!gm Why not include in the first call of dta_tsd ?76 !!gm probably associated with the use of internal damping...77 76 CALL dta_tsd_init ! Initialisation of T & S input data 78 !!gm to be moved in usrdef of C1D case 77 79 78 ! IF( lk_c1d ) CALL dta_uvd_init ! Initialization of U & V input data 80 !!gm81 79 82 rhd (:,:,: ) = 0._wp ; rhop (:,:,: ) = 0._wp ! set one for all to 0 at level jpk83 rn2b (:,:,: ) = 0._wp ; rn2 (:,:,: ) = 0._wp ! set one for all to 0 at levels 1 and jpk84 ts (:,:,:,:,Kaa) = 0._wp ! set one for all to 0 at level jpk85 rab_b(:,:,:,: ) = 0._wp ; rab_n(:,:,:,:) = 0._wp ! set one for all to 0 at level jpk80 rhd (:,:,: ) = 0._wp ; rhop (:,:,: ) = 0._wp ! set one for all to 0 at level jpk 81 rn2b (:,:,: ) = 0._wp ; rn2 (:,:,: ) = 0._wp ! set one for all to 0 at levels 1 and jpk 82 ts (:,:,:,:,Kaa) = 0._wp ! set one for all to 0 at level jpk 83 rab_b(:,:,:,: ) = 0._wp ; rab_n(:,:,:,:) = 0._wp ! set one for all to 0 at level jpk 86 84 #if defined key_agrif 87 85 uu (:,:,: ,Kaa) = 0._wp ! used in agrif_oce_sponge at initialization … … 94 92 ln_1st_euler = .true. ! Set time-step indicator at nit000 (euler forward) 95 93 CALL day_init 96 CALL agrif_istate ( Kbb, Kmm, Kaa ) ! Interp from parent94 CALL agrif_istate_oce( Kbb, Kmm, Kaa ) ! Interp from parent 97 95 ! 98 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) 99 ssh (:,:,Kmm) = ssh(:,:,Kbb) 100 uu (:,:,:,Kmm) = uu (:,:,:,Kbb) 101 vv (:,:,:,Kmm) = vv (:,:,:,Kbb) 96 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) 97 uu (:,:,: ,Kmm) = uu (:,:,: ,Kbb) 98 vv (:,:,: ,Kmm) = vv (:,:,: ,Kbb) 102 99 ELSE 103 100 #endif 104 IF( ln_rstart ) THEN ! Restart from a file 105 ! ! ------------------- 106 CALL rst_read( Kbb, Kmm ) ! Read the restart file 107 CALL day_init ! model calendar (using both namelist and restart infos) 108 ! 109 ELSE ! Start from rest 110 ! ! --------------- 111 numror = 0 ! define numror = 0 -> no restart file to read 112 l_1st_euler = .true. ! Set time-step indicator at nit000 (euler forward) 113 CALL day_init ! model calendar (using both namelist and restart infos) 114 ! ! Initialization of ocean to zero 115 ! 116 IF( ln_tsd_init ) THEN 117 CALL dta_tsd( nit000, 'ini', ts(:,:,:,:,Kbb) ) ! read 3D T and S data at nit000 101 IF( ln_rstart ) THEN ! Restart from a file 102 ! ! ------------------- 103 CALL rst_read( Kbb, Kmm ) ! Read the restart file 104 CALL day_init ! model calendar (using both namelist and restart infos) 118 105 ! 119 uu (:,:,:,Kbb) = 0._wp 120 vv (:,:,:,Kbb) = 0._wp 106 ELSE ! Start from rest 107 ! ! --------------- 108 numror = 0 ! define numror = 0 -> no restart file to read 109 l_1st_euler = .true. ! Set time-step indicator at nit000 (euler forward) 110 CALL day_init ! model calendar (using both namelist and restart infos) 111 ! ! Initialization of ocean to zero 121 112 ! 122 ELSE ! user defined initial T and S123 DO jk = 1, jpk124 zgdept(:,:,jk) = gdept(:,:,jk,Kbb)125 END DO126 CALL usr_def_istate( zgdept, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb) )127 ENDIF128 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones129 uu (:,:,:,Kmm) = uu (:,:,:,Kbb)130 vv (:,:,:,Kmm) = vv (:,:,:,Kbb)131 132 !!gm POTENTIAL BUG : 133 !!gm ISSUE : if ssh(:,:,Kbb) /= 0 then, in non linear free surface, the e3._n, e3._b should be recomputed 134 !! as well as gdept_ and gdepw_.... !!!!! 135 !! ===>>>> probably a call to domvvl initialisation here.... 136 113 IF( ln_tsd_init ) THEN 114 CALL dta_tsd( nit000, 'ini', ts(:,:,:,:,Kbb) ) ! read 3D T and S data at nit000 115 ! 116 uu (:,:,:,Kbb) = 0._wp ! set the ocean at rest 117 vv (:,:,:,Kbb) = 0._wp 118 ! 119 ELSE ! user defined initial T and S 120 DO jk = 1, jpk 121 zgdept(:,:,jk) = gdept(:,:,jk,Kbb) 122 END DO 123 CALL usr_def_istate( zgdept, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb) ) 124 ENDIF 125 ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones 126 uu (:,:,:,Kmm) = uu (:,:,:,Kbb) 127 vv (:,:,:,Kmm) = vv (:,:,:,Kbb) 137 128 138 129 ! 139 !!gm to be moved in usrdef of C1D case140 !IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000141 !ALLOCATE( zuvd(jpi,jpj,jpk,2) )142 ! CALL dta_uvd( nit000, zuvd )143 ! uu(:,:,:,Kbb) = zuvd(:,:,:,1); uu(:,:,:,Kmm) = uu(:,:,:,Kbb)144 ! vv(:,:,:,Kbb) = zuvd(:,:,:,2); vv(:,:,:,Kmm) = vv(:,:,:,Kbb)145 !DEALLOCATE( zuvd )146 !ENDIF130 !!gm ==>>> to be moved in usrdef_istate of C1D case 131 IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000 132 ALLOCATE( zuvd(jpi,jpj,jpk,2) ) 133 CALL dta_uvd( nit000, Kbb, zuvd ) 134 uu(:,:,:,Kbb) = zuvd(:,:,:,1) ; uu(:,:,:,Kmm) = uu(:,:,:,Kbb) 135 vv(:,:,:,Kbb) = zuvd(:,:,:,2) ; vv(:,:,:,Kmm) = vv(:,:,:,Kbb) 136 DEALLOCATE( zuvd ) 137 ENDIF 147 138 ! 148 !!gm This is to be changed !!!!149 ! ! - ML - ssh(:,:,Kmm) could be modified by istate_eel, so that initialization of e3t(:,:,:,Kbb) is done here150 ! IF( .NOT.ln_linssh ) THEN151 ! DO jk = 1, jpk152 ! e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm)153 ! END DO154 ! ENDIF155 !!gm156 139 ! 157 ENDIF140 ENDIF 158 141 #if defined key_agrif 159 142 ENDIF -
NEMO/trunk/tests/ISOMIP+/MY_SRC/sbcfwb.F90
r13583 r14995 24 24 ! 25 25 USE in_out_manager ! I/O manager 26 USE iom ! IOM 26 27 USE lib_mpp ! distribued memory computing library 27 28 USE timing ! Timing … … 34 35 PUBLIC sbc_fwb ! routine called by step 35 36 36 REAL(wp) :: a_fwb_b ! annual domain averaged freshwater budget37 REAL(wp) :: a_fwb ! for 2 year before (_b) and before year.38 REAL(wp) :: fwfold ! fwfold to be suppressed37 REAL(wp) :: rn_fwb0 ! initial freshwater adjustment flux [kg/m2/s] (nn_fwb = 2 only) 38 REAL(wp) :: a_fwb ! annual domain averaged freshwater budget from the 39 ! previous year 39 40 REAL(wp) :: area ! global mean ocean surface (interior domain) 40 41 … … 65 66 INTEGER, INTENT( in ) :: Kmm ! ocean time level index 66 67 ! 67 INTEGER :: i num, ikty, iyear! local integers68 INTEGER :: ios, inum, ikty ! local integers 68 69 REAL(wp) :: z_fwf, z_fwf_nsrf, zsum_fwf, zsum_erp ! local scalars 69 70 REAL(wp) :: zsurf_neg, zsurf_pos, zsurf_tospread, zcoef ! - - … … 72 73 REAL(wp) ,DIMENSION(1) :: z_fwfprv 73 74 COMPLEX(dp),DIMENSION(1) :: y_fwfnow 75 ! 76 NAMELIST/namsbc_fwb/rn_fwb0 74 77 !!---------------------------------------------------------------------- 75 78 ! 76 79 IF( kt == nit000 ) THEN 80 READ( numnam_ref, namsbc_fwb, IOSTAT = ios, ERR = 901 ) 81 901 IF( ios /= 0 ) CALL ctl_nam( ios, 'namsbc_fwb in reference namelist' ) 82 READ( numnam_cfg, namsbc_fwb, IOSTAT = ios, ERR = 902 ) 83 902 IF( ios > 0 ) CALL ctl_nam( ios, 'namsbc_fwb in configuration namelist' ) 84 IF(lwm) WRITE( numond, namsbc_fwb ) 77 85 IF(lwp) THEN 78 86 WRITE(numout,*) … … 80 88 WRITE(numout,*) '~~~~~~~' 81 89 IF( kn_fwb == 1 ) WRITE(numout,*) ' instantaneously set to zero' 82 IF( kn_fwb == 2 ) WRITE(numout,*) ' adjusted from previous year budget'90 IF( kn_fwb == 4 ) WRITE(numout,*) ' instantaneously set to zero with heat and salt flux correction (ISOMIP+)' 83 91 IF( kn_fwb == 3 ) WRITE(numout,*) ' fwf set to zero and spread out over erp area' 84 IF( kn_fwb == 4 ) WRITE(numout,*) ' instantaneously set to zero with heat and salt flux correction (ISOMIP+)' 92 IF( kn_fwb == 2 ) THEN 93 WRITE(numout,*) ' adjusted from previous year budget' 94 WRITE(numout,*) 95 WRITE(numout,*) ' Namelist namsbc_fwb' 96 WRITE(numout,*) ' Initial freshwater adjustment flux [kg/m2/s] = ', rn_fwb0 97 END IF 85 98 ENDIF 86 99 ! … … 111 124 emp(:,:) = emp(:,:) - z_fwfprv(1) * tmask(:,:,1) 112 125 qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! account for change to the heat budget due to fw correction 126 ! outputs 127 IF( iom_use('hflx_fwb_cea') ) CALL iom_put( 'hflx_fwb_cea', zcoef * sst_m(:,:) * tmask(:,:,1) ) 128 IF( iom_use('vflx_fwb_cea') ) CALL iom_put( 'vflx_fwb_cea', z_fwfprv(1) * tmask(:,:,1) ) 113 129 ENDIF 114 130 ! … … 131 147 qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! (Eq. 35 AD2015) ! use sst_m to avoid generation of any bouyancy fluxes 132 148 sfx(:,:) = sfx(:,:) + z_fwf * sss_m(:,:) * tmask(:,:,1) ! (Eq. 36 AD2015) ! use sss_m to avoid generation of any bouyancy fluxes 133 !qns(:,:) = qns(:,:) + zcoef * ( -1.9 ) * tmask(:,:,1) ! (Eq. 35 AD2015) ! could be sst_m if we don't want any bouyancy fluxes 134 !sfx(:,:) = sfx(:,:) + z_fwf * ( 33.8 ) * tmask(:,:,1) ! (Eq. 36 AD2015) ! could be sss_m if we don't want any bouyancy fluxes 135 !qns(:,:) = qns(:,:) + zcoef * ( -1.0 ) * tmask(:,:,1) ! use for ISOMIP+ coupling sanity check (keep ssh cst while playing with cpl conservation option) 136 !sfx(:,:) = sfx(:,:) + z_fwf * ( 34.2 ) * tmask(:,:,1) ! use for ISOMIP+ coupling sanity check (keep ssh cst while playing with cpl conservation option) 137 ENDIF 138 ! 139 CASE ( 2 ) !== fwf budget adjusted from the previous year ==! 140 ! 141 IF( kt == nit000 ) THEN ! initialisation 142 ! ! Read the corrective factor on precipitations (fwfold) 143 CALL ctl_opn( inum, 'EMPave_old.dat', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 144 READ ( inum, "(24X,I8,2ES24.16)" ) iyear, a_fwb_b, a_fwb 145 CLOSE( inum ) 146 fwfold = a_fwb ! current year freshwater budget correction 147 ! ! estimate from the previous year budget 149 ENDIF 150 ! 151 CASE ( 2 ) !== fw adjustment based on fw budget at the end of the previous year ==! 152 ! 153 IF( kt == nit000 ) THEN ! initialisation 154 ! ! set the fw adjustment (a_fwb) 155 IF ( ln_rstart .AND. iom_varid( numror, 'a_fwb', ldstop = .FALSE. ) > 0 ) THEN ! as read from restart file 156 IF(lwp) WRITE(numout,*) 'sbc_fwb : reading FW-budget adjustment from restart file' 157 CALL iom_get( numror, 'a_fwb', a_fwb ) 158 ELSE ! as specified in namelist 159 a_fwb = rn_fwb0 160 END IF 161 ! 148 162 IF(lwp)WRITE(numout,*) 149 IF(lwp)WRITE(numout,*)'sbc_fwb : year = ',iyear , ' freshwater budget correction = ', fwfold 150 IF(lwp)WRITE(numout,*)' year = ',iyear-1, ' freshwater budget read = ', a_fwb 151 IF(lwp)WRITE(numout,*)' year = ',iyear-2, ' freshwater budget read = ', a_fwb_b 163 IF(lwp)WRITE(numout,*)'sbc_fwb : initial freshwater-budget adjustment = ', a_fwb, 'kg/m2/s' 164 ! 152 165 ENDIF 153 ! ! Update fwfoldif new year start166 ! ! Update a_fwb if new year start 154 167 ikty = 365 * 86400 / rn_Dt !!bug use of 365 days leap year or 360d year !!!!!!! 155 168 IF( MOD( kt, ikty ) == 0 ) THEN 156 a_fwb_b = a_fwb! mean sea level taking into account the ice+snow169 ! mean sea level taking into account the ice+snow 157 170 ! sum over the global domain 158 171 a_fwb = glob_sum( 'sbcfwb', e1e2t(:,:) * ( ssh(:,:,Kmm) + snwice_mass(:,:) * r1_rho0 ) ) 159 172 a_fwb = a_fwb * 1.e+3 / ( area * rday * 365. ) ! convert in Kg/m3/s = mm/s 160 173 !!gm ! !!bug 365d year 161 fwfold = a_fwb ! current year freshwater budget correction162 ! ! estimate from the previous year budget163 174 ENDIF 164 175 ! 165 176 IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN ! correct the freshwater fluxes 166 zcoef = fwfold* rcp167 emp(:,:) = emp(:,:) + fwfold* tmask(:,:,1)177 zcoef = a_fwb * rcp 178 emp(:,:) = emp(:,:) + a_fwb * tmask(:,:,1) 168 179 qns(:,:) = qns(:,:) - zcoef * sst_m(:,:) * tmask(:,:,1) ! account for change to the heat budget due to fw correction 169 ENDIF 170 ! 171 IF( kt == nitend .AND. lwm ) THEN ! save fwfold value in a file (only one required) 172 CALL ctl_opn( inum, 'EMPave.dat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea ) 173 WRITE( inum, "(24X,I8,2ES24.16)" ) nyear, a_fwb_b, a_fwb 174 CLOSE( inum ) 175 ENDIF 180 ! outputs 181 IF( iom_use('hflx_fwb_cea') ) CALL iom_put( 'hflx_fwb_cea', -zcoef * sst_m(:,:) * tmask(:,:,1) ) 182 IF( iom_use('vflx_fwb_cea') ) CALL iom_put( 'vflx_fwb_cea', -a_fwb * tmask(:,:,1) ) 183 ENDIF 184 ! Output restart information 185 IF( lrst_oce ) THEN 186 IF(lwp) WRITE(numout,*) 187 IF(lwp) WRITE(numout,*) 'sbc_fwb : writing FW-budget adjustment to ocean restart file at it = ', kt 188 IF(lwp) WRITE(numout,*) '~~~~' 189 CALL iom_rstput( kt, nitrst, numrow, 'a_fwb', a_fwb ) 190 END IF 191 ! 192 IF( kt == nitend .AND. lwp ) WRITE(numout,*) 'sbc_fwb : final freshwater-budget adjustment = ', a_fwb, 'kg/m2/s' 176 193 ! 177 194 CASE ( 3 ) !== global fwf set to zero and spread out over erp area ==! … … 211 228 qns(:,:) = qns(:,:) - zerp_cor(:,:) * rcp * sst_m(:,:) ! account for change to the heat budget due to fw correction 212 229 erp(:,:) = erp(:,:) + zerp_cor(:,:) 230 ! outputs 231 IF( iom_use('hflx_fwb_cea') ) CALL iom_put( 'hflx_fwb_cea', -zerp_cor(:,:) * rcp * sst_m(:,:) ) 232 IF( iom_use('vflx_fwb_cea') ) CALL iom_put( 'vflx_fwb_cea', -zerp_cor(:,:) ) 213 233 ! 214 234 IF( lwp ) THEN ! control print -
NEMO/trunk/tests/ISOMIP+/MY_SRC/tradmp.F90
r13982 r14995 11 11 !! NEMO 1.0 ! 2002-08 (G. Madec, E. Durand) free form + modules 12 12 !! 3.2 ! 2009-08 (G. Madec, C. Talandier) DOCTOR norm for namelist parameter 13 !! 3.3 ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 13 !! 3.3 ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 14 14 !! 3.4 ! 2011-04 (G. Madec, C. Ethe) Merge of dtatem and dtasal + suppression of CPP keys 15 15 !! 3.6 ! 2015-06 (T. Graham) read restoring coefficient in a file … … 26 26 USE c1d ! 1D vertical configuration 27 27 USE trd_oce ! trends: ocean variables 28 USE trdtra ! trends manager: tracers 28 USE trdtra ! trends manager: tracers 29 29 USE zdf_oce ! ocean: vertical physics 30 30 USE phycst ! physical constants … … 75 75 !!---------------------------------------------------------------------- 76 76 !! *** ROUTINE tra_dmp *** 77 !! 77 !! 78 78 !! ** Purpose : Compute the tracer trend due to a newtonian damping 79 79 !! of the tracer field towards given data field and add it to the 80 80 !! general tracer trends. 81 81 !! 82 !! ** Method : Newtonian damping towards t_dta and s_dta computed 82 !! ** Method : Newtonian damping towards t_dta and s_dta computed 83 83 !! and add to the general tracer trends: 84 84 !! ta = ta + resto * (t_dta - tb) … … 101 101 IF( ln_timing ) CALL timing_start('tra_dmp') 102 102 ! 103 IF( l_trdtra ) THEN!* Save ta and sa trends104 ALLOCATE( ztrdts(jpi,jpj,jpk,jpts) ) 105 ztrdts(:,:,:,:) = pts(:,:,:,:,Krhs) 103 IF( l_trdtra .OR. iom_use('hflx_dmp_cea') .OR. iom_use('sflx_dmp_cea') ) THEN !* Save ta and sa trends 104 ALLOCATE( ztrdts(jpi,jpj,jpk,jpts) ) 105 ztrdts(:,:,:,:) = pts(:,:,:,:,Krhs) 106 106 ENDIF 107 107 ! !== input T-S data at kt ==! … … 140 140 END SELECT 141 141 ! 142 ! outputs (clem trunk) 143 IF( iom_use('hflx_dmp_cea') ) & 144 & CALL iom_put('hflx_dmp_cea', & 145 & SUM( ( pts(:,:,:,jp_tem,Krhs) - ztrdts(:,:,:,jp_tem) ) * e3t(:,:,:,Kmm), dim=3 ) * rcp * rho0 ) ! W/m2 146 IF( iom_use('sflx_dmp_cea') ) & 147 & CALL iom_put('sflx_dmp_cea', & 148 & SUM( ( pts(:,:,:,jp_sal,Krhs) - ztrdts(:,:,:,jp_sal) ) * e3t(:,:,:,Kmm), dim=3 ) * rho0 ) ! g/m2/s 149 ! 142 150 IF( l_trdtra ) THEN ! trend diagnostic 143 151 ztrdts(:,:,:,:) = pts(:,:,:,:,Krhs) - ztrdts(:,:,:,:) 144 152 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) ) 145 153 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) ) 146 DEALLOCATE( ztrdts ) 154 DEALLOCATE( ztrdts ) 147 155 ENDIF 148 156 ! ! Control print … … 158 166 !!---------------------------------------------------------------------- 159 167 !! *** ROUTINE tra_dmp_init *** 160 !! 161 !! ** Purpose : Initialization for the newtonian damping 168 !! 169 !! ** Purpose : Initialization for the newtonian damping 162 170 !! 163 171 !! ** Method : read the namtra_dmp namelist and check the parameters 164 172 !!---------------------------------------------------------------------- 165 INTEGER :: ios, imask ! local integers 173 INTEGER :: ios, imask ! local integers 166 174 ! 167 175 NAMELIST/namtra_dmp/ ln_tradmp, nn_zdmp, cn_resto
Note: See TracChangeset
for help on using the changeset viewer.
