MODULE icestp !!====================================================================== !! *** MODULE icestp *** !! Sea-Ice model : LIM Sea ice model time-stepping !!====================================================================== #if defined key_lim3 !!---------------------------------------------------------------------- !! 'key_lim3' : LIM3 sea-ice model !!---------------------------------------------------------------------- !! ice_stp : sea-ice model time-stepping !!---------------------------------------------------------------------- USE par_ice USE dom_oce USE oce ! dynamics and tracers variables USE in_out_manager USE ice_oce ! ice variables USE flx_oce ! forcings variables USE dom_ice USE cpl_oce ! nemo add USE daymod USE phycst ! Define parameters for the routines USE taumod USE ice USE iceini USE ocesbc USE lbclnk USE limdyn USE limtrp USE limthd USE limitd_th USE limitd_me USE limflx USE limdia USE limwri USE limrst USE limupdate ! update of global variables USE limvar USE prtctl ! Print control IMPLICIT NONE PRIVATE !! * Routine accessibility PUBLIC ice_stp ! called by step.F90 !! * Substitutions # include "domzgr_substitute.h90" # include "vectopt_loop_substitute.h90" !!----------------------------------------------------- !! LIM 2.0 , UCL-LODYC-IPSL (2003) !!----------------------------------------------------- CONTAINS SUBROUTINE ice_stp ( kt ) !!--------------------------------------------------------------------- !! *** ROUTINE ice_stp *** !! !! ** Purpose : Louvain la Neuve Sea Ice Model time stepping !! !! ** Action : - call the ice dynamics routine !! - call the ice advection/diffusion routine !! - call the ice thermodynamics routine !! - call the routine that computes mass and !! heat fluxes at the ice/ocean interface !! - save the outputs !! - save the outputs for restart when necessary !! !! History : !! 1.0 ! 99-11 (M. Imbard) Original code !! ! 01-03 (D. Ludicone, E. Durand, G. Madec) free surf. !! 2.0 ! 02-09 (G. Madec, C. Ethe) F90: Free form and module !! 3.0 ! 05-10 (M. Vancoppenolle) Trend terms, ice salinity, multi-layer model !!---------------------------------------------------------------------- !! * Arguments INTEGER, INTENT( in ) :: kt ! ocean time-step index !! * Local declarations INTEGER :: & ji, jj, jk, jl, & ! loop indexes index, & ! indexes for ice points numaltests, & ! number of alert tests (max 20) alert_id ! number of the current alert INTEGER, DIMENSION(20) :: & numal ! number of alerts positive REAL(wp) :: & ztmelts ! ice layer melting point REAL(wp) , DIMENSION(jpi,jpj) :: & zsss_io, zsss2_io, zsss3_io ! tempory workspaces CHARACTER (len=30), DIMENSION(20) :: & ! name of alert alname !!---------------------------------------------------------------------- !------------------------------------------------------------------------------- ! 1) First time step !------------------------------------------------------------------------------- !+++++ index = 12 jiindex = 44 jjindex = 140 WRITE(numout,*) ' The debugging point is : jiindex : ',jiindex, & ' jjindex : ',jjindex !+++++ IF( kt == nit000 ) THEN IF( lk_cpl ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'ice_stp : Louvain la Neuve Ice Model (LIM)' IF(lwp) WRITE(numout,*) '~~~~~~~ coupled case' ELSE IF(lwp) WRITE(numout,*) 'ice_stp : Louvain la Neuve Ice Model (LIM)' IF(lwp) WRITE(numout,*) '~~~~~~~ forced case using bulk formulea' ENDIF gtaux(:,:) = 0.e0 ! Set wind stress to 0 gtauy(:,:) = 0.e0 ENDIF ! Initial time step IF(lwp) WRITE(numout,*) ' ~~~ LIM-@ ~~~ ' IF(lwp) WRITE(numout,*) ' Time step : ', kt !------------------------------------------------------------------------------- ! 2) SST, SSS, wind stress, ocean velocity and seawater freezing point !------------------------------------------------------------------------------- !------------------------------------- ! Update SST, SSS and ocean velocity !------------------------------------- sst_io(:,:) = sst_io(:,:) + tn(:,:,1) + rt0 sss_io(:,:) = sss_io(:,:) + sn(:,:,1) ! vectors at F-point DO jj = 2, jpj DO ji = fs_2, jpi u_io(ji,jj) = u_io(ji,jj) + un(ji,jj,1) v_io(ji,jj) = v_io(ji,jj) + vn(ji,jj,1) END DO END DO !------------------------ ! Ice model loop starts !------------------------ IF( MOD( kt-1, nfice ) == 0 ) THEN !------------------- ! Average SST, SSS !------------------- sst_io(:,:) = sst_io(:,:) / FLOAT( nfice ) * tmask(:,:,1) sss_io(:,:) = sss_io(:,:) / FLOAT( nfice ) !---------------------------------------- ! Atmospheric stress and ocean velocity !---------------------------------------- DO jj = 2, jpj DO ji = fs_2, jpi gtaux(ji,jj) = cai / cao * tauxw(ji,jj) gtauy(ji,jj) = cai / cao * tauyw(ji,jj) u_io (ji,jj) = u_io(ji,jj) / FLOAT( nfice ) * tmu(ji,jj) v_io (ji,jj) = v_io(ji,jj) / FLOAT( nfice ) * tmv(ji,jj) END DO END DO ! Boundary conditions CALL lbc_lnk( gtaux(:,:), 'U', -1. ) ! I-point (i.e. ice U-V point) CALL lbc_lnk( gtauy(:,:), 'V', -1. ) ! I-point (i.e. ice U-V point) CALL lbc_lnk( u_io (:,:), 'U', -1. ) ! I-point (i.e. ice U-V point) CALL lbc_lnk( v_io (:,:), 'V', -1. ) ! I-point (i.e. ice U-V point) !-------------------------- ! Seawater freezing point !-------------------------- zsss_io (:,:) = SQRT( sss_io(:,:) ) zsss2_io(:,:) = sss_io(:,:) * sss_io(:,:) zsss3_io(:,:) = zsss_io(:,:) * zsss_io(:,:) * zsss_io(:,:) DO jj = 1, jpj DO ji = 1, jpi t_bo(ji,jj) = ABS ( rt0 - 0.0575 * sss_io(ji,jj) & & + 1.710523e-03 * zsss3_io(ji,jj) & & - 2.154996e-04 * zsss2_io(ji,jj) ) & * tms(ji,jj) END DO END DO !+++++ IF(ln_ctl) THEN ! print mean trends (used for debugging) CALL prt_ctl_info('Ice Forcings ') CALL prt_ctl(tab2d_1=qsr_oce ,clinfo1=' qsr_oce : ') CALL prt_ctl(tab2d_1=qnsr_oce,clinfo1=' qnsr_oce : ') CALL prt_ctl(tab2d_1=evap ,clinfo1=' evap : ') CALL prt_ctl(tab2d_1=tprecip ,clinfo1=' precip : ') CALL prt_ctl(tab2d_1=gtaux ,clinfo1=' u-stress : ', tab2d_2=gtauy , clinfo2=' v-stress : ') CALL prt_ctl(tab2d_1=sst_io ,clinfo1=' sst : ', tab2d_2=sss_io , clinfo2=' sss : ') CALL prt_ctl(tab2d_1=u_io ,clinfo1=' u_io : ', tab2d_2=v_io , clinfo2=' v_io : ') CALL prt_ctl(tab2d_1=frld ,clinfo1=' frld 1 : ') DO jl = 1, jpl CALL prt_ctl_info('* - category number ', ivar1=jl) CALL prt_ctl(tab3d_1=t_su , clinfo1=' t_su : ', kdim=jl) CALL prt_ctl(tab3d_1=qsr_ice , clinfo1=' qsr_ice : ', kdim=jl) CALL prt_ctl(tab3d_1=qnsr_ice, clinfo1=' qnsr_ice : ', kdim=jl) CALL prt_ctl(tab3d_1=dqns_ice, clinfo1=' dqns_ice : ', kdim=jl) CALL prt_ctl(tab3d_1=qla_ice , clinfo1=' qla_ice : ', kdim=jl) CALL prt_ctl(tab3d_1=dqla_ice, clinfo1=' dqla_ice : ', kdim=jl) END DO ENDIF !+++++ !---------------------------------------------------------------------------- ! 3) Store old values of ice model global variables !---------------------------------------------------------------------------- old_a_i(:,:,:) = a_i(:,:,:) ! ice area old_e_i(:,:,:,:) = e_i(:,:,:,:) ! ice thermal energy old_v_i(:,:,:) = v_i(:,:,:) ! ice volume old_v_s(:,:,:) = v_s(:,:,:) ! snow volume old_e_s(:,:,:,:) = e_s(:,:,:,:) ! snow thermal energy old_smv_i(:,:,:) = smv_i(:,:,:) ! salt content old_oa_i(:,:,:) = oa_i(:,:,:) ! areal age content d_a_i_thd(:,:,:) = 0.0 ; d_a_i_trp(:,:,:) = 0.0 d_v_i_thd(:,:,:) = 0.0 ; d_v_i_trp(:,:,:) = 0.0 d_e_i_thd(:,:,:,:) = 0.0 ; d_e_i_trp(:,:,:,:) = 0.0 d_v_s_thd(:,:,:) = 0.0 ; d_v_s_trp(:,:,:) = 0.0 d_e_s_thd(:,:,:,:) = 0.0 ; d_e_s_trp(:,:,:,:) = 0.0 d_smv_i_thd(:,:,:) = 0.0 ; d_smv_i_trp(:,:,:) = 0.0 d_oa_i_thd(:,:,:) = 0.0 ; d_oa_i_trp(:,:,:) = 0.0 fsalt(:,:) = 0.0 ; fseqv(:,:) = 0.0 fsbri(:,:) = 0.0 ; fsalt_res(:,:) = 0.0 fsalt_rpo(:,:) = 0.0 fhmec(:,:) = 0.0 ; fhbri(:,:) = 0.0; fmmec(:,:) = 0.0 ; fheat_res(:,:) = 0.0 fheat_rpo(:,:) = 0.0 ; focea2D(:,:) = 0.0 fsup2D(:,:) = 0.0 diag_sni_gr(:,:) = 0.0 ; diag_lat_gr(:,:) = 0.0 diag_bot_gr(:,:) = 0.0 ; diag_dyn_gr(:,:) = 0.0 diag_bot_me(:,:) = 0.0 ; diag_sur_me(:,:) = 0.0 ! dynamical invariants delta_i(:,:) = 0.0 divu_i(:,:) = 0.0 shear_i(:,:) = 0.0 !---------------------- ! Ice model time step !---------------------- numit = numit + nfice !---------------------------------------------------------------------------- ! 4) Ice routines, compute the trend terms of global variables !---------------------------------------------------------------------------- ! !-----------------------! CALL lim_rst_opn( kt ) ! Open Ice restart file ! ! !-----------------------! !+++++ WRITE(numout,*) ' - Beginning the time step - ' CALL lim_inst_state(jiindex,jjindex,1) WRITE(numout,*) ' ' !+++++ !--------------------------- ! 4.1) Dynamical processes !--------------------------- ! !--------------------! CALL lim_dyn ! Ice dynamics ! ( rheology/dynamics ) ! !--------------------! ! !--------------------! CALL lim_trp ! Ice transport ! ! !--------------------! CALL lim_var_agg(1) ! aggregate categories, requested CALL lim_var_glo2eqv ! equivalent variables, requested for rafting !+++++ WRITE(numout,*) ' - After ice dynamics and transport ' CALL lim_inst_state(jiindex,jjindex,1) WRITE(numout,*) ' ' WRITE(numout,*) WRITE(numout,*) ' Mechanical Check ************** ' WRITE(numout,*) ' Check what means ice divergence ' WRITE(numout,*) ' Total ice concentration ', at_i (jiindex,jjindex) WRITE(numout,*) ' Total lead fraction ', ato_i(jiindex,jjindex) WRITE(numout,*) ' Sum of both ', ato_i(jiindex,jjindex) + at_i(jiindex,jjindex) WRITE(numout,*) ' Sum of both minus 1 ', ato_i(jiindex,jjindex) + at_i(jiindex,jjindex) - 1.00 !+++++ ! !------------------------------------------------! CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) ! ! !------------------------------------------------! !------------------------- ! 4.2 Ice thermodynamics !------------------------- CALL lim_var_glo2eqv ! equivalent variables CALL lim_var_agg(1) ! aggregate ice categories CALL lim_var_bv ! bulk brine volume (diag) ! !--------------------------------! CALL lim_thd ! Heat diffusion, growth, melt ! ! !--------------------------------! ! !---------------------! ! ! Ice natural aging ! ! !---------------------! oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rdt_ice / 86400.00 !+++++ CALL lim_var_glo2eqv ! except the info message that follows, i dont ! think this one is necessary WRITE(numout,*) ' - After ice thermodynamics ' CALL lim_inst_state(jiindex,jjindex,1) !+++++ ! !-------------------------------------------! CALL lim_itd_th ! Remap ice categories, lateral accretion ! ! !-------------------------------------------! !---------------------------------------------------------------------------- ! 5) Global variables update !---------------------------------------------------------------------------- CALL lim_var_agg(1) ! requested by limupdate ! !-------------------------! CALL lim_update ! Global variables update ! ! !-------------------------! CALL lim_var_glo2eqv ! equivalent variables (outputs) CALL lim_var_agg(2) ! aggregate ice thickness categories !+++++ IF(ln_nicep) THEN WRITE(numout,*) ' - Final ice state after lim_update ' CALL lim_inst_state(jiindex,jjindex,2) WRITE(numout,*) ' ' ENDIF !+++++ !---------------------------------------------------------------------------- ! 6) Fluxes of mass and heat to the ocean !---------------------------------------------------------------------------- ! !------------------------------! CALL lim_flx ! Ice/Ocean Mass & Heat fluxes ! ! !------------------------------! !+++++ WRITE(numout,*) ' - Final ice state after lim_flx ' CALL lim_inst_state(jiindex,jjindex,3) WRITE(numout,*) ' ' !+++++ !---------------------------------------------------------------------------- ! 7) Diagnostics and outputs !---------------------------------------------------------------------------- IF( MOD( numit, ninfo ) == 0 .OR. ntmoy == 1 ) THEN ! !-------------------------------! CALL lim_dia ! Ice Diagnostics in evolu file ! ENDIF !-------------------------------! ! !----------------------------! CALL lim_wri( 1 ) ! Ice outputs in icemod file ! ! !----------------------------! ! !------------------! IF( lrst_ice ) CALL lim_rst_write( kt ) ! Ice restart file ! ! !------------------! CALL lim_var_glo2eqv ! Re-initialization of forcings qsr_oce (:,:) = 0.e0 qsr_ice (:,:,:) = 0.e0 qnsr_oce(:,:) = 0.e0 qnsr_ice(:,:,:) = 0.e0 dqns_ice(:,:,:) = 0.e0 tprecip (:,:) = 0.e0 sprecip (:,:) = 0.e0 #if defined key_coupled rrunoff (:,:) = 0.e0 calving (:,:) = 0.e0 #else qla_ice (:,:,:) = 0.e0 dqla_ice(:,:,:) = 0.e0 fr1_i0 (:,:) = 0.e0 fr2_i0 (:,:) = 0.e0 evap (:,:) = 0.e0 #endif !---------------------------------------------------------------------------- ! 7) Alerts in case of model crash !---------------------------------------------------------------------------- numaltests = 10 numal(:) = 0 ! Alert if incompatible volume and concentration alert_id = 2 ! reference number of this alert alname(alert_id) = ' Incompat vol and con ' ! name of the alert DO jl = 1, jpl DO jj = 1, jpj DO ji = 1, jpi IF ((v_i(ji,jj,jl).NE.0.0).AND.(a_i(ji,jj,jl).EQ.0.0)) THEN WRITE(numout,*) ' ALERTE 2 ' WRITE(numout,*) ' Incompatible volume and concentration ' WRITE(numout,*) ' at_i ', at_i(ji,jj) WRITE(numout,*) ' Point - category', ji, jj, jl WRITE(numout,*) WRITE(numout,*) ' a_i *** a_i_old ', a_i(ji,jj,jl), old_a_i(ji,jj,jl) WRITE(numout,*) ' v_i *** v_i_old ', v_i(ji,jj,jl), old_v_i(ji,jj,jl) WRITE(numout,*) ' d_a_i_thd/trp ', d_a_i_thd(ji,jj,jl), d_a_i_trp(ji,jj,jl) WRITE(numout,*) ' d_v_i_thd/trp ', d_v_i_thd(ji,jj,jl), d_v_i_trp(ji,jj,jl) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO END DO ! Alerte if very thick ice alert_id = 3 ! reference number of this alert alname(alert_id) = ' Very thick ice ' ! name of the alert jl = jpl DO jj = 1, jpj DO ji = 1, jpi IF (ht_i(ji,jj,jl) .GT. 50.0 ) THEN WRITE(numout,*) ' ALERTE 3 ' WRITE(numout,*) ' Very thick ice ' CALL lim_inst_state(ji,jj,2) WRITE(numout,*) ' ' numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO ! Alert if very fast ice alert_id = 4 ! reference number of this alert alname(alert_id) = ' Very fast ice ' ! name of the alert DO jj = 1, jpj DO ji = 1, jpi IF ( ( ( ABS( u_ice(ji,jj) ) .GT. 0.50) .OR. & ( ABS( v_ice(ji,jj) ) .GT. 0.50) ) .AND. & ( at_i(ji,jj) .GT. 0.0 ) ) THEN WRITE(numout,*) ' ALERTE 4 ' WRITE(numout,*) ' Very fast ice ' CALL lim_inst_state(ji,jj,1) WRITE(numout,*) ' ice strength : ', strength(ji,jj) WRITE(numout,*) ' oceanic stress ftaux : ', ftaux(ji,jj) WRITE(numout,*) ' oceanic stress ftauy : ', ftauy(ji,jj) WRITE(numout,*) ' atmosph stress gtaux : ', gtaux(ji,jj) WRITE(numout,*) ' atmosph stress gtauy : ', gtauy(ji,jj) WRITE(numout,*) ' oceanic speed u : ', u_io(ji,jj) WRITE(numout,*) ' oceanic speed v : ', v_io(ji,jj) WRITE(numout,*) ' sst : ', sst_io(ji,jj) WRITE(numout,*) ' sss : ', sss_io(ji,jj) WRITE(numout,*) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO ! Alert if there is ice on continents alert_id = 6 ! reference number of this alert alname(alert_id) = ' Ice on continents ' ! name of the alert DO jj = 1, jpj DO ji = 1, jpi IF ( ( tms(ji,jj) .LE. 0.0 ) .AND. ( at_i(ji,jj) .GT. 0.0 ) ) THEN WRITE(numout,*) ' ALERTE 6 ' WRITE(numout,*) ' Ice on continents ' CALL lim_inst_state(ji,jj,1) WRITE(numout,*) ' masks s, u, v : ', tms(ji,jj), & tmu(ji,jj), & tmv(ji,jj) WRITE(numout,*) ' sst : ', sst_io(ji,jj) WRITE(numout,*) ' sss : ', sss_io(ji,jj) WRITE(numout,*) ' at_i(ji,jj) : ', at_i(ji,jj) WRITE(numout,*) ' v_ice(ji,jj) : ', v_ice(ji,jj) WRITE(numout,*) ' v_ice(ji,jj-1) : ', v_ice(ji,jj-1) WRITE(numout,*) ' u_ice(ji-1,jj) : ', u_ice(ji-1,jj) WRITE(numout,*) ' u_ice(ji,jj) : ', v_ice(ji,jj) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO ! ! ! Alert if very fresh ice alert_id = 7 ! reference number of this alert alname(alert_id) = ' Very fresh ice ' ! name of the alert DO jl = 1, jpl DO jj = 1, jpj DO ji = 1, jpi IF ( ( ( ABS( sm_i(ji,jj,jl) ) .LT. 0.50) .OR. & ( ABS( sm_i(ji,jj,jl) ) .LT. 0.50) ) .AND. & ( a_i(ji,jj,jl) .GT. 0.0 ) ) THEN ! WRITE(numout,*) ' ALERTE 7 ' ! WRITE(numout,*) ' Very fresh ice ' ! CALL lim_inst_state(ji,jj,1) ! WRITE(numout,*) ' sst : ', sst_io(ji,jj) ! WRITE(numout,*) ' sss : ', sss_io(ji,jj) ! WRITE(numout,*) ' s_i_newice : ', s_i_newice(ji,jj,1:jpl) ! WRITE(numout,*) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO END DO ! ! ! Alert if too old ice alert_id = 9 ! reference number of this alert alname(alert_id) = ' Very old ice ' ! name of the alert DO jl = 1, jpl DO jj = 1, jpj DO ji = 1, jpi IF ( ( ( ABS( o_i(ji,jj,jl) ) .GT. rdt_ice ) .OR. & ( ABS( o_i(ji,jj,jl) ) .LT. 0.00) ) .AND. & ( a_i(ji,jj,jl) .GT. 0.0 ) ) THEN WRITE(numout,*) ' ALERTE 9 ' WRITE(numout,*) ' Wrong ice age ' CALL lim_inst_state(ji,jj,1) WRITE(numout,*) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO END DO ! Alert on salt flux alert_id = 5 ! reference number of this alert alname(alert_id) = ' High salt flux ' ! name of the alert DO jj = 1, jpj DO ji = 1, jpi IF (ABS(fsalt(ji,jj)).gt.1.0e-2) THEN WRITE(numout,*) ' ALERTE 5 ' WRITE(numout,*) ' High salt flux ' CALL lim_inst_state(ji,jj,3) WRITE(numout,*) ' ' DO jl = 1, jpl WRITE(numout,*) ' Category no: ', jl WRITE(numout,*) ' a_i : ', a_i(ji,jj,jl) , & ' old_a_i : ', old_a_i(ji,jj,jl) WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(ji,jj,jl) , & ' d_a_i_thd : ', d_a_i_thd(ji,jj,jl) WRITE(numout,*) ' v_i : ', v_i(ji,jj,jl) , & ' old_v_i : ', old_v_i(ji,jj,jl) WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(ji,jj,jl) , & ' d_v_i_thd : ', d_v_i_thd(ji,jj,jl) WRITE(numout,*) ' ' END DO numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO ! Alert if fnsolar very big alert_id = 8 ! reference number of this alert alname(alert_id) = ' fnsolar very big ' ! name of the alert DO jj = 1, jpj DO ji = 1, jpi IF ( ( ABS( fnsolar(ji,jj) ) .GT. 1500.0 ) .AND. & ( at_i(ji,jj) .GT. 0.0 ) ) THEN WRITE(numout,*) ' ALERTE 8 ' WRITE(numout,*) ' ji, jj : ', ji, jj WRITE(numout,*) ' fnsolar : ', fnsolar(ji,jj) WRITE(numout,*) ' Very high non-solar heat flux ' WRITE(numout,*) ' sst : ', sst_io(ji,jj) WRITE(numout,*) ' sss : ', sss_io(ji,jj) WRITE(numout,*) ' qcmif : ', qcmif(jiindex,jjindex) WRITE(numout,*) ' qldif : ', qldif(jiindex,jjindex) WRITE(numout,*) ' qcmif : ', qcmif(jiindex,jjindex) / rdt_ice WRITE(numout,*) ' qldif : ', qldif(jiindex,jjindex) / rdt_ice WRITE(numout,*) ' qfvbq : ', qfvbq(jiindex,jjindex) WRITE(numout,*) ' qdtcn : ', qdtcn(jiindex,jjindex) WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindex,jjindex) / rdt_ice WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindex,jjindex) / rdt_ice WRITE(numout,*) ' fdtcn : ', fdtcn(jiindex,jjindex) WRITE(numout,*) ' fhmec : ', fhmec(jiindex,jjindex) WRITE(numout,*) ' fheat_rpo : ', fheat_rpo(jiindex,jjindex) WRITE(numout,*) ' fheat_res : ', fheat_res(jiindex,jjindex) WRITE(numout,*) ' fhbri : ', fhbri(jiindex,jjindex) CALL lim_inst_state(ji,jj,2) numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO !+++++ ! Alert if very warm ice alert_id = 10 ! reference number of this alert alname(alert_id) = ' Very warm ice ' ! name of the alert numal(alert_id) = 0 DO jl = 1, jpl DO jk = 1, nlay_i DO jj = 1, jpj DO ji = 1, jpi ztmelts = -tmut*s_i(ji,jj,jk,jl) + rtt IF ( ( t_i(ji,jj,jk,jl) .GE. ztmelts) .AND. & ( v_i(ji,jj,jl) .GT. 1.0e-6) .AND. & ( a_i(ji,jj,jl) .GT. 0.0 ) ) THEN WRITE(numout,*) ' ALERTE 10 ' WRITE(numout,*) ' ji, jj, jk, jl : ', ji, jj, jk, jl WRITE(numout,*) ' Very warm ice ' WRITE(numout,*) ' t_i : ', t_i(ji,jj,jk,jl) WRITE(numout,*) ' e_i : ', e_i(ji,jj,jk,jl) WRITE(numout,*) ' s_i : ', s_i(ji,jj,jk,jl) WRITE(numout,*) ' ztmelts : ', ztmelts numal(alert_id) = numal(alert_id) + 1 ENDIF END DO END DO END DO END DO alert_id = 1 ! reference number of this alert alname(alert_id) = ' Il n''y a pas d''alerte 1 ' ! name of the alert WRITE(numout,*) WRITE(numout,*) ' All alerts at the end of ice model ' DO alert_id = 1, numaltests WRITE(numout,*) alert_id, alname(alert_id)//' : ', numal(alert_id), ' times ! ' END DO ! ENDIF ! kt EQ ice time step END SUBROUTINE ice_stp !------------------------------------------------------------------------------ SUBROUTINE lim_inst_state(i,j,n) !!----------------------------------------------------------------------- !! *** ROUTINE lim_inst_state *** !! !! ** Purpose : Writes global ice state on the (i,j) point !! in ocean.ouput !! 3 possibilities exist !! n = 1 -> simple ice state !! n = 2 -> exhaustive state !! n = 3 -> ice/ocean salt fluxes !! !! ** input : point coordinates (i,j) !! n : number of the option !! !! history : !! 9.0 ! 06-12 (M. Vancoppenolle) Original code !!------------------------------------------------------------------- INTEGER, INTENT( in ) :: i,j,n ! ocean time-step index INTEGER :: jl !!------------------------------------------------------------------- !---------------- ! Simple state !---------------- IF ( n .EQ. 1 ) THEN WRITE(numout,*) ' lim_inst_state - Point : ',i,j WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' WRITE(numout,*) ' Simple state ' WRITE(numout,*) ' masks s,u,v : ', tms(i,j), tmu(i,j), tmv(i,j) WRITE(numout,*) ' lat - long : ', gphit(i,j), glamt(i,j) WRITE(numout,*) ' Time step : ', numit WRITE(numout,*) ' - Ice drift ' WRITE(numout,*) ' ~~~~~~~~~~~ ' WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(i-1,j) WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(i,j) WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(i,j-1) WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(i,j) WRITE(numout,*) ' strength : ', strength(i,j) WRITE(numout,*) WRITE(numout,*) ' - Cell values ' WRITE(numout,*) ' ~~~~~~~~~~~ ' WRITE(numout,*) ' cell area : ', area(i,j) WRITE(numout,*) ' at_i : ', at_i(i,j) WRITE(numout,*) ' vt_i : ', vt_i(i,j) WRITE(numout,*) ' vt_s : ', vt_s(i,j) DO jl = 1, jpl WRITE(numout,*) ' - Category (', jl,')' WRITE(numout,*) ' a_i : ', a_i(i,j,jl) WRITE(numout,*) ' ht_i : ', ht_i(i,j,jl) WRITE(numout,*) ' ht_s : ', ht_s(i,j,jl) WRITE(numout,*) ' v_i : ', v_i(i,j,jl) WRITE(numout,*) ' v_s : ', v_s(i,j,jl) WRITE(numout,*) ' e_s : ', e_s(i,j,1,jl)/1.0e9 WRITE(numout,*) ' e_i : ', e_i(i,j,1:nlay_i,jl)/1.0e9 WRITE(numout,*) ' t_su : ', t_su(i,j,jl) WRITE(numout,*) ' t_snow : ', t_s(i,j,1,jl) WRITE(numout,*) ' t_i : ', t_i(i,j,1:nlay_i,jl) WRITE(numout,*) ' sm_i : ', sm_i(i,j,jl) WRITE(numout,*) ' smv_i : ', smv_i(i,j,jl) WRITE(numout,*) WRITE(numout,*) ' Pathological case : ', patho_case(i,j,jl) END DO ENDIF !-------------------- ! Exhaustive state !-------------------- IF ( n .EQ. 2 ) THEN WRITE(numout,*) ' lim_inst_state - Point : ',i,j WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' WRITE(numout,*) ' Exhaustive state ' WRITE(numout,*) ' lat - long ', gphit(i,j), glamt(i,j) WRITE(numout,*) ' Time step ', numit WRITE(numout,*) WRITE(numout,*) ' - Cell values ' WRITE(numout,*) ' ~~~~~~~~~~~ ' WRITE(numout,*) ' cell area : ', area(i,j) WRITE(numout,*) ' at_i : ', at_i(i,j) WRITE(numout,*) ' vt_i : ', vt_i(i,j) WRITE(numout,*) ' vt_s : ', vt_s(i,j) WRITE(numout,*) ' sdvt : ', sdvt(i,j) WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(i-1,j) WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(i,j) WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(i,j-1) WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(i,j) WRITE(numout,*) ' strength : ', strength(i,j) WRITE(numout,*) ' d_u_ice_dyn : ', d_u_ice_dyn(i,j), & ' d_v_ice_dyn : ', d_v_ice_dyn(i,j) WRITE(numout,*) ' old_u_ice : ', old_u_ice(i,j) , & ' old_v_ice : ', old_v_ice(i,j) WRITE(numout,*) DO jl = 1, jpl WRITE(numout,*) ' - Category (',jl,')' WRITE(numout,*) ' ~~~~~~~~ ' WRITE(numout,*) ' ht_i : ', ht_i(i,j,jl) , & ' ht_s : ', ht_s(i,j,jl) WRITE(numout,*) ' t_i : ', t_i(i,j,1:nlay_i,jl) WRITE(numout,*) ' t_su : ', t_su(i,j,jl) , & ' t_s : ', t_s(i,j,1,jl) WRITE(numout,*) ' sm_i : ', sm_i(i,j,jl) , & ' o_i : ', o_i(i,j,jl) WRITE(numout,*) ' a_i : ', a_i(i,j,jl) , & ' old_a_i : ', old_a_i(i,j,jl) WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(i,j,jl) , & ' d_a_i_thd : ', d_a_i_thd(i,j,jl) WRITE(numout,*) ' v_i : ', v_i(i,j,jl) , & ' old_v_i : ', old_v_i(i,j,jl) WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(i,j,jl) , & ' d_v_i_thd : ', d_v_i_thd(i,j,jl) WRITE(numout,*) ' v_s : ', v_s(i,j,jl) , & ' old_v_s : ', old_v_s(i,j,jl) WRITE(numout,*) ' d_v_s_trp : ', d_v_s_trp(i,j,jl) , & ' d_v_s_thd : ', d_v_s_thd(i,j,jl) WRITE(numout,*) ' e_i1 : ', e_i(i,j,1,jl)/1.0e9 , & ' old_ei1 : ', old_e_i(i,j,1,jl)/1.0e9 WRITE(numout,*) ' de_i1_trp : ', d_e_i_trp(i,j,1,jl)/1.0e9, & ' de_i1_thd : ', d_e_i_thd(i,j,1,jl)/1.0e9 WRITE(numout,*) ' e_i2 : ', e_i(i,j,2,jl)/1.0e9 , & ' old_ei2 : ', old_e_i(i,j,2,jl)/1.0e9 WRITE(numout,*) ' de_i2_trp : ', d_e_i_trp(i,j,2,jl)/1.0e9, & ' de_i2_thd : ', d_e_i_thd(i,j,2,jl)/1.0e9 WRITE(numout,*) ' e_snow : ', e_s(i,j,1,jl) , & ' old_e_snow : ', old_e_s(i,j,1,jl) WRITE(numout,*) ' d_e_s_trp : ', d_e_s_trp(i,j,1,jl) , & ' d_e_s_thd : ', d_e_s_thd(i,j,1,jl) WRITE(numout,*) ' smv_i : ', smv_i(i,j,jl) , & ' old_smv_i : ', old_smv_i(i,j,jl) WRITE(numout,*) ' d_smv_i_trp: ', d_smv_i_trp(i,j,jl) , & ' d_smv_i_thd: ', d_smv_i_thd(i,j,jl) WRITE(numout,*) ' oa_i : ', oa_i(i,j,jl) , & ' old_oa_i : ', old_oa_i(i,j,jl) WRITE(numout,*) ' d_oa_i_trp : ', d_oa_i_trp(i,j,jl) , & ' d_oa_i_thd : ', d_oa_i_thd(i,j,jl) WRITE(numout,*) ' Path. case : ', patho_case(i,j,jl) END DO !jl WRITE(numout,*) WRITE(numout,*) ' - Heat / FW fluxes ' WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' ! WRITE(numout,*) ' fsalt : ', fsalt(i,j) ! WRITE(numout,*) ' fmass : ', fmass(i,j) ! WRITE(numout,*) ' fsbri : ', fsbri(i,j) ! WRITE(numout,*) ' fseqv : ', fseqv(i,j) ! WRITE(numout,*) ' fsalt_res : ', fsalt_res(i,j) WRITE(numout,*) ' fmmec : ', fmmec(i,j) WRITE(numout,*) ' fhmec : ', fhmec(i,j) WRITE(numout,*) ' fhbri : ', fhbri(i,j) WRITE(numout,*) ' fheat_rpo : ', fheat_rpo(i,j) WRITE(numout,*) WRITE(numout,*) ' sst : ', sst_io(i,j) WRITE(numout,*) ' sss : ', sss_io(i,j) WRITE(numout,*) WRITE(numout,*) ' - Stresses ' WRITE(numout,*) ' ~~~~~~~~ ' WRITE(numout,*) ' tauxw : ', tauxw(i,j) WRITE(numout,*) ' tauyw : ', tauyw(i,j) WRITE(numout,*) ' taux : ', taux(i,j) WRITE(numout,*) ' tauy : ', tauy(i,j) WRITE(numout,*) ' ftaux : ', ftaux(i,j) WRITE(numout,*) ' ftauy : ', ftauy(i,j) WRITE(numout,*) ' gtaux : ', gtaux(i,j) WRITE(numout,*) ' gtauy : ', gtauy(i,j) WRITE(numout,*) ' oc. vel. u : ', u_io(i,j) WRITE(numout,*) ' oc. vel. v : ', v_io(i,j) ENDIF !--------------------- ! Salt / heat fluxes !--------------------- IF ( n .EQ. 3 ) THEN WRITE(numout,*) ' lim_inst_state - Point : ',i,j WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' WRITE(numout,*) ' - Salt / Heat Fluxes ' WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' WRITE(numout,*) ' lat - long ', gphit(i,j), glamt(i,j) WRITE(numout,*) ' Time step ', numit WRITE(numout,*) WRITE(numout,*) ' - Heat fluxes at bottom interface ***' WRITE(numout,*) ' fsolar : ', fsolar(i,j) WRITE(numout,*) ' fnsolar : ', fnsolar(i,j) WRITE(numout,*) WRITE(numout,*) ' - Salt fluxes at bottom interface ***' WRITE(numout,*) ' fsalt : ', fsalt(i,j) WRITE(numout,*) ' fmass : ', fmass(i,j) WRITE(numout,*) ' fsbri : ', fsbri(i,j) WRITE(numout,*) ' fseqv : ', fseqv(i,j) WRITE(numout,*) ' fsalt_res : ', fsalt_res(i,j) WRITE(numout,*) ' fsalt_rpo : ', fsalt_rpo(i,j) WRITE(numout,*) ' - Heat fluxes at bottom interface ***' WRITE(numout,*) ' fheat_res : ', fheat_res(i,j) WRITE(numout,*) WRITE(numout,*) ' - Momentum fluxes ' WRITE(numout,*) ' ftaux : ', ftaux(i,j) WRITE(numout,*) ' ftauy : ', ftauy(i,j) WRITE(numout,*) ' gtaux : ', gtaux(i,j) WRITE(numout,*) ' gtauy : ', gtauy(i,j) ENDIF END SUBROUTINE #else !!---------------------------------------------------------------------- !! Default option NO LIM sea-ice model !!---------------------------------------------------------------------- USE in_out_manager CONTAINS SUBROUTINE ice_stp ( kt ) ! Empty routine INTEGER, INTENT( in ) :: kt ! ocean time-step index IF( kt == nit000 ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'No Sea Ice Model' IF(lwp) WRITE(numout,*) '~~~~~~~' ENDIF END SUBROUTINE ice_stp SUBROUTINE lim_inst_state(i,j,n) ! Empty routine INTEGER, INTENT( in ) :: i, j, n END SUBROUTINE #endif !!====================================================================== END MODULE icestp