Changeset 4649
- Timestamp:
- 2014-05-27T11:28:12+02:00 (10 years ago)
- Location:
- branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM
- Files:
-
- 26 edited
-
CONFIG/ORCA2_LIM3/EXP00/field_def.xml (modified) (7 diffs)
-
CONFIG/ORCA2_LIM3/EXP00/iodef.xml (modified) (10 diffs)
-
CONFIG/ORCA2_LIM3/EXP00/iodef_default.xml (modified) (10 diffs)
-
CONFIG/ORCA2_LIM3/EXP00/namelist_ice (modified) (2 diffs)
-
CONFIG/ORCA2_LIM3/EXP00/namelist_ice_lim3 (modified) (2 diffs)
-
NEMO/LIM_SRC_3/ice.F90 (modified) (5 diffs)
-
NEMO/LIM_SRC_3/limcons.F90 (modified) (4 diffs)
-
NEMO/LIM_SRC_3/limdiahsb.F90 (modified) (13 diffs)
-
NEMO/LIM_SRC_3/limdyn.F90 (modified) (5 diffs)
-
NEMO/LIM_SRC_3/limistate.F90 (modified) (4 diffs)
-
NEMO/LIM_SRC_3/limitd_me.F90 (modified) (12 diffs)
-
NEMO/LIM_SRC_3/limitd_th.F90 (modified) (4 diffs)
-
NEMO/LIM_SRC_3/limsbc.F90 (modified) (5 diffs)
-
NEMO/LIM_SRC_3/limthd.F90 (modified) (19 diffs)
-
NEMO/LIM_SRC_3/limthd_dh.F90 (modified) (20 diffs)
-
NEMO/LIM_SRC_3/limthd_dif.F90 (modified) (8 diffs)
-
NEMO/LIM_SRC_3/limthd_ent.F90 (modified) (6 diffs)
-
NEMO/LIM_SRC_3/limthd_lac.F90 (modified) (21 diffs)
-
NEMO/LIM_SRC_3/limthd_sal.F90 (modified) (5 diffs)
-
NEMO/LIM_SRC_3/limtrp.F90 (modified) (9 diffs)
-
NEMO/LIM_SRC_3/limupdate1.F90 (modified) (10 diffs)
-
NEMO/LIM_SRC_3/limupdate2.F90 (modified) (12 diffs)
-
NEMO/LIM_SRC_3/limwri.F90 (modified) (6 diffs)
-
NEMO/LIM_SRC_3/thd_ice.F90 (modified) (4 diffs)
-
NEMO/OPA_SRC/SBC/sbcice_lim.F90 (modified) (5 diffs)
-
NEMO/OPA_SRC/SBC/sbcmod.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/CONFIG/ORCA2_LIM3/EXP00/field_def.xml
r4634 r4649 152 152 153 153 <!-- *_oce variables available with ln_blk_clio or ln_blk_core --> 154 <field id="qns_oce" long_name="Non solar Downward Heat Flux over open ocean" unit="W/m2" />155 154 <field id="qlw_oce" long_name="Longwave Downward Heat Flux over open ocean" unit="W/m2" /> 156 155 <field id="qsb_oce" long_name="Sensible Downward Heat Flux over open ocean" unit="W/m2" /> … … 222 221 <field id="qsr_oce" long_name="solar heat flux at ocean surface" unit="W/m2" /> 223 222 <field id="qns_oce" long_name="non-solar heat flux at ocean surface" unit="W/m2" /> 223 <field id="qt_ice" long_name="total flux at ice surface" unit="W/m2" /> 224 <field id="qsr_ice" long_name="solar heat flux at ice surface" unit="W/m2" /> 225 <field id="qns_ice" long_name="non-solar heat flux at ice surface" unit="W/m2" /> 226 <field id="qtr_ice" long_name="solar heat flux transmitted thru the ice" unit="W/m2" /> 224 227 <field id="utau_ice" long_name="Wind stress along i-axis over the ice at i-point" unit="N/m2" /> 225 228 <field id="vtau_ice" long_name="Wind stress along j-axis over the ice at i-point" unit="N/m2" /> 226 <field id="qsr_io" long_name="Ice-Oce downward solar heat flux" unit="W/m2" />227 <field id="qns_io" long_name="Ice-Oce downward non-solar heat flux" unit="W/m2" />228 229 <field id="micesalt" long_name="Mean ice salinity" unit="psu" /> 229 230 <field id="miceage" long_name="Mean ice age" unit="years" /> … … 237 238 238 239 <field id="micet" long_name="Mean ice temperature" unit="degC" /> 239 <field id="icehc" long_name="ice total heat content" unit="10^9 240 <field id="icehc" long_name="ice total heat content" unit="10^9J" /> 240 241 <field id="isnowhc" long_name="snow total heat content" unit="10^9J" /> 241 242 <field id="icest" long_name="ice surface temperature" unit="degC" /> … … 251 252 <field id="icetrp" long_name="ice volume transport" unit="m/day" /> 252 253 <field id="snwtrp" long_name="snw volume transport" unit="m/day" /> 253 <field id="deitrp" long_name="advected ice enhalpy" unit="W/ 2" />254 <field id="destrp" long_name="advected snw enhalpy" unit="W/ 2" />254 <field id="deitrp" long_name="advected ice enhalpy" unit="W/m2" /> 255 <field id="destrp" long_name="advected snw enhalpy" unit="W/m2" /> 255 256 256 257 <field id="sfxbri" long_name="brine salt flux" unit="psu*kg/m2/day" /> … … 274 275 <field id="vfxsnw" long_name="snw melt/growth" unit="m/day" /> 275 276 <field id="vfxsub" long_name="snw sublimation" unit="m/day" /> 276 277 <field id="hfxdhc1" long_name="Heat content variation in snow and ice" unit="W/m2" /> 278 <field id="hfxspr" long_name="Heat content of snow precip" unit="W/m2" /> 279 <field id="hfxqsr" long_name="solar fluxes given to ocean" unit="W/m2" /> 280 <field id="hfxqns" long_name="non solar fluxes given to ocean" unit="W/m2" /> 281 282 <field id="hfxthd" long_name="heat fluxes from ice-ocean exchange during thermo" unit="W/m2" /> 283 <field id="hfxdyn" long_name="heat fluxes from ice-ocean exchange during dynamic" unit="W/m2" /> 284 <field id="hfxres" long_name="heat fluxes from ice-ocean exchange during resultant" unit="W/m2" /> 285 <field id="hfxsnw" long_name="heat fluxes from snow-ocean exchange" unit="W/m2" /> 286 <field id="hfxsub" long_name="heat fluxes from sublimation" unit="W/m2" /> 287 <field id="hfxerr" long_name="heat fluxes error after heat diffusion" unit="W/m2" /> 288 <field id="hfxerr_rem" long_name="heat fluxes error after remapping" unit="W/m2" /> 289 <field id="hfxtot" long_name="heat fluxes total used by ice" unit="W/m2" /> 290 <field id="hfxout" long_name="non solar heat fluxes received by the ocean" unit="W/m2" /> 291 <field id="hfxin" long_name="total heat fluxes at the ice surface" unit="W/m2" /> 277 <field id="vfxspr" long_name="snw precipitation on ice" unit="m/day" /> 278 279 <field id="hfxsum" long_name="heat fluxes causing surface ice melt" unit="W/m2" /> 280 <field id="hfxbom" long_name="heat fluxes causing bottom ice melt" unit="W/m2" /> 281 <field id="hfxbog" long_name="heat fluxes causing bottom ice growth" unit="W/m2" /> 282 <field id="hfxdif" long_name="heat fluxes causing ice temperature change" unit="W/m2" /> 283 <field id="hfxopw" long_name="heat fluxes causing open water ice formation" unit="W/m2" /> 284 <field id="hfxsnw" long_name="heat fluxes causing snow melt" unit="W/m2" /> 285 <field id="hfxerr" long_name="heat fluxes error after heat diffusion" unit="W/m2" /> 286 <field id="hfxerr_rem" long_name="heat fluxes error after remapping" unit="W/m2" /> 287 <field id="hfxout" long_name="non solar heat fluxes received by the ocean" unit="W/m2" /> 288 <field id="hfxin" long_name="total heat fluxes at the ice surface" unit="W/m2" /> 289 290 <!-- heat flux associated with mass exchange --> 291 <field id="hfxthd" long_name="heat fluxes from ice-ocean mass exchange during thermo" unit="W/m2" /> 292 <field id="hfxdyn" long_name="heat fluxes from ice-ocean mass exchange during dynamic" unit="W/m2" /> 293 <field id="hfxres" long_name="heat fluxes from ice-ocean mass exchange during resultant" unit="W/m2" /> 294 <field id="hfxsub" long_name="heat fluxes from ice-atm. mass exchange during sublimation" unit="W/m2" /> 295 <field id="hfxspr" long_name="heat fluxes from ice-atm. mass exchange during snow precip" unit="W/m2" /> 296 297 <!-- diags --> 298 <field id="hfxdhc" long_name="Heat content variation in snow and ice" unit="W/m2" /> 299 <field id="hfxtur" long_name="turbulent heat flux at the ice base" unit="W/m2" /> 292 300 293 301 </field_group> … … 403 411 <field id="ibgvfxsum" long_name="global mean volume flux (surface melt)" unit="m/day" /> 404 412 <field id="ibgvfxres" long_name="global mean volume flux (resultant)" unit="m/day" /> 413 <field id="ibgvfxspr" long_name="global mean volume flux (snow precip)" unit="m/day" /> 414 <field id="ibgvfxsnw" long_name="global mean volume flux (snow melt)" unit="m/day" /> 415 <field id="ibgvfxsub" long_name="global mean volume flux (snow sublimation)" unit="m/day" /> 405 416 406 417 <field id="ibgsfx" long_name="global mean salt flux (total)" unit="psu*m/day" /> … … 415 426 416 427 417 <field id="ibghfxdhc 1"long_name="Heat content variation in snow and ice" unit="W" />428 <field id="ibghfxdhc" long_name="Heat content variation in snow and ice" unit="W" /> 418 429 <field id="ibghfxspr" long_name="Heat content of snow precip" unit="W" /> 419 <field id="ibghfxqsr" long_name="solar fluxes given to ocean" unit="W" />420 <field id="ibghfxqns" long_name="non solar fluxes given to ocean" unit="W" />421 430 422 431 <field id="ibghfxthd" long_name="heat fluxes from ice-ocean exchange during thermo" unit="W" /> 432 <field id="ibghfxsum" long_name="heat fluxes causing surface ice melt" unit="W" /> 433 <field id="ibghfxbom" long_name="heat fluxes causing bottom ice melt" unit="W" /> 434 <field id="ibghfxbog" long_name="heat fluxes causing bottom ice growth" unit="W" /> 435 <field id="ibghfxdif" long_name="heat fluxes causing ice temperature change" unit="W" /> 436 <field id="ibghfxopw" long_name="heat fluxes causing open water ice formation" unit="W" /> 423 437 <field id="ibghfxdyn" long_name="heat fluxes from ice-ocean exchange during dynamic" unit="W" /> 424 438 <field id="ibghfxres" long_name="heat fluxes from ice-ocean exchange during resultant" unit="W" /> 425 439 <field id="ibghfxsub" long_name="heat fluxes from sublimation" unit="W" /> 426 440 <field id="ibghfxsnw" long_name="heat fluxes from snow-ocean exchange" unit="W" /> 427 <field id="ibghfxtot" long_name="heat fluxes total used by ice" unit="W" />428 441 <field id="ibghfxout" long_name="non solar heat fluxes received by the ocean" unit="W" /> 429 442 <field id="ibghfxin" long_name="total heat fluxes at the ice surface" unit="W" /> -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/CONFIG/ORCA2_LIM3/EXP00/iodef.xml
r4634 r4649 53 53 <file id="file2" name_suffix="_SBC" description="surface fluxes variables" > <!-- time step automaticaly defined based on nn_fsbc --> 54 54 <field field_ref="empmr" name="wfo" long_name="water_flux_into_sea_water" /> 55 <field field_ref="qsr" name="hfxqsr" long_name="surface_net_downward_shortwave_flux" /> 56 <field field_ref="qt" name="hfxtot" long_name="surface_net_downward_total_heat_flux" /> 55 <field field_ref="qsr_oce" name="qsr_oce" long_name="downward shortwave flux at ocean surface" /> 56 <field field_ref="qns_oce" name="qns_oce" long_name="downward non solar flux at ocean surface" /> 57 <field field_ref="qt_oce" name="qt_oce" long_name="downward total flux at ocean surface" /> 58 <field field_ref="qsr_ice" name="qsr_ice" long_name="downward shortwave flux at ice surface" /> 59 <field field_ref="qns_ice" name="qns_ice" long_name="downward non solar flux at ice surface" /> 60 <field field_ref="qtr_ice" name="qtr_ice" long_name="shortwave flux transmitted thru the ice" /> 61 <field field_ref="qt_ice" name="qt_ice" long_name="downward total flux at ice surface" /> 57 62 <field field_ref="saltflx" name="sfx" /> 58 63 <field field_ref="taum" name="taum" /> … … 63 68 <field field_ref="utau_ice" name="utau_ice" /> 64 69 <field field_ref="vtau_ice" name="vtau_ice" /> 65 <!-- clem66 <field field_ref="qsr_io" name="iicesflx" />67 <field field_ref="qns_io" name="iicenflx" />68 -->69 70 70 71 </file> … … 97 98 98 99 <file id="file6" name_suffix="_icemod" description="ice variables" enabled=".true." > 99 <field field_ref="snowthic_cea" name="sn dept" long_name="surface_snow_thickness" />100 <field field_ref="snowthic_cea" name="snthic" long_name="surface_snow_thickness" /> 100 101 <field field_ref="icethic_cea" name="sithic" long_name="sea_ice_thickness" /> 101 102 <field field_ref="icevolu" name="sivolu" /> … … 113 114 <field field_ref="vfxsnw" name="vfxsnw" /> 114 115 <field field_ref="vfxsub" name="vfxsub" /> 116 <field field_ref="vfxspr" name="vfxspr" /> 115 117 116 118 <field field_ref="icetrp" name="sivtrp" /> … … 129 131 <field field_ref="sfx" name="sfx" /> 130 132 131 <field field_ref="hfxdhc1" name="hfxdhc1" /> 132 <field field_ref="hfxspr" name="hfxspr" /> 133 <field field_ref="hfxqsr" name="hfxqsr" /> 134 <field field_ref="hfxqns" name="hfxqns" /> 135 136 <field field_ref="hfxthd" name="hfxthd" /> 133 <field field_ref="hfxsum" name="hfxsum" /> 134 <field field_ref="hfxbom" name="hfxbom" /> 135 <field field_ref="hfxbog" name="hfxbog" /> 136 <field field_ref="hfxdif" name="hfxdif" /> 137 <field field_ref="hfxopw" name="hfxopw" /> 138 <field field_ref="hfxout" name="hfxout" /> 139 <field field_ref="hfxin" name="hfxin" /> 140 <field field_ref="hfxsnw" name="hfxsnw" /> 141 <field field_ref="hfxerr" name="hfxerr" /> 142 <field field_ref="hfxerr_rem" name="hfxerr_rem" /> 143 144 <!-- ice-ocean heat flux from mass exchange --> 137 145 <field field_ref="hfxdyn" name="hfxdyn" /> 138 146 <field field_ref="hfxres" name="hfxres" /> 139 <field field_ref="hfxout" name="hfxout" /> 140 <field field_ref="hfxin" name="hfxin" /> 141 <field field_ref="hfxtot" name="hfxtot" /> 142 <field field_ref="hfxsnw" name="hfxsnw" /> 147 <field field_ref="hfxthd" name="hfxthd" /> 148 <!-- ice-atm. heat flux from mass exchange --> 143 149 <field field_ref="hfxsub" name="hfxsub" /> 144 <field field_ref="hfxerr" name="hfxerr" /> 145 <field field_ref="hfxerr_rem" name="hfxerr_rem" /> 150 <field field_ref="hfxspr" name="hfxspr" /> 151 152 <!-- diags --> 153 <field field_ref="hfxdhc" name="hfxdhc" /> 154 <field field_ref="hfxtur" name="hfxtur" /> 146 155 147 156 <field field_ref="isst" name="sst" /> … … 167 176 <field field_ref="iceconc_cat" name="siconcat"/> 168 177 <field field_ref="icethic_cat" name="sithicat"/> 169 <field field_ref="snowthic_cat" name="sn deptcat"/>178 <field field_ref="snowthic_cat" name="snthicat"/> 170 179 <field field_ref="salinity_cat" name="salincat"/> 171 <field field_ref="brinevol_cat" name="sibri ncat"/>180 <field field_ref="brinevol_cat" name="sibricat"/> 172 181 173 182 </file> … … 191 200 <field field_ref="bgfrctem" name="bgfrctem" /> 192 201 <field field_ref="bgfrcsal" name="bgfrcsal" /> 193 <!-- 194 <field field_ref="bgmistem" name="bgmistem" /> 195 <field field_ref="bgmissal" name="bgmissal" /> 196 --> 202 197 203 </file> 198 204 … … 215 221 <field field_ref="ibgvfxsum" name="ibgvfxsum" /> 216 222 <field field_ref="ibgvfxres" name="ibgvfxres" /> 223 <field field_ref="ibgvfxspr" name="ibgvfxspr" /> 224 <field field_ref="ibgvfxsnw" name="ibgvfxsnw" /> 225 <field field_ref="ibgvfxsub" name="ibgvfxsub" /> 217 226 218 227 <field field_ref="ibgsfx" name="ibgsfx" /> … … 226 235 <field field_ref="ibgsfxsum" name="ibgsfxsum" /> 227 236 228 <field field_ref="ibghfxdhc 1" name="ibghfxdhc1" />237 <field field_ref="ibghfxdhc" name="ibghfxdhc" /> 229 238 <field field_ref="ibghfxspr" name="ibghfxspr" /> 230 <field field_ref="ibghfxqsr" name="ibghfxqsr" />231 <field field_ref="ibghfxqns" name="ibghfxqns" />232 239 233 240 <field field_ref="ibghfxres" name="ibghfxres" /> … … 235 242 <field field_ref="ibghfxdyn" name="ibghfxdyn" /> 236 243 <field field_ref="ibghfxthd" name="ibghfxthd" /> 244 <field field_ref="ibghfxsum" name="ibghfxsum" /> 245 <field field_ref="ibghfxbom" name="ibghfxbom" /> 246 <field field_ref="ibghfxbog" name="ibghfxbog" /> 247 <field field_ref="ibghfxdif" name="ibghfxdif" /> 248 <field field_ref="ibghfxopw" name="ibghfxopw" /> 237 249 <field field_ref="ibghfxout" name="ibghfxout" /> 238 250 <field field_ref="ibghfxin" name="ibghfxin" /> 239 <field field_ref="ibghfxtot" name="ibghfxtot" />240 251 <field field_ref="ibghfxsnw" name="ibghfxsnw" /> 241 252 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/CONFIG/ORCA2_LIM3/EXP00/iodef_default.xml
r4634 r4649 53 53 <file id="file2" name_suffix="_SBC" description="surface fluxes variables" > <!-- time step automaticaly defined based on nn_fsbc --> 54 54 <field field_ref="empmr" name="wfo" long_name="water_flux_into_sea_water" /> 55 <field field_ref="qsr" name="hfxqsr" long_name="surface_net_downward_shortwave_flux" /> 56 <field field_ref="qt" name="hfxtot" long_name="surface_net_downward_total_heat_flux" /> 55 <field field_ref="qsr_oce" name="qsr_oce" long_name="downward shortwave flux at ocean surface" /> 56 <field field_ref="qns_oce" name="qns_oce" long_name="downward non solar flux at ocean surface" /> 57 <field field_ref="qt_oce" name="qt_oce" long_name="downward total flux at ocean surface" /> 58 <field field_ref="qsr_ice" name="qsr_ice" long_name="downward shortwave flux at ice surface" /> 59 <field field_ref="qns_ice" name="qns_ice" long_name="downward non solar flux at ice surface" /> 60 <field field_ref="qtr_ice" name="qtr_ice" long_name="shortwave flux transmitted thru the ice" /> 61 <field field_ref="qt_ice" name="qt_ice" long_name="downward total flux at ice surface" /> 57 62 <field field_ref="saltflx" name="sfx" /> 58 63 <field field_ref="taum" name="taum" /> … … 63 68 <field field_ref="utau_ice" name="utau_ice" /> 64 69 <field field_ref="vtau_ice" name="vtau_ice" /> 65 <!-- clem66 <field field_ref="qsr_io" name="iicesflx" />67 <field field_ref="qns_io" name="iicenflx" />68 -->69 70 70 71 </file> … … 97 98 98 99 <file id="file6" name_suffix="_icemod" description="ice variables" enabled=".true." > 99 <field field_ref="snowthic_cea" name="sn dept" long_name="surface_snow_thickness" />100 <field field_ref="snowthic_cea" name="snthic" long_name="surface_snow_thickness" /> 100 101 <field field_ref="icethic_cea" name="sithic" long_name="sea_ice_thickness" /> 101 102 <field field_ref="icevolu" name="sivolu" /> … … 113 114 <field field_ref="vfxsnw" name="vfxsnw" /> 114 115 <field field_ref="vfxsub" name="vfxsub" /> 116 <field field_ref="vfxspr" name="vfxspr" /> 115 117 116 118 <field field_ref="icetrp" name="sivtrp" /> … … 129 131 <field field_ref="sfx" name="sfx" /> 130 132 131 <field field_ref="hfxdhc1" name="hfxdhc1" /> 132 <field field_ref="hfxspr" name="hfxspr" /> 133 <field field_ref="hfxqsr" name="hfxqsr" /> 134 <field field_ref="hfxqns" name="hfxqns" /> 135 136 <field field_ref="hfxthd" name="hfxthd" /> 133 <field field_ref="hfxsum" name="hfxsum" /> 134 <field field_ref="hfxbom" name="hfxbom" /> 135 <field field_ref="hfxbog" name="hfxbog" /> 136 <field field_ref="hfxdif" name="hfxdif" /> 137 <field field_ref="hfxopw" name="hfxopw" /> 138 <field field_ref="hfxout" name="hfxout" /> 139 <field field_ref="hfxin" name="hfxin" /> 140 <field field_ref="hfxsnw" name="hfxsnw" /> 141 <field field_ref="hfxerr" name="hfxerr" /> 142 <field field_ref="hfxerr_rem" name="hfxerr_rem" /> 143 144 <!-- ice-ocean heat flux from mass exchange --> 137 145 <field field_ref="hfxdyn" name="hfxdyn" /> 138 146 <field field_ref="hfxres" name="hfxres" /> 139 <field field_ref="hfxout" name="hfxout" /> 140 <field field_ref="hfxin" name="hfxin" /> 141 <field field_ref="hfxtot" name="hfxtot" /> 142 <field field_ref="hfxsnw" name="hfxsnw" /> 147 <field field_ref="hfxthd" name="hfxthd" /> 148 <!-- ice-atm. heat flux from mass exchange --> 143 149 <field field_ref="hfxsub" name="hfxsub" /> 144 <field field_ref="hfxerr" name="hfxerr" /> 145 <field field_ref="hfxerr_rem" name="hfxerr_rem" /> 150 <field field_ref="hfxspr" name="hfxspr" /> 151 152 <!-- diags --> 153 <field field_ref="hfxdhc" name="hfxdhc" /> 154 <field field_ref="hfxtur" name="hfxtur" /> 146 155 147 156 <field field_ref="isst" name="sst" /> … … 167 176 <field field_ref="iceconc_cat" name="siconcat"/> 168 177 <field field_ref="icethic_cat" name="sithicat"/> 169 <field field_ref="snowthic_cat" name="sn deptcat"/>178 <field field_ref="snowthic_cat" name="snthicat"/> 170 179 <field field_ref="salinity_cat" name="salincat"/> 171 <field field_ref="brinevol_cat" name="sibri ncat"/>180 <field field_ref="brinevol_cat" name="sibricat"/> 172 181 173 182 </file> … … 191 200 <field field_ref="bgfrctem" name="bgfrctem" /> 192 201 <field field_ref="bgfrcsal" name="bgfrcsal" /> 193 <!-- 194 <field field_ref="bgmistem" name="bgmistem" /> 195 <field field_ref="bgmissal" name="bgmissal" /> 196 --> 202 197 203 </file> 198 204 … … 215 221 <field field_ref="ibgvfxsum" name="ibgvfxsum" /> 216 222 <field field_ref="ibgvfxres" name="ibgvfxres" /> 223 <field field_ref="ibgvfxspr" name="ibgvfxspr" /> 224 <field field_ref="ibgvfxsnw" name="ibgvfxsnw" /> 225 <field field_ref="ibgvfxsub" name="ibgvfxsub" /> 217 226 218 227 <field field_ref="ibgsfx" name="ibgsfx" /> … … 226 235 <field field_ref="ibgsfxsum" name="ibgsfxsum" /> 227 236 228 <field field_ref="ibghfxdhc 1" name="ibghfxdhc1" />237 <field field_ref="ibghfxdhc" name="ibghfxdhc" /> 229 238 <field field_ref="ibghfxspr" name="ibghfxspr" /> 230 <field field_ref="ibghfxqsr" name="ibghfxqsr" />231 <field field_ref="ibghfxqns" name="ibghfxqns" />232 239 233 240 <field field_ref="ibghfxres" name="ibghfxres" /> … … 235 242 <field field_ref="ibghfxdyn" name="ibghfxdyn" /> 236 243 <field field_ref="ibghfxthd" name="ibghfxthd" /> 244 <field field_ref="ibghfxsum" name="ibghfxsum" /> 245 <field field_ref="ibghfxbom" name="ibghfxbom" /> 246 <field field_ref="ibghfxbog" name="ibghfxbog" /> 247 <field field_ref="ibghfxdif" name="ibghfxdif" /> 248 <field field_ref="ibghfxopw" name="ibghfxopw" /> 237 249 <field field_ref="ibghfxout" name="ibghfxout" /> 238 250 <field field_ref="ibghfxin" name="ibghfxin" /> 239 <field field_ref="ibghfxtot" name="ibghfxtot" />240 251 <field field_ref="ibghfxsnw" name="ibghfxsnw" /> 241 252 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/CONFIG/ORCA2_LIM3/EXP00/namelist_ice
r4634 r4649 15 15 !----------------------------------------------------------------------- 16 16 cn_icerst_in = "restart_ice" ! suffix of ice restart name (input) 17 cn_icerst_out = "restart_ice" 18 ln_limdyn = .true.! ice dynamics (T) or thermodynamics only (F)19 amax = 0.999! maximum ice concentration20 cai = 1.40e-3 ! atmospheric drag over sea ice21 cao = 1.00e-3 ! atmospheric drag over ocean22 ln_nicep = .false.! Ice points output for debug (yes or no)23 ln_limdiahsb = .false. ! check the heat and salt budgets (T) or not (F)24 ln_limdiaout = .true. ! output the heat and salt budgets (T) or not (F)17 cn_icerst_out = "restart_ice" ! suffix of ice restart name (output) 18 ln_limdyn = .true. ! ice dynamics (T) or thermodynamics only (F) 19 amax = 0.999 ! maximum ice concentration 20 cai = 1.40e-3 ! atmospheric drag over sea ice (clio) 21 cao = 1.00e-3 ! atmospheric drag over ocean (clio) 22 ln_nicep = .false. ! Ice points output for debug (yes or no) 23 ln_limdiahsb = .false. ! check the heat and salt budgets (T) or not (F) 24 ln_limdiaout = .true. ! output the heat and salt budgets (T) or not (F) 25 25 / 26 26 !----------------------------------------------------------------------- … … 147 147 !&namicehsb ! Heat and salt budgets 148 148 !!----------------------------------------------------------------------- 149 !150 149 !/ 151 !-----------------------------------------------------------------------152 &namiceout ! parameters for outputs153 !-----------------------------------------------------------------------154 ! noumef = 43 ! number of fields155 !SF add_diag_swi= 1 ! 1 -> diagnose distribution in thickness space156 !SF ! 0 -> only simple diagnostics157 !158 ! ! title of the field ! name ! units ! save ! multipl. ! additive !159 ! ! ! ! ! or not ! factor ! factor !160 /161 150 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/CONFIG/ORCA2_LIM3/EXP00/namelist_ice_lim3
r4634 r4649 15 15 !----------------------------------------------------------------------- 16 16 cn_icerst_in = "restart_ice" ! suffix of ice restart name (input) 17 cn_icerst_out = "restart_ice" 18 ln_limdyn = .true.! ice dynamics (T) or thermodynamics only (F)19 amax = 0.999! maximum ice concentration20 cai = 1.40e-3 ! atmospheric drag over sea ice21 cao = 1.00e-3 ! atmospheric drag over ocean22 ln_nicep = .false.! Ice points output for debug (yes or no)23 ln_limdiahsb = .false. ! check the heat and salt budgets (T) or not (F)24 ln_limdiaout = .true. ! output the heat and salt budgets (T) or not (F)17 cn_icerst_out = "restart_ice" ! suffix of ice restart name (output) 18 ln_limdyn = .true. ! ice dynamics (T) or thermodynamics only (F) 19 amax = 0.999 ! maximum ice concentration 20 cai = 1.40e-3 ! atmospheric drag over sea ice (clio) 21 cao = 1.00e-3 ! atmospheric drag over ocean (clio) 22 ln_nicep = .false. ! Ice points output for debug (yes or no) 23 ln_limdiahsb = .false. ! check the heat and salt budgets (T) or not (F) 24 ln_limdiaout = .true. ! output the heat and salt budgets (T) or not (F) 25 25 / 26 26 !----------------------------------------------------------------------- … … 147 147 !&namicehsb ! Heat and salt budgets 148 148 !!----------------------------------------------------------------------- 149 !150 149 !/ 151 !-----------------------------------------------------------------------152 &namiceout ! parameters for outputs153 !-----------------------------------------------------------------------154 ! noumef = 43 ! number of fields155 !SF add_diag_swi= 1 ! 1 -> diagnose distribution in thickness space156 !SF ! 0 -> only simple diagnostics157 !158 ! ! title of the field ! name ! units ! save ! multipl. ! additive !159 ! ! ! ! ! or not ! factor ! factor !160 /161 150 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
r4634 r4649 270 270 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_sum ! vertical surface melt 271 271 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_res ! production (growth+melt) due to limupdate 272 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_spr ! snow precipitation on ice 272 273 273 274 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_bog !: salt flux due to ice growth/melt [PSU/m2/s] … … 280 281 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx_res !: residual salt flux due to correction of ice thickness [PSU/m2/s] 281 282 282 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_ thd !: ice-ocean heat flux from thermo processes (limthd_dh)283 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_ dyn !: ice-ocean heat flux from mecanical processes (limitd_me)284 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_ tot !: total heat flux lost/gained by ice285 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_ spr !: heat flux of the snow precipitation286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_ res !: residual heat flux due to correction of ice thickness283 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_bog !: total heat flux causing bottom ice growth 284 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_bom !: total heat flux causing bottom ice melt 285 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_sum !: total heat flux causing surface ice melt 286 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_opw !: total heat flux causing open water ice formation 287 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_dif !: total heat flux causing Temp change in the ice 287 288 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_snw !: heat flux for snow melt 288 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_sub !: heat flux for sublimation289 289 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_err !: heat flux error after heat diffusion 290 290 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_err_rem !: heat flux error after heat remapping 291 291 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_in !: heat flux available for thermo transformations 292 292 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_out !: heat flux remaining at the end of thermo transformations 293 294 ! heat flux associated with ice-atmosphere mass exchange 295 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_sub !: heat flux for sublimation 296 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_spr !: heat flux of the snow precipitation 297 298 ! heat flux associated with ice-ocean mass exchange 299 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_thd !: ice-ocean heat flux from thermo processes (limthd_dh) 300 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_dyn !: ice-ocean heat flux from mecanical processes (limitd_me) 301 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hfx_res !: residual heat flux due to correction of ice thickness 293 302 294 303 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ftr_ice !: transmitted solar radiation under ice … … 416 425 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_trp_es ! transport of snw enthalpy (W/m2) 417 426 ! 418 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_heat_dhc 1! snw/ice heat content variation [W/m2]427 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_heat_dhc ! snw/ice heat content variation [W/m2] 419 428 ! 420 429 INTEGER , PUBLIC :: jiindx, jjindx !: indexes of the debugging point … … 450 459 451 460 ii = ii + 1 452 ALLOCATE( sist (jpi,jpj) , icethi (jpi,jpj) , t_bo (jpi,jpj) , & 453 & frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , & 454 & wfx_snw (jpi,jpj) , wfx_ice(jpi,jpj) , wfx_sub(jpi,jpj) , & 455 & wfx_bog(jpi,jpj) , wfx_dyn(jpi,jpj) , wfx_bom(jpi,jpj) , wfx_sum(jpi,jpj) , & 456 & wfx_res(jpi,jpj) , wfx_sni(jpi,jpj) , wfx_opw(jpi,jpj) , qlead (jpi,jpj) , & 457 & fhtur (jpi,jpj) , ftr_ice(jpi,jpj,jpl) , & 458 & sfx_res (jpi,jpj) , sfx_bri(jpi,jpj) , sfx_dyn(jpi,jpj) , & 459 & sfx_bog (jpi,jpj) , sfx_bom (jpi,jpj) , sfx_sum (jpi,jpj) , sfx_sni (jpi,jpj) , sfx_opw (jpi,jpj) , & 460 & hfx_res (jpi,jpj) , hfx_snw (jpi,jpj) , hfx_sub(jpi,jpj) , hfx_err(jpi,jpj), hfx_err_rem(jpi,jpj), & 461 & hfx_in (jpi,jpj) , hfx_out(jpi,jpj) , fhld(jpi,jpj) , & 462 & hfx_tot (jpi,jpj) , hfx_thd (jpi,jpj) , hfx_dyn(jpi,jpj) , hfx_spr(jpi,jpj), STAT=ierr(ii) ) 461 ALLOCATE( sist (jpi,jpj) , icethi (jpi,jpj) , t_bo (jpi,jpj) , & 462 & frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , & 463 & wfx_snw(jpi,jpj) , wfx_ice(jpi,jpj) , wfx_sub(jpi,jpj) , & 464 & wfx_bog(jpi,jpj) , wfx_dyn(jpi,jpj) , wfx_bom(jpi,jpj) , wfx_sum(jpi,jpj) , & 465 & wfx_res(jpi,jpj) , wfx_sni(jpi,jpj) , wfx_opw(jpi,jpj) , wfx_spr(jpi,jpj) , qlead (jpi,jpj) , & 466 & fhtur (jpi,jpj) , ftr_ice(jpi,jpj,jpl) , & 467 & sfx_res(jpi,jpj) , sfx_bri(jpi,jpj) , sfx_dyn(jpi,jpj) , & 468 & sfx_bog(jpi,jpj) , sfx_bom(jpi,jpj) , sfx_sum(jpi,jpj) , sfx_sni(jpi,jpj) , sfx_opw(jpi,jpj) , & 469 & hfx_res(jpi,jpj) , hfx_snw(jpi,jpj) , hfx_sub(jpi,jpj) , hfx_err(jpi,jpj) , hfx_err_rem(jpi,jpj), & 470 & hfx_in (jpi,jpj) , hfx_out(jpi,jpj) , fhld(jpi,jpj) , & 471 & hfx_sum(jpi,jpj) , hfx_bom(jpi,jpj) , hfx_bog(jpi,jpj) , hfx_dif(jpi,jpj) , hfx_opw(jpi,jpj) , & 472 & hfx_thd(jpi,jpj) , hfx_dyn(jpi,jpj) , hfx_spr(jpi,jpj) , STAT=ierr(ii) ) 463 473 464 474 ii = ii + 1 … … 526 536 ALLOCATE( dv_dt_thd(jpi,jpj,jpl) , & 527 537 & izero (jpi,jpj,jpl) , diag_trp_vi(jpi,jpj) , diag_trp_vs(jpi,jpj), diag_trp_ei(jpi,jpj), diag_trp_es(jpi,jpj), & 528 & diag_heat_dhc 1(jpi,jpj) , STAT=ierr(ii) )538 & diag_heat_dhc(jpi,jpj) , STAT=ierr(ii) ) 529 539 530 540 ice_alloc = MAXVAL( ierr(:) ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90
r4045 r4649 7 7 !! 3.0 ! 2004-06 (M. Vancoppenolle) Energy Conservation 8 8 !! 4.0 ! 2011-02 (G. Madec) add mpp considerations 9 !! - ! 2014-05 (C. Rousset) add lim_cons_hsm 9 10 !!---------------------------------------------------------------------- 10 11 #if defined key_lim3 … … 14 15 !! lim_cons : checks whether energy, mass and salt are conserved 15 16 !!---------------------------------------------------------------------- 17 USE phycst ! physical constants 16 18 USE par_ice ! LIM-3 parameter 17 19 USE ice ! LIM-3 variables … … 28 30 PUBLIC lim_column_sum_energy 29 31 PUBLIC lim_cons_check 32 PUBLIC lim_cons_hsm 30 33 31 34 !!---------------------------------------------------------------------- … … 151 154 END SUBROUTINE lim_cons_check 152 155 156 157 SUBROUTINE lim_cons_hsm( icount, cd_routine, zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b ) 158 !!------------------------------------------------------------------- 159 !! *** ROUTINE lim_cons_hsm *** 160 !! 161 !! ** Purpose : Test the conservation of heat, salt and mass for each routine 162 !! 163 !! ** Method : 164 !!--------------------------------------------------------------------- 165 INTEGER , INTENT(in) :: icount ! determine wether this is the beggining of the routine (0) or the end (1) 166 CHARACTER(len=*), INTENT(in) :: cd_routine ! name of the routine 167 REAL(wp) , INTENT(inout) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 168 REAL(wp) :: zvi, zsmv, zei, zfs, zfw, zft 169 REAL(wp) :: zvmin, zamin, zamax 170 171 IF( icount == 0 ) THEN 172 173 zvi_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 174 zsmv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 175 zei_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 176 zfw_b = glob_sum( - ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + & 177 & wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) + wfx_sub(:,:) + wfx_spr(:,:) ) * area(:,:) * tms(:,:) ) 178 zfs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + & 179 & sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 180 zft_b = glob_sum( ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:) & 181 & - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 182 183 ELSEIF( icount == 1 ) THEN 184 185 zfs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + & 186 & sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zfs_b 187 zfw = glob_sum( - ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + & 188 & wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) + wfx_sub(:,:) + wfx_spr(:,:) ) * area(:,:) * tms(:,:) ) - zfw_b 189 zft = glob_sum( ( hfx_sum(:,:) + hfx_bom(:,:) + hfx_bog(:,:) + hfx_dif(:,:) + hfx_opw(:,:) + hfx_snw(:,:) & 190 & - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) - hfx_sub(:,:) - hfx_spr(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zft_b 191 192 zvi = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zvi_b ) * r1_rdtice - zfw 193 zsmv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zsmv_b ) * r1_rdtice + ( zfs / rhoic ) 194 zei = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zei_b * r1_rdtice + zft 195 196 zvmin = glob_min(v_i) 197 zamax = glob_max(SUM(a_i,dim=3)) 198 zamin = glob_min(a_i) 199 200 IF(lwp) THEN 201 IF ( ABS( zvi ) > 1.e-4 ) WRITE(numout,*) 'violation volume [kg/day] (',cd_routine,') = ',(zvi * rday) 202 IF ( ABS( zsmv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (',cd_routine,') = ',(zsmv * rday) 203 IF ( ABS( zei ) > 1. ) WRITE(numout,*) 'violation enthalpy [1e9 J] (',cd_routine,') = ',(zei) 204 IF ( zvmin < 0. ) WRITE(numout,*) 'violation v_i<0 [m] (',cd_routine,') = ',(zvmin) 205 IF( cd_routine /= 'limtrp' .AND. cd_routine /= 'limitd_me' .AND. zamax > amax+1.e-10 ) THEN 206 WRITE(numout,*) 'violation a_i>amax (',cd_routine,') = ',zamax 207 ENDIF 208 IF ( zamin < 0. ) WRITE(numout,*) 'violation a_i<0 (',cd_routine,') = ',zamin 209 ENDIF 210 211 ENDIF 212 213 END SUBROUTINE lim_cons_hsm 214 153 215 #else 154 216 !!---------------------------------------------------------------------- -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limdiahsb.F90
r4634 r4649 35 35 !!PUBLIC lim_diahsb_rst ! routine called by ice_init.F90 36 36 37 REAL( wp) :: frc_sal, frc_vol ! global forcing trends38 REAL( wp) :: bg_grme ! global ice growth+melt trends37 REAL(dp) :: frc_sal, frc_vol ! global forcing trends 38 REAL(dp) :: bg_grme ! global ice growth+melt trends 39 39 REAL(wp) :: epsi06 = 1.e-6_wp ! small number 40 40 … … 47 47 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 48 48 !!---------------------------------------------------------------------- 49 49 50 CONTAINS 50 51 … … 57 58 !!--------------------------------------------------------------------------- 58 59 !! 59 REAL(wp) :: zbg_ivo, zbg_svo, zbg_are, zbg_sal ,zbg_tem ,zbg_ihc ,zbg_shc 60 REAL(wp) :: zbg_sfx, zbg_sfx_bri, zbg_sfx_bog, zbg_sfx_bom, zbg_sfx_sum, zbg_sfx_sni, zbg_sfx_opw, zbg_sfx_res, zbg_sfx_dyn 61 REAL(wp) :: zbg_vfx, zbg_vfx_bog, zbg_vfx_opw, zbg_vfx_sni, zbg_vfx_dyn, zbg_vfx_bom, zbg_vfx_sum, zbg_vfx_res 62 REAL(wp) :: zbg_hfx_dhc1, zbg_hfx_spr, zbg_hfx_qsr, zbg_hfx_qns 63 REAL(wp) :: zbg_hfx_res, zbg_hfx_sub, zbg_hfx_dyn, zbg_hfx_thd, zbg_hfx_snw, zbg_hfx_tot, zbg_hfx_out, zbg_hfx_in 64 REAL(wp) :: z_frc_vol, z_frc_sal, z_bg_grme 65 REAL(wp) :: z1_area, zcoef 66 REAL(wp) :: zinda, zindb 60 REAL(dp) :: zbg_ivo, zbg_svo, zbg_are, zbg_sal ,zbg_tem ,zbg_ihc ,zbg_shc 61 REAL(dp) :: zbg_sfx, zbg_sfx_bri, zbg_sfx_bog, zbg_sfx_bom, zbg_sfx_sum, zbg_sfx_sni, zbg_sfx_opw, zbg_sfx_res, zbg_sfx_dyn 62 REAL(dp) :: zbg_vfx, zbg_vfx_bog, zbg_vfx_opw, zbg_vfx_sni, zbg_vfx_dyn 63 REAL(dp) :: zbg_vfx_bom, zbg_vfx_sum, zbg_vfx_res, zbg_vfx_spr, zbg_vfx_snw, zbg_vfx_sub 64 REAL(dp) :: zbg_hfx_dhc, zbg_hfx_spr 65 REAL(dp) :: zbg_hfx_res, zbg_hfx_sub, zbg_hfx_dyn, zbg_hfx_thd, zbg_hfx_snw, zbg_hfx_out, zbg_hfx_in 66 REAL(dp) :: zbg_hfx_sum, zbg_hfx_bom, zbg_hfx_bog, zbg_hfx_dif, zbg_hfx_opw 67 REAL(dp) :: z_frc_vol, z_frc_sal, z_bg_grme 68 REAL(dp) :: z1_area ! - - 69 REAL(dp) :: zinda, zindb 67 70 !!--------------------------------------------------------------------------- 68 71 IF( nn_timing == 1 ) CALL timing_start('lim_diahsb') … … 70 73 IF( numit == nstart ) CALL lim_diahsb_init 71 74 72 z1_area = 1._wp / MAX( glob_sum( area(:,:) * tms(:,:) ), epsi06 ) 73 74 zinda = MAX( 0._wp , SIGN( 1._wp , glob_sum( area(:,:) * tms(:,:) ) - epsi06 ) ) 75 ! 1/area 76 z1_area = 1.d0 / MAX( glob_sum( area(:,:) * tms(:,:) ), epsi06 ) 77 78 zinda = MAX( 0.d0 , SIGN( 1.d0 , glob_sum( area(:,:) * tms(:,:) ) - epsi06 ) ) 75 79 ! ----------------------- ! 76 80 ! 1 - Content variations ! … … 79 83 zbg_svo = glob_sum( vt_s(:,:) * area(:,:) * tms(:,:) ) ! volume snow 80 84 zbg_are = glob_sum( at_i(:,:) * area(:,:) * tms(:,:) ) ! area 81 zbg_sal = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) )! mean salt content85 zbg_sal = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) ! mean salt content 82 86 zbg_tem = glob_sum( ( tm_i(:,:) - rtt ) * vt_i(:,:) * area(:,:) * tms(:,:) ) ! mean temp content 83 87 84 zcoef = zinda * z1_area * r1_rau0 * rday 88 !zbg_ihc = glob_sum( et_i(:,:) * area(:,:) * tms(:,:) ) / MAX( zbg_ivo,epsi06 ) ! ice heat content 89 !zbg_shc = glob_sum( et_s(:,:) * area(:,:) * tms(:,:) ) / MAX( zbg_svo,epsi06 ) ! snow heat content 90 85 91 ! Volume 86 zbg_vfx = glob_sum( emp(:,:) * area(:,:) * tms(:,:) ) * zcoef 87 zbg_vfx_bog = glob_sum( wfx_bog(:,:) * area(:,:) * tms(:,:) ) * zcoef 88 zbg_vfx_opw = glob_sum( wfx_opw(:,:) * area(:,:) * tms(:,:) ) * zcoef 89 zbg_vfx_sni = glob_sum( wfx_sni(:,:) * area(:,:) * tms(:,:) ) * zcoef 90 zbg_vfx_dyn = glob_sum( wfx_dyn(:,:) * area(:,:) * tms(:,:) ) * zcoef 91 zbg_vfx_bom = glob_sum( wfx_bom(:,:) * area(:,:) * tms(:,:) ) * zcoef 92 zbg_vfx_sum = glob_sum( wfx_sum(:,:) * area(:,:) * tms(:,:) ) * zcoef 93 zbg_vfx_res = glob_sum( wfx_res(:,:) * area(:,:) * tms(:,:) ) * zcoef 92 zbg_vfx = zinda * glob_sum( emp(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 93 zbg_vfx_bog = zinda * glob_sum( wfx_bog(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 94 zbg_vfx_opw = zinda * glob_sum( wfx_opw(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 95 zbg_vfx_sni = zinda * glob_sum( wfx_sni(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 96 zbg_vfx_dyn = zinda * glob_sum( wfx_dyn(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 97 zbg_vfx_bom = zinda * glob_sum( wfx_bom(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 98 zbg_vfx_sum = zinda * glob_sum( wfx_sum(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 99 zbg_vfx_res = zinda * glob_sum( wfx_res(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 100 zbg_vfx_spr = zinda * glob_sum( wfx_spr(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 101 zbg_vfx_snw = zinda * glob_sum( wfx_snw(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 102 zbg_vfx_sub = zinda * glob_sum( wfx_sub(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 94 103 95 104 ! Salt 96 zbg_sfx = glob_sum( sfx(:,:) * area(:,:) * tms(:,:) ) * zcoef97 zbg_sfx_bri = glob_sum( sfx_bri(:,:) * area(:,:) * tms(:,:) ) * zcoef98 zbg_sfx_res = glob_sum( sfx_res(:,:) * area(:,:) * tms(:,:) ) * zcoef99 zbg_sfx_dyn = glob_sum( sfx_dyn(:,:) * area(:,:) * tms(:,:) ) * zcoef100 101 zbg_sfx_bog = glob_sum( sfx_bog(:,:) * area(:,:) * tms(:,:) ) * zcoef102 zbg_sfx_opw = glob_sum( sfx_opw(:,:) * area(:,:) * tms(:,:) ) * zcoef103 zbg_sfx_sni = glob_sum( sfx_sni(:,:) * area(:,:) * tms(:,:) ) * zcoef104 zbg_sfx_bom = glob_sum( sfx_bom(:,:) * area(:,:) * tms(:,:) ) * zcoef105 zbg_sfx_sum = glob_sum( sfx_sum(:,:) * area(:,:) * tms(:,:) ) * zcoef105 zbg_sfx = zinda * glob_sum( sfx(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 106 zbg_sfx_bri = zinda * glob_sum( sfx_bri(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 107 zbg_sfx_res = zinda * glob_sum( sfx_res(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 108 zbg_sfx_dyn = zinda * glob_sum( sfx_dyn(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 109 110 zbg_sfx_bog = zinda * glob_sum( sfx_bog(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 111 zbg_sfx_opw = zinda * glob_sum( sfx_opw(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 112 zbg_sfx_sni = zinda * glob_sum( sfx_sni(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 113 zbg_sfx_bom = zinda * glob_sum( sfx_bom(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 114 zbg_sfx_sum = zinda * glob_sum( sfx_sum(:,:) * area(:,:) * tms(:,:) ) * z1_area * r1_rau0 * rday 106 115 107 116 ! Heat budget 108 117 zbg_ihc = glob_sum( et_i(:,:) * 1.e-20 ) ! ice heat content [1.e-20 J] 109 118 zbg_shc = glob_sum( et_s(:,:) * 1.e-20 ) ! snow heat content [1.e-20 J] 110 zbg_hfx_dhc 1 = glob_sum( diag_heat_dhc1(:,:) * area(:,:) * tms(:,:) ) ! [in W]119 zbg_hfx_dhc = glob_sum( diag_heat_dhc(:,:) * area(:,:) * tms(:,:) ) ! [in W] 111 120 zbg_hfx_spr = glob_sum( hfx_spr(:,:) * area(:,:) * tms(:,:) ) ! [in W] 112 zbg_hfx_qsr = glob_sum( qsr(:,:) * area(:,:) * tms(:,:) ) ! [in W]113 zbg_hfx_qns = glob_sum( qns(:,:) * area(:,:) * tms(:,:) ) ! [in W]114 121 115 122 zbg_hfx_thd = glob_sum( hfx_thd(:,:) * area(:,:) * tms(:,:) ) ! [in W] … … 118 125 zbg_hfx_sub = glob_sum( hfx_sub(:,:) * area(:,:) * tms(:,:) ) ! [in W] 119 126 zbg_hfx_snw = glob_sum( hfx_snw(:,:) * area(:,:) * tms(:,:) ) ! [in W] 120 zbg_hfx_tot = glob_sum( hfx_tot(:,:) * area(:,:) * tms(:,:) ) ! [in W] 127 zbg_hfx_sum = glob_sum( hfx_sum(:,:) * area(:,:) * tms(:,:) ) ! [in W] 128 zbg_hfx_bom = glob_sum( hfx_bom(:,:) * area(:,:) * tms(:,:) ) ! [in W] 129 zbg_hfx_bog = glob_sum( hfx_bog(:,:) * area(:,:) * tms(:,:) ) ! [in W] 130 zbg_hfx_dif = glob_sum( hfx_dif(:,:) * area(:,:) * tms(:,:) ) ! [in W] 131 zbg_hfx_opw = glob_sum( hfx_opw(:,:) * area(:,:) * tms(:,:) ) ! [in W] 121 132 zbg_hfx_out = glob_sum( hfx_out(:,:) * area(:,:) * tms(:,:) ) ! [in W] 122 133 zbg_hfx_in = glob_sum( hfx_in(:,:) * area(:,:) * tms(:,:) ) ! [in W] … … 127 138 z_frc_vol = r1_rau0 * glob_sum( - emp(:,:) * area(:,:) * tms(:,:) ) ! volume fluxes 128 139 z_frc_sal = r1_rau0 * glob_sum( sfx(:,:) * area(:,:) * tms(:,:) ) ! salt fluxes 129 z_bg_grme = glob_sum( ( wfx_bog(:,:) + wfx_opw(:,:) + wfx_sni(:,:) + wfx_dyn(:,:) + & 130 & wfx_bom(:,:) + wfx_sum(:,:) + wfx_res(:,:) ) / rhoic * area(:,:) * tms(:,:) ) ! volume fluxes 131 ! 132 frc_vol = frc_vol + z_frc_vol * rdt_ice 133 frc_sal = frc_sal + z_frc_sal * rdt_ice 134 bg_grme = bg_grme + z_bg_grme * rdt_ice 135 140 z_bg_grme = glob_sum( - ( wfx_bog(:,:) + wfx_opw(:,:) + wfx_sni(:,:) + wfx_dyn(:,:) + & 141 & wfx_bom(:,:) + wfx_sum(:,:) + wfx_res(:,:) + wfx_snw(:,:) + wfx_sub(:,:) ) * area(:,:) * tms(:,:) ) ! volume fluxes 142 ! 143 frc_vol = frc_vol + z_frc_vol * rdt_ice 144 frc_sal = frc_sal + z_frc_sal * rdt_ice 145 bg_grme = bg_grme + z_bg_grme * rdt_ice 146 147 ! difference 148 !frc_vol = zbg_ivo - frc_vol 149 !frc_sal = zbg_sal - frc_sal 150 136 151 ! ----------------------- ! 137 152 ! 3 - Diagnostics writing ! 138 153 ! ----------------------- ! 139 zindb = MAX( 0. _wp , SIGN( 1._wp, zbg_ivo - epsi06 ) )154 zindb = MAX( 0.d0 , SIGN( 1.d0 , zbg_ivo - epsi06 ) ) 140 155 ! 141 156 CALL iom_put( 'ibgvoltot' , zbg_ivo * rhoic * r1_rau0 * 1.e-9 ) ! ice volume (km3 equivalent liquid) … … 156 171 CALL iom_put( 'ibgvfxsum' , zbg_vfx_sum ) ! volume flux surface melt - 157 172 CALL iom_put( 'ibgvfxres' , zbg_vfx_res ) ! volume flux resultant - 173 CALL iom_put( 'ibgvfxspr' , zbg_vfx_spr ) ! volume flux from snow precip - 174 CALL iom_put( 'ibgvfxsnw' , zbg_vfx_snw ) ! volume flux from snow melt - 175 CALL iom_put( 'ibgvfxsub' , zbg_vfx_sub ) ! volume flux from sublimation - 158 176 159 177 CALL iom_put( 'ibgsfx' , zbg_sfx ) ! salt flux -(psu*m/day equivalent liquid) … … 167 185 CALL iom_put( 'ibgsfxsum' , zbg_sfx_sum ) ! salt flux surface melt - 168 186 169 CALL iom_put( 'ibghfxdhc 1', zbg_hfx_dhc1) ! Heat content variation in snow and ice [W]187 CALL iom_put( 'ibghfxdhc' , zbg_hfx_dhc ) ! Heat content variation in snow and ice [W] 170 188 CALL iom_put( 'ibghfxspr' , zbg_hfx_spr ) ! Heat content of snow precip [W] 171 CALL iom_put( 'ibghfxqsr' , zbg_hfx_qsr ) ! solar fluxes used by snw/ice [W]172 CALL iom_put( 'ibghfxqns' , zbg_hfx_qns ) ! non solar fluxes used by snw/ice [W]173 189 174 190 CALL iom_put( 'ibghfxres' , zbg_hfx_res ) ! … … 177 193 CALL iom_put( 'ibghfxthd' , zbg_hfx_thd ) ! 178 194 CALL iom_put( 'ibghfxsnw' , zbg_hfx_snw ) ! 179 CALL iom_put( 'ibghfxtot' , zbg_hfx_tot ) ! 195 CALL iom_put( 'ibghfxsum' , zbg_hfx_sum ) ! 196 CALL iom_put( 'ibghfxbom' , zbg_hfx_bom ) ! 197 CALL iom_put( 'ibghfxbog' , zbg_hfx_bog ) ! 198 CALL iom_put( 'ibghfxdif' , zbg_hfx_dif ) ! 199 CALL iom_put( 'ibghfxopw' , zbg_hfx_opw ) ! 180 200 CALL iom_put( 'ibghfxout' , zbg_hfx_out ) ! 181 CALL iom_put( 'ibghfxin' , zbg_hfx_in ) !201 CALL iom_put( 'ibghfxin' , zbg_hfx_in ) ! 182 202 183 203 CALL iom_put( 'ibgfrcvol' , frc_vol * 1.e-9 ) ! vol - forcing (km3 equivalent liquid) 184 204 CALL iom_put( 'ibgfrcsfx' , frc_sal * 1.e-9 ) ! sal - forcing (psu*km3 equivalent liquid) 185 CALL iom_put( 'ibgvolgrm' , bg_grme * r hoic * r1_rau0 * 1.e-9) ! vol growth + melt (km3 equivalent liquid)205 CALL iom_put( 'ibgvolgrm' , bg_grme * r1_rau0 * 1.e-9 ) ! vol growth + melt (km3 equivalent liquid) 186 206 187 207 ! … … 189 209 ! 190 210 IF( nn_timing == 1 ) CALL timing_stop('lim_diahsb') 191 211 ! 192 212 END SUBROUTINE lim_diahsb 193 213 … … 223 243 ! 2 - initial conservation variables ! 224 244 ! ---------------------------------- ! 245 !frc_vol = 0.d0 ! volume trend due to forcing 246 !frc_sal = 0.d0 ! salt content - - - - 247 !bg_grme = 0.d0 ! ice growth + melt volume trend 225 248 ! 226 249 CALL lim_diahsb_rst( nstart, 'READ' ) !* read or initialize all required files … … 256 279 IF(lwp) WRITE(numout,*) ' lim_diahsb at initial state ' 257 280 IF(lwp) WRITE(numout,*) '~~~~~~~' 258 frc_vol = 0. _wp259 frc_sal = 0. _wp260 bg_grme = 0. _wp281 frc_vol = 0.d0 282 frc_sal = 0.d0 283 bg_grme = 0.d0 261 284 ENDIF 262 285 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90
r4634 r4649 30 30 USE lib_fortran ! glob_sum 31 31 USE timing ! Timing 32 USE limcons ! conservation tests 32 33 33 34 IMPLICIT NONE … … 66 67 REAL(wp), POINTER, DIMENSION(:) :: zmsk ! i-averaged of tmask 67 68 REAL(wp), POINTER, DIMENSION(:,:) :: zu_io, zv_io ! ice-ocean velocity 68 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset)69 REAL(wp) :: z chk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset)69 ! 70 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 70 71 !!--------------------------------------------------------------------- 71 72 … … 75 76 CALL wrk_alloc( jpj, zind, zmsk ) 76 77 77 ! -------------------------------78 !- check conservation (C Rousset)79 IF (ln_limdiahsb) THEN80 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) )81 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) )82 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) ) * area(:,:) * tms(:,:) )83 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) )84 ENDIF85 !- check conservation (C Rousset)86 ! -------------------------------87 88 78 IF( kt == nit000 ) CALL lim_dyn_init ! Initialization (first time-step only) 89 79 90 80 IF( ln_limdyn ) THEN 91 81 ! 82 ! conservation test 83 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limdyn', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 84 92 85 old_u_ice(:,:) = u_ice(:,:) * tmu(:,:) 93 86 old_v_ice(:,:) = v_ice(:,:) * tmv(:,:) … … 171 164 END DO 172 165 END DO 166 ! 167 ! conservation test 168 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limdyn', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 173 169 ! 174 170 ELSE ! no ice dynamics : transmit directly the atmospheric stress to the ocean … … 224 220 ENDIF 225 221 ! 226 ! -------------------------------227 !- check conservation (C Rousset)228 IF (ln_limdiahsb) THEN229 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b230 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b231 232 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - ( zchk_fw / rhoic )233 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic )234 235 zchk_vmin = glob_min(v_i)236 zchk_amax = glob_max(SUM(a_i,dim=3))237 zchk_amin = glob_min(a_i)238 239 IF(lwp) THEN240 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limdyn) = ',(zchk_v_i * rday)241 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limdyn) = ',(zchk_smv * rday)242 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limdyn) = ',(zchk_vmin * 1.e-3)243 !IF ( zchk_amax > amax+1.e-10 ) WRITE(numout,*) 'violation a_i>amax (limdyn) = ',zchk_amax244 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limdyn) = ',zchk_amin245 ENDIF246 ENDIF247 !- check conservation (C Rousset)248 ! -------------------------------249 250 222 CALL wrk_dealloc( jpi, jpj, zu_io, zv_io ) 251 223 CALL wrk_dealloc( jpj, zind, zmsk ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90
r4634 r4649 26 26 USE dom_ice ! sea-ice domain 27 27 USE in_out_manager ! I/O manager 28 USE lbclnk ! lateral boundary condition - MPP exchanges29 28 USE lib_mpp ! MPP library 30 29 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) … … 390 389 tn_ice (:,:,:) = t_su (:,:,:) 391 390 392 !--------------------------------------------------------------------393 ! 4) Lateral boundary conditions |394 !--------------------------------------------------------------------395 DO jl = 1, jpl396 397 CALL lbc_lnk( a_i(:,:,jl) , 'T', 1. )398 CALL lbc_lnk( v_i(:,:,jl) , 'T', 1. )399 CALL lbc_lnk( v_s(:,:,jl) , 'T', 1. )400 CALL lbc_lnk( smv_i(:,:,jl), 'T', 1. )401 CALL lbc_lnk( oa_i(:,:,jl) , 'T', 1. )402 403 CALL lbc_lnk( ht_i(:,:,jl) , 'T', 1. )404 CALL lbc_lnk( ht_s(:,:,jl) , 'T', 1. )405 CALL lbc_lnk( sm_i(:,:,jl) , 'T', 1. )406 CALL lbc_lnk( o_i(:,:,jl) , 'T', 1. )407 CALL lbc_lnk( t_su(:,:,jl) , 'T', 1. )408 CALL lbc_lnk( tn_ice(:,:,jl) , 'T', 1. )409 DO jk = 1, nlay_s410 CALL lbc_lnk(t_s(:,:,jk,jl), 'T', 1. )411 CALL lbc_lnk(e_s(:,:,jk,jl), 'T', 1. )412 END DO413 DO jk = 1, nlay_i414 CALL lbc_lnk(t_i(:,:,jk,jl), 'T', 1. )415 CALL lbc_lnk(e_i(:,:,jk,jl), 'T', 1. )416 END DO417 !418 ! a_i(:,:,jl) = tms(:,:) * a_i(:,:,jl)419 END DO420 421 391 ELSE 422 392 ! if ln_limini=false … … 449 419 at_i (:,:) = at_i (:,:) + a_i (:,:,jl) 450 420 END DO 451 CALL lbc_lnk( at_i , 'T', 1. )452 ! at_i(:,:) = tms(:,:) * at_i(:,:)453 421 ! 454 422 !-------------------------------------------------------------------- 455 ! 5) Global ice variables for output diagnostics |423 ! 4) Global ice variables for output diagnostics | 456 424 !-------------------------------------------------------------------- 457 425 u_ice (:,:) = 0._wp … … 462 430 463 431 !-------------------------------------------------------------------- 464 ! 6) Moments for advection432 ! 5) Moments for advection 465 433 !-------------------------------------------------------------------- 466 434 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r4634 r4649 22 22 USE limthd_lac ! LIM 23 23 USE limvar ! LIM 24 USE limcons ! LIM25 24 USE in_out_manager ! I/O manager 26 25 USE lbclnk ! lateral boundary condition - MPP exchanges … … 33 32 USE limdiahsb 34 33 USE timing ! Timing 34 USE limcons ! conservation tests 35 35 36 36 IMPLICIT NONE … … 143 143 REAL(wp), POINTER, DIMENSION(:,:) :: esnow_mlt ! energy needed to melt snow in ocean (J m-2) 144 144 REAL(wp), POINTER, DIMENSION(:,:) :: vt_i_init, vt_i_final ! ice volume summed over categories 145 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b146 REAL(wp) :: z chk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset)145 ! 146 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 147 147 !!----------------------------------------------------------------------------- 148 148 IF( nn_timing == 1 ) CALL timing_start('limitd_me') … … 158 158 159 159 IF( ln_limdyn ) THEN ! Start ridging and rafting ! 160 ! ------------------------------- 161 !- check conservation (C Rousset) 162 IF (ln_limdiahsb) THEN 163 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 164 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 165 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 166 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) 167 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 168 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 169 ENDIF 170 !- check conservation (C Rousset) 171 ! ------------------------------- 160 161 ! conservation test 162 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limitd_me', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 172 163 173 164 !-----------------------------------------------------------------------------! … … 355 346 ! 5) Heat, salt and freshwater fluxes 356 347 !-----------------------------------------------------------------------------! 357 wfx_snw(ji,jj) = wfx_snw(ji,jj) -msnow_mlt(ji,jj) * r1_rdtice ! fresh water source for ocean348 wfx_snw(ji,jj) = wfx_snw(ji,jj) + msnow_mlt(ji,jj) * r1_rdtice ! fresh water source for ocean 358 349 hfx_dyn(ji,jj) = hfx_dyn(ji,jj) + esnow_mlt(ji,jj) * unit_fac / area(ji,jj) * r1_rdtice ! heat sink for ocean (<0, W.m-2) 359 350 … … 425 416 ENDIF 426 417 427 ! ------------------------------- 428 !- check conservation (C Rousset) 429 IF (ln_limdiahsb) THEN 430 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 431 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b 432 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b 433 434 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw 435 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic ) 436 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft 437 438 zchk_vmin = glob_min(v_i) 439 zchk_amax = glob_max(SUM(a_i,dim=3)) 440 zchk_amin = glob_min(a_i) 441 442 IF(lwp) THEN 443 IF ( ABS( zchk_v_i ) > 1.e-2 ) WRITE(numout,*) 'violation volume [kg/day] (limitd_me) = ',(zchk_v_i * rday) 444 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limitd_me) = ',(zchk_smv * rday) 445 IF ( ABS( zchk_e_i ) > 1.e-4 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limitd_me) = ',(zchk_e_i) 446 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_me) = ',(zchk_vmin * 1.e-3) 447 IF ( zchk_amax > kamax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_me) = ',zchk_amax 448 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limitd_me) = ',zchk_amin 449 ENDIF 450 ENDIF 451 !- check conservation (C Rousset) 452 ! ------------------------------- 418 ! conservation test 419 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limitd_me', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 453 420 454 421 ENDIF ! ln_limdyn=.true. … … 1108 1075 1109 1076 sfx_dyn(ji,jj) = sfx_dyn(ji,jj) - smsw(ji,jj) * rhoic * r1_rdtice 1110 wfx_dyn(ji,jj) = wfx_dyn(ji,jj) + vsw (ji,jj) * rhoic * r1_rdtice ! gurvan: increase in ice volume du to seawater frozen in voids 1111 ! MV HC 2014 this previous line seems ok, i'm not sure at this moment of the sign convention 1077 wfx_dyn(ji,jj) = wfx_dyn(ji,jj) - vsw (ji,jj) * rhoic * r1_rdtice ! gurvan: increase in ice volume du to seawater frozen in voids 1112 1078 1113 1079 !------------------------------------ … … 1440 1406 CALL wrk_alloc( jpi, jpj, zmask ) 1441 1407 1408 ! to be sure that at_i is the sum of a_i(jl) 1409 at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) 1410 1442 1411 DO jl = 1, jpl 1443 1444 1412 !----------------------------------------------------------------- 1445 1413 ! Count categories to be zapped. 1446 ! Abort model in case of negative area.1447 1414 !----------------------------------------------------------------- 1448 1415 icells = 0 … … 1450 1417 DO jj = 1, jpj 1451 1418 DO ji = 1, jpi 1452 ! IF( ( a_i(ji,jj,jl) >= -epsi10 .AND. a_i(ji,jj,jl) < 0._wp ) .OR. & 1453 ! & ( a_i(ji,jj,jl) > 0._wp .AND. a_i(ji,jj,jl) <= epsi10 ) .OR. & 1454 ! & ( v_i(ji,jj,jl) == 0._wp .AND. a_i(ji,jj,jl) > 0._wp ) .OR. & 1455 ! & ( v_i(ji,jj,jl) > 0._wp .AND. v_i(ji,jj,jl) <= epsi10 ) ) zmask(ji,jj) = 1._wp 1456 IF( ( a_i(ji,jj,jl) >= -epsi10 .AND. a_i(ji,jj,jl) <= epsi10 ) .OR. & 1457 & ( v_i(ji,jj,jl) >= 0._wp .AND. v_i(ji,jj,jl) <= epsi10 ) ) zmask(ji,jj) = 1._wp 1419 IF( a_i(ji,jj,jl) <= epsi10 .OR. v_i(ji,jj,jl) <= epsi10 .OR. at_i(ji,jj) <= epsi10 ) THEN 1420 zmask(ji,jj) = 1._wp 1421 ENDIF 1458 1422 END DO 1459 1423 END DO … … 1470 1434 zei = e_i(ji,jj,jk,jl) 1471 1435 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * ( 1._wp - zmask(ji,jj) ) 1436 t_i(ji,jj,jk,jl) = t_i(ji,jj,jk,jl) * ( 1._wp - zmask(ji,jj) ) + rtt * zmask(ji,jj) 1472 1437 ! update exchanges with ocean 1473 1438 hfx_res(ji,jj) = hfx_res(ji,jj) + ( e_i(ji,jj,jk,jl) - zei ) * unit_fac / area(ji,jj) * r1_rdtice ! W.m-2 <0 … … 1497 1462 ! zap ice and snow volume, add water and salt to ocean 1498 1463 !----------------------------------------------------------------- 1499 1500 ! xtmp = (rhoi*vicen(i,j,n) + rhos*vsnon(i,j,n)) / dt 1501 ! sfx_res(ji,jj) = sfx_res(ji,jj) + ( sss_m(ji,jj) ) & 1502 ! * rhosn * v_s(ji,jj,jl) * r1_rdtice 1503 ! sfx_res(ji,jj) = sfx_res(ji,jj) + ( sss_m(ji,jj) - sm_i(ji,jj,jl) ) & 1504 ! * rhoic * v_i(ji,jj,jl) * r1_rdtice 1505 ! sfx (i,j) = sfx (i,j) + xtmp 1506 1507 ato_i(ji,jj) = a_i (ji,jj,jl) * zmask(ji,jj) + ato_i(ji,jj) 1464 ato_i(ji,jj) = a_i (ji,jj,jl) * zmask(ji,jj) + ato_i(ji,jj) 1508 1465 a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 1509 1466 v_i (ji,jj,jl) = v_i (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) … … 1520 1477 ! update exchanges with ocean 1521 1478 sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsal ) * rhoic * r1_rdtice 1522 wfx_res(ji,jj) = wfx_res(ji,jj) +( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice1523 wfx_snw(ji,jj) = wfx_snw(ji,jj) +( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice1479 wfx_res(ji,jj) = wfx_res(ji,jj) - ( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice 1480 wfx_snw(ji,jj) = wfx_snw(ji,jj) - ( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice 1524 1481 hfx_res(ji,jj) = hfx_res(ji,jj) + ( e_s(ji,jj,1,jl) - zes ) * unit_fac / area(ji,jj) * r1_rdtice ! W.m-2 <0 1525 1482 END DO 1526 1483 END DO 1527 ! 1528 END DO ! jl 1484 END DO ! jl 1485 1486 ! to be sure that at_i is the sum of a_i(jl) 1487 at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) 1529 1488 ! 1530 1489 CALL wrk_dealloc( jpi, jpj, zmask ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90
r4634 r4649 35 35 USE lib_fortran ! to use key_nosignedzero 36 36 USE timing ! Timing 37 USE limcons ! conservation tests 37 38 38 39 IMPLICIT NONE … … 66 67 ! 67 68 INTEGER :: ji,jj, jk, jl, ja, jm, jbnd1, jbnd2 ! ice types dummy loop index 68 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b ! Check conservation (C Rousset)69 REAL(wp) :: z chk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset)69 ! 70 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 70 71 !!------------------------------------------------------------------ 71 72 IF( nn_timing == 1 ) CALL timing_start('limitd_th') 72 73 73 ! ------------------------------- 74 !- check conservation (C Rousset) 75 IF (ln_limdiahsb) THEN 76 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 77 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 78 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 79 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) 80 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 81 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 82 ENDIF 83 !- check conservation (C Rousset) 84 ! ------------------------------- 74 ! conservation test 75 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limitd_th', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 85 76 86 77 IF( kt == nit000 .AND. lwp ) THEN … … 107 98 ! 3) Add frazil ice growing in leads. 108 99 !------------------------------------------------------------------------------| 109 110 100 CALL lim_thd_lac 111 101 CALL lim_var_glo2eqv ! only for info … … 143 133 ENDIF 144 134 ! 145 ! ------------------------------- 146 !- check conservation (C Rousset) 147 IF( ln_limdiahsb ) THEN 148 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 149 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b 150 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b 151 152 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw 153 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic ) 154 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft 155 156 zchk_vmin = glob_min(v_i) 157 zchk_amax = glob_max(SUM(a_i,dim=3)) 158 zchk_amin = glob_min(a_i) 159 160 IF(lwp) THEN 161 IF ( ABS( zchk_v_i ) > 1.e-4 ) WRITE(numout,*) 'violation volume [kg/day] (limitd_th) = ',(zchk_v_i * rday) 162 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limitd_th) = ',(zchk_smv * rday) 163 IF ( ABS( zchk_e_i ) > 1.e-2 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limitd_th) = ',(zchk_e_i) 164 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_th) = ',(zchk_vmin * 1.e-3) 165 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_th) = ',zchk_amax 166 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limitd_th) = ',zchk_amin 167 ENDIF 168 ENDIF 169 !- check conservation (C Rousset) 170 ! ------------------------------- 135 ! conservation test 136 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limitd_th', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 171 137 ! 172 138 IF( nn_timing == 1 ) CALL timing_stop('limitd_th') -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90
r4634 r4649 42 42 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 43 43 USE traqsr ! clem: add penetration of solar flux into the calculation of heat budget 44 USE iom 44 45 45 46 IMPLICIT NONE … … 104 105 INTEGER :: ji, jj, jl, jk ! dummy loop indices 105 106 REAL(wp) :: zinda, zemp ! local scalars 106 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace107 REAL(wp) :: ztmelts ! clem 2014: for HC diags108 109 107 REAL(wp) :: zf_mass ! Heat flux associated with mass exchange ice->ocean (W.m-2) 110 108 REAL(wp) :: zfcm1 ! New solar flux received by the ocean 109 REAL(wp), POINTER, DIMENSION(:,:,:) :: zalb, zalbp ! 2D/3D workspace 111 110 !!--------------------------------------------------------------------- 112 111 113 112 IF( lk_cpl ) CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp ) 113 114 ! make calls for heat fluxes before it is modified 115 CALL iom_put( "qsr_oce" , qsr(:,:) * pfrld(:,:) ) ! solar flux at ocean surface 116 CALL iom_put( "qns_oce" , qns(:,:) * pfrld(:,:) ) ! non-solar flux at ocean surface 117 CALL iom_put( "qsr_ice" , SUM( qsr_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) ) ! solar flux at ice surface 118 CALL iom_put( "qns_ice" , SUM( qns_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) ) ! non-solar flux at ice surface 119 CALL iom_put( "qtr_ice" , SUM( ftr_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) ) ! solar flux transmitted thru ice 120 CALL iom_put( "qt_oce" , ( qsr(:,:) + qns(:,:) ) * pfrld(:,:) ) 121 CALL iom_put( "qt_ice" , SUM( ( qns_ice(:,:,:) + qsr_ice(:,:,:) ) * old_a_i(:,:,:), dim=3 ) ) 114 122 115 123 ! pfrld is the lead fraction at the previous time step (actually between TRP and THD) … … 165 173 IF( lk_cpl ) THEN 166 174 zemp = - emp_tot(ji,jj) + emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! 167 & -wfx_snw(ji,jj)175 & + wfx_snw(ji,jj) 168 176 ELSE 169 177 zemp = emp(ji,jj) * pfrld(ji,jj) & ! evaporation over oceanic fraction … … 176 184 177 185 ! mass flux at the ocean/ice interface 178 fmmflx(ji,jj) = wfx_ice(ji,jj) * rdt_ice ! F/M mass flux save at least for biogeochemical model179 emp(ji,jj) = zemp + wfx_ice(ji,jj) + wfx_snw(ji,jj)! mass flux + F/M mass flux (always ice/ocean mass exchange)186 fmmflx(ji,jj) = - wfx_ice(ji,jj) * rdt_ice ! F/M mass flux save at least for biogeochemical model 187 emp(ji,jj) = zemp - wfx_ice(ji,jj) - wfx_snw(ji,jj) - wfx_sub(ji,jj) ! mass flux + F/M mass flux (always ice/ocean mass exchange) 180 188 181 189 END DO … … 213 221 ENDIF 214 222 215 ! -------------------------------------------------216 ! C. Rousset Begin Diagnostics for heat in W/m2217 ! -------------------------------------------------218 DO jj = 1, jpj219 DO ji = 1, jpi220 diag_heat_dhc1(ji,jj) = ( SUM( d_e_i_trp(ji,jj,1:nlay_i,:) + d_e_i_thd(ji,jj,1:nlay_i,:) ) + &221 & SUM( d_e_s_trp(ji,jj,1:nlay_s,:) + d_e_s_thd(ji,jj,1:nlay_s,:) ) ) * unit_fac * r1_rdtice / area(ji,jj)222 END DO223 END DO224 ! -------------------------------------------------225 ! C. Rousset End Diagnostics226 ! -------------------------------------------------227 223 228 224 IF(ln_ctl) THEN -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r4634 r4649 44 44 USE timing ! Timing 45 45 USE cpl_oasis3, ONLY : lk_cpl 46 USE limcons ! conservation tests 46 47 47 48 IMPLICIT NONE … … 86 87 INTEGER :: nbpb ! nb of icy pts for thermo. cal. 87 88 INTEGER :: ii, ij ! temporary dummy loop index 88 REAL(wp) :: zfric_umin = 5e-03_wp ! lower bound for the friction velocity 89 REAL(wp) :: zfric_umax = 2e-02_wp ! upper bound for the friction velocity 90 REAL(wp) :: zinda, zindb, zfric_u ! local scalar 91 REAL(wp) :: zareamin ! - - 92 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b 93 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 94 REAL(wp) :: zqld, zqfr 89 REAL(wp) :: zfric_umin = 0._wp ! lower bound for the friction velocity (cice value=5.e-04) 90 REAL(wp) :: zch = 0.0057_wp ! heat transfer coefficient 91 REAL(wp) :: zinda, zindb, zareamin 92 REAL(wp) :: zfric_u, zqld, zqfr 95 93 REAL(wp), POINTER, DIMENSION(:) :: zdq, zq_ini, zhfx, zqfx 96 94 REAL(wp) :: zhfx_err, ztest 95 ! 96 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 97 97 !!------------------------------------------------------------------- 98 98 IF( nn_timing == 1 ) CALL timing_start('limthd') … … 103 103 zdq(:) = 0._wp ; zq_ini(:) = 0._wp ; zhfx(:) = 0._wp ; zqfx(:) = 0._wp 104 104 105 ! ------------------------------- 106 !- check conservation (C Rousset) 107 IF (ln_limdiahsb) THEN 108 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 109 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 110 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 111 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) 112 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 113 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 114 ENDIF 115 !- check conservation (C Rousset) 116 ! ------------------------------- 105 ! conservation test 106 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 117 107 118 108 !------------------------------------------------------------------------------! … … 158 148 !CDIR NOVERRCHK 159 149 DO ji = 1, jpi 160 zinda = tms(ji,jj) * ( 1. 0- MAX( 0._wp , SIGN( 1._wp , - at_i(ji,jj) + epsi10 ) ) ) ! 0 if no ice150 zinda = tms(ji,jj) * ( 1._wp - MAX( 0._wp , SIGN( 1._wp , - at_i(ji,jj) + epsi10 ) ) ) ! 0 if no ice 161 151 ! 162 152 ! ! solar irradiance transmission at the mixed layer bottom and used in the lead heat budget … … 165 155 ! ! net downward heat flux from the ice to the ocean, expressed as a function of ocean 166 156 ! ! temperature and turbulent mixing (McPhee, 1992) 167 168 ! friction velocity 169 zfric_u = MAX ( MIN( SQRT( ust2s(ji,jj) ) , zfric_umax ) , zfric_umin ) 170 171 !-- Energy from the turbulent oceanic heat flux. here the drag will depend on ice thickness and type (0.006) 172 fhtur(ji,jj) = zinda * rau0 * rcp * 0.006 * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ! W.m-2 173 ! clem: why not the following? 174 !fhtur(ji,jj) = zinda * rau0 * rcp * 0.006 * SQRT( ust2s(ji,jj) ) * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) 175 176 !-- Energy received in the lead, zqld is defined everywhere (J.m-2) 177 ! It includes turbulent ocean heat flux (only in the leads, the rest is used for bottom melting) 178 zqld = tms(ji,jj) * rdt_ice * & 179 & ( pfrld(ji,jj) * ( qsr(ji,jj) * oatte(ji,jj) & ! solar heat + clem modif 180 & + qns(ji,jj) & ! non solar heat 181 & + fhtur(ji,jj) ) & ! turbulent ice-ocean heat (0 if no ice) 157 ! 158 ! --- Energy received in the lead, zqld is defined everywhere (J.m-2) --- ! 159 zqld = tms(ji,jj) * rdt_ice * & 160 & ( pfrld(ji,jj) * ( qsr(ji,jj) * oatte(ji,jj) & ! solar heat + clem modif 161 & + qns(ji,jj) ) & ! non solar heat 182 162 ! latent heat of precip (note that precip is included in qns but not in qns_ice) 183 163 & + ( pfrld(ji,jj)**betas - pfrld(ji,jj) ) * sprecip(ji,jj) * ( cpic * ( MIN( tatm_ice(ji,jj), rt0_snow ) - rtt ) - lfus ) & 184 164 & + ( 1._wp - pfrld(ji,jj) ) * ( tprecip(ji,jj) - sprecip(ji,jj) ) * rcp * ( tatm_ice(ji,jj) - rtt ) ) 185 165 186 !-- Energy needed to bring ocean surface layer until its freezing (<0, J.m-2) 166 !-- Energy needed to bring ocean surface layer until its freezing (<0, J.m-2) --- ! 187 167 zqfr = tms(ji,jj) * rau0 * rcp * fse3t_m(ji,jj,1) * ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ) 188 168 … … 192 172 ! If there is ice and leads are warming, then transfer energy from the lead budget and use it for bottom melting 193 173 IF( at_i(ji,jj) > epsi10 .AND. zqld > 0._wp ) THEN 194 fhld (ji,jj) = zqld * r1_rdtice / at_i(ji,jj) ! divided by a _i since this is (re)multiplied by a_i in limthd_dh.F90174 fhld (ji,jj) = zqld * r1_rdtice / at_i(ji,jj) ! divided by at_i since this is (re)multiplied by a_i in limthd_dh.F90 195 175 qlead(ji,jj) = 0._wp 196 176 ENDIF 197 177 ! 198 IF( qlead(ji,jj) == 0._wp ) zqld = 0._wp ; zqfr = 0._wp 199 ! 178 !-- Energy from the turbulent oceanic heat flux --- ! 179 !clem zfric_u = MAX ( MIN( SQRT( ust2s(ji,jj) ) , zfric_umax ) , zfric_umin ) 180 zfric_u = MAX( SQRT( ust2s(ji,jj) ), zfric_umin ) 181 fhtur(ji,jj) = MAX( 0._wp, zinda * rau0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ) ! W.m-2 182 ! upper bound for fhtur: we do not want SST to drop below Tfreeze. 183 ! So we say that the heat retrieved from the ocean (fhtur+fhld) must be < to the heat necessary to reach Tfreeze (zqfr) 184 ! This is not a clean budget, so that should be corrected at some point 185 fhtur(ji,jj) = zinda * MIN( fhtur(ji,jj), - fhld(ji,jj) - zqfr * r1_rdtice / MAX( at_i(ji,jj), epsi10 ) ) 186 200 187 ! ----------------------------------------- 201 188 ! Net heat flux on top of ice-ocean [W.m-2] … … 317 304 CALL tab_2d_1d( nbpb, wfx_sum_1d (1:nbpb), wfx_sum , jpi, jpj, npb(1:nbpb) ) 318 305 CALL tab_2d_1d( nbpb, wfx_sni_1d (1:nbpb), wfx_sni , jpi, jpj, npb(1:nbpb) ) 306 CALL tab_2d_1d( nbpb, wfx_res_1d (1:nbpb), wfx_res , jpi, jpj, npb(1:nbpb) ) 307 CALL tab_2d_1d( nbpb, wfx_spr_1d (1:nbpb), wfx_spr , jpi, jpj, npb(1:nbpb) ) 319 308 320 309 CALL tab_2d_1d( nbpb, sfx_bog_1d (1:nbpb), sfx_bog , jpi, jpj, npb(1:nbpb) ) … … 323 312 CALL tab_2d_1d( nbpb, sfx_sni_1d (1:nbpb), sfx_sni , jpi, jpj, npb(1:nbpb) ) 324 313 CALL tab_2d_1d( nbpb, sfx_bri_1d (1:nbpb), sfx_bri , jpi, jpj, npb(1:nbpb) ) 314 CALL tab_2d_1d( nbpb, sfx_res_1d (1:nbpb), sfx_res , jpi, jpj, npb(1:nbpb) ) 325 315 326 316 CALL tab_2d_1d( nbpb, iatte_1d (1:nbpb), iatte , jpi, jpj, npb(1:nbpb) ) … … 329 319 CALL tab_2d_1d( nbpb, hfx_thd_1d (1:nbpb), hfx_thd , jpi, jpj, npb(1:nbpb) ) 330 320 CALL tab_2d_1d( nbpb, hfx_spr_1d (1:nbpb), hfx_spr , jpi, jpj, npb(1:nbpb) ) 331 CALL tab_2d_1d( nbpb, hfx_tot_1d (1:nbpb), hfx_tot , jpi, jpj, npb(1:nbpb) ) 321 CALL tab_2d_1d( nbpb, hfx_sum_1d (1:nbpb), hfx_sum , jpi, jpj, npb(1:nbpb) ) 322 CALL tab_2d_1d( nbpb, hfx_bom_1d (1:nbpb), hfx_bom , jpi, jpj, npb(1:nbpb) ) 323 CALL tab_2d_1d( nbpb, hfx_bog_1d (1:nbpb), hfx_bog , jpi, jpj, npb(1:nbpb) ) 324 CALL tab_2d_1d( nbpb, hfx_dif_1d (1:nbpb), hfx_dif , jpi, jpj, npb(1:nbpb) ) 325 CALL tab_2d_1d( nbpb, hfx_opw_1d (1:nbpb), hfx_opw , jpi, jpj, npb(1:nbpb) ) 332 326 CALL tab_2d_1d( nbpb, hfx_snw_1d (1:nbpb), hfx_snw , jpi, jpj, npb(1:nbpb) ) 333 327 CALL tab_2d_1d( nbpb, hfx_sub_1d (1:nbpb), hfx_sub , jpi, jpj, npb(1:nbpb) ) … … 365 359 ii = MOD( npb(ji) - 1, jpi ) + 1 ; ij = ( npb(ji) - 1 ) / jpi + 1 366 360 hfx_in (ii,ij) = hfx_in (ii,ij) + a_i_b(ji) * ( qsr_ice_1d(ji) + qns_ice_1d(ji) ) 367 368 361 END DO 369 362 … … 379 372 CALL lim_thd_dh( 1, nbpb, jl ) 380 373 381 ! --- Ice /Snowenthalpy remapping --- !382 CALL lim_thd_ent( 1, nbpb, jl )374 ! --- Ice enthalpy remapping --- ! 375 CALL lim_thd_ent( 1, nbpb, jl, q_i_b(1:nbpb,:) ) 383 376 ! 384 377 ! --- diag error on heat remapping - PART 2 --- ! … … 391 384 392 385 !---------------------------------! 393 ! Ice salinity!386 ! --- Ice salinity --- ! 394 387 !---------------------------------! 395 388 CALL lim_thd_sal( 1, nbpb ) 396 389 397 ! CALL lim_thd_enmelt(1,nbpb) ! computes sea ice energy of melting 390 !---------------------------------! 391 ! --- temperature update --- ! 392 !---------------------------------! 393 CALL lim_thd_temp( 1, nbpb ) 394 398 395 !-------------------------------- 399 396 ! 4.4) Move 1D to 2D vectors … … 424 421 CALL tab_1d_2d( nbpb, wfx_sum , npb, wfx_sum_1d(1:nbpb) , jpi, jpj ) 425 422 CALL tab_1d_2d( nbpb, wfx_sni , npb, wfx_sni_1d(1:nbpb) , jpi, jpj ) 423 CALL tab_1d_2d( nbpb, wfx_res , npb, wfx_res_1d(1:nbpb) , jpi, jpj ) 424 CALL tab_1d_2d( nbpb, wfx_spr , npb, wfx_spr_1d(1:nbpb) , jpi, jpj ) 426 425 427 426 CALL tab_1d_2d( nbpb, sfx_bog , npb, sfx_bog_1d(1:nbpb) , jpi, jpj ) … … 429 428 CALL tab_1d_2d( nbpb, sfx_sum , npb, sfx_sum_1d(1:nbpb) , jpi, jpj ) 430 429 CALL tab_1d_2d( nbpb, sfx_sni , npb, sfx_sni_1d(1:nbpb) , jpi, jpj ) 430 CALL tab_1d_2d( nbpb, sfx_res , npb, sfx_res_1d(1:nbpb) , jpi, jpj ) 431 431 ! 432 432 IF( num_sal == 2 ) THEN … … 436 436 CALL tab_1d_2d( nbpb, hfx_thd , npb, hfx_thd_1d(1:nbpb) , jpi, jpj ) 437 437 CALL tab_1d_2d( nbpb, hfx_spr , npb, hfx_spr_1d(1:nbpb) , jpi, jpj ) 438 CALL tab_1d_2d( nbpb, hfx_tot , npb, hfx_tot_1d(1:nbpb) , jpi, jpj ) 438 CALL tab_1d_2d( nbpb, hfx_sum , npb, hfx_sum_1d(1:nbpb) , jpi, jpj ) 439 CALL tab_1d_2d( nbpb, hfx_bom , npb, hfx_bom_1d(1:nbpb) , jpi, jpj ) 440 CALL tab_1d_2d( nbpb, hfx_bog , npb, hfx_bog_1d(1:nbpb) , jpi, jpj ) 441 CALL tab_1d_2d( nbpb, hfx_dif , npb, hfx_dif_1d(1:nbpb) , jpi, jpj ) 442 CALL tab_1d_2d( nbpb, hfx_opw , npb, hfx_opw_1d(1:nbpb) , jpi, jpj ) 439 443 CALL tab_1d_2d( nbpb, hfx_snw , npb, hfx_snw_1d(1:nbpb) , jpi, jpj ) 440 444 CALL tab_1d_2d( nbpb, hfx_sub , npb, hfx_sub_1d(1:nbpb) , jpi, jpj ) … … 521 525 ENDIF 522 526 ! 523 ! ------------------------------- 524 !- check conservation (C Rousset) 525 IF (ln_limdiahsb) THEN 526 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 527 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b 528 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b 529 530 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw 531 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic ) 532 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft 533 534 zchk_vmin = glob_min(v_i) 535 zchk_amax = glob_max(SUM(a_i,dim=3)) 536 zchk_amin = glob_min(a_i) 537 538 IF(lwp) THEN 539 IF ( ABS( zchk_v_i ) > 1.e-4 ) WRITE(numout,*) 'violation volume [kg/day] (limthd) = ',(zchk_v_i * rday) 540 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limthd) = ',(zchk_smv * rday) 541 IF ( ABS( zchk_e_i ) > 1.e-2 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limthd) = ',(zchk_e_i) 542 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limthd) = ',(zchk_vmin * 1.e-3) 543 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limthd) = ',zchk_amax 544 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limthd) = ',zchk_amin 545 ENDIF 546 ENDIF 547 !- check conservation (C Rousset) 548 ! ------------------------------- 527 ! conservation test 528 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 549 529 ! 550 530 CALL wrk_dealloc( jpij, zdq, zq_ini, zhfx, zqfx ) … … 558 538 !! *** ROUTINE lim_thd_enmelt *** 559 539 !! 560 !! ** Purpose : Computes sea ice energy of melting q_i (J.m-3) 540 !! ** Purpose : Computes sea ice energy of melting q_i (J.m-3) from temperature 561 541 !! 562 542 !! ** Method : Formula (Bitz and Lipscomb, 1999) … … 565 545 !! 566 546 INTEGER :: ji, jk ! dummy loop indices 567 REAL(wp) :: ztmelts ! local scalar547 REAL(wp) :: ztmelts, zindb ! local scalar 568 548 !!------------------------------------------------------------------- 569 549 ! 570 550 DO jk = 1, nlay_i ! Sea ice energy of melting 571 551 DO ji = kideb, kiut 572 ztmelts = - tmut * s_i_b(ji,jk) + rtt 573 q_i_b(ji,jk) = rhoic * ( cpic * ( ztmelts - t_i_b(ji,jk) ) & 574 & + lfus * ( 1.0 - (ztmelts-rtt) / MIN( t_i_b(ji,jk)-rtt, -epsi10 ) ) & 575 & - rcp * ( ztmelts-rtt ) ) 552 ztmelts = - tmut * s_i_b(ji,jk) + rtt 553 zindb = MAX( 0._wp , SIGN( 1._wp , -(t_i_b(ji,jk) - rtt) - epsi10 ) ) 554 q_i_b(ji,jk) = rhoic * ( cpic * ( ztmelts - t_i_b(ji,jk) ) & 555 & + lfus * ( 1.0 - zindb * ( ztmelts-rtt ) / MIN( t_i_b(ji,jk)-rtt, -epsi10 ) ) & 556 & - rcp * ( ztmelts-rtt ) ) 576 557 END DO 577 558 END DO … … 584 565 END SUBROUTINE lim_thd_enmelt 585 566 567 SUBROUTINE lim_thd_temp( kideb, kiut ) 568 !!----------------------------------------------------------------------- 569 !! *** ROUTINE lim_thd_temp *** 570 !! 571 !! ** Purpose : Computes sea ice temperature (Kelvin) from enthalpy 572 !! 573 !! ** Method : Formula (Bitz and Lipscomb, 1999) 574 !!------------------------------------------------------------------- 575 INTEGER, INTENT(in) :: kideb, kiut ! bounds for the spatial loop 576 !! 577 INTEGER :: ji, jk ! dummy loop indices 578 REAL(wp) :: ztmelts, zswitch, zaaa, zbbb, zccc, zdiscrim ! local scalar 579 !!------------------------------------------------------------------- 580 ! Recover ice temperature 581 DO jk = 1, nlay_i 582 DO ji = kideb, kiut 583 ztmelts = -tmut * s_i_b(ji,jk) + rtt 584 ! Conversion q(S,T) -> T (second order equation) 585 zaaa = cpic 586 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk) / rhoic - lfus 587 zccc = lfus * ( ztmelts - rtt ) 588 zdiscrim = SQRT( MAX( zbbb * zbbb - 4._wp * zaaa * zccc, 0._wp ) ) 589 t_i_b(ji,jk) = rtt - ( zbbb + zdiscrim ) / ( 2._wp * zaaa ) 590 591 ! mask temperature 592 zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_i_b(ji) ) ) 593 t_i_b(ji,jk) = zswitch * t_i_b(ji,jk) + ( 1._wp - zswitch ) * rtt 594 END DO 595 END DO 596 597 END SUBROUTINE lim_thd_temp 586 598 587 599 SUBROUTINE lim_thd_init -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r4634 r4649 72 72 INTEGER :: ji , jk ! dummy loop indices 73 73 INTEGER :: ii, ij ! 2D corresponding indices to ji 74 INTEGER :: i_ice_switch ! ice thickness above a certain treshold or not75 74 INTEGER :: iter 76 75 … … 114 113 115 114 ! mass and salt flux (clem) 116 REAL(wp) :: zdvres, zswitch_sal 115 REAL(wp) :: zdvres, zswitch_sal, zswitch 117 116 118 117 ! Heat conservation … … 180 179 zq_su (ji) = MAX( 0._wp, zfdum * rdt_ice ) * MAX( 0._wp , SIGN( 1._wp, t_su_b(ji) - ztmelts ) ) 181 180 zq_bo (ji) = MAX( 0._wp, zf_tt(ji) * rdt_ice ) 182 END DO ! ji181 END DO 183 182 184 183 ! … … 192 191 hfx_res_1d(ji) = hfx_res_1d(ji) + q_s_b(ji,1) * ht_s_b(ji) * a_i_b(ji) * r1_rdtice 193 192 ! Contribution to mass flux 194 wfx_snw_1d(ji) = wfx_snw_1d(ji) -rhosn * ht_s_b(ji) * a_i_b(ji) * r1_rdtice193 wfx_snw_1d(ji) = wfx_snw_1d(ji) + rhosn * ht_s_b(ji) * a_i_b(ji) * r1_rdtice 195 194 ! updates 196 195 ht_s_b(ji) = 0._wp … … 247 246 zdh_s_pre(ji) = zcoeff * sprecip_1d(ji) * rdt_ice / rhosn 248 247 ! enthalpy of the precip (>0, J.m-3) (tatm_ice is now in K) 249 zqprec 248 zqprec (ji) = rhosn * ( cpic * ( rtt - MIN( tatm_ice_1d(ji), rt0_snow) ) + lfus ) 250 249 IF( sprecip_1d(ji) == 0._wp ) zqprec(ji) = 0._wp 251 250 ! heat flux from snow precip (>0, W.m-2) 252 251 hfx_spr_1d(ji) = hfx_spr_1d(ji) + zdh_s_pre(ji) * a_i_b(ji) * zqprec(ji) * r1_rdtice 252 ! mass flux, <0 253 wfx_spr_1d(ji) = wfx_spr_1d(ji) - rhosn * a_i_b(ji) * zdh_s_pre(ji) * r1_rdtice 253 254 ! update thickness 254 255 ht_s_b (ji) = MAX( 0._wp , ht_s_b(ji) + zdh_s_pre(ji) ) … … 258 259 !--------------------- 259 260 ! thickness change 261 IF( zdh_s_pre(ji) > 0._wp ) THEN 260 262 zindq = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zqprec(ji) + epsi20 ) ) 261 263 zdh_s_mel (ji) = - zindq * zq_su(ji) / MAX( zqprec(ji) , epsi20 ) 262 264 zdh_s_mel (ji) = MAX( - zdh_s_pre(ji), zdh_s_mel(ji) ) ! bound melting 263 ! Heat flux associated with snow melt of falling snow (W.m-2, <0)264 hfx_snw_1d(ji) = hfx_snw_1d(ji) + zdh_s_mel(ji) * a_i_b(ji) * zqprec(ji) * r1_rdtice265 265 ! heat used to melt snow (W.m-2, >0) 266 hfx_ tot_1d(ji) = hfx_tot_1d(ji) - zdh_s_mel(ji) * a_i_b(ji) * zqprec(ji) * r1_rdtice267 ! snow melting only = water into the ocean (then without snow precip) 268 wfx_snw_1d(ji) = wfx_snw_1d(ji) +rhosn * a_i_b(ji) * zdh_s_mel(ji) * r1_rdtice266 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdh_s_mel(ji) * a_i_b(ji) * zqprec(ji) * r1_rdtice 267 ! snow melting only = water into the ocean (then without snow precip), >0 268 wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_b(ji) * zdh_s_mel(ji) * r1_rdtice 269 269 270 270 ! updates available heat + thickness … … 275 275 ! clem debug: variation of enthalpy (J.m-2) 276 276 dq_s(ji) = dq_s(ji) + ( zdh_s_pre(ji) + zdh_s_mel(ji) ) * zqprec(ji) 277 277 ENDIF 278 278 END DO 279 279 … … 288 288 zdeltah (ji,jk) = MAX( zdeltah(ji,jk) , - zh_s(ji) ) ! bound melting 289 289 zdh_s_mel(ji) = zdh_s_mel(ji) + zdeltah(ji,jk) 290 ! heat flux associated with snow melt(W.m-2, <0)291 hfx_snw_1d(ji) = hfx_snw_1d(ji) + zdeltah(ji,jk) * a_i_b(ji) * q_s_b(ji,jk) * r1_rdtice292 290 ! heat used to melt snow(W.m-2, >0) 293 hfx_ tot_1d(ji) = hfx_tot_1d(ji) - zdeltah(ji,jk) * a_i_b(ji) * q_s_b(ji,jk) * r1_rdtice291 hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,jk) * a_i_b(ji) * q_s_b(ji,jk) * r1_rdtice 294 292 ! snow melting only = water into the ocean (then without snow precip) 295 wfx_snw_1d(ji) = wfx_snw_1d(ji) +rhosn * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice293 wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice 296 294 297 295 ! updates available heat + thickness … … 308 306 !---------------------- 309 307 ! qla_ice is always >=0 (upwards), heat goes to the atmosphere, therefore snow sublimates 308 ! clem comment: not counted in mass exchange in limsbc since this is an exchange with atm. (not ocean) 309 ! clem comment: ice should also sublimate 310 310 IF( lk_cpl ) THEN 311 311 ! coupled mode: sublimation already included in emp_ice (to do in limsbc_ice) … … 320 320 & ( zdh_s_sub(ji) - MAX( zdh_s_sub(ji), - MAX( 0._wp, zdh_s_pre(ji) + zdh_s_mel(ji) ) ) ) * q_s_b(ji,1) ) & 321 321 & * a_i_b(ji) * r1_rdtice 322 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zcoeff ! diag only (to close heat budget) 323 ! heat used for sublimation (>0, W.m-2) 324 !!? hfx_tot_1d(ji) = hfx_tot_1d(ji) - zcoeff 322 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zcoeff 325 323 ! Mass flux by sublimation 326 wfx_sub_1d(ji) = wfx_sub_1d(ji) + rhosn * a_i_b(ji) * zdh_s_sub(ji) * r1_rdtice ! diag only 327 wfx_snw_1d(ji) = wfx_snw_1d(ji) + rhosn * a_i_b(ji) * zdh_s_sub(ji) * r1_rdtice 324 wfx_sub_1d(ji) = wfx_sub_1d(ji) - rhosn * a_i_b(ji) * zdh_s_sub(ji) * r1_rdtice 328 325 ! new snow thickness 329 326 ht_s_b(ji) = MAX( 0._wp , ht_s_b(ji) + zdh_s_sub(ji) ) … … 388 385 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_b(ji) * zEw * r1_rdtice 389 386 390 ! Total heat flux used in this process [W.m-2], <0391 hfx_ tot_1d(ji) = hfx_tot_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice387 ! Total heat flux used in this process [W.m-2], > 0 388 hfx_sum_1d(ji) = hfx_sum_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice 392 389 393 390 ! Contribution to mass flux 394 wfx_sum_1d(ji) = wfx_sum_1d(ji) +rhoic * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice391 wfx_sum_1d(ji) = wfx_sum_1d(ji) - rhoic * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice 395 392 396 393 ! record which layers have disappeared (for bottom melting) … … 424 421 ! Basal growth is driven by heat imbalance at the ice-ocean interface, 425 422 ! between the inner conductive flux (fc_bo_i), from the open water heat flux 426 ! (fhld b) and the turbulent ocean flux (fhtur).423 ! (fhld) and the turbulent ocean flux (fhtur). 427 424 ! fc_bo_i is positive downwards. fhtur and fhld are positive to the ice 428 425 … … 490 487 ! New ice growth 491 488 492 zfmdt = - rhoic * dh_i_bott(ji) ! Mass flux x time step (kg/m2, < 0) 489 zfmdt = - rhoic * dh_i_bott(ji) ! Mass flux x time step (kg/m2, < 0) 490 491 ztmelts = - tmut * s_i_new(ji) + rtt ! New ice melting point (K) 492 493 zt_i_new = zswitch_sal * t_bo_b(ji) + ( 1. - zswitch_sal) * t_i_b(ji, nlay_i) 494 495 zEi = cpic * ( zt_i_new - ztmelts ) & ! Specific enthalpy of forming ice (J/kg, <0) 496 & - lfus * ( 1.0 - ( ztmelts - rtt ) / ( zt_i_new - rtt ) ) & 497 & + rcp * ( ztmelts-rtt ) 498 499 zEw = rcp * ( t_bo_b(ji) - rt0 ) ! Specific enthalpy of seawater (J/kg, < 0) 500 501 zdE = zEi - zEw ! Specific enthalpy difference (J/kg, <0) 493 502 494 503 ! Contribution to heat flux to the ocean [W.m-2], >0 495 504 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_b(ji) * zEw * r1_rdtice 496 ! Total heat flux used in this process [W.m-2] 497 hfx_tot_1d(ji) = hfx_tot_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice 505 506 ! Total heat flux used in this process [W.m-2], <0 507 hfx_bog_1d(ji) = hfx_bog_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice 498 508 499 ! Contribution to salt flux ()509 ! Contribution to salt flux, <0 500 510 sfx_bog_1d(ji) = sfx_bog_1d(ji) + s_i_new(ji) * a_i_b(ji) * zfmdt * r1_rdtice 501 511 502 ! Contribution to mass flux 503 wfx_bog_1d(ji) = wfx_bog_1d(ji) +rhoic * a_i_b(ji) * dh_i_bott(ji) * r1_rdtice512 ! Contribution to mass flux, <0 513 wfx_bog_1d(ji) = wfx_bog_1d(ji) - rhoic * a_i_b(ji) * dh_i_bott(ji) * r1_rdtice 504 514 505 515 ! clem debug: variation of enthalpy (J.m-2) … … 542 552 hfx_res_1d(ji) = hfx_res_1d(ji) + zfmdt * a_i_b(ji) * zEi * r1_rdtice 543 553 554 ! Contribution to salt flux (clem: using sm_i_b and not s_i_b(jk) is ok) 555 sfx_res_1d(ji) = sfx_res_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice 556 557 ! Contribution to mass flux 558 wfx_res_1d(ji) = wfx_res_1d(ji) - rhoic * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice 559 544 560 ! clem debug: variation of enthalpy (J.m-2) 545 561 dq_i(ji) = dq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) … … 573 589 ! Contribution to heat flux to the ocean [W.m-2], <0 574 590 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * a_i_b(ji) * zEw * r1_rdtice 591 592 ! Contribution to salt flux (clem: using sm_i_b and not s_i_b(jk) is ok) 593 sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice 594 595 ! Total heat flux used in this process [W.m-2], >0 596 hfx_bom_1d(ji) = hfx_bom_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice 597 598 ! Contribution to mass flux 599 wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoic * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice 575 600 576 601 ! clem debug: variation of enthalpy (J.m-2) … … 581 606 h_i_old (ji,jk) = h_i_old (ji,jk) + zdeltah(ji,jk) 582 607 ENDIF 583 584 ! Contribution to salt flux (clem: using sm_i_b and not s_i_b(jk) is ok)585 sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdeltah(ji,jk) * rhoic * r1_rdtice586 587 ! Total heat flux used in this process [W.m-2]588 hfx_tot_1d(ji) = hfx_tot_1d(ji) - zfmdt * a_i_b(ji) * zdE * r1_rdtice589 590 ! Contribution to mass flux591 wfx_bom_1d(ji) = wfx_bom_1d(ji) + rhoic * a_i_b(ji) * zdeltah(ji,jk) * r1_rdtice592 608 593 609 ENDIF … … 602 618 ! clem bug: I think this should be included above, so we would not have to 603 619 ! track heat/salt/mass fluxes backwards 604 IF( jpl == 1 ) THEN605 DO ji = kideb, kiut606 IF( zf_tt(ji) >= 0._wp ) THEN607 zdh = MAX( hmelt , dh_i_bott(ji) )608 zdvres = zdh - dh_i_bott(ji) ! >=0609 dh_i_bott(ji) = zdh610 611 ! excessive energy is sent to lateral ablation612 zinda = MAX( 0._wp, SIGN( 1._wp , 1._wp - at_i_b(ji) - epsi20 ) )613 zq_1cat(ji) = zinda * rhoic * lfus * at_i_b(ji) / MAX( 1._wp - at_i_b(ji) , epsi20 ) * zdvres ! J.m-2 >=0614 615 ! correct salt and mass fluxes616 sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdvres * rhoic * r1_rdtice ! this is only a raw approximation617 wfx_bom_1d(ji) = wfx_bom_1d(ji) +rhoic * a_i_b(ji) * zdvres * r1_rdtice618 ENDIF619 END DO620 ENDIF620 ! IF( jpl == 1 ) THEN 621 ! DO ji = kideb, kiut 622 ! IF( zf_tt(ji) >= 0._wp ) THEN 623 ! zdh = MAX( hmelt , dh_i_bott(ji) ) 624 ! zdvres = zdh - dh_i_bott(ji) ! >=0 625 ! dh_i_bott(ji) = zdh 626 ! 627 ! ! excessive energy is sent to lateral ablation 628 ! zinda = MAX( 0._wp, SIGN( 1._wp , 1._wp - at_i_b(ji) - epsi20 ) ) 629 ! zq_1cat(ji) = zinda * rhoic * lfus * at_i_b(ji) / MAX( 1._wp - at_i_b(ji) , epsi20 ) * zdvres ! J.m-2 >=0 630 ! 631 ! ! correct salt and mass fluxes 632 ! sfx_bom_1d(ji) = sfx_bom_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdvres * rhoic * r1_rdtice ! this is only a raw approximation 633 ! wfx_bom_1d(ji) = wfx_bom_1d(ji) - rhoic * a_i_b(ji) * zdvres * r1_rdtice 634 ! ENDIF 635 ! END DO 636 ! ENDIF 621 637 622 638 !------------------------------------------- … … 644 660 ! 645 661 ! zq_rema(ji) = zq_rema(ji) + zdeltah(ji,1) * zq_s(ji) ! update available heat (J.m-2) 646 ! ! Heat flux associated with snow melt647 ! hfx_snw_1d(ji) = hfx_snw_1d(ji) + zdeltah(ji,1) * a_i_b(ji) * zq_s(ji) * r1_rdtice ! W.m-2 (<0)648 662 ! ! heat used to melt snow 649 ! hfx_ tot_1d(ji) = hfx_tot_1d(ji) - zdeltah(ji,1) * a_i_b(ji) * zq_s(ji) * r1_rdtice ! W.m-2 (>0)663 ! hfx_snw_1d(ji) = hfx_snw_1d(ji) - zdeltah(ji,1) * a_i_b(ji) * zq_s(ji) * r1_rdtice ! W.m-2 (>0) 650 664 ! ! Contribution to mass flux 651 ! wfx_snw_1d(ji) = wfx_snw_1d(ji) +rhosn * a_i_b(ji) * zdeltah(ji,1) * r1_rdtice665 ! wfx_snw_1d(ji) = wfx_snw_1d(ji) - rhosn * a_i_b(ji) * zdeltah(ji,1) * r1_rdtice 652 666 ! ! clem debug: variation of enthalpy (J.m-2) 653 667 ! dq_s(ji) = dq_s(ji) + zdeltah(ji,1) * q_s_b(ji,1) … … 680 694 ! new salinity difference stored (to be used in limthd_ent.F90) 681 695 IF ( num_sal == 2 ) THEN 682 i_ice_switch = MAX( 0._wp , SIGN( 1._wp , ht_i_b(ji) - epsi10 ) )696 zswitch = MAX( 0._wp , SIGN( 1._wp , ht_i_b(ji) - epsi10 ) ) 683 697 ! salinity dif due to snow-ice formation 684 dsm_i_si_1d(ji) = ( zs_snic - sm_i_b(ji) ) * dh_snowice(ji) / MAX( ht_i_b(ji), epsi10 ) * i_ice_switch698 dsm_i_si_1d(ji) = ( zs_snic - sm_i_b(ji) ) * dh_snowice(ji) / MAX( ht_i_b(ji), epsi10 ) * zswitch 685 699 ! salinity dif due to bottom growth 686 700 IF ( zf_tt(ji) < 0._wp ) THEN 687 dsm_i_se_1d(ji) = ( s_i_new(ji) - sm_i_b(ji) ) * dh_i_bott(ji) / MAX( ht_i_b(ji), epsi10 ) * i_ice_switch701 dsm_i_se_1d(ji) = ( s_i_new(ji) - sm_i_b(ji) ) * dh_i_bott(ji) / MAX( ht_i_b(ji), epsi10 ) * zswitch 688 702 ENDIF 689 703 ENDIF … … 704 718 ! Contribution to mass flux 705 719 ! All snow is thrown in the ocean, and seawater is taken to replace the volume 706 wfx_sni_1d(ji) = wfx_sni_1d(ji) +a_i_b(ji) * dh_snowice(ji) * rhoic * r1_rdtice707 wfx_snw_1d(ji) = wfx_snw_1d(ji) -a_i_b(ji) * dh_snowice(ji) * rhosn * r1_rdtice720 wfx_sni_1d(ji) = wfx_sni_1d(ji) - a_i_b(ji) * dh_snowice(ji) * rhoic * r1_rdtice 721 wfx_snw_1d(ji) = wfx_snw_1d(ji) + a_i_b(ji) * dh_snowice(ji) * rhosn * r1_rdtice 708 722 709 723 ! clem debug: variation of enthalpy (J.m-2) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_dif.F90
r4634 r4649 32 32 PUBLIC lim_thd_dif ! called by lim_thd 33 33 34 REAL(wp) :: epsi10 =1.e-10_wp !34 REAL(wp) :: epsi10 = 1.e-10_wp ! 35 35 !!---------------------------------------------------------------------- 36 36 !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011) … … 92 92 !! (04-2007) Energy conservation tested by M. Vancoppenolle 93 93 !!------------------------------------------------------------------ 94 INTEGER , INTENT (in) :: kideb ! Start point on which the the computation is applied 95 INTEGER , INTENT (in) :: kiut ! End point on which the the computation is applied 96 INTEGER , INTENT (in) :: jl ! Category number 94 INTEGER , INTENT(in) :: kideb, kiut ! Start/End point on which the the computation is applied 95 INTEGER , INTENT(in) :: jl ! Thickness cateogry number 97 96 98 97 !! * Local variables … … 103 102 INTEGER :: nconv ! number of iterations in iterative procedure 104 103 INTEGER :: minnumeqmin, maxnumeqmax 105 INTEGER, DIMENSION(kiut) :: numeqmin ! reference number of top equation106 INTEGER, DIMENSION(kiut) :: numeqmax ! reference number of bottom equation107 INTEGER, DIMENSION(kiut) :: isnow ! switch for presence (1) or absence (0) of snow104 INTEGER, POINTER, DIMENSION(:) :: numeqmin ! reference number of top equation 105 INTEGER, POINTER, DIMENSION(:) :: numeqmax ! reference number of bottom equation 106 INTEGER, POINTER, DIMENSION(:) :: isnow ! switch for presence (1) or absence (0) of snow 108 107 REAL(wp) :: zg1s = 2._wp ! for the tridiagonal system 109 108 REAL(wp) :: zg1 = 2._wp ! … … 115 114 REAL(wp) :: ztmelt_i ! ice melting temperature 116 115 REAL(wp) :: zerritmax ! current maximal error on temperature 117 REAL(wp), DIMENSION(kiut) :: ztfs ! ice melting point 118 REAL(wp), DIMENSION(kiut) :: ztsuold ! old surface temperature (before the iterative procedure ) 119 REAL(wp), DIMENSION(kiut) :: ztsuoldit ! surface temperature at previous iteration 120 REAL(wp), DIMENSION(kiut) :: zh_i ! ice layer thickness 121 REAL(wp), DIMENSION(kiut) :: zh_s ! snow layer thickness 122 REAL(wp), DIMENSION(kiut) :: zfsw ! solar radiation absorbed at the surface 123 REAL(wp), DIMENSION(kiut) :: zf ! surface flux function 124 REAL(wp), DIMENSION(kiut) :: dzf ! derivative of the surface flux function 125 REAL(wp), DIMENSION(kiut) :: zerrit ! current error on temperature 126 REAL(wp), DIMENSION(kiut) :: zdifcase ! case of the equation resolution (1->4) 127 REAL(wp), DIMENSION(kiut) :: zftrice ! solar radiation transmitted through the ice 128 REAL(wp), DIMENSION(kiut) :: zihic, zhsu 129 REAL(wp), DIMENSION(kiut,0:nlay_i) :: ztcond_i ! Ice thermal conductivity 130 REAL(wp), DIMENSION(kiut,0:nlay_i) :: zradtr_i ! Radiation transmitted through the ice 131 REAL(wp), DIMENSION(kiut,0:nlay_i) :: zradab_i ! Radiation absorbed in the ice 132 REAL(wp), DIMENSION(kiut,0:nlay_i) :: zkappa_i ! Kappa factor in the ice 133 REAL(wp), DIMENSION(kiut,0:nlay_i) :: ztiold ! Old temperature in the ice 134 REAL(wp), DIMENSION(kiut,0:nlay_i) :: zeta_i ! Eta factor in the ice 135 REAL(wp), DIMENSION(kiut,0:nlay_i) :: ztitemp ! Temporary temperature in the ice to check the convergence 136 REAL(wp), DIMENSION(kiut,0:nlay_i) :: zspeche_i ! Ice specific heat 137 REAL(wp), DIMENSION(kiut,0:nlay_i) :: z_i ! Vertical cotes of the layers in the ice 138 REAL(wp), DIMENSION(kiut,0:nlay_s) :: zradtr_s ! Radiation transmited through the snow 139 REAL(wp), DIMENSION(kiut,0:nlay_s) :: zradab_s ! Radiation absorbed in the snow 140 REAL(wp), DIMENSION(kiut,0:nlay_s) :: zkappa_s ! Kappa factor in the snow 141 REAL(wp), DIMENSION(kiut,0:nlay_s) :: zeta_s ! Eta factor in the snow 142 REAL(wp), DIMENSION(kiut,0:nlay_s) :: ztstemp ! Temporary temperature in the snow to check the convergence 143 REAL(wp), DIMENSION(kiut,0:nlay_s) :: ztsold ! Temporary temperature in the snow 144 REAL(wp), DIMENSION(kiut,0:nlay_s) :: z_s ! Vertical cotes of the layers in the snow 145 REAL(wp), DIMENSION(kiut,jkmax+2) :: zindterm ! Independent term 146 REAL(wp), DIMENSION(kiut,jkmax+2) :: zindtbis ! temporary independent term 147 REAL(wp), DIMENSION(kiut,jkmax+2) :: zdiagbis 148 REAL(wp), DIMENSION(kiut,jkmax+2,3) :: ztrid ! tridiagonal system terms 149 REAL(wp) :: ztemp ! local scalar 116 REAL(wp), POINTER, DIMENSION(:) :: ztfs ! ice melting point 117 REAL(wp), POINTER, DIMENSION(:) :: ztsuold ! old surface temperature (before the iterative procedure ) 118 REAL(wp), POINTER, DIMENSION(:) :: ztsuoldit ! surface temperature at previous iteration 119 REAL(wp), POINTER, DIMENSION(:) :: zh_i ! ice layer thickness 120 REAL(wp), POINTER, DIMENSION(:) :: zh_s ! snow layer thickness 121 REAL(wp), POINTER, DIMENSION(:) :: zfsw ! solar radiation absorbed at the surface 122 REAL(wp), POINTER, DIMENSION(:) :: zf ! surface flux function 123 REAL(wp), POINTER, DIMENSION(:) :: dzf ! derivative of the surface flux function 124 REAL(wp), POINTER, DIMENSION(:) :: zerrit ! current error on temperature 125 REAL(wp), POINTER, DIMENSION(:) :: zdifcase ! case of the equation resolution (1->4) 126 REAL(wp), POINTER, DIMENSION(:) :: zftrice ! solar radiation transmitted through the ice 127 REAL(wp), POINTER, DIMENSION(:) :: zihic, zhsu 128 REAL(wp), POINTER, DIMENSION(:,:) :: ztcond_i ! Ice thermal conductivity 129 REAL(wp), POINTER, DIMENSION(:,:) :: zradtr_i ! Radiation transmitted through the ice 130 REAL(wp), POINTER, DIMENSION(:,:) :: zradab_i ! Radiation absorbed in the ice 131 REAL(wp), POINTER, DIMENSION(:,:) :: zkappa_i ! Kappa factor in the ice 132 REAL(wp), POINTER, DIMENSION(:,:) :: ztiold ! Old temperature in the ice 133 REAL(wp), POINTER, DIMENSION(:,:) :: zeta_i ! Eta factor in the ice 134 REAL(wp), POINTER, DIMENSION(:,:) :: ztitemp ! Temporary temperature in the ice to check the convergence 135 REAL(wp), POINTER, DIMENSION(:,:) :: zspeche_i ! Ice specific heat 136 REAL(wp), POINTER, DIMENSION(:,:) :: z_i ! Vertical cotes of the layers in the ice 137 REAL(wp), POINTER, DIMENSION(:,:) :: zradtr_s ! Radiation transmited through the snow 138 REAL(wp), POINTER, DIMENSION(:,:) :: zradab_s ! Radiation absorbed in the snow 139 REAL(wp), POINTER, DIMENSION(:,:) :: zkappa_s ! Kappa factor in the snow 140 REAL(wp), POINTER, DIMENSION(:,:) :: zeta_s ! Eta factor in the snow 141 REAL(wp), POINTER, DIMENSION(:,:) :: ztstemp ! Temporary temperature in the snow to check the convergence 142 REAL(wp), POINTER, DIMENSION(:,:) :: ztsold ! Temporary temperature in the snow 143 REAL(wp), POINTER, DIMENSION(:,:) :: z_s ! Vertical cotes of the layers in the snow 144 REAL(wp), POINTER, DIMENSION(:,:) :: zindterm ! Independent term 145 REAL(wp), POINTER, DIMENSION(:,:) :: zindtbis ! temporary independent term 146 REAL(wp), POINTER, DIMENSION(:,:) :: zdiagbis 147 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrid ! tridiagonal system terms 150 148 !!------------------------------------------------------------------ 151 149 ! 150 CALL wrk_alloc( jpij, numeqmin, numeqmax, isnow ) 151 CALL wrk_alloc( jpij, ztfs, ztsuold, ztsuoldit, zh_i, zh_s, zfsw ) 152 CALL wrk_alloc( jpij, zf, dzf, zerrit, zdifcase, zftrice, zihic, zhsu ) 153 CALL wrk_alloc( jpij, nlay_i+1, ztcond_i, zradtr_i, zradab_i, zkappa_i, ztiold, zeta_i, ztitemp, z_i, zspeche_i, kjstart=0) 154 CALL wrk_alloc( jpij, nlay_s+1, zradtr_s, zradab_s, zkappa_s, ztsold, zeta_s, ztstemp, z_s, kjstart=0) 155 CALL wrk_alloc( jpij, jkmax+2, zindterm, zindtbis, zdiagbis ) 156 CALL wrk_alloc( jpij, jkmax+2, 3, ztrid ) 157 152 158 !------------------------------------------------------------------------------! 153 159 ! 1) Initialization ! … … 318 324 DO layer = 1, nlay_i-1 319 325 DO ji = kideb , kiut 320 ztemp = t_i_b(ji,layer) + t_i_b(ji,layer+1) - 2._wp * rtt 321 ztcond_i(ji,layer) = rcdic + 0.0900_wp * ( s_i_b(ji,layer) + s_i_b(ji,layer+1) ) & 322 & / MIN( -2.0_wp * epsi10, ztemp ) & 323 & - 0.0055_wp * ztemp 326 ztcond_i(ji,layer) = rcdic + 0.090_wp * ( s_i_b(ji,layer) + s_i_b(ji,layer+1) ) & 327 & / MIN(-2.0_wp * epsi10, t_i_b(ji,layer)+t_i_b(ji,layer+1) - 2.0_wp * rtt) & 328 & - 0.0055_wp* ( t_i_b(ji,layer) + t_i_b(ji,layer+1) - 2.0*rtt ) 324 329 ztcond_i(ji,layer) = MAX( ztcond_i(ji,layer), zkimin ) 325 330 END DO 326 331 END DO 327 332 DO ji = kideb , kiut 328 ztemp = t_bo_b(ji) - rtt 329 ztcond_i(ji,nlay_i) = rcdic + 0.090_wp * s_i_b(ji,nlay_i) / MIN( -epsi10, ztemp ) & 330 & - 0.011_wp * ztemp 333 ztcond_i(ji,nlay_i) = rcdic + 0.090_wp * s_i_b(ji,nlay_i) / MIN(-epsi10,t_bo_b(ji)-rtt) & 334 & - 0.011_wp * ( t_bo_b(ji) - rtt ) 331 335 ztcond_i(ji,nlay_i) = MAX( ztcond_i(ji,nlay_i), zkimin ) 332 336 END DO … … 395 399 ! 396 400 DO ji = kideb , kiut 397 398 401 ! update of the non solar flux according to the update in T_su 399 402 qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( t_su_b(ji) - ztsuoldit(ji) ) … … 403 406 + qns_ice_1d(ji) ! non solar total flux 404 407 ! (LWup, LWdw, SH, LH) 405 406 ! heat flux used to change surface temperature407 !hfx_tot_1d(ji) = hfx_tot_1d(ji) + dqns_ice_1d(ji) * ( t_su_b(ji) - ztsuoldit(ji) ) * a_i_b(ji)408 408 END DO 409 409 … … 721 721 DO ji = kideb, kiut 722 722 IF( t_su_b(ji) < rtt ) THEN ! case T_su < 0degC 723 hfx_ tot_1d(ji) = hfx_tot_1d(ji) + ( qns_ice_1d(ji) + qsr_ice_1d(ji) - zradtr_i(ji,nlay_i) - fc_bo_i(ji) ) * a_i_b(ji)723 hfx_dif_1d(ji) = hfx_dif_1d(ji) + ( qns_ice_1d(ji) + qsr_ice_1d(ji) - zradtr_i(ji,nlay_i) - fc_bo_i(ji) ) * a_i_b(ji) 724 724 ELSE ! case T_su = 0degC 725 hfx_ tot_1d(ji) = hfx_tot_1d(ji) + ( fc_su(ji) + i0(ji) * qsr_ice_1d(ji) - zradtr_i(ji,nlay_i) - fc_bo_i(ji) ) * a_i_b(ji)725 hfx_dif_1d(ji) = hfx_dif_1d(ji) + ( fc_su(ji) + i0(ji) * qsr_ice_1d(ji) - zradtr_i(ji,nlay_i) - fc_bo_i(ji) ) * a_i_b(ji) 726 726 ENDIF 727 727 END DO 728 728 729 729 ! 730 CALL wrk_dealloc( jpij, numeqmin, numeqmax, isnow ) 731 CALL wrk_dealloc( jpij, ztfs, ztsuold, ztsuoldit, zh_i, zh_s, zfsw ) 732 CALL wrk_dealloc( jpij, zf, dzf, zerrit, zdifcase, zftrice, zihic, zhsu ) 733 CALL wrk_dealloc( jpij, nlay_i+1, ztcond_i, zradtr_i, zradab_i, zkappa_i, ztiold, zeta_i, ztitemp, z_i, zspeche_i, kjstart = 0 ) 734 CALL wrk_dealloc( jpij, nlay_s+1, zradtr_s, zradab_s, zkappa_s, ztsold, zeta_s, ztstemp, z_s, kjstart = 0 ) 735 CALL wrk_dealloc( jpij, jkmax+2, zindterm, zindtbis, zdiagbis ) 736 CALL wrk_dealloc( jpij, jkmax+2, 3, ztrid ) 737 730 738 END SUBROUTINE lim_thd_dif 731 739 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_ent.F90
r4634 r4649 48 48 CONTAINS 49 49 50 SUBROUTINE lim_thd_ent( kideb, kiut, jl )50 SUBROUTINE lim_thd_ent( kideb, kiut, jl, qnew ) 51 51 !!------------------------------------------------------------------- 52 52 !! *** ROUTINE lim_thd_ent *** … … 62 62 !! ** Steps : 1) cumulative integrals of old enthalpies/thicknesses 63 63 !! 2) linear remapping on the new layers 64 !! 3) Ice salinity update + recover temperature from enthalpies 64 !! 65 !! ------------ cum0(0) ------------- cum1(0) 66 !! NEW ------------- 67 !! ------------ cum0(1) ==> ------------- 68 !! ... ------------- 69 !! ------------ ------------- 70 !! ------------ cum0(nlay_i+2) ------------- cum1(nlay_i) 71 !! 65 72 !! 66 73 !! References : Bitz & Lipscomb, JGR 99; Vancoppenolle et al., GRL, 2005 … … 68 75 INTEGER , INTENT(in) :: kideb, kiut ! Start/End point on which the the computation is applied 69 76 INTEGER , INTENT(in) :: jl ! Thickness cateogry number 77 78 REAL(wp), INTENT(inout), DIMENSION(:,:) :: qnew ! new enthlapies (remapped) 70 79 71 80 INTEGER :: ji,ii,ij ! dummy loop indices … … 76 85 REAL(wp), POINTER, DIMENSION(:,:) :: zqh_cum0, zh_cum0 ! old cumulative enthlapies and layers interfaces 77 86 REAL(wp), POINTER, DIMENSION(:,:) :: zqh_cum1, zh_cum1 ! new cumulative enthlapies and layers interfaces 87 REAL(wp), POINTER, DIMENSION(:) :: zhnew ! new layers thicknesses 78 88 !!------------------------------------------------------------------- 79 89 80 90 CALL wrk_alloc( jpij, nlay_i+3, zqh_cum0, zh_cum0, kjstart = 0 ) 81 91 CALL wrk_alloc( jpij, nlay_i+1, zqh_cum1, zh_cum1, kjstart = 0 ) 92 CALL wrk_alloc( jpij, zhnew ) 82 93 83 94 !-------------------------------------------------------------------------- … … 96 107 ! 2) Interpolation on the new layers 97 108 !------------------------------------ 109 ! new layer thickesses 110 DO ji = kideb, kiut 111 zhnew(ji) = SUM( h_i_old(ji,0:nlay_i+1) ) / REAL( nlay_i ) 112 ENDDO 113 98 114 ! new layers interfaces 99 115 zh_cum1(:,0:nlay_i) = 0._wp 100 116 DO jk1 = 1, nlay_i 101 117 DO ji = kideb, kiut 102 zh_cum1(ji,jk1) = zh_cum1(ji,jk1-1) + ht_i_b(ji) / REAL( nlay_i)118 zh_cum1(ji,jk1) = zh_cum1(ji,jk1-1) + zhnew(ji) 103 119 ENDDO 104 120 ENDDO … … 123 139 DO jk1 = 1, nlay_i 124 140 DO ji = kideb, kiut 125 zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_i_b(ji) + epsi20 ) )126 q _i_b(ji,jk1) = zswitch * ( zqh_cum1(ji,jk1) - zqh_cum1(ji,jk1-1) ) * REAL( nlay_i ) / MAX( ht_i_b(ji), epsi20 )141 zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zhnew(ji) + epsi10 ) ) 142 qnew(ji,jk1) = zswitch * ( zqh_cum1(ji,jk1) - zqh_cum1(ji,jk1-1) ) / MAX( zhnew(ji), epsi10 ) 127 143 ENDDO 128 144 ENDDO 129 130 !---------------------------------------------------------131 ! 3) Update ice salinity and recover ice temperature132 !---------------------------------------------------------133 ! Update salinity (basal entrapment, snow ice formation)134 DO ji = kideb, kiut135 sm_i_b(ji) = sm_i_b(ji) + dsm_i_se_1d(ji) + dsm_i_si_1d(ji)136 END DO !ji137 138 ! Recover ice temperature139 DO jk1 = 1, nlay_i140 DO ji = kideb, kiut141 ztmelts = -tmut * s_i_b(ji,jk1) + rtt142 ! Conversion q(S,T) -> T (second order equation)143 zaaa = cpic144 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk1) / rhoic - lfus145 zccc = lfus * ( ztmelts - rtt )146 zdiscrim = SQRT( MAX( zbbb * zbbb - 4._wp * zaaa * zccc, 0._wp ) )147 t_i_b(ji,jk1) = rtt - ( zbbb + zdiscrim ) / ( 2._wp * zaaa )148 149 ! mask temperature150 zswitch = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_i_b(ji) ) )151 t_i_b(ji,jk1) = zswitch * t_i_b(ji,jk1) + ( 1._wp - zswitch ) * rtt152 END DO153 END DO154 155 145 ! 156 146 CALL wrk_dealloc( jpij, nlay_i+3, zqh_cum0, zh_cum0, kjstart = 0 ) 157 147 CALL wrk_dealloc( jpij, nlay_i+1, zqh_cum1, zh_cum1, kjstart = 0 ) 148 CALL wrk_dealloc( jpij, zhnew ) 158 149 ! 159 150 END SUBROUTINE lim_thd_ent -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90
r4634 r4649 30 30 USE wrk_nemo ! work arrays 31 31 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 32 USE limthd_ent 32 33 33 34 IMPLICIT NONE … … 85 86 REAL(wp) :: zEw ! seawater specific enthalpy (J/kg) 86 87 REAL(wp) :: zfmdt ! mass flux x time step (kg/m2, >0 towards ocean) 87 88 INTEGER , POINTER, DIMENSION(:) :: zcatac ! indexes of categories where new ice grows 88 89 REAL(wp) :: zv_newfra 90 91 INTEGER , POINTER, DIMENSION(:) :: jcat ! indexes of categories where new ice grows 89 92 REAL(wp), POINTER, DIMENSION(:) :: zswinew ! switch for new ice or not 90 93 … … 97 100 REAL(wp), POINTER, DIMENSION(:) :: zdv_res ! residual volume in case of excessive heat budget 98 101 REAL(wp), POINTER, DIMENSION(:) :: zda_res ! residual area in case of excessive heat budget 99 REAL(wp), POINTER, DIMENSION(:) :: zat_i_ ac! total ice fraction102 REAL(wp), POINTER, DIMENSION(:) :: zat_i_1d ! total ice fraction 100 103 REAL(wp), POINTER, DIMENSION(:) :: zat_i_lev ! total ice fraction for level ice only (type 1) 101 104 REAL(wp), POINTER, DIMENSION(:) :: zv_frazb ! accretion of frazil ice at the ice bottom 102 REAL(wp), POINTER, DIMENSION(:) :: zvrel_ ac! relative ice / frazil velocity (1D vector)105 REAL(wp), POINTER, DIMENSION(:) :: zvrel_1d ! relative ice / frazil velocity (1D vector) 103 106 104 107 REAL(wp), POINTER, DIMENSION(:,:) :: zv_old ! old volume of ice in category jl 105 108 REAL(wp), POINTER, DIMENSION(:,:) :: za_old ! old area of ice in category jl 106 REAL(wp), POINTER, DIMENSION(:,:) :: za_i_ac ! 1-D version of a_i 107 REAL(wp), POINTER, DIMENSION(:,:) :: zv_i_ac ! 1-D version of v_i 108 REAL(wp), POINTER, DIMENSION(:,:) :: zoa_i_ac ! 1-D version of oa_i 109 REAL(wp), POINTER, DIMENSION(:,:) :: zsmv_i_ac ! 1-D version of smv_i 110 111 REAL(wp), POINTER, DIMENSION(:,:,:) :: ze_i_ac !: 1-D version of e_i 112 113 REAL(wp), POINTER, DIMENSION(:,:,:) :: zqm0 ! old layer-system heat content 114 REAL(wp), POINTER, DIMENSION(:,:,:) :: zthick0 ! old ice thickness 115 116 REAL(wp), POINTER, DIMENSION(:,:) :: vt_i_init, vt_i_final ! ice volume summed over categories 117 REAL(wp), POINTER, DIMENSION(:,:) :: vt_s_init, vt_s_final ! snow volume summed over categories 118 REAL(wp), POINTER, DIMENSION(:,:) :: et_i_init, et_i_final ! ice energy summed over categories 119 REAL(wp), POINTER, DIMENSION(:,:) :: et_s_init ! snow energy summed over categories 109 REAL(wp), POINTER, DIMENSION(:,:) :: za_i_1d ! 1-D version of a_i 110 REAL(wp), POINTER, DIMENSION(:,:) :: zv_i_1d ! 1-D version of v_i 111 REAL(wp), POINTER, DIMENSION(:,:) :: zoa_i_1d ! 1-D version of oa_i 112 REAL(wp), POINTER, DIMENSION(:,:) :: zsmv_i_1d ! 1-D version of smv_i 113 114 REAL(wp), POINTER, DIMENSION(:,:,:) :: ze_i_1d !: 1-D version of e_i 115 120 116 REAL(wp), POINTER, DIMENSION(:,:) :: zvrel ! relative ice / frazil velocity 121 117 !!-----------------------------------------------------------------------! 122 118 123 CALL wrk_alloc( jpij, zcatac) ! integer119 CALL wrk_alloc( jpij, jcat ) ! integer 124 120 CALL wrk_alloc( jpij, zswinew, zv_newice, za_newice, zh_newice, ze_newice, zs_newice, zo_newice ) 125 CALL wrk_alloc( jpij, zdv_res, zda_res, zat_i_ac, zat_i_lev, zv_frazb, zvrel_ac ) 126 CALL wrk_alloc( jpij,jpl, zv_old, za_old, za_i_ac, zv_i_ac, zoa_i_ac, zsmv_i_ac ) 127 CALL wrk_alloc( jpij,jkmax,jpl, ze_i_ac ) 128 CALL wrk_alloc( jpij,jkmax+1,jpl, zqm0, zthick0 ) 129 CALL wrk_alloc( jpi,jpj, vt_i_init, vt_i_final, vt_s_init, vt_s_final, et_i_init, et_i_final, et_s_init, zvrel ) 130 131 et_i_init(:,:) = 0._wp 132 et_s_init(:,:) = 0._wp 133 vt_i_init(:,:) = 0._wp 134 vt_s_init(:,:) = 0._wp 135 136 !------------------------------------------------------------------------------! 137 ! 1) Conservation check and changes in each ice category 138 !------------------------------------------------------------------------------! 139 IF( con_i ) THEN 140 CALL lim_column_sum ( jpl, v_i , vt_i_init) 141 CALL lim_column_sum ( jpl, v_s , vt_s_init) 142 CALL lim_column_sum_energy ( jpl, nlay_i , e_i , et_i_init) 143 CALL lim_column_sum ( jpl, e_s(:,:,1,:) , et_s_init) 144 ENDIF 121 CALL wrk_alloc( jpij, zdv_res, zda_res, zat_i_1d, zat_i_lev, zv_frazb, zvrel_1d ) 122 CALL wrk_alloc( jpij,jpl, zv_old, za_old, za_i_1d, zv_i_1d, zoa_i_1d, zsmv_i_1d ) 123 CALL wrk_alloc( jpij,jkmax,jpl, ze_i_1d ) 124 CALL wrk_alloc( jpi,jpj, zvrel ) 145 125 146 126 !------------------------------------------------------------------------------| … … 204 184 & + vtau_ice(ji ,jj ) * tmv(ji ,jj ) ) * 0.5_wp 205 185 ! Square root of wind stress 206 ztenagm = SQRT( SQRT( ztaux * ztaux + ztauy * ztauy) )186 ztenagm = SQRT( SQRT( ztaux**2 + ztauy**2 ) ) 207 187 208 188 !--------------------- 209 189 ! Frazil ice velocity 210 190 !--------------------- 211 zvfrx = zgamafr * zsqcd * ztaux / MAX(ztenagm,epsi10) 212 zvfry = zgamafr * zsqcd * ztauy / MAX(ztenagm,epsi10) 191 zindb = MAX( 0._wp, SIGN( 1._wp , ztenagm - epsi10 ) ) 192 zvfrx = zindb * zgamafr * zsqcd * ztaux / MAX( ztenagm, epsi10 ) 193 zvfry = zindb * zgamafr * zsqcd * ztauy / MAX( ztenagm, epsi10 ) 213 194 214 195 !------------------- … … 305 286 IF ( nbpac > 0 ) THEN 306 287 307 CALL tab_2d_1d( nbpac, zat_i_ ac(1:nbpac) , at_i , jpi, jpj, npac(1:nbpac) )288 CALL tab_2d_1d( nbpac, zat_i_1d (1:nbpac) , at_i , jpi, jpj, npac(1:nbpac) ) 308 289 DO jl = 1, jpl 309 CALL tab_2d_1d( nbpac, za_i_ ac(1:nbpac,jl), a_i (:,:,jl), jpi, jpj, npac(1:nbpac) )310 CALL tab_2d_1d( nbpac, zv_i_ ac(1:nbpac,jl), v_i (:,:,jl), jpi, jpj, npac(1:nbpac) )311 CALL tab_2d_1d( nbpac, zoa_i_ ac(1:nbpac,jl), oa_i (:,:,jl), jpi, jpj, npac(1:nbpac) )312 CALL tab_2d_1d( nbpac, zsmv_i_ ac(1:nbpac,jl), smv_i(:,:,jl), jpi, jpj, npac(1:nbpac) )290 CALL tab_2d_1d( nbpac, za_i_1d (1:nbpac,jl), a_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 291 CALL tab_2d_1d( nbpac, zv_i_1d (1:nbpac,jl), v_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 292 CALL tab_2d_1d( nbpac, zoa_i_1d (1:nbpac,jl), oa_i (:,:,jl), jpi, jpj, npac(1:nbpac) ) 293 CALL tab_2d_1d( nbpac, zsmv_i_1d(1:nbpac,jl), smv_i(:,:,jl), jpi, jpj, npac(1:nbpac) ) 313 294 DO jk = 1, nlay_i 314 CALL tab_2d_1d( nbpac, ze_i_ ac(1:nbpac,jk,jl), e_i(:,:,jk,jl) , jpi, jpj, npac(1:nbpac) )295 CALL tab_2d_1d( nbpac, ze_i_1d(1:nbpac,jk,jl), e_i(:,:,jk,jl) , jpi, jpj, npac(1:nbpac) ) 315 296 END DO ! jk 316 297 END DO ! jl … … 322 303 CALL tab_2d_1d( nbpac, wfx_opw_1d(1:nbpac) , wfx_opw, jpi, jpj, npac(1:nbpac) ) 323 304 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , jpi, jpj, npac(1:nbpac) ) 324 CALL tab_2d_1d( nbpac, zvrel_ ac(1:nbpac) , zvrel , jpi, jpj, npac(1:nbpac) )305 CALL tab_2d_1d( nbpac, zvrel_1d (1:nbpac) , zvrel , jpi, jpj, npac(1:nbpac) ) 325 306 326 307 CALL tab_2d_1d( nbpac, hfx_thd_1d(1:nbpac) , hfx_thd, jpi, jpj, npac(1:nbpac) ) 327 CALL tab_2d_1d( nbpac, hfx_ tot_1d(1:nbpac) , hfx_tot, jpi, jpj, npac(1:nbpac) )308 CALL tab_2d_1d( nbpac, hfx_opw_1d(1:nbpac) , hfx_opw, jpi, jpj, npac(1:nbpac) ) 328 309 329 310 !------------------------------------------------------------------------------! … … 334 315 ! Keep old ice areas and volume in memory 335 316 !----------------------------------------- 336 zv_old(:,:) = zv_i_ ac(:,:)337 za_old(:,:) = za_i_ ac(:,:)317 zv_old(:,:) = zv_i_1d(:,:) 318 za_old(:,:) = za_i_1d(:,:) 338 319 339 320 !---------------------- … … 343 324 zh_newice(ji) = hiccrit(1) 344 325 END DO 345 IF( fraz_swi == 1.0 ) 326 IF( fraz_swi == 1.0 ) zh_newice(:) = hicol_b(:) 346 327 347 328 !---------------------- … … 370 351 ztmelts = - tmut * zs_newice(ji) + rtt ! Melting point (K) 371 352 ze_newice(ji) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 372 & + lfus * ( 1.0 - ( ztmelts - rtt ) / ( t_bo_b(ji) - rtt) ) &353 & + lfus * ( 1.0 - ( ztmelts - rtt ) / MIN( t_bo_b(ji) - rtt, -epsi10 ) ) & 373 354 & - rcp * ( ztmelts - rtt ) ) 374 ! MV HC 2014 comment I dont see why this line below is here... ?375 ! This implies that ze_newice gets to rhoic*Lfus if it was negative, but this should never happen376 ze_newice(ji) = MAX( ze_newice(ji) , 0._wp ) &377 & + MAX( 0.0 , SIGN( 1.0 , - ze_newice(ji) ) ) * rhoic * lfus378 355 END DO ! ji 379 356 !---------------- … … 405 382 hfx_thd_1d(ji) = hfx_thd_1d(ji) + zfmdt * zEw * r1_rdtice 406 383 ! Total heat flux used in this process [W.m-2] 407 hfx_ tot_1d(ji) = hfx_tot_1d(ji) - zfmdt * zdE * r1_rdtice384 hfx_opw_1d(ji) = hfx_opw_1d(ji) - zfmdt * zdE * r1_rdtice 408 385 ! mass flux 409 wfx_opw_1d(ji) = wfx_opw_1d(ji) +zv_newice(ji) * rhoic * r1_rdtice386 wfx_opw_1d(ji) = wfx_opw_1d(ji) - zv_newice(ji) * rhoic * r1_rdtice 410 387 ! salt flux 411 388 sfx_opw_1d(ji) = sfx_opw_1d(ji) - zv_newice(ji) * rhoic * zs_newice(ji) * r1_rdtice 412 389 413 390 ! A fraction zfrazb of frazil ice is accreted at the ice bottom 414 zfrazb = ( TANH ( Cfrazb * ( zvrel_ac(ji) - vfrazb ) ) + 1.0 ) * 0.5 * maxfrazb 391 zinda = 1._wp - MAX( 0._wp, SIGN( 1._wp , - zat_i_1d(ji) ) ) 392 zfrazb = zinda * ( TANH ( Cfrazb * ( zvrel_1d(ji) - vfrazb ) ) + 1.0 ) * 0.5 * maxfrazb 415 393 zv_frazb(ji) = zfrazb * zv_newice(ji) 416 394 zv_newice(ji) = ( 1.0 - zfrazb ) * zv_newice(ji) 417 395 END DO 418 396 419 !------------------------------------420 ! Diags for energy conservation test421 !------------------------------------422 DO ji = 1, nbpac423 ii = MOD( npac(ji) - 1 , jpi ) + 1424 ij = ( npac(ji) - 1 ) / jpi + 1425 !426 zde = ze_newice(ji) / unit_fac * area(ii,ij) * zv_newice(ji)427 !zde = ze_newice(ji) * area(ii,ij) * zv_newice(ji)428 !429 ! clem: change that?430 vt_i_init(ii,ij) = vt_i_init(ii,ij) + zv_newice(ji) ! volume431 et_i_init(ii,ij) = et_i_init(ii,ij) + zde ! Energy432 433 END DO434 397 435 398 !----------------- … … 438 401 DO ji = 1, nbpac 439 402 za_newice(ji) = zv_newice(ji) / zh_newice(ji) 440 END DO !ji403 END DO 441 404 442 405 !------------------------------------------------------------------------------! … … 450 413 ! we keep the excessive volume in memory and attribute it later to bottom accretion 451 414 DO ji = 1, nbpac 452 IF ( za_newice(ji) > ( amax - zat_i_ ac(ji) ) ) THEN453 zda_res(ji) = za_newice(ji) - ( amax - zat_i_ ac(ji) )415 IF ( za_newice(ji) > ( amax - zat_i_1d(ji) ) ) THEN 416 zda_res(ji) = za_newice(ji) - ( amax - zat_i_1d(ji) ) 454 417 zdv_res(ji) = zda_res (ji) * zh_newice(ji) 455 418 za_newice(ji) = za_newice(ji) - zda_res (ji) … … 464 427 ! Laterally redistribute new ice volume and area 465 428 !------------------------------------------------ 466 zat_i_ ac(:) = 0._wp429 zat_i_1d(:) = 0._wp 467 430 DO jl = 1, jpl 468 431 DO ji = 1, nbpac 469 IF( hi_max (jl-1) < zh_newice(ji) .AND. & 470 & zh_newice(ji) <= hi_max (jl) ) THEN 471 za_i_ac (ji,jl) = za_i_ac (ji,jl) + za_newice(ji) 472 zv_i_ac (ji,jl) = zv_i_ac (ji,jl) + zv_newice(ji) 473 zcatac (ji) = jl 432 IF( zh_newice(ji) > hi_max(jl-1) .AND. zh_newice(ji) <= hi_max(jl) ) THEN 433 za_i_1d (ji,jl) = za_i_1d (ji,jl) + za_newice(ji) 434 zv_i_1d (ji,jl) = zv_i_1d (ji,jl) + zv_newice(ji) 435 jcat (ji) = jl 474 436 ENDIF 475 zat_i_ ac(ji) = zat_i_ac(ji) + za_i_ac(ji,jl)437 zat_i_1d(ji) = zat_i_1d(ji) + za_i_1d (ji,jl) 476 438 END DO 477 439 END DO … … 481 443 !---------------------------------- 482 444 DO ji = 1, nbpac 483 jl = zcatac(ji)! categroy in which new ice is put484 zswinew (ji) = MAX( 0._wp , SIGN( 1._wp , - za_old(ji,jl) + epsi10) ) ! 0 if old ice445 jl = jcat(ji) ! categroy in which new ice is put 446 zswinew (ji) = MAX( 0._wp , SIGN( 1._wp , - za_old(ji,jl) ) ) ! 0 if old ice 485 447 END DO 486 448 487 449 DO jk = 1, nlay_i 488 450 DO ji = 1, nbpac 489 jl = zcatac(ji) 490 ze_i_ac(ji,jk,jl) = zswinew(ji) * ze_newice(ji) & 491 & + ( 1.0 - zswinew(ji) ) * ( ze_newice(ji) * zv_newice(ji) + ze_i_ac(ji,jk,jl) * zv_old(ji,jl) ) / zv_i_ac(ji,jl) 451 jl = jcat(ji) 452 zinda = MAX( 0._wp, SIGN( 1._wp , zv_i_1d(ji,jl) - epsi20 ) ) 453 ze_i_1d(ji,jk,jl) = zswinew(ji) * ze_newice(ji) + & 454 & ( 1.0 - zswinew(ji) ) * ( ze_newice(ji) * zv_newice(ji) + ze_i_1d(ji,jk,jl) * zv_old(ji,jl) ) & 455 & * zinda / MAX( zv_i_1d(ji,jl), epsi20 ) 492 456 END DO 493 457 END DO … … 496 460 ! Add excessive volume of new ice at the bottom 497 461 !----------------------------------------------- 498 jm = 1 499 500 ! --- Redistributing energy on the new grid (energy is equally distributed in every layer) --- ! 501 ! DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 502 ! DO jk = 1, nlay_i 503 ! DO ji = 1, nbpac 504 ! ze_i_ac(ji,jk,jl) = ( ze_i_ac(ji,jk,jl) * zv_i_ac(ji,jl) + ze_newice(ji) * ( zdv_res(ji) + zv_frazb(ji) ) ) / & 505 ! & ( zv_i_ac(ji,jl) + ( zdv_res(ji) + zv_frazb(ji) ) ) 506 ! END DO 507 ! END DO 508 ! END DO 509 510 ! --- Redistributing energy on the new grid (energy is sent to the bottom) PART 1 --- ! 511 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 462 DO jl = 1, jpl 463 h_i_old (1:nbpac,0:nlay_i+1) = 0._wp 464 qh_i_old(1:nbpac,0:nlay_i+1) = 0._wp 465 512 466 DO jk = 1, nlay_i 513 467 DO ji = 1, nbpac 514 zthick0(ji,jk,jl) = zv_i_ac(ji,jl) / REAL( nlay_i )515 zqm0 (ji,jk,jl) = ze_i_ac(ji,jk,jl) * zthick0(ji,jk,jl)468 h_i_old (ji,jk) = zv_i_1d(ji,jl) / REAL( nlay_i ) 469 qh_i_old(ji,jk) = ze_i_1d(ji,jk,jl) * h_i_old(ji,jk) 516 470 END DO 517 471 END DO 518 END DO 519 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 472 520 473 DO ji = 1, nbpac 521 zinda = MAX( 0._wp, SIGN( 1._wp , zat_i_ac(ji) - epsi10 ) ) 522 zthick0(ji,nlay_i+1,jl) = zinda * ( zdv_res(ji) + zv_frazb(ji) ) * za_i_ac(ji,jl) / MAX( zat_i_ac(ji) , epsi10 ) 523 zqm0 (ji,nlay_i+1,jl) = ze_newice(ji) * zthick0(ji,nlay_i+1,jl) 524 END DO ! ji 525 END DO ! jl 526 527 ze_i_ac(:,:,:) = 0._wp 528 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 529 DO jk = 1, nlay_i 530 DO layer = 1, nlay_i + 1 531 DO ji = 1, nbpac 532 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zthick0(ji,layer,jl) + epsi10 ) ) 533 zweight = zindb * MAX( 0._wp, & 534 & MIN( zv_i_ac(ji,jl) * REAL( layer ), ( zv_i_ac(ji,jl) + zthick0(ji,nlay_i+1,jl) ) * REAL( jk ) ) & 535 & - MAX( zv_i_ac(ji,jl) * REAL( layer - 1 ), ( zv_i_ac(ji,jl) + zthick0(ji,nlay_i+1,jl) ) * REAL( jk - 1 ) ) ) & 536 & / ( REAL( nlay_i ) * MAX( zthick0(ji,layer,jl), epsi10 ) ) 537 ze_i_ac(ji,jk,jl) = ze_i_ac(ji,jk,jl) + zweight * zqm0(ji,layer,jl) 538 END DO 539 END DO 540 END DO 541 END DO 542 543 ! --- new volumes and layer thickness --- 544 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 545 DO ji = 1, nbpac 546 zinda = MAX( 0._wp, SIGN( 1._wp , zat_i_ac(ji) - epsi10 ) ) 547 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + zinda * ( zdv_res(ji) + zv_frazb(ji) ) * za_i_ac(ji,jl) / MAX( zat_i_ac(ji) , epsi10 ) 548 END DO 549 END DO 550 551 ! --- Redistributing energy on the new grid (energy is sent to the bottom) PART 2 --- ! 552 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 553 DO jk = 1, nlay_i 554 DO ji = 1, nbpac 555 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zv_i_ac(ji,jl) + epsi10 ) ) 556 ze_i_ac(ji,jk,jl) = zindb * ze_i_ac(ji,jk,jl) / MAX( zv_i_ac(ji,jl), epsi10 ) * REAL( nlay_i ) 557 END DO 558 END DO 559 END DO 560 474 zinda = MAX( 0._wp, SIGN( 1._wp , zat_i_1d(ji) - epsi20 ) ) 475 zv_newfra = zinda * ( zdv_res(ji) + zv_frazb(ji) ) * za_i_1d(ji,jl) / MAX( zat_i_1d(ji) , epsi20 ) 476 za_i_1d(ji,jl) = zinda * za_i_1d(ji,jl) 477 478 h_i_old (ji,nlay_i+1) = zv_newfra 479 qh_i_old(ji,nlay_i+1) = ze_newice(ji) * zv_newfra 480 481 zv_i_1d(ji,jl) = zv_i_1d(ji,jl) + zv_newfra 482 ENDDO 483 484 ! --- Ice enthalpy remapping --- ! 485 CALL lim_thd_ent( 1, nbpac, jl, ze_i_1d(1:nbpac,:,jl) ) 486 487 ENDDO 561 488 562 489 !------------ … … 565 492 DO jl = 1, jpl 566 493 DO ji = 1, nbpac 567 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - za_i_ ac(ji,jl) + epsi10 ) ) ! 0 if no ice and 1 if yes568 zoa_i_ ac(ji,jl) = za_old(ji,jl) * zoa_i_ac(ji,jl) / MAX( za_i_ac(ji,jl) , epsi10 ) * zindb494 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - za_i_1d(ji,jl) + epsi20 ) ) ! 0 if no ice and 1 if yes 495 zoa_i_1d(ji,jl) = za_old(ji,jl) * zoa_i_1d(ji,jl) / MAX( za_i_1d(ji,jl) , epsi20 ) * zindb 569 496 END DO 570 497 END DO … … 576 503 DO jl = 1, jpl 577 504 DO ji = 1, nbpac 578 zdv = zv_i_ ac(ji,jl) - zv_old(ji,jl)579 zsmv_i_ ac(ji,jl) = zsmv_i_ac(ji,jl) + zdv * zs_newice(ji)505 zdv = zv_i_1d(ji,jl) - zv_old(ji,jl) 506 zsmv_i_1d(ji,jl) = zsmv_i_1d(ji,jl) + zdv * zs_newice(ji) 580 507 END DO 581 508 END DO … … 586 513 !------------------------------------------------------------------------------! 587 514 DO jl = 1, jpl 588 CALL tab_1d_2d( nbpac, a_i (:,:,jl), npac(1:nbpac), za_i_ ac(1:nbpac,jl), jpi, jpj )589 CALL tab_1d_2d( nbpac, v_i (:,:,jl), npac(1:nbpac), zv_i_ ac(1:nbpac,jl), jpi, jpj )590 CALL tab_1d_2d( nbpac, oa_i(:,:,jl), npac(1:nbpac), zoa_i_ ac(1:nbpac,jl), jpi, jpj )591 CALL tab_1d_2d( nbpac, smv_i (:,:,jl), npac(1:nbpac), zsmv_i_ ac(1:nbpac,jl) , jpi, jpj )515 CALL tab_1d_2d( nbpac, a_i (:,:,jl), npac(1:nbpac), za_i_1d (1:nbpac,jl), jpi, jpj ) 516 CALL tab_1d_2d( nbpac, v_i (:,:,jl), npac(1:nbpac), zv_i_1d (1:nbpac,jl), jpi, jpj ) 517 CALL tab_1d_2d( nbpac, oa_i(:,:,jl), npac(1:nbpac), zoa_i_1d(1:nbpac,jl), jpi, jpj ) 518 CALL tab_1d_2d( nbpac, smv_i (:,:,jl), npac(1:nbpac), zsmv_i_1d(1:nbpac,jl) , jpi, jpj ) 592 519 DO jk = 1, nlay_i 593 CALL tab_1d_2d( nbpac, e_i(:,:,jk,jl), npac(1:nbpac), ze_i_ ac(1:nbpac,jk,jl), jpi, jpj )520 CALL tab_1d_2d( nbpac, e_i(:,:,jk,jl), npac(1:nbpac), ze_i_1d(1:nbpac,jk,jl), jpi, jpj ) 594 521 END DO 595 522 END DO … … 599 526 600 527 CALL tab_1d_2d( nbpac, hfx_thd, npac(1:nbpac), hfx_thd_1d(1:nbpac), jpi, jpj ) 601 CALL tab_1d_2d( nbpac, hfx_ tot, npac(1:nbpac), hfx_tot_1d(1:nbpac), jpi, jpj )528 CALL tab_1d_2d( nbpac, hfx_opw, npac(1:nbpac), hfx_opw_1d(1:nbpac), jpi, jpj ) 602 529 ! 603 530 ENDIF ! nbpac > 0 … … 612 539 ! heat content in Joules 613 540 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * area(ji,jj) * v_i(ji,jj,jl) / ( REAL( nlay_i ) * unit_fac ) 614 !e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * area(ji,jj) * v_i(ji,jj,jl) / ( REAL( nlay_i ) )615 541 END DO 616 542 END DO … … 618 544 END DO 619 545 620 !------------------------------------------------------------------------------|621 ! 10) Conservation check and changes in each ice category622 !------------------------------------------------------------------------------|623 IF( con_i ) THEN624 CALL lim_column_sum (jpl, v_i, vt_i_final)625 fieldid = 'v_i, limthd_lac'626 CALL lim_cons_check (vt_i_init, vt_i_final, 1.0e-6, fieldid)627 !628 CALL lim_column_sum_energy(jpl, nlay_i, e_i, et_i_final)629 fieldid = 'e_i, limthd_lac'630 CALL lim_cons_check (et_i_final, et_i_final, 1.0e-3, fieldid)631 !632 CALL lim_column_sum (jpl, v_s, vt_s_final)633 fieldid = 'v_s, limthd_lac'634 CALL lim_cons_check (vt_s_init, vt_s_final, 1.0e-6, fieldid)635 !636 ! CALL lim_column_sum (jpl, e_s(:,:,1,:) , et_s_init)637 ! fieldid = 'e_s, limthd_lac'638 ! CALL lim_cons_check (et_s_init, et_s_final, 1.0e-3, fieldid)639 IF( ln_nicep ) THEN640 DO ji = mi0(jiindx), mi1(jiindx)641 DO jj = mj0(jjindx), mj1(jjindx)642 WRITE(numout,*) ' vt_i_init : ', vt_i_init (ji,jj)643 WRITE(numout,*) ' vt_i_final: ', vt_i_final(ji,jj)644 WRITE(numout,*) ' et_i_init : ', et_i_init (ji,jj)645 WRITE(numout,*) ' et_i_final: ', et_i_final(ji,jj)646 END DO647 END DO648 ENDIF649 !650 ENDIF651 546 ! 652 CALL wrk_dealloc( jpij, zcatac) ! integer547 CALL wrk_dealloc( jpij, jcat ) ! integer 653 548 CALL wrk_dealloc( jpij, zswinew, zv_newice, za_newice, zh_newice, ze_newice, zs_newice, zo_newice ) 654 CALL wrk_dealloc( jpij, zdv_res, zda_res, zat_i_ac, zat_i_lev, zv_frazb, zvrel_ac ) 655 CALL wrk_dealloc( jpij,jpl, zv_old, za_old, za_i_ac, zv_i_ac, zoa_i_ac, zsmv_i_ac ) 656 CALL wrk_dealloc( jpij,jkmax,jpl, ze_i_ac ) 657 CALL wrk_dealloc( jpij,jkmax+1,jpl, zqm0, zthick0 ) 658 CALL wrk_dealloc( jpi,jpj, vt_i_init, vt_i_final, vt_s_init, vt_s_final, et_i_init, et_i_final, et_s_init, zvrel ) 549 CALL wrk_dealloc( jpij, zdv_res, zda_res, zat_i_1d, zat_i_lev, zv_frazb, zvrel_1d ) 550 CALL wrk_dealloc( jpij,jpl, zv_old, za_old, za_i_1d, zv_i_1d, zoa_i_1d, zsmv_i_1d ) 551 CALL wrk_dealloc( jpij,jkmax,jpl, ze_i_1d ) 552 CALL wrk_dealloc( jpi,jpj, zvrel ) 659 553 ! 660 554 END SUBROUTINE lim_thd_lac -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90
r4634 r4649 53 53 ! 54 54 INTEGER :: ji, jk ! dummy loop indices 55 REAL(wp) :: iflush, igravdr, ztmelts ! local scalars 56 REAL(wp) :: zaaa, zbbb, zccc, zdiscrim ! local scalars 55 REAL(wp) :: iflush, igravdr ! local scalars 57 56 !!--------------------------------------------------------------------- 58 57 58 !--------------------------------------------------------- 59 ! 0) Update ice salinity from snow-ice and bottom growth 60 !--------------------------------------------------------- 61 DO ji = kideb, kiut 62 sm_i_b(ji) = sm_i_b(ji) + dsm_i_se_1d(ji) + dsm_i_si_1d(ji) 63 END DO 64 59 65 !------------------------------------------------------------------------------| 60 66 ! 1) Constant salinity, constant in time | … … 71 77 ! Module 2 : Constant salinity varying in time | 72 78 !------------------------------------------------------------------------------| 73 74 79 IF( num_sal == 2 ) THEN 75 80 … … 78 83 ! Switches 79 84 !---------- 80 iflush =MAX( 0._wp , SIGN( 1._wp , t_su_b(ji) - rtt ) ) ! =1 if summer81 igravdr =MAX( 0._wp , SIGN( 1._wp , t_bo_b(ji) - t_su_b(ji) ) ) ! =1 if t_su < t_bo85 iflush = MAX( 0._wp , SIGN( 1._wp , t_su_b(ji) - rtt ) ) ! =1 if summer 86 igravdr = MAX( 0._wp , SIGN( 1._wp , t_bo_b(ji) - t_su_b(ji) ) ) ! =1 if t_su < t_bo 82 87 83 88 !--------------------- 84 89 ! Salinity tendencies 85 90 !--------------------- 86 ! 87 dsm_i_gd_1d(ji) = - igravdr * MAX( sm_i_b(ji) - sal_G , 0._wp ) / time_G * rdt_ice ! drainage by gravity 88 dsm_i_fl_1d(ji) = - iflush * MAX( sm_i_b(ji) - sal_F , 0._wp ) / time_F * rdt_ice ! drainage by flushing 89 ! 91 ! drainage by gravity drainage 92 dsm_i_gd_1d(ji) = - igravdr * MAX( sm_i_b(ji) - sal_G , 0._wp ) / time_G * rdt_ice 93 ! drainage by flushing 94 dsm_i_fl_1d(ji) = - iflush * MAX( sm_i_b(ji) - sal_F , 0._wp ) / time_F * rdt_ice 95 90 96 !----------------- 91 97 ! Update salinity … … 104 110 ! Salinity profile 105 111 CALL lim_var_salprof1d( kideb, kiut ) 106 107 108 ! Only necessary for conservation check since salinity is modified109 !--------------------110 ! Temperature update111 !--------------------112 DO jk = 1, nlay_i113 DO ji = kideb, kiut114 ztmelts = -tmut*s_i_b(ji,jk) + rtt115 !Conversion q(S,T) -> T (second order equation)116 zaaa = cpic117 zbbb = ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk) / rhoic - lfus118 zccc = lfus * ( ztmelts - rtt )119 zdiscrim = SQRT( MAX( zbbb*zbbb - 4.0*zaaa*zccc, 0._wp ) )120 t_i_b(ji,jk) = rtt - ( zbbb + zdiscrim ) / ( 2.0 *zaaa )121 END DO122 END DO123 112 ! 124 113 ENDIF … … 127 116 ! Module 3 : Profile of salinity, constant in time | 128 117 !------------------------------------------------------------------------------| 129 130 118 IF( num_sal == 3 ) CALL lim_var_salprof1d( kideb, kiut ) 131 119 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90
r4634 r4649 30 30 USE limvar ! clem for ice thickness correction 31 31 USE timing ! Timing 32 USE limcons ! conservation tests 32 33 33 34 IMPLICIT NONE … … 72 73 REAL(wp), POINTER, DIMENSION(:,:,:) :: zs0ice, zs0sn, zs0a, zs0c0 , zs0sm , zs0oi 73 74 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zs0e 74 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b ! Check conservation (C Rousset)75 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax, zchk_umax ! Check errors (C Rousset)76 75 ! mass and salt flux (clem) 77 76 REAL(wp), POINTER, DIMENSION(:,:,:) :: zviold, zvsold ! old ice volume... … … 79 78 REAL(wp), POINTER, DIMENSION(:,:) :: zeiold, zesold ! old enthalpies 80 79 REAL(wp) :: zdv, zda, zvi, zvs, zsmv, zes, zei 80 ! 81 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 81 82 !!--------------------------------------------------------------------- 82 83 IF( nn_timing == 1 ) CALL timing_start('limtrp') … … 87 88 88 89 CALL wrk_alloc( jpi, jpj, jpl, zaiold, zhimax, zviold, zvsold ) ! clem 89 90 ! -------------------------------91 !- check conservation (C Rousset)92 IF( ln_limdiahsb ) THEN93 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) )94 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) )95 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) )96 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) )97 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) )98 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) )99 ENDIF100 !- check conservation (C Rousset)101 ! -------------------------------102 90 103 91 IF( numit == nstart .AND. lwp ) THEN … … 114 102 IF( ln_limdyn ) THEN ! Advection of sea ice properties ! 115 103 ! !-------------------------------------! 104 105 ! conservation test 106 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limtrp', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 107 116 108 ! mass and salt flux init (clem) 117 109 zviold(:,:,:) = v_i(:,:,:) … … 368 360 369 361 ! Update fluxes 370 wfx_res(ji,jj) = wfx_res(ji,jj) +( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice371 wfx_snw(ji,jj) = wfx_snw(ji,jj) +( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice362 wfx_res(ji,jj) = wfx_res(ji,jj) - ( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice 363 wfx_snw(ji,jj) = wfx_snw(ji,jj) - ( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice 372 364 sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsmv ) * rhoic * r1_rdtice 373 365 hfx_res(ji,jj) = hfx_res(ji,jj) + ( e_s(ji,jj,1,jl) - zes ) * unit_fac / area(ji,jj) * r1_rdtice ! W.m-2 <0 … … 425 417 426 418 ! Update mass fluxes 427 wfx_res(ji,jj) = wfx_res(ji,jj) +( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice428 wfx_snw(ji,jj) = wfx_snw(ji,jj) +( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice419 wfx_res(ji,jj) = wfx_res(ji,jj) - ( v_i(ji,jj,jl) - zvi ) * rhoic * r1_rdtice 420 wfx_snw(ji,jj) = wfx_snw(ji,jj) - ( v_s(ji,jj,jl) - zvs ) * rhosn * r1_rdtice 429 421 sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsmv ) * rhoic * r1_rdtice 430 422 hfx_res(ji,jj) = hfx_res(ji,jj) + ( e_s(ji,jj,1,jl) - zes ) * unit_fac / area(ji,jj) * r1_rdtice ! W.m-2 <0 … … 474 466 END DO 475 467 END DO 468 469 ! conservation test 470 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limtrp', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 476 471 477 472 ENDIF … … 508 503 END DO 509 504 ENDIF 510 ! -------------------------------511 !- check conservation (C Rousset)512 IF( ln_limdiahsb ) THEN513 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b514 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b515 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b516 517 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw518 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic )519 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft520 521 zchk_vmin = glob_min(v_i)522 zchk_amax = glob_max(SUM(a_i,dim=3))523 zchk_amin = glob_min(a_i)524 zchk_umax = glob_max(SQRT(u_ice**2 + v_ice**2))525 526 IF(lwp) THEN527 IF ( ABS( zchk_v_i ) > 1.e-4 ) THEN528 WRITE(numout,*) 'violation volume [kg/day] (limtrp) = ',(zchk_v_i * rday)529 WRITE(numout,*) 'u_ice max [m/s] (limtrp) = ',zchk_umax530 WRITE(numout,*) 'number of time steps (limtrp) =',kt531 ENDIF532 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limtrp) = ',(zchk_smv * rday)533 IF ( ABS( zchk_e_i ) > 1.e-2 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limtrp) = ',(zchk_e_i)534 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limtrp) = ',(zchk_vmin * 1.e-3)535 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limtrp) = ',zchk_amin536 ENDIF537 ENDIF538 !- check conservation (C Rousset)539 ! -------------------------------540 505 ! 541 506 CALL wrk_dealloc( jpi, jpj, zui_u, zvi_v, zsm, zs0at, zs0ow, zeiold, zesold ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90
r4634 r4649 38 38 USE wrk_nemo ! work arrays 39 39 USE lib_fortran ! glob_sum 40 ! Check budget (Rousset)41 40 USE in_out_manager ! I/O manager 42 41 USE iom ! I/O manager 43 42 USE lib_mpp ! MPP library 44 43 USE timing ! Timing 44 USE limcons ! conservation tests 45 45 46 46 IMPLICIT NONE … … 82 82 REAL(wp) :: zh, zdvres, zsal 83 83 REAL(wp) :: zat_i_old, ztmelts 84 85 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b ! Check conservation (C Rousset) 86 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 84 ! 85 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 87 86 !!------------------------------------------------------------------- 88 87 IF( nn_timing == 1 ) CALL timing_start('limupdate1') 89 88 90 !------------------------------------------------------------------------------ 91 ! 1. Update of Global variables | 92 !------------------------------------------------------------------------------ 93 94 !----------------- 95 ! Trend terms 96 !----------------- 97 98 ! ------------------------------- 99 !- check conservation (C Rousset) 100 IF (ln_limdiahsb) THEN 101 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 102 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 103 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 104 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) 105 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 106 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 107 ENDIF 108 !- check conservation (C Rousset) 109 ! ------------------------------- 110 89 IF( ln_limdyn ) THEN 90 91 ! conservation test 92 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limupdate1', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 93 94 !----------------- 111 95 ! zap small values 112 96 !----------------- … … 115 99 CALL lim_var_glo2eqv 116 100 117 at_i(:,:) = 0._wp118 DO jl = 1, jpl119 at_i(:,:) = a_i(:,:,jl) + at_i(:,:)120 END DO121 122 101 !---------------------------------------------------- 123 ! 2.2)Rebin categories with thickness out of bounds102 ! Rebin categories with thickness out of bounds 124 103 !---------------------------------------------------- 125 104 DO jm = 1, jpm … … 134 113 END DO 135 114 136 !--- 2.13 ice concentration should not exceed amax 115 !---------------------------------------------------- 116 ! ice concentration should not exceed amax 137 117 !----------------------------------------------------- 138 118 DO jl = 1, jpl … … 147 127 END DO 148 128 149 at_i(:,:) = 0. 0129 at_i(:,:) = 0._wp 150 130 DO jl = 1, jpl 151 131 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 152 132 END DO 153 154 133 134 ! -------------------------------------- 155 135 ! Final thickness distribution rebinning 156 136 ! -------------------------------------- … … 163 143 END DO 164 144 165 145 !----------------- 166 146 ! zap small values 167 147 !----------------- … … 169 149 170 150 !--------------------- 171 ! 2.11) Ice salinity151 ! Ice salinity bounds 172 152 !--------------------- 173 153 IF ( num_sal == 2 ) THEN … … 179 159 ! salinity stays in bounds 180 160 i_ice_switch = 1._wp - MAX( 0._wp, SIGN( 1._wp, - v_i(ji,jj,jl) ) ) 181 smv_i(ji,jj,jl) = i_ice_switch * MAX( MIN( s_i_max * v_i(ji,jj,jl), smv_i(ji,jj,jl) ), s_i_min * v_i(ji,jj,jl) ) !+ s_i_min * ( 1._wp - i_ice_switch ) * v_i(ji,jj,jl)161 smv_i(ji,jj,jl) = i_ice_switch * MAX( MIN( s_i_max * v_i(ji,jj,jl), smv_i(ji,jj,jl) ), s_i_min * v_i(ji,jj,jl) ) 182 162 ! associated salt flux 183 163 sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsal ) * rhoic * r1_rdtice 184 END DO ! ji185 END DO ! jj186 END DO !jl164 END DO 165 END DO 166 END DO 187 167 ENDIF 188 189 ! -------------------190 at_i(:,:) = a_i(:,:,1)191 DO jl = 2, jpl192 at_i(:,:) = a_i(:,:,jl) + at_i(:,:)193 END DO194 195 168 196 169 ! ------------------------------------------------- … … 206 179 d_oa_i_trp (:,:,:) = oa_i (:,:,:) - old_oa_i (:,:,:) 207 180 d_smv_i_trp(:,:,:) = 0._wp 208 IF( num_sal == 2 ) d_smv_i_trp(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:) 209 210 ! ------------------------------- 211 !- check conservation (C Rousset) 212 IF( ln_limdiahsb ) THEN 213 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 214 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b 215 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b 216 217 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw 218 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic ) 219 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft 220 221 zchk_vmin = glob_min(v_i) 222 zchk_amax = glob_max(SUM(a_i,dim=3)) 223 zchk_amin = glob_min(a_i) 224 225 IF(lwp) THEN 226 IF ( ABS( zchk_v_i ) > 1.e-4 ) WRITE(numout,*) 'violation volume [kg/day] (limupdate1) = ',(zchk_v_i * rday) 227 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limupdate1) = ',(zchk_smv * rday) 228 IF ( ABS( zchk_e_i ) > 1.e-2 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limupdate1) = ',(zchk_e_i) 229 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate1) = ',(zchk_vmin * 1.e-3) 230 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate1) = ',zchk_amax 231 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate1) = ',zchk_amin 232 ENDIF 233 ENDIF 234 !- check conservation (C Rousset) 235 ! ------------------------------- 181 IF( num_sal == 2 ) d_smv_i_trp(:,:,:) = smv_i(:,:,:) - old_smv_i(:,:,:) 182 183 ! conservation test 184 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limupdate1', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 236 185 237 186 IF(ln_ctl) THEN ! Control print … … 302 251 ENDIF 303 252 253 ENDIF ! ln_limdyn 254 304 255 IF( nn_timing == 1 ) CALL timing_stop('limupdate1') 305 256 END SUBROUTINE lim_update1 -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90
r4634 r4649 39 39 USE lib_fortran ! glob_sum 40 40 USE timing ! Timing 41 USE limcons ! conservation tests 41 42 42 43 IMPLICIT NONE … … 84 85 REAL(wp) :: zdE ! specific enthalpy difference (J/kg) 85 86 REAL(wp) :: zfmdt ! exchange mass flux x time step (J/m2), >0 towards the ocean 86 87 REAL(wp) :: zchk_v_i, zchk_smv, zchk_e_i, zchk_fs, zchk_fw, zchk_ft, zchk_v_i_b, zchk_smv_b, zchk_e_i_b, zchk_fs_b, zchk_fw_b, zchk_ft_b ! Check conservation (C Rousset) 88 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 87 ! 88 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 89 89 !!------------------------------------------------------------------- 90 90 IF( nn_timing == 1 ) CALL timing_start('limupdate2') 91 91 92 !---------------------------------------------------------------------------------------- 93 ! 1. Computation of trend terms 94 !---------------------------------------------------------------------------------------- 95 96 ! ------------------------------- 97 !- check conservation (C Rousset) 98 IF (ln_limdiahsb) THEN 99 zchk_v_i_b = glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) 100 zchk_smv_b = glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) 101 zchk_e_i_b = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) 102 zchk_fw_b = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) 103 zchk_fs_b = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) 104 zchk_ft_b = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) 105 ENDIF 106 !- check conservation (C Rousset) 107 ! ------------------------------- 108 92 ! conservation test 93 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limupdate2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 94 95 !----------------- 109 96 ! zap small values 110 97 !----------------- … … 113 100 CALL lim_var_glo2eqv 114 101 115 !--------------------------------------116 ! 2. Review of all pathological cases117 !--------------------------------------118 at_i(:,:) = 0._wp119 DO jl = 1, jpl120 at_i(:,:) = a_i(:,:,jl) + at_i(:,:)121 END DO122 123 102 !---------------------------------------------------- 124 ! 2.2)Rebin categories with thickness out of bounds103 ! Rebin categories with thickness out of bounds 125 104 !---------------------------------------------------- 126 105 DO jm = 1, jpm … … 130 109 END DO 131 110 132 133 !clem debug: it is done in limthd_dh now 134 ! ! Melt of snow 135 ! !-------------- 136 ! DO jl = 1, jpl 137 ! DO jj = 1, jpj 138 ! DO ji = 1, jpi 139 ! IF( v_s(ji,jj,jl) >= epsi20 ) THEN 140 ! ! If snow energy of melting smaller then Lf 141 ! ! Then all snow melts and heat go to the ocean 142 ! !IF ( zEs <= lfus ) THEN 143 ! IF( t_s(ji,jj,1,jl) >= rtt ) THEN 144 ! zdvres = - v_s(ji,jj,jl) 145 ! zEs = - e_s(ji,jj,1,jl) * unit_fac / ( area(ji,jj) * MAX( v_s(ji,jj,jl), epsi20 ) ) ! snow energy of melting (J.m-3) 146 ! ! Contribution to heat flux to the ocean [W.m-2], < 0 147 ! hfx_res(ji,jj) = hfx_res(ji,jj) - zEs * zdvres * r1_rdtice 148 ! ! Contribution to mass flux 149 ! wfx_snw(ji,jj) = wfx_snw(ji,jj) + rhosn * zdvres * r1_rdtice 150 ! ! updates 151 ! v_s (ji,jj,jl) = 0._wp 152 ! ht_s(ji,jj,jl) = 0._wp 153 ! e_s (ji,jj,1,jl) = 0._wp 154 ! t_s (ji,jj,1,jl) = rtt 155 ! ENDIF 156 ! ENDIF 157 ! END DO 158 ! END DO 159 ! END DO 160 !clem debug 161 162 !--- 2.12 Constrain the thickness of the smallest category above 10 cm 111 !---------------------------------------------------------------------- 112 ! Constrain the thickness of the smallest category above hiclim 163 113 !---------------------------------------------------------------------- 164 114 DO jm = 1, jpm … … 176 126 END DO !jm 177 127 178 !--- 2.13 ice concentration should not exceed amax179 128 !----------------------------------------------------- 180 at_i(:,:) = 0.0 129 ! ice concentration should not exceed amax 130 !----------------------------------------------------- 131 at_i(:,:) = 0._wp 181 132 DO jl = 1, jpl 182 133 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) … … 199 150 END DO 200 151 152 ! -------------------------------------- 201 153 ! Final thickness distribution rebinning 202 154 ! -------------------------------------- … … 209 161 END DO 210 162 163 !----------------- 211 164 ! zap small values 212 165 !----------------- … … 232 185 ENDIF 233 186 234 235 ! -------------------236 at_i(:,:) = a_i(:,:,1)237 DO jl = 2, jpl238 at_i(:,:) = a_i(:,:,jl) + at_i(:,:)239 END DO240 241 187 !------------------------------------------------------------------------------ 242 188 ! 2) Corrections to avoid wrong values | … … 246 192 DO jj = 2, jpjm1 247 193 DO ji = 2, jpim1 248 IF ( at_i(ji,jj) .EQ. 0.0) THEN ! what to do if there is no ice249 IF ( at_i(ji+1,jj) .EQ. 0.0 ) u_ice(ji,jj) = 0.0! right side250 IF ( at_i(ji-1,jj) .EQ. 0.0 ) u_ice(ji-1,jj) = 0.0! left side251 IF ( at_i(ji,jj+1) .EQ. 0.0 ) v_ice(ji,jj) = 0.0! upper side252 IF ( at_i(ji,jj-1) .EQ. 0.0 ) v_ice(ji,jj-1) = 0.0! bottom side194 IF ( at_i(ji,jj) == 0._wp ) THEN ! what to do if there is no ice 195 IF ( at_i(ji+1,jj) == 0._wp ) u_ice(ji,jj) = 0._wp ! right side 196 IF ( at_i(ji-1,jj) == 0._wp ) u_ice(ji-1,jj) = 0._wp ! left side 197 IF ( at_i(ji,jj+1) == 0._wp ) v_ice(ji,jj) = 0._wp ! upper side 198 IF ( at_i(ji,jj-1) == 0._wp ) v_ice(ji,jj-1) = 0._wp ! bottom side 253 199 ENDIF 254 200 END DO … … 261 207 v_ice(:,:) = v_ice(:,:) * tmv(:,:) 262 208 263 264 209 ! ------------------------------------------------- 265 210 ! Diagnostics … … 276 221 dv_dt_thd(:,:,:) = d_v_i_thd(:,:,:) * r1_rdtice * rday 277 222 278 ! ------------------------------- 279 !- check conservation (C Rousset) 280 IF( ln_limdiahsb ) THEN 281 zchk_fs = glob_sum( ( sfx_bri(:,:) + sfx_bog(:,:) + sfx_bom(:,:) + sfx_sum(:,:) + sfx_sni(:,:) + sfx_opw(:,:) + sfx_res(:,:) + sfx_dyn(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 282 zchk_fw = glob_sum( ( wfx_bog(:,:) + wfx_bom(:,:) + wfx_sum(:,:) + wfx_sni(:,:) + wfx_opw(:,:) + wfx_res(:,:) + wfx_dyn(:,:) + wfx_snw(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fw_b 283 zchk_ft = glob_sum( ( hfx_tot(:,:) - hfx_thd(:,:) - hfx_dyn(:,:) - hfx_res(:,:) ) * area(:,:) / unit_fac * tms(:,:) ) - zchk_ft_b 284 285 zchk_v_i = ( glob_sum( SUM( v_i(:,:,:)*rhoic + v_s(:,:,:)*rhosn, dim=3 ) * area(:,:) * tms(:,:) ) - zchk_v_i_b ) * r1_rdtice - zchk_fw 286 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) * r1_rdtice + ( zchk_fs / rhoic ) 287 zchk_e_i = glob_sum( SUM( e_i(:,:,1:nlay_i,:), dim=3 ) + SUM( e_s(:,:,1:nlay_s,:), dim=3 ) ) * r1_rdtice - zchk_e_i_b * r1_rdtice + zchk_ft 288 289 zchk_vmin = glob_min(v_i) 290 zchk_amax = glob_max(SUM(a_i,dim=3)) 291 zchk_amin = glob_min(a_i) 292 293 IF(lwp) THEN 294 IF ( ABS( zchk_v_i ) > 1.e-4 ) WRITE(numout,*) 'violation volume [kg/day] (limupdate2) = ',(zchk_v_i * rday) 295 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limupdate2) = ',(zchk_smv * rday) 296 IF ( ABS( zchk_e_i ) > 1.e-2 ) WRITE(numout,*) 'violation enthalpy [1e9 J] (limupdate2) = ',(zchk_e_i) 297 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate2) = ',(zchk_vmin * 1.e-3) 298 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate2) = ',zchk_amax 299 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate2) = ',zchk_amin 300 ENDIF 301 ENDIF 302 !- check conservation (C Rousset) 303 ! ------------------------------- 223 ! heat content variation (W.m-2) 224 DO jj = 1, jpj 225 DO ji = 1, jpi 226 diag_heat_dhc(ji,jj) = ( SUM( d_e_i_trp(ji,jj,1:nlay_i,:) + d_e_i_thd(ji,jj,1:nlay_i,:) ) + & 227 & SUM( d_e_s_trp(ji,jj,1:nlay_s,:) + d_e_s_thd(ji,jj,1:nlay_s,:) ) ) * unit_fac * r1_rdtice / area(ji,jj) 228 END DO 229 END DO 230 231 ! conservation test 232 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limupdate2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 304 233 305 234 IF(ln_ctl) THEN ! Control print … … 370 299 371 300 IF( nn_timing == 1 ) CALL timing_stop('limupdate2') 301 372 302 END SUBROUTINE lim_update2 373 303 #else -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90
r4635 r4649 102 102 CALL iom_put( "isst" , sst_m ) ! sea surface temperature 103 103 CALL iom_put( "isss" , sss_m ) ! sea surface salinity 104 CALL iom_put( "qt_oce" , qns + qsr ) ! total flux at ocean surface105 104 ! 106 105 DO jj = 2 , jpjm1 … … 120 119 CALL iom_put( "vice_ipa" , z2db ) ! ice velocity v component 121 120 CALL iom_put( "icevel" , z2d ) ! ice velocity module 122 !!SF BE CAREFUL : qsr_oce qnd qns_oce are after penetration over ice123 CALL iom_put( "qsr_oce" , qsr ) ! solar flux at ocean surface124 CALL iom_put( "qns_oce" , qns ) ! non-solar flux at ocean surface125 !!SF end be careful126 121 CALL iom_put( "utau_ice" , utau_ice ) ! wind stress over ice along i-axis at I-point 127 122 CALL iom_put( "vtau_ice" , vtau_ice ) ! wind stress over ice along j-axis at I-point 128 !!SF commented because this computation is not ok129 !SF because qsr is not qsr_ocean but it contains already qsr_ice130 !SF131 !SF DO jj = 1 , jpj132 !SF DO ji = 1 , jpi133 !SF z2d(ji,jj) = ( 1._wp - at_i(ji,jj) ) * qsr(ji,jj)134 !SF END DO135 !SF END DO136 !SF CALL iom_put( "qsr_io" , z2d ) ! solar flux at ice/ocean surface137 !SF DO jj = 1 , jpj138 !SF DO ji = 1 , jpi139 !SF z2d(ji,jj) = ( 1._wp - at_i(ji,jj) ) * qns(ji,jj)140 !SF END DO141 !SF END DO142 !SF CALL iom_put( "qns_io" , z2d ) ! non-solar flux at ice/ocean surface143 123 CALL iom_put( "snowpre" , sprecip ) ! snow precipitation 144 124 CALL iom_put( "micesalt" , smt_i ) ! mean ice salinity … … 211 191 CALL iom_put( "vfxsnw" , wfx_snw * rday / rhoic ) ! total snw growth/melt 212 192 CALL iom_put( "vfxsub" , wfx_sub * rday / rhoic ) ! sublimation (snow) 213 214 CALL iom_put ('hfxdhc1', diag_heat_dhc1(:,:) ) ! Heat content variation in snow and ice 215 CALL iom_put ('hfxspr', hfx_spr(:,:) ) ! Heat content of snow precip 216 CALL iom_put ('hfxqsr', qsr(:,:) ) ! solar fluxes used by snw/ice 217 CALL iom_put ('hfxqns', qns(:,:) ) ! non solar fluxes used by snw/ice 193 CALL iom_put( "vfxspr" , wfx_spr * rday / rhoic ) ! precip (snow) 218 194 219 195 CALL iom_put ('hfxthd', hfx_thd(:,:) ) ! … … 222 198 CALL iom_put ('hfxout', hfx_out(:,:) ) ! 223 199 CALL iom_put ('hfxin' , hfx_in(:,:) ) ! 224 CALL iom_put ('hfxtot', hfx_tot(:,:) ) !225 200 CALL iom_put ('hfxsnw', hfx_snw(:,:) ) ! 226 201 CALL iom_put ('hfxsub', hfx_sub(:,:) ) ! 227 202 CALL iom_put ('hfxerr', hfx_err(:,:) ) ! 228 203 CALL iom_put ('hfxerr_rem', hfx_err_rem(:,:) ) ! 204 205 CALL iom_put ('hfxsum', hfx_sum(:,:) ) ! 206 CALL iom_put ('hfxbom', hfx_bom(:,:) ) ! 207 CALL iom_put ('hfxbog', hfx_bog(:,:) ) ! 208 CALL iom_put ('hfxdif', hfx_dif(:,:) ) ! 209 CALL iom_put ('hfxopw', hfx_opw(:,:) ) ! 210 CALL iom_put ('hfxtur', fhtur(:,:) * at_i(:,:) ) ! turbulent heat flux at ice base 211 CALL iom_put ('hfxdhc', diag_heat_dhc(:,:) ) ! Heat content variation in snow and ice 212 CALL iom_put ('hfxspr', hfx_spr(:,:) ) ! Heat content of snow precip 229 213 230 214 !-------------------------------- … … 327 311 CALL histdef( kid, "iicebome", "Ice bottom melt" , "m/s" , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 328 312 CALL histdef( kid, "iicesume", "Ice surface melt" , "m/s" , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 329 CALL histdef( kid, "iisfxthd", "Salt flux from thermo" , "" , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )330 313 CALL histdef( kid, "iisfxdyn", "Salt flux from dynmics" , "" , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) 331 314 CALL histdef( kid, "iisfxres", "Salt flux from limupdate", "" , jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt ) … … 354 337 CALL histwrite( kid, "iicebome", kt, wfx_bom , jpi*jpj, (/1/) ) 355 338 CALL histwrite( kid, "iicesume", kt, wfx_sum , jpi*jpj, (/1/) ) 356 !CALL histwrite( kid, "iisfxthd", kt, sfx_thd , jpi*jpj, (/1/) )357 339 CALL histwrite( kid, "iisfxdyn", kt, sfx_dyn , jpi*jpj, (/1/) ) 358 340 CALL histwrite( kid, "iisfxres", kt, sfx_res , jpi*jpj, (/1/) ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90
r4634 r4649 57 57 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: t_bo_b !: <==> the 2D t_bo 58 58 59 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_tot_1d 60 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_thd_1d 61 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_spr_1d 59 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_sum_1d 60 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_bom_1d 61 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_bog_1d 62 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_dif_1d 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_opw_1d 62 64 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_snw_1d 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_sub_1d64 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_res_1d65 65 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_err_1d 66 66 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_err_rem_1d 67 68 ! heat flux associated with ice-atmosphere mass exchange 69 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_sub_1d 70 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_spr_1d 71 72 ! heat flux associated with ice-ocean mass exchange 73 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_thd_1d 74 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: hfx_res_1d 67 75 68 76 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_ice_1d !: <==> the 2D wfx_ice … … 75 83 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_sni_1d !: <==> the 2D wfx_ice 76 84 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_opw_1d !: <==> the 2D wfx_ice 85 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_res_1d !: <==> the 2D wfx_ice 86 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_spr_1d !: <==> the 2D wfx_ice 77 87 78 88 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: sfx_bri_1d !: <==> the 2D sfx_bri … … 82 92 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: sfx_sni_1d !: 83 93 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: sfx_opw_1d !: 94 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: sfx_res_1d !: 84 95 85 96 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: sprecip_1d !: <==> the 2D sprecip 86 97 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: frld_1d !: <==> the 2D frld 87 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: at_i_b !: <==> the 2D frld98 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: at_i_b !: <==> the 2D at_i 88 99 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: fhtur_1d !: <==> the 2D fhtur 89 100 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: fhld_1d !: <==> the 2D fhld … … 150 161 & qsr_ice_1d (jpij) , & 151 162 & fr1_i0_1d(jpij) , fr2_i0_1d(jpij) , qns_ice_1d(jpij) , & 152 & t_bo_b (jpij) , iatte_1d (jpij) , & 153 & oatte_1d (jpij) , hfx_tot_1d(jpij), hfx_thd_1d(jpij) , hfx_spr_1d(jpij) , & 154 & hfx_snw_1d(jpij), hfx_sub_1d(jpij), hfx_err_1d(jpij) , hfx_res_1d(jpij) , hfx_err_rem_1d(jpij), STAT=ierr(1) ) 163 & t_bo_b (jpij) , iatte_1d (jpij) , oatte_1d (jpij) , & 164 & hfx_sum_1d(jpij) , hfx_bom_1d(jpij) ,hfx_bog_1d(jpij) ,hfx_dif_1d(jpij) ,hfx_opw_1d(jpij) , & 165 & hfx_thd_1d(jpij) , hfx_spr_1d(jpij) , & 166 & hfx_snw_1d(jpij) , hfx_sub_1d(jpij) , hfx_err_1d(jpij) , hfx_res_1d(jpij) , hfx_err_rem_1d(jpij), STAT=ierr(1) ) 155 167 ! 156 168 ALLOCATE( sprecip_1d (jpij) , frld_1d (jpij) , at_i_b (jpij) , & 157 & fhtur_1d (jpij) , wfx_ice_1d (jpij) , wfx_snw_1d (jpij) , &158 & fhld_1d (jpij) , wfx_sub_1d (jpij) , wfx_bog_1d(jpij) , wfx_bom_1d(jpij) , wfx_sum_1d(jpij) , wfx_sni_1d (jpij) , wfx_opw_1d (jpij) , &169 & fhtur_1d (jpij) , wfx_ice_1d (jpij) , wfx_snw_1d (jpij) , wfx_spr_1d (jpij) , & 170 & fhld_1d (jpij) , wfx_sub_1d (jpij) , wfx_bog_1d(jpij) , wfx_bom_1d(jpij) , wfx_sum_1d(jpij) , wfx_sni_1d (jpij) , wfx_opw_1d (jpij) , wfx_res_1d (jpij) , & 159 171 & dqns_ice_1d(jpij) , qla_ice_1d (jpij) , dqla_ice_1d(jpij) , & 160 172 & tatm_ice_1d(jpij) , & 161 173 & i0 (jpij) , & 162 & sfx_bri_1d (jpij) , sfx_bog_1d (jpij) , sfx_bom_1d (jpij) ,sfx_sum_1d (jpij) ,sfx_sni_1d (jpij) , sfx_opw_1d (jpij) , &174 & sfx_bri_1d (jpij) , sfx_bog_1d (jpij) , sfx_bom_1d (jpij) ,sfx_sum_1d (jpij) ,sfx_sni_1d (jpij) , sfx_opw_1d (jpij) , sfx_res_1d (jpij) , & 163 175 & dsm_i_fl_1d(jpij) , dsm_i_gd_1d(jpij) , dsm_i_se_1d(jpij) , & 164 176 & dsm_i_si_1d(jpij) , hicol_b (jpij) , STAT=ierr(2) ) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r4634 r4649 69 69 70 70 PUBLIC sbc_ice_lim ! routine called by sbcmod.F90 71 PUBLIC lim_prt_state 71 72 72 73 !! * Substitutions … … 170 171 ! 171 172 IF( ln_nicep ) THEN ! control print at a given point 172 jiindx = 3 ; jjindx = 49173 jiindx = 15 ; jjindx = 44 173 174 IF(lwp) WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx 174 175 ENDIF … … 183 184 u_oce(:,:) = ssu_m(:,:) ! mean surface ocean current at ice velocity point 184 185 v_oce(:,:) = ssv_m(:,:) ! (C-grid dynamics : U- & V-points as the ocean) 185 ! 186 t_bo(:,:) = tfreez( sss_m ) + rt0 ! masked sea surface freezing temperature [Kelvin] 187 ! ! (set to rt0 over land) 186 187 ! masked sea surface freezing temperature [Kelvin] 188 t_bo(:,:) = ( tfreez( sss_m ) + rt0 ) * tmask(:,:,1) + rt0 * ( 1. - tmask(:,:,1) ) 189 188 190 CALL albedo_ice( t_su, ht_i, ht_s, zalb_ice_cs, zalb_ice_os ) ! ... ice albedo 189 191 … … 310 312 ! salt, heat and mass fluxes 311 313 sfx (:,:) = 0._wp ; 312 sfx_bri(:,:) = 0._wp ; sfx_dyn(:,:) = 0._wp ; sfx_res(:,:) = 0._wp314 sfx_bri(:,:) = 0._wp ; sfx_dyn(:,:) = 0._wp 313 315 sfx_sni(:,:) = 0._wp ; sfx_opw(:,:) = 0._wp 314 316 sfx_bog(:,:) = 0._wp ; sfx_dyn(:,:) = 0._wp 315 317 sfx_bom(:,:) = 0._wp ; sfx_sum(:,:) = 0._wp 316 317 hfx_thd(:,:) = 0._wp ; hfx_dyn(:,:) = 0._wp ; hfx_snw(:,:) = 0._wp 318 hfx_tot(:,:) = 0._wp ; hfx_spr(:,:) = 0._wp ; hfx_res(:,:) = 0._wp 319 hfx_sub(:,:) = 0._wp ; hfx_err(:,:) = 0._wp ; hfx_in (:,:) = 0._wp ; hfx_out(:,:) = 0._wp 320 hfx_err_rem(:,:) = 0._wp 321 322 wfx_snw(:,:) = 0._wp ; wfx_ice(:,:) = 0._wp ; wfx_sub(:,:) = 0._wp 318 sfx_res(:,:) = 0._wp 319 320 wfx_snw(:,:) = 0._wp ; wfx_ice(:,:) = 0._wp 323 321 wfx_sni(:,:) = 0._wp ; wfx_opw(:,:) = 0._wp 324 322 wfx_bog(:,:) = 0._wp ; wfx_dyn(:,:) = 0._wp 325 323 wfx_bom(:,:) = 0._wp ; wfx_sum(:,:) = 0._wp 326 wfx_res(:,:) = 0._wp ; 327 ! 328 fhld (:,:) = 0._wp 329 fmmflx(:,:) = 0._wp 330 ftr_ice(:,:,:) = 0._wp ! part of solar radiation transmitted through the ice 324 wfx_res(:,:) = 0._wp ; wfx_sub(:,:) = 0._wp 325 wfx_spr(:,:) = 0._wp ; 326 327 hfx_in (:,:) = 0._wp ; hfx_out(:,:) = 0._wp 328 hfx_thd(:,:) = 0._wp ; 329 hfx_snw(:,:) = 0._wp ; hfx_opw(:,:) = 0._wp 330 hfx_bog(:,:) = 0._wp ; hfx_dyn(:,:) = 0._wp 331 hfx_bom(:,:) = 0._wp ; hfx_sum(:,:) = 0._wp 332 hfx_res(:,:) = 0._wp ; hfx_sub(:,:) = 0._wp 333 hfx_spr(:,:) = 0._wp ; hfx_dif(:,:) = 0._wp 334 hfx_err(:,:) = 0._wp ; hfx_err_rem(:,:) = 0._wp 335 336 ! 337 fhld (:,:) = 0._wp 338 fmmflx(:,:) = 0._wp 339 ! part of solar radiation transmitted through the ice 340 ftr_ice(:,:,:) = 0._wp 331 341 332 342 ! diags 333 diag_trp_vi (:,:) = 0._wp ; diag_trp_vs(:,:) = 0._wp ; diag_trp_ei(:,:) = 0._wp ; diag_trp_es(:,:) = 0._wp ;334 diag_heat_dhc 1(:,:) = 0._wp ;343 diag_trp_vi (:,:) = 0._wp ; diag_trp_vs(:,:) = 0._wp ; diag_trp_ei(:,:) = 0._wp ; diag_trp_es(:,:) = 0._wp 344 diag_heat_dhc(:,:) = 0._wp 335 345 336 346 ! dynamical invariants … … 820 830 WRITE(numout,*) ' hfx_in : ', hfx_in(ji,jj) 821 831 WRITE(numout,*) ' hfx_out : ', hfx_out(ji,jj) 822 WRITE(numout,*) ' hfx_tot : ', hfx_tot(ji,jj) 823 WRITE(numout,*) ' dhc : ', diag_heat_dhc1(ji,jj) 832 WRITE(numout,*) ' dhc : ', diag_heat_dhc(ji,jj) 824 833 WRITE(numout,*) 825 834 WRITE(numout,*) ' hfx_dyn : ', hfx_dyn(ji,jj) -
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r4634 r4649 391 391 CALL iom_put( "qt" , qns + qsr ) ! total heat flux 392 392 CALL iom_put( "qns" , qns ) ! solar heat flux 393 CALL iom_put( "qsr" , qsr) ! solar heat flux393 CALL iom_put( "qsr" , qsr ) ! solar heat flux 394 394 IF( nn_ice > 0 ) CALL iom_put( "ice_cover", fr_i ) ! ice fraction 395 395 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
