Changeset 1857
- Timestamp:
- 2010-05-03T13:59:46+02:00 (14 years ago)
- Location:
- branches/DEV_r1837_mass_heat_salt_fluxes/NEMO
- Files:
-
- 21 edited
- 1 copied
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/dom_ice_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE dom_ice *** 4 !! LIM -2Sea Ice : Domain variables4 !! LIM 2.0 Sea Ice : Domain variables 5 5 !!====================================================================== 6 !! History : 2.0 ! 2003-08 (C. Ethe) Free form and module6 !! History : 2.0 ! 03-08 (C. Ethe) Free form and module 7 7 !!---------------------------------------------------------------------- 8 8 #if defined key_lim2 9 9 !!---------------------------------------------------------------------- 10 !! 'key_lim2' : LIM 2.0 sea-ice model 10 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 11 !! $Id$ 12 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 11 13 !!---------------------------------------------------------------------- 12 14 USE par_ice_2 … … 16 18 17 19 LOGICAL, PUBLIC :: l_jeq = .TRUE. !: Equator inside the domain flag 20 18 21 INTEGER, PUBLIC :: njeq , njeqm1 !: j-index of the equator if it is inside the domain 22 ! ! (otherwise = jpj+10 (SH) or -10 (SH) ) 19 23 20 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fs2cor , fcor!: coriolis factor and coeficient21 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: covrai!: sine of geographic latitude22 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: area!: surface of grid cell23 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: tms , tmu!: temperature and velocity points masks24 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2) :: wght 25 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2) :: akappa, bkappa!: first and third group of metric coefficients26 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2,2,2) :: alambd 24 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fs2cor , fcor, & !: coriolis factor and coeficient 25 & covrai , & !: sine of geographic latitude 26 & area , & !: surface of grid cell 27 & tms , tmu !: temperature and velocity points masks 28 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2) :: wght , & !: weight of the 4 neighbours to compute averages 29 & akappa , bkappa !: first and third group of metric coefficients 30 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2,2,2) :: alambd !: second group of metric coefficients 27 31 28 #else 29 !!---------------------------------------------------------------------- 30 !! Default option Dummy module NO LIM 2.0 sea-ice model 31 !!---------------------------------------------------------------------- 32 !!====================================================================== 32 33 #endif 33 34 !!-------- -------------------------------------------------------------35 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)36 !! $Id$37 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)38 !!======================================================================39 34 END MODULE dom_ice_2 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/ice_2.F90
r1855 r1857 4 4 !! Sea Ice physics: diagnostics variables of ice defined in memory 5 5 !!===================================================================== 6 !! History : 2.0 ! 2003-08 (C. Ethe) F90: Free form and module 7 !! - ! 2010-05 (G. Madec) suppression of thd_ice_2 6 !! History : 2.0 ! 03-08 (C. Ethe) F90: Free form and module 8 7 !!---------------------------------------------------------------------- 9 8 #if defined key_lim2 … … 11 10 !! 'key_lim2' : LIM 2.0 sea-ice model 12 11 !!---------------------------------------------------------------------- 12 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 13 !! $Id$ 14 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 15 !!---------------------------------------------------------------------- 16 !! * Modules used 13 17 USE par_ice_2 ! LIM sea-ice parameters 14 18 … … 16 20 PRIVATE 17 21 18 !---------------------------------------------------! 19 ! Sea-Ice namelist ! 20 !---------------------------------------------------! 21 ! !!! ** Namelist namicerun ** 22 CHARACTER(len=32), PUBLIC :: cn_icerst_in = "restart_ice_in" !: suffix of ice restart name (input) 23 CHARACTER(len=32), PUBLIC :: cn_icerst_out = "restart_ice" !: suffix of ice restart name (output) 24 LOGICAL , PUBLIC :: ln_limdyn = .TRUE. !: flag for ice dynamics (T) or not (F) 25 LOGICAL , PUBLIC :: ln_limdmp = .FALSE. !: Ice damping 26 REAL(wp) , PUBLIC :: hsndif = 0.e0 !: computation of temp. in snow (0) or not (9999) 27 REAL(wp) , PUBLIC :: hicdif = 0.e0 !: computation of temp. in ice (0) or not (9999) 28 REAL(wp) , PUBLIC :: acrit(2) = (/ 1.e-06 , 1.e-06 /) !: min lead fraction in north and south hemispheres 29 30 ! !!! ** Namelist namicedyn ** 31 INTEGER , PUBLIC :: nbiter = 1 !: number of sub-time steps for relaxation 32 INTEGER , PUBLIC :: nbitdr = 250 !: maximum number of iterations for relaxation 33 REAL(wp), PUBLIC :: rdt_ice !: ice time step 34 REAL(wp), PUBLIC :: epsd = 1.0e-20 !: tolerance parameter for dynamic 35 REAL(wp), PUBLIC :: alpha = 0.5 !: coefficient for semi-implicit coriolis 36 REAL(wp), PUBLIC :: dm = 0.6e+03 !: diffusion constant for dynamics 37 REAL(wp), PUBLIC :: om = 0.5 !: relaxation constant 38 REAL(wp), PUBLIC :: resl = 5.0e-05 !: maximum value for the residual of relaxation 39 REAL(wp), PUBLIC :: cw = 5.0e-03 !: drag coefficient for oceanic stress 40 REAL(wp), PUBLIC :: angvg = 0.e0 !: turning angle for oceanic stress 41 REAL(wp), PUBLIC :: pstar = 1.0e+04 !: first bulk-rheology parameter 42 REAL(wp), PUBLIC :: c_rhg = 20.e0 !: second bulk-rhelogy parameter 43 REAL(wp), PUBLIC :: etamn = 0.e+07 !: minimun value for viscosity 44 REAL(wp), PUBLIC :: creepl = 2.e-08 !: creep limit 45 REAL(wp), PUBLIC :: ecc = 2.e0 !: eccentricity of the elliptical yield curve 46 REAL(wp), PUBLIC :: ahi0 = 350.e0 !: sea-ice hor. eddy diffusivity coeff. (m2/s) 22 !!* Share parameters namelist (namicerun read in iceini) * 23 CHARACTER(len=32) , PUBLIC :: cn_icerst_in = "restart_ice_in" !: suffix of ice restart name (input) 24 CHARACTER(len=32) , PUBLIC :: cn_icerst_out = "restart_ice" !: suffix of ice restart name (output) 25 LOGICAL , PUBLIC :: ln_limdyn = .TRUE. !: flag for ice dynamics (T) or not (F) 26 LOGICAL , PUBLIC :: ln_limdmp = .FALSE. !: Ice damping 27 REAL(wp) , PUBLIC :: hsndif = 0.e0 !: computation of temp. in snow (0) or not (9999) 28 REAL(wp) , PUBLIC :: hicdif = 0.e0 !: computation of temp. in ice (0) or not (9999) 29 REAL(wp), DIMENSION(2), PUBLIC :: acrit = (/ 1.e-06 , 1.e-06 /) !: minimum fraction for leads in 30 ! !: north and south hemisphere 31 !!* ice-dynamic namelist (namicedyn) * 32 INTEGER , PUBLIC :: nbiter = 1 !: number of sub-time steps for relaxation 33 INTEGER , PUBLIC :: nbitdr = 250 !: maximum number of iterations for relaxation 34 REAL(wp), PUBLIC :: rdt_ice !: ice time step 35 REAL(wp), PUBLIC :: epsd = 1.0e-20 !: tolerance parameter for dynamic 36 REAL(wp), PUBLIC :: alpha = 0.5 !: coefficient for semi-implicit coriolis 37 REAL(wp), PUBLIC :: dm = 0.6e+03 !: diffusion constant for dynamics 38 REAL(wp), PUBLIC :: om = 0.5 !: relaxation constant 39 REAL(wp), PUBLIC :: resl = 5.0e-05 !: maximum value for the residual of relaxation 40 REAL(wp), PUBLIC :: cw = 5.0e-03 !: drag coefficient for oceanic stress 41 REAL(wp), PUBLIC :: angvg = 0.e0 !: turning angle for oceanic stress 42 REAL(wp), PUBLIC :: pstar = 1.0e+04 !: first bulk-rheology parameter 43 REAL(wp), PUBLIC :: c_rhg = 20.e0 !: second bulk-rhelogy parameter 44 REAL(wp), PUBLIC :: etamn = 0.e+07 !: minimun value for viscosity 45 REAL(wp), PUBLIC :: creepl = 2.e-08 !: creep limit 46 REAL(wp), PUBLIC :: ecc = 2.e0 !: eccentricity of the elliptical yield curve 47 REAL(wp), PUBLIC :: ahi0 = 350.e0 !: sea-ice hor. eddy diffusivity coeff. (m2/s) 47 48 48 ! !!! ** Namelist namicethd ** 49 REAL(wp), PUBLIC :: hmelt = -0.15 !: maximum melting at the bottom 50 REAL(wp), PUBLIC :: hicmin = 0.2 !: ice th. corr. to max. ener. in brine pocket 51 REAL(wp), PUBLIC :: hiclim = 0.05 !: minimum ice thickness 52 REAL(wp), PUBLIC :: amax = 0.999 !: maximum lead fraction 53 REAL(wp), PUBLIC :: swiqst = 1.0 !: energy stored in brine pocket (1) or not (0) 54 REAL(wp), PUBLIC :: sbeta = 1.0 !: numerical scheme for diffusion in ice 55 REAL(wp), PUBLIC :: parlat = 0.0 !: percent. of energy used for lateral ablation 56 REAL(wp), PUBLIC :: hakspl = 0.5 !: slope of distr. for Hakkinen-Mellro's lat. melt 57 REAL(wp), PUBLIC :: hibspl = 0.5 !: slope of distribution for Hibler's lat. melt 58 REAL(wp), PUBLIC :: exld = 2.0 !: exponent for leads-closure rate 59 REAL(wp), PUBLIC :: hakdif = 1.0 !: coefficient for diffusions of ice and snow 60 REAL(wp), PUBLIC :: thth = 0.2 !: thick. for comp. of eq. thermal conduct 61 REAL(wp), PUBLIC :: hnzst = 0.1 !: thick. of the surf. layer in temp. comp. 62 REAL(wp), PUBLIC :: parsub = 1.0 !: switch for snow sublimation or not 63 REAL(wp), PUBLIC :: alphs = 1.0 !: coef. for snow density when snow-ice formation 64 REAL(wp), PUBLIC :: hiccrit(2) = (/0.3,0.3/) !: ice th. for lateral accretion in the NH (SH) (m) 65 66 !---------------------------------------------------! 67 ! Sea-Ice variables ! 68 !---------------------------------------------------! 69 70 REAL(wp), PUBLIC :: usecc2 !: = 1.0 / ( ecc * ecc ) 71 REAL(wp), PUBLIC :: rhoco !: = rau0 * cw 72 REAL(wp), PUBLIC :: sangvg, cangvg !: sin and cos of the turning angle for ocean stress 73 REAL(wp), PUBLIC :: pstarh !: pstar / 2.0 74 REAL(wp), PUBLIC :: uscomi !: inverse of minimum lead fraction !!gm to be suppress 75 REAL(wp), PUBLIC :: cnscg !: ratio rcpsn/rcpic !!gm to be suppress 49 REAL(wp), PUBLIC :: usecc2 !: = 1.0 / ( ecc * ecc ) 50 REAL(wp), PUBLIC :: rhoco !: = rau0 * cw 51 REAL(wp), PUBLIC :: sangvg, cangvg !: sin and cos of the turning angle for ocean stress 52 REAL(wp), PUBLIC :: pstarh !: pstar / 2.0 76 53 77 54 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: ahiu , ahiv !: hor. diffusivity coeff. at ocean U- and V-points (m2/s) 78 55 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: pahu , pahv !: ice hor. eddy diffusivity coef. at ocean U- and V-points 79 56 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: hsnm , hicm !: mean snow and ice thicknesses 80 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: ust2s !: friction velocity57 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: ust2s !: friction velocity 81 58 82 INTEGER , PUBLIC, DIMENSION(jpij) :: npb !: number of points where computations has to be done 83 INTEGER , PUBLIC, DIMENSION(jpij) :: npac !: correspondance between the points 59 !!* diagnostic quantities 60 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: sist !: Sea-Ice Surface Temperature (Kelvin) 61 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: tfu !: Freezing/Melting point temperature of sea water at SSS 62 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: hicif !: Ice thickness 63 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: hsnif !: Snow thickness 64 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: hicifp !: Ice production/melting 65 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: frld !: Leads fraction = 1-a/totalarea 66 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: phicif !: ice thickness at previous time 67 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: pfrld !: Leads fraction at previous time 68 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qstoif !: Energy stored in the brine pockets 69 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fbif !: Heat flux at the ice base 70 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: rdmsnif !: Variation of snow mass 71 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: rdmicif !: Variation of ice mass 72 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qldif !: heat balance of the lead (or of the open ocean) 73 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qcmif !: Energy needed to bring the ocean surface layer until its freezing 74 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fdtcn !: net downward heat flux from the ice to the ocean 75 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qdtcn !: energy from the ice to the ocean point (at a factor 2) 76 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: thcm !: part of the solar energy used in the lead heat budget 77 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fstric !: Solar flux transmitted trough the ice 78 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: ffltbif !: Array linked with the max heat contained in brine pockets (?) 79 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fscmbq !: Linked with the solar flux below the ice (?) 80 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fsbbq !: Also linked with the solar flux below the ice (?) 81 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qfvbq !: Array used to store energy in case of toral lateral ablation (?) 82 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: dmgwi !: Variation of the mass of snow ice 84 83 85 ! !!* sea-ice variables * 86 REAL(wp), PUBLIC :: sist (jpi,jpj), sist_1d (jpij) !: Sea-Ice Surface Temperature (Kelvin) 87 REAL(wp), PUBLIC :: tfu (jpi,jpj), tfu_1d (jpij) !: Freezing/Melting point temperature of sea water at SSS 88 REAL(wp), PUBLIC :: hicif (jpi,jpj), h_ice_1d (jpij) !: ice thickness 89 REAL(wp), PUBLIC :: phicif (jpi,jpj) !: ice thickness at previous time 90 REAL(wp), PUBLIC :: hsnif (jpi,jpj), h_snow_1d (jpij) !: Snow thickness 91 REAL(wp), PUBLIC :: hicifp (jpi,jpj) !: Ice production/melting 92 REAL(wp), PUBLIC :: frld (jpi,jpj), frld_1d (jpij) !: Leads fraction = 1-a/totalarea 93 REAL(wp), PUBLIC :: pfrld (jpi,jpj) !: Leads fraction at previous time 94 REAL(wp), PUBLIC :: qstoif (jpi,jpj), qstbif_1d (jpij) !: Energy stored in the brine pockets !!gm 95 REAL(wp), PUBLIC :: fbif (jpi,jpj), fbif_1d (jpij) !: Heat flux at the ice base 96 REAL(wp), PUBLIC :: rdmsnif(jpi,jpj), rdmsnif_1d(jpij) !: Variation of snow mass 97 REAL(wp), PUBLIC :: rdmicif(jpi,jpj), rdmicif_1d(jpij) !: Variation of ice mass 98 REAL(wp), PUBLIC :: qldif (jpi,jpj), qldif_1d (jpij) !: heat balance of the lead (or of the open ocean) 99 REAL(wp), PUBLIC :: qcmif (jpi,jpj), qcmif_1d (jpij) !: Energy needed to freeze the ocean surface layer 100 REAL(wp), PUBLIC :: fdtcn (jpi,jpj) !: net downward heat flux from the ice to the ocean 101 REAL(wp), PUBLIC :: qdtcn (jpi,jpj) !: energy from the ice to the ocean point (at a factor 2) 102 REAL(wp), PUBLIC :: thcm (jpi,jpj), thcm_1d (jpij) !: part of the solar energy used in the lead heat budget 103 REAL(wp), PUBLIC :: fstric (jpi,jpj), fstric_1d (jpij) !: Solar flux transmitted trough the ice 104 REAL(wp), PUBLIC :: ffltbif(jpi,jpj), fltbif_1d (jpij) !: Array linked with the max heat contained in brine pockets (?) 105 REAL(wp), PUBLIC :: fscmbq (jpi,jpj), fscbq_1d (jpij) !: Linked with the solar flux below the ice (?) 106 REAL(wp), PUBLIC :: fsbbq (jpi,jpj) !: Also linked with the solar flux below the ice (?) 107 REAL(wp), PUBLIC :: qfvbq (jpi,jpj), qfvbq_1d (jpij) !: Array used to store energy in case of total lateral ablation (?) 108 REAL(wp), PUBLIC :: dmgwi(jpi,jpj), dmgwi_1d(jpij) !: Variation of the mass of snow ice 84 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: u_ice, v_ice !: two components of the ice velocity at I-point (m/s) 85 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: u_oce, v_oce !: two components of the ocean velocity at I-point (m/s) 109 86 110 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: u_ice, v_ice !: two components of the ice velocity at I-point (m/s) 111 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: u_oce, v_oce !: two components of the ocean velocity at I-point (m/s) 87 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jplayersp1) :: tbif !: Temperature inside the ice/snow layer 112 88 113 REAL(wp), PUBLIC :: tbif(jpi,jpj,jplayersp1), tbif_1d(jpij,jplayersp1) !: Temperature inside the ice+snow layer 114 115 ! Surface forcing transford into 1D field (2D field defined in /OPA_SRC/SBC/sbc_ice.F90) 116 REAL(wp), PUBLIC, DIMENSION(jpij) :: qns_ice_1d !: qns_ice 117 REAL(wp), PUBLIC, DIMENSION(jpij) :: qsr_ice_1d !: qsr_ice 118 REAL(wp), PUBLIC, DIMENSION(jpij) :: qla_ice_1d !: qla_ice 119 REAL(wp), PUBLIC, DIMENSION(jpij) :: dqla_ice_1d !: dqla_ice 120 REAL(wp), PUBLIC, DIMENSION(jpij) :: dqns_ice_1d !: dqns_ice 121 REAL(wp), PUBLIC, DIMENSION(jpij) :: fr1_i0_1d !: fr1_i0 122 REAL(wp), PUBLIC, DIMENSION(jpij) :: fr2_i0_1d !: fr2_i0 123 124 ! Surface forcing transformed into 1D field (2D field defined in /OPA_SRC/SBC/sbc_oce.F90) 125 REAL(wp), PUBLIC, DIMENSION(jpij) :: sprecip_1d !: sprecip 126 127 ! local variable transformed into 1D field (2D field defined in /LIM_SRC_2/limthd_2.F90) 128 REAL(wp), PUBLIC, DIMENSION(jpij) :: qlbbq_1d !: zqlbsbq 129 REAL(wp), PUBLIC, DIMENSION(jpij) :: dvsbq_1d !: zdvosif 130 REAL(wp), PUBLIC, DIMENSION(jpij) :: dvbbq_1d !: zdvobif 131 REAL(wp), PUBLIC, DIMENSION(jpij) :: dvlbq_1d !: zdvolif 132 REAL(wp), PUBLIC, DIMENSION(jpij) :: dvnbq_1d !: zdvonif 133 134 135 REAL(wp), PUBLIC, DIMENSION(jpij) :: rdvomif_1d !: rdvomif <<<== !gm should be suppressed ! 136 137 !---------------------------------------------------! 138 ! Prather' advection scheme : 1st & 2nd moments ! 139 !---------------------------------------------------! 89 !!* moment used in the advection scheme 140 90 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: sxice, syice, sxxice, syyice, sxyice !: for ice volume 141 91 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: sxsn, sysn, sxxsn, syysn, sxysn !: for snow volume … … 152 102 #endif 153 103 154 !!-------- -------------------------------------------------------------155 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)156 !! $Id$157 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)158 104 !!====================================================================== 159 105 END MODULE ice_2 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/iceini_2.F90
r1855 r1857 4 4 !! Sea-ice model : LIM 2.0 Sea ice model Initialization 5 5 !!====================================================================== 6 !! History : 1.0 ! 2002-08 (G. Madec) F90: Free form and modules7 !! 2.0 ! 2003-08 (C. Ethe) add ice_run6 !! History : 1.0 ! 02-08 (G. Madec) F90: Free form and modules 7 !! 2.0 ! 03-08 (C. Ethe) add ice_run 8 8 !!---------------------------------------------------------------------- 9 9 #if defined key_lim2 … … 11 11 !! 'key_lim2' : LIM 2.0 sea-ice model 12 12 !!---------------------------------------------------------------------- 13 !! ice_init_2 : sea-ice model initialization14 13 !!---------------------------------------------------------------------- 15 USE dom_oce ! ocean domain 16 USE dom_ice_2 ! LIM-2 ice domain 17 USE sbc_oce ! surface boundary condition: ocean 18 USE sbc_ice ! surface boundary condition: sea ice 19 USE phycst ! physical constant 20 USE ice_2 ! LIM-2 ice varibles 21 USE limmsh_2 ! LIM-2 mesh 22 USE limistate_2 ! LIM-2 initial state 23 USE limrst_2 ! LIM-2 restart 24 USE in_out_manager ! I/O manager 14 !! ice_init_2 : sea-ice model initialization 15 !! ice_run_2 : Definition some run parameter for ice model 16 !!---------------------------------------------------------------------- 17 USE dom_oce 18 USE dom_ice_2 19 USE sbc_oce ! surface boundary condition: ocean 20 USE sbc_ice ! surface boundary condition: ice 21 USE phycst ! Define parameters for the routines 22 USE ice_2 23 USE limmsh_2 24 USE limistate_2 25 USE limrst_2 26 USE in_out_manager 25 27 26 28 IMPLICIT NONE 27 29 PRIVATE 28 30 29 PUBLIC ice_init_2 ! called by sbcice_lim_2.F9031 PUBLIC ice_init_2 ! called by sbcice_lim_2.F90 30 32 31 33 !!---------------------------------------------------------------------- 32 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)34 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 33 35 !! $Id$ 34 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)36 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 35 37 !!---------------------------------------------------------------------- 36 38 … … 43 45 !! ** purpose : 44 46 !!---------------------------------------------------------------------- 47 ! 48 ! Open the namelist file 49 CALL ctl_opn( numnam_ice, 'namelist_ice', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 50 CALL ice_run_2 ! read in namelist some run parameters 51 52 ! Louvain la Neuve Ice model 53 rdt_ice = nn_fsbc * rdttra(1) 54 55 CALL lim_msh_2 ! ice mesh initialization 56 57 ! Initial sea-ice state 58 IF( .NOT.ln_rstart ) THEN 59 CALL lim_istate_2 ! start from rest: sea-ice deduced from sst 60 ELSE 61 CALL lim_rst_read_2 ! start from a restart file 62 ENDIF 63 64 tn_ice(:,:,1) = sist(:,:) ! initialisation of ice temperature 65 fr_i (:,:) = 1.0 - frld(:,:) ! initialisation of sea-ice fraction 66 ! 67 END SUBROUTINE ice_init_2 68 69 70 SUBROUTINE ice_run_2 71 !!------------------------------------------------------------------- 72 !! *** ROUTINE ice_run_2 *** 73 !! 74 !! ** Purpose : Definition some run parameter for ice model 75 !! 76 !! ** Method : Read the namicerun namelist and check the parameter 77 !! values called at the first timestep (nit000) 78 !! 79 !! ** input : Namelist namicerun 80 !!------------------------------------------------------------------- 45 81 NAMELIST/namicerun/ cn_icerst_in, cn_icerst_out, ln_limdyn, ln_limdmp, acrit, hsndif, hicdif 46 82 !!------------------------------------------------------------------- 47 ! 48 ! ! Open the ice namelist file 49 CALL ctl_opn( numnam_ice, 'namelist_ice', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 50 51 REWIND ( numnam_ice ) ! Read namicerun namelist 83 ! 84 REWIND ( numnam_ice ) ! Read Namelist namicerun 52 85 READ ( numnam_ice , namicerun ) 53 ! 54 IF(lwp) THEN ! control print86 87 IF(lwp) THEN 55 88 WRITE(numout,*) 56 89 WRITE(numout,*) 'ice_run : ice share parameters for dynamics/advection/thermo of sea-ice' … … 62 95 WRITE(numout,*) ' computation of temp. in ice (=0) or not (=9999) hicdif = ', hicdif 63 96 ENDIF 64 65 rdt_ice = nn_fsbc * rdttra(1) ! ice time step66 67 CALL lim_msh_2 ! ice mesh initialization68 69 ! ! Initial sea-ice state70 IF( .NOT.ln_rstart ) THEN ; CALL lim_istate_2 ! start from rest: sea-ice deduced from sst71 ELSE ; CALL lim_rst_read_2 ! restart from a file72 ENDIF73 74 ! ! ice variables see by the ocean75 tn_ice(:,:,1) = sist(:,:) ! surface ice temperature76 fr_i (:,:) = 1.0 - frld(:,:) ! sea-ice fraction77 97 ! 78 END SUBROUTINE ice_ init_298 END SUBROUTINE ice_run_2 79 99 80 100 #else … … 82 102 !! Default option : Empty module NO LIM 2.0 sea-ice model 83 103 !!---------------------------------------------------------------------- 104 CONTAINS 105 SUBROUTINE ice_init_2 ! Empty routine 106 END SUBROUTINE ice_init_2 84 107 #endif 85 108 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limadv_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limadv_2 *** 4 !! LIM -2sea-ice model : sea-ice advection4 !! LIM 2.0 sea-ice model : sea-ice advection 5 5 !!====================================================================== 6 6 !! History : OPA ! 2000-01 (LIM) Original code -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdia_2.F90
r1855 r1857 12 12 !! 'key_lim2' : LIM 2.0 sea-ice model 13 13 !!---------------------------------------------------------------------- 14 !!---------------------------------------------------------------------- 14 15 !! lim_dia_2 : computation of the time evolution of keys var. 15 16 !! lim_dia_init_2 : initialization and namelist read … … 27 28 PRIVATE 28 29 29 PUBLIC lim_dia_2 ! called by sbc_ice_lim_2 30 31 INTEGER, PUBLIC :: ntmoy = 1 !: instantaneous values of ice evolution or averaging ntmoy 32 INTEGER, PUBLIC :: ninfo = 1 !: frequency of ouputs on file ice_evolu in case of averaging 33 34 ! !!! Parameters for outputs to files "evolu" 35 INTEGER, PARAMETER :: jpinfmx = 100 ! maximum number of key variables 36 INTEGER, PARAMETER :: jpchinf = 5 ! ??? 37 INTEGER, PARAMETER :: jpchsep = jpchinf + 2 ! ??? 38 39 INTEGER :: nfrinf = 4 ! number of variables written in one line 40 INTEGER :: nferme ! last time step at which the var. are written on file 41 INTEGER :: nvinfo ! number of total variables 42 INTEGER :: nbvt ! number of time variables 43 INTEGER :: naveg ! number of step for accumulation before averaging 44 45 CHARACTER(len= 8) :: fmtinf = '1PE13.5 ' ! format of the output values 46 CHARACTER(len=30) :: fmtw, fmtr, fmtitr ! various formats 47 CHARACTER(len=jpchsep), DIMENSION(jpinfmx) :: titvar ! title of key variables 30 PUBLIC lim_dia_2 ! called by sbc_ice_lim_2 31 INTEGER, PUBLIC :: ntmoy = 1 , & !: instantaneous values of ice evolution or averaging ntmoy 32 & ninfo = 1 !: frequency of ouputs on file ice_evolu in case of averaging 33 34 INTEGER, PARAMETER :: & ! Parameters for outputs to files "evolu" 35 jpinfmx = 100 , & ! maximum number of key variables 36 jpchinf = 5 , & ! ??? 37 jpchsep = jpchinf + 2 ! ??? 38 39 INTEGER :: & 40 nfrinf = 4 , & ! number of variables written in one line 41 nferme , & ! last time step at which the var. are written on file 42 nvinfo , & ! number of total variables 43 nbvt , & ! number of time variables 44 naveg ! number of step for accumulation before averaging 45 46 CHARACTER(len= 8) :: fmtinf = '1PE13.5 ' ! format of the output values 47 CHARACTER(len=30) :: fmtw , & ! formats 48 & fmtr , & ! ??? 49 & fmtitr ! ??? 50 CHARACTER(len=jpchsep), DIMENSION(jpinfmx) :: titvar ! title of key variables 48 51 49 52 REAL(wp) :: epsi06 = 1.e-06 ! ??? … … 54 57 # include "vectopt_loop_substitute.h90" 55 58 !!---------------------------------------------------------------------- 56 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)59 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 57 60 !! $Id$ 58 61 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 68 71 !! the temporal evolution of some key variables 69 72 !!------------------------------------------------------------------- 70 INTEGER, INTENT(in) :: kt ! number of iteration73 INTEGER, INTENT(in) :: kt ! number of iteration 71 74 !! 72 75 INTEGER :: jv,ji, jj ! dummy loop indices 73 76 INTEGER :: nv ! indice of variable 74 REAL(wp) :: zarea , zldarea 75 REAL(wp) :: zextent15, zextent85! sea-ice extent (15% and 85%)76 REAL(wp) :: zicevol , zsnwvol! sea-ice and snow volume volume77 REAL(wp) ::zicespd ! sea-ice velocity77 REAL(wp) :: zarea , zldarea , & ! sea-ice and leads area 78 & zextent15, zextent85, & ! sea-ice extent (15% and 85%) 79 & zicevol , zsnwvol , & ! sea-ice and snow volume volume 80 & zicespd ! sea-ice velocity 78 81 REAL(wp), DIMENSION(jpinfmx) :: vinfor ! temporary working space 79 82 !!------------------------------------------------------------------- … … 117 120 vinfor(nv+13) = SQRT( vinfor(nv+13) / MAX( vinfor(nv+9) , epsi06 ) ) 118 121 122 119 123 ! variables in southern Hemis 120 124 nv = nv + 1 … … 173 177 INTEGER :: nv ! indice of variable 174 178 REAL(wp) :: zxx0, zxx1 ! temporary scalars 175 !! 179 176 180 NAMELIST/namicedia/fmtinf, nfrinf, ninfo, ntmoy 177 181 !!------------------------------------------------------------------- 178 182 179 REWIND ( numnam_ice ) ! Read Namelist namicedia 183 ! Read Namelist namicedia 184 REWIND ( numnam_ice ) 180 185 READ ( numnam_ice , namicedia ) 181 186 182 IF(lwp) THEN ! control print187 IF(lwp) THEN 183 188 WRITE(numout,*) 184 189 WRITE(numout,*) 'lim_dia_init_2 : ice parameters for ice diagnostics ' … … 269 274 !! Default option : NO LIM 2.0 sea-ice model 270 275 !!---------------------------------------------------------------------- 276 CONTAINS 277 SUBROUTINE lim_dia_2 ! Empty routine 278 END SUBROUTINE lim_dia_2 271 279 #endif 272 280 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdmp_2.F90
r1855 r1857 4 4 !! Ice model : restoring Ice thickness and Fraction leads 5 5 !!====================================================================== 6 !! History : 2.0 ! 2004-04 (S. Theetten) Original code6 !! History : 2.0 ! 04-04 (S. Theetten) Original code 7 7 !!---------------------------------------------------------------------- 8 8 #if defined key_lim2 && defined key_tradmp … … 11 11 !! 'key_tradmp' Damping 12 12 !!---------------------------------------------------------------------- 13 !! lim_dmp_2 : ice model damping14 13 !!---------------------------------------------------------------------- 15 USE oce ! ocean variables 16 USE dom_oce ! ocean domain 17 USE phycst ! physical constants 18 USE ice_2 ! LIM-2 variables 19 USE tradmp ! traceur damping 20 USE in_out_manager ! I/O manager 21 USE iom ! 14 !! lim_dmp_2 : ice model damping 15 !!---------------------------------------------------------------------- 16 USE in_out_manager ! I/O manager 17 USE phycst ! physical constants 18 USE ice_2 19 USE tradmp 20 USE dom_oce 21 USE oce 22 USE iom 22 23 23 24 IMPLICIT NONE … … 26 27 PUBLIC lim_dmp_2 ! called by ice_step_2 27 28 28 INTEGER :: nice1, nice2 29 INTEGER :: inumice_dmp! logical unit for ice variables (damping)30 REAL(wp), DIMENSION(jpi,jpj) :: hicif_dta 31 REAL(wp), DIMENSION(jpi,jpj) :: frld_dta! fraction lead at a given time32 REAL(wp), DIMENSION(jpi,jpj,2) :: hicif_data 33 REAL(wp), DIMENSION(jpi,jpj,2) :: frld_data! fraction lead data at two consecutive times29 INTEGER :: nice1, nice2, & ! first and second record used 30 & inumice_dmp ! logical unit for ice variables (damping) 31 REAL(wp), DIMENSION(jpi,jpj) :: hicif_dta , & ! ice thickness at a given time 32 & frld_dta ! fraction lead at a given time 33 REAL(wp), DIMENSION(jpi,jpj,2) :: hicif_data , & ! ice thickness data at two consecutive times 34 & frld_data ! fraction lead data at two consecutive times 34 35 35 36 !! * Substitution 36 37 # include "vectopt_loop_substitute.h90" 37 38 !!---------------------------------------------------------------------- 38 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)39 !! LIM 2.0 , UCL-LOCEAN-IPSL (2006) 39 40 !! $Id$ 40 41 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 43 44 CONTAINS 44 45 45 SUBROUTINE lim_dmp_2( kt)46 SUBROUTINE lim_dmp_2(kt) 46 47 !!------------------------------------------------------------------- 47 48 !! *** ROUTINE lim_dmp_2 *** … … 52 53 !! ** method : the key_tradmp must be used to compute resto(:,:) coef. 53 54 !!--------------------------------------------------------------------- 54 INTEGER, INTENT(in) :: kt ! ocean time-step55 ! !56 INTEGER :: ji, jj! dummy loop indices55 INTEGER, INTENT(in) :: kt ! ocean time-step 56 ! 57 INTEGER :: ji, jj ! dummy loop indices 57 58 !!--------------------------------------------------------------------- 58 59 ! 59 60 CALL dta_lim_2( kt ) 60 ! 61 61 62 DO jj = 2, jpjm1 62 63 DO ji = fs_2, fs_jpim1 ! vector opt. 63 64 hicif(ji,jj) = hicif(ji,jj) - rdt_ice * resto(ji,jj,1) * ( hicif(ji,jj) - hicif_dta(ji,jj) ) 64 frld(ji,jj) = frld (ji,jj) - rdt_ice * resto(ji,jj,1) * ( frld (ji,jj)- frld_dta (ji,jj) )65 frld(ji,jj) = frld (ji,jj) - rdt_ice * resto(ji,jj,1) * ( frld(ji,jj) - frld_dta (ji,jj) ) 65 66 END DO 66 67 END DO -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdyn_2.F90
r1855 r1857 17 17 !!---------------------------------------------------------------------- 18 18 USE dom_oce ! ocean space and time domain 19 USE sbc_oce ! surface boundary condition : ocean20 USE phycst ! physical constants21 USE ice_2 ! LIM-2 sea-ice variables22 USE dom_ice_2 ! LIM-2 domain23 USE limistate_2 ! LIM-2 initial state24 USE limrhg_2 ! LIM-2rheology19 USE sbc_oce ! 20 USE phycst ! 21 USE ice_2 ! 22 USE dom_ice_2 ! 23 USE limistate_2 ! 24 USE limrhg_2 ! ice rheology 25 25 26 26 USE lbclnk ! … … 34 34 PUBLIC lim_dyn_2 ! routine called by sbc_ice_lim 35 35 36 !! * Module variables 36 37 REAL(wp) :: rone = 1.e0 ! constant value 37 38 38 39 # include "vectopt_loop_substitute.h90" 39 40 !!---------------------------------------------------------------------- 40 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)41 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 41 42 !! $Id$ 42 43 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 58 59 !! - treatment of the case if no ice dynamic 59 60 !!--------------------------------------------------------------------- 60 INTEGER, INTENT(in) :: kt ! number of iteration61 !! 62 INTEGER :: ji, jj ! dummy loop indices63 INTEGER :: i_j1, i_jpj ! Starting/ending j-indices for rheology64 REAL(wp) :: zcoef ! temporary scalar61 INTEGER, INTENT(in) :: kt ! number of iteration 62 !! 63 INTEGER :: ji, jj ! dummy loop indices 64 INTEGER :: i_j1, i_jpj ! Starting/ending j-indices for rheology 65 REAL(wp) :: zcoef ! temporary scalar 65 66 REAL(wp), DIMENSION(jpj) :: zind ! i-averaged indicator of sea-ice 66 67 REAL(wp), DIMENSION(jpj) :: zmsk ! i-averaged of tmask … … 143 144 ENDIF 144 145 145 IF(ln_ctl) CALL prt_ctl( tab2d_1=u_ice, clinfo1=' lim_dyn : u_ice :', tab2d_2=v_ice, clinfo2=' v_ice :')146 IF(ln_ctl) CALL prt_ctl(tab2d_1=u_ice , clinfo1=' lim_dyn : u_ice :', tab2d_2=v_ice , clinfo2=' v_ice :') 146 147 147 148 ! computation of friction velocity … … 178 179 CALL lbc_lnk( ust2s, 'T', 1. ) ! T-point 179 180 ! 180 IF(ln_ctl) CALL prt_ctl( tab2d_1=ust2s, clinfo1=' lim_dyn : ust2s :')181 IF(ln_ctl) CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn : ust2s :') 181 182 182 183 END SUBROUTINE lim_dyn_2 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limhdf_2.F90
r1855 r1857 4 4 !! LIM 2.0 ice model : horizontal diffusion of sea-ice quantities 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (UCL) Original code7 !! 1.0 ! 2001-05 (G. Madec, R. Hordoir) opa norm8 !! 2.0 ! 2002-08 (C. Ethe) F90, free form9 !!----------------------------------------------------------------------10 6 #if defined key_lim2 11 7 !!---------------------------------------------------------------------- … … 14 10 !! lim_hdf_2 : diffusion trend on sea-ice variable 15 11 !!---------------------------------------------------------------------- 16 USE dom_oce ! ocean domain 17 USE ice_2 ! LIM-2 variables 18 USE prtctl ! Print control 19 USE lbclnk ! 20 USE lib_mpp ! 21 USE in_out_manager ! I/O manager 12 !! * Modules used 13 USE dom_oce 14 USE in_out_manager 15 USE ice_2 16 USE lbclnk 17 USE lib_mpp 18 USE prtctl ! Print control 22 19 23 20 IMPLICIT NONE 24 21 PRIVATE 25 22 26 PUBLIC lim_hdf_2 ! called by lim_tra_2 23 !! * Routine accessibility 24 PUBLIC lim_hdf_2 ! called by lim_tra_2 27 25 28 LOGICAL :: linit = .TRUE. ! indictor of initialisation 26 !! * Module variables 27 LOGICAL :: linit = .TRUE. ! ??? 29 28 REAL(wp) :: epsi04 = 1e-04 ! constant 30 REAL(wp), DIMENSION(jpi,jpj) :: zfact ! metric term29 REAL(wp), DIMENSION(jpi,jpj) :: zfact ! ??? 31 30 32 31 !! * Substitution 33 32 # include "vectopt_loop_substitute.h90" 34 33 !!---------------------------------------------------------------------- 35 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)34 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 36 35 !! $Id$ 37 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)36 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 38 37 !!---------------------------------------------------------------------- 39 38 … … 44 43 !! *** ROUTINE lim_hdf_2 *** 45 44 !! 46 !! ** purpose : Compute and add the diffusive trend on sea-ice variables 45 !! ** purpose : Compute and add the diffusive trend on sea-ice 46 !! variables 47 47 !! 48 48 !! ** method : Second order diffusive operator evaluated using a 49 !! 49 !! Cranck-Nicholson time Scheme. 50 50 !! 51 !! ** Action : update ptab with the diffusive contribution 51 !! ** Action : update ptab with the diffusive contribution 52 !! 53 !! History : 54 !! ! 00-01 (LIM) Original code 55 !! ! 01-05 (G. Madec, R. Hordoir) opa norm 56 !! ! 02-08 (C. Ethe) F90, free form 52 57 !!------------------------------------------------------------------- 53 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj) :: ptab ! Field on which the diffusion is applied 54 !! 58 ! * Arguments 59 REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) :: & 60 ptab ! Field on which the diffusion is applied 61 REAL(wp), DIMENSION(jpi,jpj) :: & 62 ptab0 ! ??? 63 64 ! * Local variables 55 65 INTEGER :: ji, jj ! dummy loop indices 56 INTEGER :: its, iter ! temporary integers 66 INTEGER :: & 67 its, iter ! temporary integers 57 68 CHARACTER (len=55) :: charout 58 REAL(wp) :: zalfa, zrlxint, zconv, zeps ! temporary scalars 59 REAL(wp), DIMENSION(jpi,jpj) :: zrlx, zflu, zflv, zdiv0, zdiv ! temporary workspaces 60 REAL(wp), DIMENSION(jpi,jpj) :: ztab0 ! ??? 69 REAL(wp) :: & 70 zalfa, zrlxint, zconv, zeps ! temporary scalars 71 REAL(wp), DIMENSION(jpi,jpj) :: & 72 zrlx, zflu, zflv, & ! temporary workspaces 73 zdiv0, zdiv ! " " 61 74 !!------------------------------------------------------------------- 62 75 … … 69 82 70 83 ! Arrays initialization 71 ztab0(:, : ) = ptab(:,:) 84 ptab0 (:, : ) = ptab(:,:) 85 !bug zflu (:,jpj) = 0.e0 86 !bug zflv (:,jpj) = 0.e0 72 87 zdiv0(:, 1 ) = 0.e0 73 88 zdiv0(:,jpj) = 0.e0 … … 83 98 DO jj = 2, jpjm1 84 99 DO ji = fs_2 , fs_jpim1 ! vector opt. 85 zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj) + e1v(ji,jj) + e1v(ji,jj-1) ) / ( e1t(ji,jj) * e2t(ji,jj) ) 100 zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj ) + e1v(ji,jj) + e1v(ji,jj-1) ) & 101 & / ( e1t(ji,jj) * e2t(ji,jj) ) 86 102 END DO 87 103 END DO … … 94 110 iter = 0 95 111 96 ! !======================! 97 DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) ) ! Sub-time step loop ! 98 ! !======================! 99 iter = iter + 1 ! incrementation of the sub-time step number 112 ! !=================== 113 DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) ) ! Sub-time step loop 114 ! !=================== 115 ! incrementation of the sub-time step number 116 iter = iter + 1 100 117 101 DO jj = 1, jpjm1 ! diffusive fluxes in U- and V- direction 118 ! diffusive fluxes in U- and V- direction 119 DO jj = 1, jpjm1 102 120 DO ji = 1 , fs_jpim1 ! vector opt. 103 121 zflu(ji,jj) = pahu(ji,jj) * e2u(ji,jj) / e1u(ji,jj) * ( ptab(ji+1,jj) - ptab(ji,jj) ) … … 105 123 END DO 106 124 END DO 107 DO jj= 2, jpjm1 ! diffusive trend : divergence of the fluxes 125 126 ! diffusive trend : divergence of the fluxes 127 DO jj= 2, jpjm1 108 128 DO ji = fs_2 , fs_jpim1 ! vector opt. 109 129 zdiv (ji,jj) = ( zflu(ji,jj) - zflu(ji-1,jj ) & … … 111 131 END DO 112 132 END DO 113 ! ! save the first evaluation of the diffusive trend in zdiv0 133 134 ! save the first evaluation of the diffusive trend in zdiv0 114 135 IF( iter == 1 ) zdiv0(:,:) = zdiv(:,:) 115 136 116 DO jj = 2, jpjm1 ! XXXX iterative evaluation????? 137 ! XXXX iterative evaluation????? 138 DO jj = 2, jpjm1 117 139 DO ji = fs_2 , fs_jpim1 ! vector opt. 118 zrlxint = ( ztab0(ji,jj) &140 zrlxint = ( ptab0(ji,jj) & 119 141 & + rdt_ice * ( zalfa * ( zdiv(ji,jj) + zfact(ji,jj) * ptab(ji,jj) ) & 120 142 & + ( 1.0 - zalfa ) * zdiv0(ji,jj) ) ) & … … 123 145 END DO 124 146 END DO 147 148 ! lateral boundary condition on ptab 125 149 CALL lbc_lnk( zrlx, 'T', 1. ) 126 150 127 zconv = 0.e0 ! convergence test 151 ! convergence test 152 zconv = 0.e0 128 153 DO jj = 2, jpjm1 129 154 DO ji = 2, jpim1 … … 132 157 END DO 133 158 IF( lk_mpp ) CALL mpp_max( zconv ) ! max over the global domain 134 ! 135 ptab(:,:) = zrlx(:,:) ! update value 136 ! !=============================! 137 END DO ! end of sub-time step loop ! 138 ! !=============================! 159 160 ptab(:,:) = zrlx(:,:) 161 162 ! !========================== 163 END DO ! end of sub-time step loop 164 ! !========================== 139 165 140 166 IF(ln_ctl) THEN 141 zrlx(:,:) = ptab(:,:) - ztab0(:,:)167 zrlx(:,:) = ptab(:,:) - ptab0(:,:) 142 168 WRITE(charout,FMT="(' lim_hdf : zconv =',D23.16, ' iter =',I4,2X)") zconv, iter 143 CALL prt_ctl( tab2d_1=zrlx, clinfo1=charout)169 CALL prt_ctl(tab2d_1=zrlx, clinfo1=charout) 144 170 ENDIF 145 ! 171 146 172 END SUBROUTINE lim_hdf_2 147 173 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limistate_2.F90
r1855 r1857 4 4 !! Initialisation of diagnostics ice variables 5 5 !!====================================================================== 6 !! History : 1.0 ! 2001-04 (C. Ethe, G. Madec) Original code7 !! 2.0 ! 2003-08 (G. Madec) add lim_istate_init8 !! - ! 2004-04 (S. Theetten) initialization from a file9 !! - ! 2006-07 (S. Masson) IOM to read the restart10 !! - ! 2007-10 (G. Madec) surface module6 !! History : 1.0 ! 01-04 (C. Ethe, G. Madec) Original code 7 !! 2.0 ! 03-08 (G. Madec) add lim_istate_init 8 !! ! 04-04 (S. Theetten) initialization from a file 9 !! ! 06-07 (S. Masson) IOM to read the restart 10 !! ! 07-10 (G. Madec) surface module 11 11 !!-------------------------------------------------------------------- 12 12 #if defined key_lim2 … … 14 14 !! 'key_lim2' : LIM 2.0 sea-ice model 15 15 !!---------------------------------------------------------------------- 16 !!---------------------------------------------------------------------- 16 17 !! lim_istate_2 : Initialisation of diagnostics ice variables 17 18 !! lim_istate_init_2 : initialization of ice state and namelist read 18 19 !!---------------------------------------------------------------------- 19 USE oce ! ocean variables20 USE ice_2 ! LIM-2 variables21 USE par_ice_2 ! LIM-2 ice parameters22 USE dom_ice_2 ! LIM-2 domain23 USE phycst ! physical constants24 USE eosbn2 ! equation of state25 USE lbclnk !26 USE iom !27 USE in_out_manager !20 USE phycst 21 USE par_ice_2 ! ice parameters 22 USE dom_ice_2 23 USE eosbn2 ! equation of state 24 USE lbclnk 25 USE oce 26 USE ice_2 27 USE iom 28 USE in_out_manager 28 29 29 30 IMPLICIT NONE 30 31 PRIVATE 31 32 32 PUBLIC lim_istate_2! routine called by lim_init_2.F9033 34 !!!! ** init namelist (namiceini) **35 LOGICAL :: ln_limini = .FALSE. ! Ice initialization state33 PUBLIC lim_istate_2 ! routine called by lim_init_2.F90 34 35 !!! ** init namelist (namiceini) ** 36 LOGICAL :: ln_limini = .FALSE. !: Ice initialization state 36 37 REAL(wp) :: ttest = 2.0 ! threshold water temperature for initial sea ice 37 38 REAL(wp) :: hninn = 0.5 ! initial snow thickness in the north … … 44 45 REAL(wp) :: zero = 0.e0 ! constant value = 0 45 46 REAL(wp) :: zone = 1.e0 ! constant value = 1 46 47 !!---------------------------------------------------------------------- 48 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010) 47 !!---------------------------------------------------------------------- 48 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 49 49 !! $Id$ 50 50 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 62 62 !! or from arbitrary sea-ice conditions 63 63 !!-------------------------------------------------------------------- 64 INTEGER :: ji, jj, jk ! dummy loop indices65 REAL(wp) :: zidto ! temporary scalar64 INTEGER :: ji, jj, jk ! dummy loop indices 65 REAL(wp) :: zidto ! temporary scalar 66 66 !-------------------------------------------------------------------- 67 67 … … 69 69 70 70 IF( .NOT. ln_limini ) THEN 71 !71 72 72 tfu(:,:) = tfreez( sn(:,:,1) ) * tmask(:,:,1) ! freezing/melting point of sea water [Celcius] 73 ! 73 74 74 DO jj = 1, jpj 75 75 DO ji = 1, jpi … … 126 126 127 127 !-- lateral boundary conditions 128 CALL lbc_lnk( hicif, 'T', 1. ) ; CALL lbc_lnk( frld , 'T', 1. ) 128 CALL lbc_lnk( hicif, 'T', 1. ) 129 CALL lbc_lnk( frld , 'T', 1. ) 129 130 130 131 ! C A U T I O N frld = 1 over land and lbc_lnk put zero along … … 139 140 CALL lbc_lnk( fsbbq , 'T', 1. ) 140 141 CALL lbc_lnk( qstoif , 'T', 1. ) 141 ! 142 142 143 END SUBROUTINE lim_istate_2 143 144 … … 150 151 !! 151 152 !! ** Method : Read the namiceini namelist and check the parameter 152 !! 153 !! values called at the first timestep (nit000) 153 154 !! 154 155 !! ** input : Namelist namiceini 155 156 !!------------------------------------------------------------------- 156 INTEGER :: ji,jj ! dummy loop indices157 INTEGER :: inum_ice ! temporary integer158 !! 157 INTEGER :: inum_ice 158 INTEGER :: ji,jj 159 159 160 NAMELIST/namiceini/ ln_limini, ttest, hninn, hginn, alinn, & 160 & 161 & hnins, hgins, alins 161 162 !!------------------------------------------------------------------- 162 163 ! … … 164 165 READ ( numnam_ice , namiceini ) 165 166 ! 166 IF(lwp) THEN ! control print167 IF(lwp) THEN 167 168 WRITE(numout,*) 168 169 WRITE(numout,*) 'lim_istate_init_2 : ice parameters inititialisation ' … … 178 179 ENDIF 179 180 180 IF( ln_limini ) THEN ! Ice initialization using input file181 IF( ln_limini ) THEN ! Ice initialization using input file 181 182 ! 182 183 CALL iom_open( 'Ice_initialization.nc', inum_ice ) … … 185 186 IF(lwp) WRITE(numout,*) 186 187 IF(lwp) WRITE(numout,*) ' ice state initialization with : Ice_initialization.nc' 187 !188 188 189 CALL iom_get( inum_ice, jpdom_data, 'hicif', hicif ) 189 190 CALL iom_get( inum_ice, jpdom_data, 'hsnif', hsnif ) … … 191 192 CALL iom_get( inum_ice, jpdom_data, 'ts' , sist ) 192 193 CALL iom_get( inum_ice, jpdom_unknown, 'tbif', tbif(1:nlci,1:nlcj,:), & 193 &kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) )194 & kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) ) 194 195 ! put some values in the extra-halo... 195 196 DO jj = nlcj+1, jpj ; tbif(1:nlci,jj,:) = tbif(1:nlci,nlej,:) ; END DO 196 197 DO ji = nlci+1, jpi ; tbif(ji ,: ,:) = tbif(nlei ,: ,:) ; END DO 197 ! 198 198 199 CALL iom_close( inum_ice) 199 200 ! … … 207 208 !! Default option : Empty module NO LIM 2.0 sea-ice model 208 209 !!---------------------------------------------------------------------- 210 CONTAINS 211 SUBROUTINE lim_istate_2 ! Empty routine 212 END SUBROUTINE lim_istate_2 209 213 #endif 210 214 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limmsh_2.F90
r1855 r1857 4 4 !! LIM 2.0 ice model : definition of the ice mesh parameters 5 5 !!====================================================================== 6 !! ** History : LIM ! 2001-04 (UCL) original code7 !! 2.0 ! 2002-08 (C. Ethe, G. Madec) F90, module8 !!----------------------------------------------------------------------9 6 #if defined key_lim2 10 7 !!---------------------------------------------------------------------- … … 13 10 !! lim_msh_2 : definition of the ice mesh 14 11 !!---------------------------------------------------------------------- 12 !! * Modules used 15 13 USE phycst 16 14 USE dom_oce … … 22 20 PRIVATE 23 21 24 PUBLIC lim_msh_2 ! routine called by ice_ini_2.F90 25 26 !!---------------------------------------------------------------------- 27 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010) 22 !! * Accessibility 23 PUBLIC lim_msh_2 ! routine called by ice_ini_2.F90 24 25 !!---------------------------------------------------------------------- 26 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 28 27 !! $Id$ 29 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)28 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 30 29 !!---------------------------------------------------------------------- 31 30 … … 43 42 !! - Initialization of the ice masks (tmsk, umsk) 44 43 !! 45 !! Reference : Deleersnijder et al. Ocean Modelling 100, 7-10 44 !! ** Refer. : Deleersnijder et al. Ocean Modelling 100, 7-10 45 !! 46 !! ** History : 47 !! original : 01-04 (LIM) 48 !! addition : 02-08 (C. Ethe, G. Madec) 46 49 !!--------------------------------------------------------------------- 47 INTEGER :: ji, jj ! dummy loop indices 48 REAL(wp) :: zh1p, zd1d2p, zusden ! temporary scalars 49 REAL(wp) :: zh2p, zd2d1p, zusden2 ! - - 50 REAL(wp), DIMENSION(jpi,jpj) :: zd2d1 , zd1d2 ! 2D workspace 50 !! * Local variables 51 INTEGER :: ji, jj ! dummy loop indices 52 53 REAL(wp), DIMENSION(jpi,jpj) :: & 54 zd2d1 , zd1d2 ! Derivative of zh2 (resp. zh1) in the x direction 55 ! ! (resp. y direction) (defined at the center) 56 REAL(wp) :: & 57 zh1p , zh2p , & ! Idem zh1, zh2 for the bottom left corner of the grid 58 zd2d1p, zd1d2p , & ! Idem zd2d1, zd1d2 for the bottom left corner of the grid 59 zusden, zusden2 ! temporary scalars 51 60 !!--------------------------------------------------------------------- 52 61 … … 103 112 !------------------- 104 113 !!ibug ??? 105 akappa(:,:,:,:) 106 wght (:,:,:,:)= 0.e0114 akappa(:,:,:,:) = 0.e0 115 wght(:,:,:,:) = 0.e0 107 116 alambd(:,:,:,:,:,:) = 0.e0 108 tmu (:,:)= 0.e0117 tmu(:,:) = 0.e0 109 118 !!i 110 119 … … 229 238 END DO 230 239 END DO 231 CALL lbc_lnk( tmu(:,:), 'I', 1. ) !--lateral boundary conditions 232 233 area(:,:) = e1t(:,:) * e2t(:,:) ! unmasked and masked area of T-grid cell 240 241 !--lateral boundary conditions 242 CALL lbc_lnk( tmu(:,:), 'I', 1. ) 243 244 ! unmasked and masked area of T-grid cell 245 area(:,:) = e1t(:,:) * e2t(:,:) 234 246 235 247 END SUBROUTINE lim_msh_2 … … 239 251 !! Default option Dummy Module NO LIM sea-ice model 240 252 !!---------------------------------------------------------------------- 253 CONTAINS 254 SUBROUTINE lim_msh_2 ! Dummy routine 255 END SUBROUTINE lim_msh_2 241 256 #endif 242 257 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limrhg_2.F90
r1855 r1857 4 4 !! Ice rheology : performs sea ice rheology 5 5 !!====================================================================== 6 !! History : 0.0 ! 1993-12 (M.A. Morales Maqueda.) Original code7 !! 1.0 ! 1994-12 (H. Goosse)8 !! 2.0 ! 1903-12 (C. Ethe, G. Madec) F90, mpp9 !! - ! 2006-08 (G. Madec) surface module, ice-stress at I-point10 !! - ! 2009-09 (G. Madec) Huge verctor optimisation6 !! History : 0.0 ! 93-12 (M.A. Morales Maqueda.) Original code 7 !! 1.0 ! 94-12 (H. Goosse) 8 !! 2.0 ! 03-12 (C. Ethe, G. Madec) F90, mpp 9 !! " " ! 06-08 (G. Madec) surface module, ice-stress at I-point 10 !! " " ! 09-09 (G. Madec) Huge verctor optimisation 11 11 !!---------------------------------------------------------------------- 12 12 #if defined key_lim2 … … 14 14 !! 'key_lim2' LIM 2.0 sea-ice model 15 15 !!---------------------------------------------------------------------- 16 !! lim_rhg_2 : computes ice velocities 16 !!---------------------------------------------------------------------- 17 !! lim_rhg_2 : computes ice velocities 17 18 !!---------------------------------------------------------------------- 18 19 USE par_oce ! ocean parameter … … 33 34 PUBLIC lim_rhg_2 ! routine called by lim_dyn 34 35 35 REAL(wp) :: rzero = 0.e0 ! constant value: zero36 REAL(wp) :: rone = 1.e0 ! and one36 REAL(wp) :: rzero = 0.e0 ! constant value: zero 37 REAL(wp) :: rone = 1.e0 ! and one 37 38 38 39 !! * Substitutions 39 40 # include "vectopt_loop_substitute.h90" 40 41 !!---------------------------------------------------------------------- 41 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)42 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 42 43 !! $Id$ 43 44 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 55 56 !! viscous-plastic law including shear strength and a bulk rheology. 56 57 !! 57 !! ** Action : - compute u_ice, v_ice the sea-ice velocity defined at I-point 58 !! ** Action : - compute u_ice, v_ice the sea-ice velocity defined 59 !! at I-point 58 60 !!------------------------------------------------------------------- 59 61 INTEGER, INTENT(in) :: k_j1 ! southern j-index for ice computation … … 580 582 !! Default option Dummy module NO 2.0 LIM sea-ice model 581 583 !!---------------------------------------------------------------------- 584 CONTAINS 585 SUBROUTINE lim_rhg_2( k1 , k2 ) ! Dummy routine 586 WRITE(*,*) 'lim_rhg_2: You should not have seen this print! error?', k1, k2 587 END SUBROUTINE lim_rhg_2 582 588 #endif 583 589 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limrst_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limrst_2 *** 4 !! LIM-2 ice model : open/read/write the restart file4 !! Ice restart : write the ice restart file 5 5 !!====================================================================== 6 !! History : 2.0 ! 2001-04 (C. Ethe, G. Madec) Original code7 !! ! 2006-07 (S. Masson) use IOM for restart read/write6 !! History : 2.0 ! 01-04 (C. Ethe, G. Madec) Original code 7 !! ! 06-07 (S. Masson) use IOM for restart read/write 8 8 !!---------------------------------------------------------------------- 9 9 #if defined key_lim2 10 10 !!---------------------------------------------------------------------- 11 11 !! 'key_lim2' : LIM 2.0 sea-ice model 12 !!---------------------------------------------------------------------- 12 13 !!---------------------------------------------------------------------- 13 14 !! lim_rst_opn_2 : open ice restart file … … 19 20 USE sbc_oce 20 21 USE sbc_ice 22 21 23 USE in_out_manager 22 24 USE iom … … 33 35 34 36 !!---------------------------------------------------------------------- 35 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)37 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 36 38 !! $Id$ 37 39 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 82 84 END SUBROUTINE lim_rst_opn_2 83 85 84 85 86 SUBROUTINE lim_rst_write_2( kt ) 86 87 !!---------------------------------------------------------------------- … … 90 91 !!---------------------------------------------------------------------- 91 92 INTEGER, INTENT(in) :: kt ! number of iteration 92 ! !93 ! 93 94 INTEGER :: iter ! kt + nn_fsbc -1 94 95 !!---------------------------------------------------------------------- … … 180 181 ENDIF 181 182 182 IF ( jprstlib == jprstdimg ) THEN183 IF ( jprstlib == jprstdimg ) THEN 183 184 ! eventually read netcdf file (monobloc) for restarting on different number of processors 184 185 ! if {cn_icerst_in}.nc exists, then set jlibalt to jpnf90 … … 187 188 ENDIF 188 189 189 CALL iom_open ( cn_icerst_in, numrir, kiolib = jlibalt )190 CALL iom_open ( cn_icerst_in, numrir, kiolib = jlibalt ) 190 191 191 192 CALL iom_get( numrir, 'kt_ice' , ziter ) -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limsbc_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limsbc_2 *** 4 !! LIM 2 ice model : heat, salt, mass and momentum fluxesat the sea ice/ocean interface4 !! computation of the flux at the sea ice/ocean interface 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (H. Goosse) Original code7 !! 2.0 ! 2002-07 (C. Ethe, G. Madec) re-writing F908 !! - ! 2006-07 (G. Madec) surface module6 !! History : 00-01 (H. Goosse) Original code 7 !! 02-07 (C. Ethe, G. Madec) re-writing F90 8 !! 06-07 (G. Madec) surface module 9 9 !!---------------------------------------------------------------------- 10 10 #if defined key_lim2 11 11 !!---------------------------------------------------------------------- 12 12 !! 'key_lim2' LIM 2.0 sea-ice model 13 !!---------------------------------------------------------------------- 13 14 !!---------------------------------------------------------------------- 14 15 !! lim_sbc_2 : flux at the ice / ocean interface … … 19 20 USE sbc_oce ! surface boundary condition 20 21 USE phycst ! physical constants 21 USE ice_2 ! LIM 2sea-ice variables22 USE ice_2 ! LIM sea-ice variables 22 23 23 24 USE lbclnk ! ocean lateral boundary condition … … 32 33 PRIVATE 33 34 34 PUBLIC lim_sbc_2! called by sbc_ice_lim_235 PUBLIC lim_sbc_2 ! called by sbc_ice_lim_2 35 36 36 37 REAL(wp) :: epsi16 = 1.e-16 ! constant values … … 43 44 # include "vectopt_loop_substitute.h90" 44 45 !!---------------------------------------------------------------------- 45 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)46 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 46 47 !! $Id$ 47 48 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 108 109 sice_r(:,:) = sice 109 110 ! 110 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA_R2 configuration 111 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN 112 ! ! ======================= 113 ! ! ORCA_R2 configuration 114 ! ! ======================= 111 115 ii0 = 145 ; ii1 = 180 ! Baltic Sea 112 116 ij0 = 113 ; ij1 = 130 ; soce_r(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 … … 190 194 END DO 191 195 192 CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) 193 CALL iom_put( 'qns_io_cea' 194 CALL iom_put( 'qsr_io_cea' , fstric(:,:) * (1. - pfrld(:,:)))196 CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) ) 197 CALL iom_put( 'qns_io_cea', qns(:,:) - zqnsoce(:,:) * pfrld(:,:) ) 198 CALL iom_put( 'qsr_io_cea', fstric(:,:) * (1. - pfrld(:,:)) ) 195 199 196 200 !------------------------------------------! … … 298 302 !-----------------------------------------------! 299 303 300 IF( lk_cpl ) THEN 304 IF ( lk_cpl ) THEN 305 ! Ice surface temperature 301 306 tn_ice(:,:,1) = sist(:,:) ! sea-ice surface temperature 302 ! ! snow/ice and ocean albedos307 ! Computation of snow/ice and ocean albedo 303 308 CALL albedo_ice( tn_ice, reshape( hicif, (/jpi,jpj,1/) ), reshape( hsnif, (/jpi,jpj,1/) ), zalbp, zalb ) 304 309 alb_ice(:,:,1) = 0.5 * ( zalbp(:,:,1) + zalb (:,:,1) ) ! Ice albedo (mean clear and overcast skys) … … 307 312 308 313 IF(ln_ctl) THEN 309 CALL prt_ctl( tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ')310 CALL prt_ctl( tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=emps , clinfo2=' emps : ')311 CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau : ', mask1=umask, &312 & tab2d_2=vtau, clinfo2=' vtau : ' , mask2=vmask )313 CALL prt_ctl( tab2d_1=fr_i, clinfo1=' lim_sbc: fr_i : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice : ')314 CALL prt_ctl(tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ') 315 CALL prt_ctl(tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=emps , clinfo2=' emps : ') 316 CALL prt_ctl(tab2d_1=utau , clinfo1=' lim_sbc: utau : ', mask1=umask, & 317 & tab2d_2=vtau , clinfo2=' vtau : ' , mask2=vmask ) 318 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice : ') 314 319 ENDIF 315 ! 320 321 END SUBROUTINE lim_sbc_2 322 323 #else 324 !!---------------------------------------------------------------------- 325 !! Default option : Dummy module NO LIM 2.0 sea-ice model 326 !!---------------------------------------------------------------------- 327 CONTAINS 328 SUBROUTINE lim_sbc_2 ! Dummy routine 316 329 END SUBROUTINE lim_sbc_2 317 318 #else319 !!----------------------------------------------------------------------320 !! Default option : Dummy module NO LIM 2.0 sea-ice model321 !!----------------------------------------------------------------------322 330 #endif 323 331 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limtab_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limtab_2 *** 4 !! LIM 2 ice model : transform 1D (2D) array to a 2D (1D) array4 !! transform 1D (2D) array to a 2D (1D) table 5 5 !!====================================================================== 6 6 #if defined key_lim2 7 7 !!---------------------------------------------------------------------- 8 !! 'key_lim2' : LIM 2.0 sea-ice model 8 !! tab_2d_1d : 2-D to 1-D 9 !! tab_1d_2d : 1-D to 2-D 9 10 !!---------------------------------------------------------------------- 10 !! tab_2d_1d_2 : 2-D to 1-D 11 !! tab_1d_2d_2 : 1-D to 2-D 12 !!---------------------------------------------------------------------- 11 !! * Modules used 13 12 USE par_kind 14 13 … … 16 15 PRIVATE 17 16 18 PUBLIC tab_2d_1d_2 ! called by lim_thd_2 19 PUBLIC tab_1d_2d_2 ! called by lim_thd_2 17 !! * Routine accessibility 18 PUBLIC tab_2d_1d_2 ! called by lim_ther 19 PUBLIC tab_1d_2d_2 ! called by lim_ther 20 20 21 21 !!---------------------------------------------------------------------- 22 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)22 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 23 23 !! $Id$ 24 24 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 26 26 CONTAINS 27 27 28 SUBROUTINE tab_2d_1d_2( kdim1d, ptab1d, ptab2d, kdim2d_i, kdim2d_j, ptab_ind ) 29 !!------------------------------------------------------------------- 30 INTEGER , INTENT(in ) :: kdim1d 31 INTEGER , INTENT(in ) :: kdim2d_i, kdim2d_j 32 REAL(wp), INTENT(in ), DIMENSION(kdim2d_i, kdim2d_j) :: ptab2d 33 INTEGER , INTENT(in ), DIMENSION(kdim1d) :: ptab_ind 34 REAL(wp), INTENT( out), DIMENSION(kdim1d) :: ptab1d 35 !! 36 INTEGER :: jn , jid, jjd ! dummy loop indices 37 !!------------------------------------------------------------------- 38 ! 39 DO jn = 1, kdim1d 40 jid = MOD( ptab_ind(jn) - 1, kdim2d_i ) + 1 41 jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 42 ptab1d(jn) = ptab2d(jid,jjd) 28 SUBROUTINE tab_2d_1d_2 ( ndim1d, tab1d, tab2d, ndim2d_x, ndim2d_y, tab_ind ) 29 30 INTEGER, INTENT(in) :: & 31 ndim1d, ndim2d_x, ndim2d_y 32 33 REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT(in) :: & 34 tab2d 35 36 INTEGER, DIMENSION ( ndim1d), INTENT ( in) :: & 37 tab_ind 38 39 REAL(wp), DIMENSION(ndim1d), INTENT ( out) :: & 40 tab1d 41 42 INTEGER :: & 43 jn , jid, jjd 44 45 DO jn = 1, ndim1d 46 jid = MOD( tab_ind(jn) - 1, ndim2d_x ) + 1 47 jjd = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 48 tab1d( jn) = tab2d( jid, jjd) 43 49 END DO 44 ! 50 45 51 END SUBROUTINE tab_2d_1d_2 46 52 47 53 48 SUBROUTINE tab_1d_2d_2( kdim1d, ptab2d, ptab_ind, ptab1d, kdim2d_i, kdim2d_j ) 49 !!------------------------------------------------------------------- 50 INTEGER , INTENT(in ) :: kdim1d 51 INTEGER , INTENT(in ) :: kdim2d_i, kdim2d_j 52 INTEGER , INTENT(in ), DIMENSION(kdim1d) :: ptab_ind 53 REAL(wp), INTENT(in ), DIMENSION(kdim1d) :: ptab1d 54 REAL(wp), INTENT( out), DIMENSION(kdim2d_i, kdim2d_j) :: ptab2d 55 !! 56 INTEGER :: jn, jid, jjd ! dummy loop indices 57 !!------------------------------------------------------------------- 58 ! 59 DO jn = 1, kdim1d 60 jid = MOD( ptab_ind(jn) - 1, kdim2d_i) + 1 61 jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 62 ptab2d(jid,jjd) = ptab1d(jn) 54 SUBROUTINE tab_1d_2d_2 ( ndim1d, tab2d, tab_ind, tab1d, ndim2d_x, ndim2d_y ) 55 56 INTEGER, INTENT ( in) :: & 57 ndim1d, ndim2d_x, ndim2d_y 58 59 INTEGER, DIMENSION (ndim1d) , INTENT (in) :: & 60 tab_ind 61 62 REAL(wp), DIMENSION(ndim1d), INTENT (in) :: & 63 tab1d 64 65 REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT ( out) :: & 66 tab2d 67 68 INTEGER :: & 69 jn, jid, jjd 70 71 DO jn = 1, ndim1d 72 jid = MOD( tab_ind(jn) - 1, ndim2d_x) + 1 73 jjd = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 74 tab2d(jid, jjd) = tab1d( jn) 63 75 END DO 64 ! 76 65 77 END SUBROUTINE tab_1d_2d_2 66 78 67 #else68 !!----------------------------------------------------------------------69 !! Default option Dummy module NO LIM 2.0 sea-ice model70 !!----------------------------------------------------------------------71 79 #endif 72 73 !!======================================================================74 80 END MODULE limtab_2 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_2.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limthd_2 *** 4 !! LIM 2 ice model : thermodynamics4 !! LIM thermo ice model : ice thermodynamic 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (UCL) Original code6 !! History : 1.0 ! 2000-01 (LIM) 7 7 !! 2.0 ! 2002-07 (C. Ethe, G. Madec) F90 8 !! -! 2003-08 (C. Ethe) add lim_thd_init9 !! - ! 2008- 08 (A. Caubel, G. Madec, E. Maisonnave, S. Masson ) generic coupled interface8 !! 2.0 ! 2003-08 (C. Ethe) add lim_thd_init 9 !! - ! 2008-2008 (A. Caubel, G. Madec, E. Maisonnave, S. Masson ) generic coupled interface 10 10 !!--------------------------------------------------------------------- 11 11 #if defined key_lim2 … … 13 13 !! 'key_lim2' : LIM 2.0 sea-ice model 14 14 !!---------------------------------------------------------------------- 15 !! lim_thd_2 : thermodynamic sof sea ice15 !! lim_thd_2 : thermodynamic of sea ice 16 16 !! lim_thd_init_2 : initialisation of sea-ice thermodynamic 17 17 !!---------------------------------------------------------------------- … … 23 23 USE lib_mpp 24 24 USE iom ! IOM library 25 USE ice_2 ! LIM 2sea-ice variables25 USE ice_2 ! LIM sea-ice variables 26 26 USE sbc_oce ! 27 27 USE sbc_ice ! 28 USE dom_ice_2 ! LIM 2 sea-ice domain 28 USE thd_ice_2 ! LIM thermodynamic sea-ice variables 29 USE dom_ice_2 ! LIM sea-ice domain 29 30 USE limthd_zdf_2 30 31 USE limthd_lac_2 … … 48 49 # include "domzgr_substitute.h90" 49 50 # include "vectopt_loop_substitute.h90" 50 !!-------- --------------------------------------------------------------51 !! NEMO/LIM 3. 3, UCL-LOCEAN-IPSL (2010)51 !!-------- ------------------------------------------------------------- 52 !! NEMO/LIM 3.2, UCL-LOCEAN-IPSL (2009) 52 53 !! $Id$ 53 54 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 565 566 !! Default option Dummy module NO LIM 2.0 sea-ice model 566 567 !!---------------------------------------------------------------------- 568 CONTAINS 569 SUBROUTINE lim_thd_2 ! Dummy routine 570 END SUBROUTINE lim_thd_2 567 571 #endif 568 572 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_lac_2.F90
r1855 r1857 1 1 MODULE limthd_lac_2 2 #if defined key_lim2 2 3 !!====================================================================== 3 4 !! *** MODULE limthd_lac_2 *** 4 !! LIM 2 ice model : thermodynamics -- lateral accretion of the ice 5 !!====================================================================== 6 !! History : LIM ! 2001-04 (UCL) original code 7 !! 2.0 ! 2002-08 (C. Ethe, G. Madec) F90, mpp 5 !! lateral thermodynamic growth of the ice 6 !!====================================================================== 7 8 8 !!---------------------------------------------------------------------- 9 #if defined key_lim2 10 !!---------------------------------------------------------------------- 11 !! 'key_lim2' : LIM 2.0 sea-ice model 12 !!---------------------------------------------------------------------- 13 !! lim_lat_acr_2 : lateral accretion of ice 14 !!---------------------------------------------------------------------- 9 !! lim_lat_acr_2 : lateral accretion of ice 10 !! * Modules used 15 11 USE par_oce ! ocean parameters 16 USE phycst ! physical constants 17 USE ice_2 ! LIM 2 sea-ice variables 18 USE limistate_2 ! LIM 2 initial state 12 USE phycst 13 USE thd_ice_2 14 USE ice_2 15 USE limistate_2 19 16 20 17 IMPLICIT NONE 21 18 PRIVATE 22 19 23 PUBLIC lim_thd_lac_2 ! called by lim_thd_2 24 25 REAL(wp) :: epsi20 = 1.e-20 ! constant values 26 REAL(wp) :: epsi13 = 1.e-13 ! 27 REAL(wp) :: rzero = 0.e0 ! 28 REAL(wp) :: rone = 1.e0 ! 29 20 !! * Routine accessibility 21 PUBLIC lim_thd_lac_2 ! called by lim_thd_2 22 23 !! * Module variables 24 REAL(wp) :: & ! constant values 25 epsi20 = 1.e-20 , & 26 epsi13 = 1.e-13 , & 27 zzero = 0.e0 , & 28 zone = 1.e0 30 29 !!---------------------------------------------------------------------- 31 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)30 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 32 31 !! $Id$ 33 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)32 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 34 33 !!---------------------------------------------------------------------- 35 34 CONTAINS … … 61 60 !! update h_snow_1d, h_ice_1d and tbif_1d(:,:) 62 61 !! 63 !! References : M. Maqueda, 1995, PhD Thesis, Univesidad Complutense Madrid 64 !! Fichefet T. and M. Maqueda 1997, J. Geo. Res., 102(C6), 12609 -12646 62 !! ** References : 63 !! M. Maqueda, 1995, PhD Thesis, Univesidad Complutense Madrid 64 !! Fichefet T. and M. Maqueda 1997, J. Geo. Res., 102(C6), 65 !! 12609 -12646 66 !! History : 67 !! 1.0 ! 01-04 (LIM) original code 68 !! 2.0 ! 02-08 (C. Ethe, G. Madec) F90, mpp 65 69 !!------------------------------------------------------------------- 66 INTEGER, INTENT(in) :: kideb ! start point on which the the computation is applied 67 INTEGER, INTENT(in) :: kiut ! end point on which the the computation is applied 68 !! 69 INTEGER :: ji ! dummy loop indices 70 INTEGER :: iicefr ! 1 = existing ice ; 0 = no ice 71 INTEGER :: iiceform ! 1 = ice formed ; 0 = no ice formed 72 INTEGER :: ihemis ! hemisphere indicator 73 REAL(wp), DIMENSION(jpij) :: zqbgow ! heat budget of the open water (negative) 74 REAL(wp), DIMENSION(jpij) :: zfrl_old ! previous sea/ice fraction 75 REAL(wp), DIMENSION(jpij) :: zhice_old ! previous ice thickness 76 REAL(wp), DIMENSION(jpij) :: zhice0 ! thickness of newly formed ice in leads 77 REAL(wp), DIMENSION(jpij) :: zfrlmin ! minimum fraction for leads 78 REAL(wp), DIMENSION(jpij) :: zdhicbot ! part of thickness of newly formed ice in leads which 79 ! ! has been already used in transport for example 80 REAL(wp) :: zhemis ! hemisphere (0 = North, 1 = South) 81 REAL(wp) :: zhicenew ! new ice thickness 82 REAL(wp) :: zholds2 ! ratio of previous ice thickness and 2 83 REAL(wp) :: zhnews2 ! ratio of new ice thickness and 2 84 REAL(wp) :: zfrlnew ! new sea/ice fraction 85 REAL(wp) :: zfrld ! ratio of sea/ice fraction and minimum fraction for leads 86 REAL(wp) :: zfrrate ! leads-closure rate 87 REAL(wp) :: zdfrl ! sea-ice fraction increment 88 REAL(wp) :: zdh1 , zdh2, zdh3, zdh4, zdh5 ! tempory scalars 89 REAL(wp) :: ztint, zta1, zta2, zta3, zta4 ! - - 90 REAL(wp) :: zah, zalpha, zbeta ! - - 70 !! * Arguments 71 INTEGER , INTENT(IN):: & 72 kideb , & ! start point on which the the computation is applied 73 kiut ! end point on which the the computation is applied 74 75 !! * Local variables 76 INTEGER :: & 77 ji , & ! dummy loop indices 78 iicefr , & ! 1 = existing ice ; 0 = no ice 79 iiceform , & ! 1 = ice formed ; 0 = no ice formed 80 ihemis ! dummy indice 81 REAL(wp), DIMENSION(jpij) :: & 82 zqbgow , & ! heat budget of the open water (negative) 83 zfrl_old , & ! previous sea/ice fraction 84 zhice_old , & ! previous ice thickness 85 zhice0 , & ! thickness of newly formed ice in leads 86 zfrlmin , & ! minimum fraction for leads 87 zdhicbot ! part of thickness of newly formed ice in leads which 88 ! has been already used in transport for example 89 REAL(wp) :: & 90 zhemis , & ! hemisphere (0 = North, 1 = South) 91 zhicenew , & ! new ice thickness 92 zholds2 , & ! ratio of previous ice thickness and 2 93 zhnews2 , & ! ratio of new ice thickness and 2 94 zfrlnew , & ! new sea/ice fraction 95 zfrld , & ! ratio of sea/ice fraction and minimum fraction for leads 96 zfrrate , & ! leads-closure rate 97 zdfrl ! sea-ice fraction increment 98 REAL(wp) :: & 99 zdh1 , zdh2 , zdh3 , zdh4, zdh5 , & ! tempory scalars 100 ztint , zta1 , zta2 , zta3 , zta4 , & 101 zah, zalpha , zbeta 91 102 !!--------------------------------------------------------------------- 92 103 93 104 !-------------------------------------------------------------- 94 105 ! Computation of the heat budget of the open water (negative) 95 !-------------------------------------------------------------- 106 !-------------------------------------------------------------- 107 96 108 DO ji = kideb , kiut 97 109 zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji) … … 119 131 ! Morales Maqueda, 1995 - Fichefet and Morales Maqueda, 1997 120 132 !--------------------------------------------------------------------- 133 121 134 !CDIR NOVERRCHK 122 135 DO ji = kideb , kiut … … 132 145 !--computation of the remaining part of ice thickness which has been already used 133 146 zdhicbot(ji) = ( frld_1d(ji) - zfrlnew ) * zhice0(ji) / ( 1.0 - zfrlmin(ji) ) & 134 &- ( ( 1.0 - zfrrate ) / ( 1.0 - frld_1d(ji) ) ) * ( zqbgow(ji) / xlic )147 - ( ( 1.0 - zfrrate ) / ( 1.0 - frld_1d(ji) ) ) * ( zqbgow(ji) / xlic ) 135 148 END DO 136 149 … … 142 155 ! (1-Anew) * hsnew = (1-Aold) * hsold 143 156 !-------------------------------------------------------------------------------------------- 157 144 158 DO ji = kideb , kiut 145 159 iicefr = 1 - MAX( 0, INT( SIGN( 1.5 * zone , zfrl_old(ji) - 1.0 + epsi13 ) ) ) … … 155 169 ! Ajustement of ice internal temperatures 156 170 !------------------------------------------------------- 171 157 172 DO ji = kideb , kiut 158 173 iicefr = 1 - MAX( 0, INT( SIGN( 1.5 * zone , zfrl_old(ji) - 1.0 + epsi13 ) ) ) … … 198 213 ! dV = Vnew - Vold 199 214 !------------------------------------------------------------- 215 200 216 DO ji = kideb , kiut 201 217 dvlbq_1d (ji) = ( 1. - frld_1d(ji) ) * h_ice_1d(ji) - ( 1. - zfrl_old(ji) ) * zhice_old(ji) 202 218 rdmicif_1d(ji) = rdmicif_1d(ji) + rhoic * dvlbq_1d(ji) 203 219 END DO 204 !220 205 221 END SUBROUTINE lim_thd_lac_2 206 207 222 #else 208 !!---------------------------------------------------------------------- 209 !! Default option Dummy module NO LIM 2.0 sea-ice model 210 !!---------------------------------------------------------------------- 223 !!====================================================================== 224 !! *** MODULE limthd_lac_2 *** 225 !! no sea ice model 226 !!====================================================================== 227 CONTAINS 228 SUBROUTINE lim_thd_lac_2 ! Empty routine 229 END SUBROUTINE lim_thd_lac_2 211 230 #endif 212 213 !!======================================================================214 231 END MODULE limthd_lac_2 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_zdf_2.F90
r1855 r1857 4 4 !! thermodynamic growth and decay of the ice 5 5 !!====================================================================== 6 !! History : 1.0 ! 2001-04 (LIM) Original code7 !! 2.0 ! 2002-08 (C. Ethe, G. Madec) F906 !! History : 1.0 ! 01-04 (LIM) Original code 7 !! 2.0 ! 02-08 (C. Ethe, G. Madec) F90 8 8 !!---------------------------------------------------------------------- 9 9 #if defined key_lim2 … … 11 11 !! 'key_lim2' LIM 2.0 sea-ice model 12 12 !!---------------------------------------------------------------------- 13 !!---------------------------------------------------------------------- 13 14 !! lim_thd_zdf_2 : vertical accr./abl. and lateral ablation of sea ice 14 15 !!---------------------------------------------------------------------- 16 !! * Modules used 15 17 USE par_oce ! ocean parameters 16 USE phycst ! physical constants17 USE ice_2 ! LIM 2 sea-ice variables18 USE limistate_2 ! LIM 2 initial state19 USE in_out_manager ! I/O manager20 !!18 USE phycst ! ??? 19 USE thd_ice_2 20 USE ice_2 21 USE limistate_2 22 USE in_out_manager 21 23 USE cpl_oasis3, ONLY : lk_cpl 22 24 … … 24 26 PRIVATE 25 27 26 PUBLIC lim_thd_zdf_2 ! called by lim_thd_2 27 28 REAL(wp) :: epsi20 = 1.e-20 ! constant values 29 REAL(wp) :: epsi13 = 1.e-13 ! 30 REAL(wp) :: rzero = 0.e0 ! 31 REAL(wp) :: rone = 1.e0 ! 32 33 !!---------------------------------------------------------------------- 34 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010) 28 PUBLIC lim_thd_zdf_2 ! called by lim_thd_2 29 30 REAL(wp) :: epsi20 = 1.e-20 , & ! constant values 31 & epsi13 = 1.e-13 , & 32 & zzero = 0.e0 , & 33 & zone = 1.e0 34 !!---------------------------------------------------------------------- 35 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 35 36 !! $Id$ 36 37 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 72 73 !! 73 74 INTEGER :: ji ! dummy loop indices 74 REAL(wp), DIMENSION(jpij,2) :: zqcmlt ! energy due to surface( /1 ) and bottom melting( /2 ) 75 REAL(wp), DIMENSION(jpij) :: ztsmlt ! snow/ice surface melting temperature 76 REAL(wp), DIMENSION(jpij) :: ztbif ! int. temp. at the mid-point of the 1st layer of the snow/ice sys. 77 REAL(wp), DIMENSION(jpij) :: zksn ! effective conductivity of snow 78 REAL(wp), DIMENSION(jpij) :: zkic ! effective conductivity of ice 79 REAL(wp), DIMENSION(jpij) :: zksndh ! thermal cond. at the mid-point of the 1st layer of the snow/ice sys. 80 REAL(wp), DIMENSION(jpij) :: zfcsu ! conductive heat flux at the surface of the snow/ice system 81 REAL(wp), DIMENSION(jpij) :: zfcsudt ! = zfcsu * dt 82 REAL(wp), DIMENSION(jpij) :: zi0 ! frac. of the net SW rad. which is not absorbed at the surface 83 REAL(wp), DIMENSION(jpij) :: z1mi0 ! fraction of the net SW radiation absorbed at the surface 84 REAL(wp), DIMENSION(jpij) :: zqmax ! maximum energy stored in brine pockets 85 REAL(wp), DIMENSION(jpij) :: zrcpdt ! h_su*rho_su*cp_su/dt(h_su being the thick. of surf. layer) 86 REAL(wp), DIMENSION(jpij) :: zts_old ! previous surface temperature 87 REAL(wp), DIMENSION(jpij) :: zidsn , z1midsn , zidsnic ! tempory variables 88 !! 89 REAL(wp), DIMENSION(jpij) :: zfnet ! net heat flux at the top surface( incl. conductive heat flux) 90 REAL(wp), DIMENSION(jpij) :: zsprecip ! snow accumulation 91 REAL(wp), DIMENSION(jpij) :: zhsnw_old ! previous snow thickness 92 REAL(wp), DIMENSION(jpij) :: zdhictop ! change in ice thickness due to top surf ablation/accretion 93 REAL(wp), DIMENSION(jpij) :: zdhicbot ! change in ice thickness due to bottom surf abl/acc 94 REAL(wp), DIMENSION(jpij) :: zqsup ! energy transmitted to ocean (coming from top surf abl/acc) 95 REAL(wp), DIMENSION(jpij) :: zqocea ! energy transmitted to ocean (coming from bottom sur abl/acc) 96 REAL(wp), DIMENSION(jpij) :: zfrl_old ! previous sea/ice fraction 97 REAL(wp), DIMENSION(jpij) :: zfrld_1d ! new sea/ice fraction 98 REAL(wp), DIMENSION(jpij) :: zep ! internal temperature of the 2nd layer of the snow/ice system 99 !! 100 REAL(wp), DIMENSION(3) :: zplediag ! principle diagonal, subdiag. and supdiag. of the 101 REAL(wp), DIMENSION(3) :: zsubdiag ! tri-diagonal matrix coming from the computation 102 REAL(wp), DIMENSION(3) :: zsupdiag ! of the temperatures inside the snow-ice system 103 REAL(wp), DIMENSION(3) :: zsmbr ! second member 104 REAL(wp) :: zhsu ! thickness of surface layer 105 REAL(wp) :: zhe ! effective thickness for compu. of equ. thermal conductivity 106 REAL(wp) :: zheshth ! = zhe / thth 107 REAL(wp) :: zghe ! correction factor of the thermal conductivity 108 REAL(wp) :: zumsb ! parameter for numerical method to solve heat-diffusion eq. 109 REAL(wp) :: zkhsn ! conductivity at the snow layer 110 REAL(wp) :: zkhic ! conductivity at the ice layers 111 REAL(wp) :: zkint ! equivalent conductivity at the snow-ice interface 112 REAL(wp) :: zkhsnint ! = zkint*dt / (hsn*rhosn*cpsn) 113 REAL(wp) :: zkhicint ! = 2*zkint*dt / (hic*rhoic*cpic) 114 REAL(wp) :: zpiv1 , zpiv2 ! tempory scalars used to solve the tri-diagonal system 115 REAL(wp) :: zb2, zd2, zb3, zd3 116 REAL(wp) :: ztint ! equivalent temperature at the snow-ice interface 117 REAL(wp) :: zexp ! exponential function of the ice thickness 118 REAL(wp) :: zfsab ! part of solar radiation stored in brine pockets 119 REAL(wp) :: zfts ! value of energy balance function when the temp. equal surf. temp. 120 REAL(wp) :: zdfts ! value of derivative of ztfs when the temp. equal surf. temp. 121 REAL(wp) :: zdts ! surface temperature increment 122 REAL(wp) :: zqsnw_mlt ! energy needed to melt snow 123 REAL(wp) :: zdhsmlt ! change in snow thickness due to melt 124 REAL(wp) :: zhsn ! snow thickness (previous+accumulation-melt) 125 REAL(wp) :: zqsn_mlt_rem ! remaining heat coming from snow melting 126 REAL(wp) :: zqice_top_mlt ! energy used to melt ice at top surface 127 REAL(wp) :: zdhssub ! change in snow thick. due to sublimation or evaporation 128 REAL(wp) :: zdhisub ! change in ice thick. due to sublimation or evaporation 129 REAL(wp) :: zdhsn ! snow ice thickness increment 130 REAL(wp) :: zdtsn ! snow internal temp. increment 131 REAL(wp) :: zdtic ! ice internal temp. increment 132 REAL(wp) :: zqnes ! conductive energy due to ice melting in the first ice layer 133 !! 134 REAL(wp) :: ztbot ! temperature at the bottom surface 135 REAL(wp) :: zfcbot ! conductive heat flux at bottom surface 136 REAL(wp) :: zqice_bot ! energy used for bottom melting/growing 137 REAL(wp) :: zqice_bot_mlt ! energy used for bottom melting 138 REAL(wp) :: zqstbif_bot ! part of energy stored in brine pockets used for bottom melting 139 REAL(wp) :: zqstbif_old ! tempory var. for zqstbif_bot 140 REAL(wp) :: zdhicmlt ! change in ice thickness due to bottom melting 141 REAL(wp) :: zdhicm ! change in ice thickness var. 142 REAL(wp) :: zdhsnm ! change in snow thickness var. 143 REAL(wp) :: zhsnfi ! snow thickness var. 144 REAL(wp) :: zc1, zpc1, zc2, zpc2, zp1, zp2 ! temporary scalars 145 REAL(wp) :: ztb2, ztb3 146 !! 147 REAL(wp) :: zdrmh ! change in snow/ice thick. after snow-ice formation 148 REAL(wp) :: zhicnew ! new ice thickness 149 REAL(wp) :: zhsnnew ! new snow thickness 150 REAL(wp) :: zquot , ztneq ! tempory temp. variables 151 REAL(wp) :: zqice, zqicetot & ! total heat inside the snow/ice system 152 REAL(wp) :: zdfrl ! change in ice concentration 153 REAL(wp) :: zdvsnvol ! change in snow volume 154 REAL(wp) :: zdrfrl1, zdrfrl2 ! tempory scalars 155 REAL(wp) :: zihsn, zidhb, zihic, zihe, zihq, ziexp, ziqf, zihnf, zibmlt, ziqr, zihgnew, zind 156 !!---------------------------------------------------------------------- 157 158 !----------------------------------------------------------------------- 159 ! 1. Boundaries conditions for snow/ice system internal temperature 160 ! - If tbif_1d(ji,1) > rt0_snow, tbif_1d(ji,1) = rt0_snow 161 ! - If tbif_1d(ji,2/3) > rt0_ice, tbif_1d(ji,2/3) = rt0_ice 162 ! Computation of energies due to surface and bottom melting 163 !----------------------------------------------------------------------- 164 DO ji = kideb , kiut 165 zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 166 zihic = MAX( zzero , SIGN( zone , hicdif - h_ice_1d(ji) ) ) 167 !--computation of energy due to surface melting 168 zqcmlt(ji,1) = ( MAX ( zzero , & 169 & rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 170 !--computation of energy due to bottom melting 171 zqcmlt(ji,2) = ( MAX( zzero , & 172 & rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 173 & + MAX( zzero , & 174 & rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 175 & ) * ( 1.0 - zihic ) 176 !--limitation of snow/ice system internal temperature 177 tbif_1d(ji,1) = MIN( rt0_snow, tbif_1d(ji,1) ) 178 tbif_1d(ji,2) = MIN( rt0_ice , tbif_1d(ji,2) ) 179 tbif_1d(ji,3) = MIN( rt0_ice , tbif_1d(ji,3) ) 180 END DO 181 182 !------------------------------------------- 183 ! 2. Calculate some intermediate variables. 184 !------------------------------------------- 185 zhsu = hnzst ! initialisation of the thickness of surface layer 186 ! 187 DO ji = kideb , kiut 188 zind = MAX( zzero , SIGN( zone , zhsu - h_snow_1d(ji) ) ) 189 zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 190 zihsn = MAX( zihsn , zind ) 191 zihic = MAX( zzero , sign( zone , hicdif - h_ice_1d (ji) ) ) 192 ! 193 ! 2.1. Computation of surface melting temperature 194 !---------------------------------------------------- 195 zind = MAX( zzero , SIGN( zone , -h_snow_1d(ji) ) ) 196 ztsmlt(ji) = ( 1.0 - zind ) * rt0_snow + zind * rt0_ice 197 ! 198 ! 2.2. Effective conductivity of snow and ice 199 !----------------------------------------------- 200 201 !---computation of the correction factor on the thermal conductivity 202 !-- (Morales Maqueda, 1995 ; Fichefet and Morales Maqueda, 1997) 203 zhe = ( rcdsn / ( rcdsn + rcdic ) ) * h_ice_1d(ji) & 204 & + ( rcdic / ( rcdsn + rcdic ) ) * h_snow_1d(ji) 205 zihe = MAX( zzero , SIGN( zone , 2.0 * zhe - thth ) ) 206 zheshth = zhe / thth 207 zghe = ( 1.0 - zihe ) * zheshth * ( 2.0 - zheshth ) & 208 & + zihe * 0.5 * ( 1.5 + LOG( 2.0 * zheshth ) ) 209 210 !---effective conductivities 211 zksn(ji) = zghe * rcdsn 212 zkic(ji) = zghe * rcdic 213 ! 214 ! 2.3. Computation of the conductive heat flux from the snow/ice 215 ! system interior toward the top surface 216 !------------------------------------------------------------------ 217 218 !---Thermal conductivity at the mid-point of the first snow/ice system layer 219 zksndh(ji) = ( ( 1.0 - zihsn ) * 2.0 * zksn(ji) + zihsn * 4.0 * zkic(ji) ) & 220 & / ( ( 1.0 - zihsn ) * h_snow_1d(ji) & 221 & + zihsn * ( ( 1.0 + 3.0 * zihic ) * h_ice_1d(ji) & 222 & + 4.0 * zkic(ji)/zksn(ji) * h_snow_1d(ji) ) ) 223 224 !---internal temperature at the mid-point of the first snow/ice system layer 225 ztbif(ji) = ( 1.0 - zihsn ) * tbif_1d(ji,1) & 226 & + zihsn * ( ( 1.0 - zihic ) * tbif_1d(ji,2) & 227 & + zihic * tfu_1d(ji) ) 228 !---conductive heat flux 229 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 230 231 END DO 232 233 !-------------------------------------------------------------------- 234 ! 3. Calculate : 235 ! - fstbif_1d, part of solar radiation absorbing inside the ice 236 ! assuming an exponential absorption (Grenfell and Maykut, 1977) 237 ! - zqmax, maximum energy stored in brine pockets 238 ! - qstbif_1d, total energy stored in brine pockets (updating) 239 !------------------------------------------------------------------- 240 241 DO ji = kideb , kiut 242 zihsn = MAX( zzero , SIGN (zone , -h_snow_1d(ji) ) ) 243 zihic = MAX( zzero , 1.0 - ( h_ice_1d(ji) / zhsu ) ) 244 zind = MAX( zzero , SIGN (zone , hicdif - h_ice_1d(ji) ) ) 245 !--Computation of the fraction of the net shortwave radiation which 246 !--penetrates inside the ice cover ( See Forcat) 247 zi0(ji) = zihsn * ( fr1_i0_1d(ji) + zihic * fr2_i0_1d(ji) ) 248 zexp = MIN( zone , EXP( -1.5 * ( h_ice_1d(ji) - zhsu ) ) ) 249 fstbif_1d(ji) = zi0(ji) * qsr_ice_1d(ji) * zexp 250 !--Computation of maximum energy stored in brine pockets zqmax and update 251 !--the total energy stored in brine pockets, if less than zqmax 252 zqmax(ji) = MAX( zzero , 0.5 * xlic * ( h_ice_1d(ji) - hicmin ) ) 253 zfsab = zi0(ji) * qsr_ice_1d(ji) * ( 1.0 - zexp ) 254 zihq = ( 1.0 - zind ) * MAX(zzero, SIGN( zone , qstbif_1d(ji) - zqmax(ji) ) ) & 255 & + zind * zone 256 qstbif_1d(ji) = ( qstbif_1d(ji) + ( 1.0 - zihq ) * zfsab * rdt_ice ) * swiqst 257 !--fraction of shortwave radiation absorbed at surface 258 ziexp = zihq * zexp + ( 1.0 - zihq ) * ( swiqst + ( 1.0 - swiqst ) * zexp ) 259 z1mi0(ji) = 1.0 - zi0(ji) * ziexp 260 END DO 261 262 !-------------------------------------------------------------------------------- 263 ! 4. Computation of the surface temperature : determined by considering the 264 ! budget of a thin layer of thick. zhsu at the top surface (H. Grenier, 1995) 265 ! and based on a surface energy balance : 266 ! hsu * rcp * dT/dt = Fsr + Fnsr(T) + Fcs(T), 267 ! where - Fsr is the net absorbed solar radiation, 268 ! - Fnsr is the total non solar radiation (incoming and outgoing long-wave, 269 ! sensible and latent heat fluxes) 270 ! - Fcs the conductive heat flux at the top of surface 271 !------------------------------------------------------------------------------ 272 273 ! 4.1. Computation of intermediate values 274 !--------------------------------------------- 275 DO ji = kideb, kiut 276 zrcpdt(ji) = ( rcpsn * MIN( h_snow_1d(ji) , zhsu ) & 277 & + rcpic * MAX( zhsu - h_snow_1d(ji) , zzero ) ) / rdt_ice 278 zts_old(ji) = sist_1d(ji) 279 END DO 280 281 ! 4.2. Computation of surface temperature by expanding the eq. of energy balance 282 ! with Ts = Tp + DT. One obtain , F(Tp) + DT * DF(Tp) = 0 283 ! where - F(Tp) = Fsr + Fnsr(Tp) + Fcs(Tp) 284 ! - DF(Tp)= (dFnsr(Tp)/dT) + (dFcs(Tp)/dT) - hsu*rcp/dt 285 !--------------------------------------------------------------------------------- 286 287 DO ji = kideb, kiut 288 !---computation of the derivative of energy balance function 289 zdfts = zksndh(ji) & ! contribution of the conductive heat flux 290 & + zrcpdt(ji) & ! contribution of hsu * rcp / dt 291 & - dqns_ice_1d (ji) ! contribution of the total non solar radiation 292 !---computation of the energy balance function 293 zfts = - z1mi0 (ji) * qsr_ice_1d(ji) & ! net absorbed solar radiation 294 & - qns_ice_1d(ji) & ! total non solar radiation 295 & - zfcsu (ji) ! conductive heat flux from the surface 296 !---computation of surface temperature increment 297 zdts = -zfts / zdfts 298 !---computation of the new surface temperature 299 sist_1d(ji) = sist_1d(ji) + zdts 300 END DO 301 302 !---------------------------------------------------------------------------- 303 ! 5. Boundary condition at the top surface 304 !-- IF Tsb < Tmelt, Fnet = Fcs (the net heat flux equal the conductive heat flux) 305 ! Otherwise Tsb = Tmelt and Qnet(Tmelt) > 0 306 ! Fnet(Tmelt) is therefore the net surface flux needed for melting 307 !---------------------------------------------------------------------------- 75 REAL(wp), DIMENSION(jpij,2) :: zqcmlt ! energy due to surface( /1 ) and bottom melting( /2 ) 76 REAL(wp), DIMENSION(jpij) :: & 77 ztsmlt & ! snow/ice surface melting temperature 78 ,ztbif & ! int. temp. at the mid-point of the 1st layer of the snow/ice sys. 79 ,zksn & ! effective conductivity of snow 80 ,zkic & ! effective conductivity of ice 81 ,zksndh & ! thermal cond. at the mid-point of the 1st layer of the snow/ice sys. 82 , zfcsu & ! conductive heat flux at the surface of the snow/ice system 83 , zfcsudt & ! = zfcsu * dt 84 , zi0 & ! frac. of the net SW rad. which is not absorbed at the surface 85 , z1mi0 & ! fraction of the net SW radiation absorbed at the surface 86 , zqmax & ! maximum energy stored in brine pockets 87 , zrcpdt & ! h_su*rho_su*cp_su/dt(h_su being the thick. of surf. layer) 88 , zts_old & ! previous surface temperature 89 , zidsn , z1midsn , zidsnic ! tempory variables 90 REAL(wp), DIMENSION(jpij) :: & 91 zfnet & ! net heat flux at the top surface( incl. conductive heat flux) 92 , zsprecip & ! snow accumulation 93 , zhsnw_old & ! previous snow thickness 94 , zdhictop & ! change in ice thickness due to top surf ablation/accretion 95 , zdhicbot & ! change in ice thickness due to bottom surf abl/acc 96 , zqsup & ! energy transmitted to ocean (coming from top surf abl/acc) 97 , zqocea & ! energy transmitted to ocean (coming from bottom sur abl/acc) 98 , zfrl_old & ! previous sea/ice fraction 99 , zfrld_1d & ! new sea/ice fraction 100 , zep ! internal temperature of the 2nd layer of the snow/ice system 101 REAL(wp), DIMENSION(3) :: & 102 zplediag & ! principle diagonal, subdiag. and supdiag. of the 103 , zsubdiag & ! tri-diagonal matrix coming from the computation 104 , zsupdiag & ! of the temperatures inside the snow-ice system 105 , zsmbr ! second member 106 REAL(wp) :: & 107 zhsu & ! thickness of surface layer 108 , zhe & ! effective thickness for compu. of equ. thermal conductivity 109 , zheshth & ! = zhe / thth 110 , zghe & ! correction factor of the thermal conductivity 111 , zumsb & ! parameter for numerical method to solve heat-diffusion eq. 112 , zkhsn & ! conductivity at the snow layer 113 , zkhic & ! conductivity at the ice layers 114 , zkint & ! equivalent conductivity at the snow-ice interface 115 , zkhsnint & ! = zkint*dt / (hsn*rhosn*cpsn) 116 , zkhicint & ! = 2*zkint*dt / (hic*rhoic*cpic) 117 , zpiv1 , zpiv2 & ! tempory scalars used to solve the tri-diagonal system 118 , zb2 , zd2 , zb3 , zd3 & 119 , ztint ! equivalent temperature at the snow-ice interface 120 REAL(wp) :: & 121 zexp & ! exponential function of the ice thickness 122 , zfsab & ! part of solar radiation stored in brine pockets 123 , zfts & ! value of energy balance function when the temp. equal surf. temp. 124 , zdfts & ! value of derivative of ztfs when the temp. equal surf. temp. 125 , zdts & ! surface temperature increment 126 , zqsnw_mlt & ! energy needed to melt snow 127 , zdhsmlt & ! change in snow thickness due to melt 128 , zhsn & ! snow thickness (previous+accumulation-melt) 129 , zqsn_mlt_rem & ! remaining heat coming from snow melting 130 , zqice_top_mlt & ! energy used to melt ice at top surface 131 , zdhssub & ! change in snow thick. due to sublimation or evaporation 132 , zdhisub & ! change in ice thick. due to sublimation or evaporation 133 , zdhsn & ! snow ice thickness increment 134 , zdtsn & ! snow internal temp. increment 135 , zdtic & ! ice internal temp. increment 136 , zqnes ! conductive energy due to ice melting in the first ice layer 137 REAL(wp) :: & 138 ztbot & ! temperature at the bottom surface 139 , zfcbot & ! conductive heat flux at bottom surface 140 , zqice_bot & ! energy used for bottom melting/growing 141 , zqice_bot_mlt & ! energy used for bottom melting 142 , zqstbif_bot & ! part of energy stored in brine pockets used for bottom melting 143 , zqstbif_old & ! tempory var. for zqstbif_bot 144 , zdhicmlt & ! change in ice thickness due to bottom melting 145 , zdhicm & ! change in ice thickness var. 146 , zdhsnm & ! change in snow thickness var. 147 , zhsnfi & ! snow thickness var. 148 , zc1, zpc1, zc2, zpc2, zp1, zp2 & ! tempory variables 149 , ztb2, ztb3 150 REAL(wp) :: & 151 zdrmh & ! change in snow/ice thick. after snow-ice formation 152 , zhicnew & ! new ice thickness 153 , zhsnnew & ! new snow thickness 154 , zquot , ztneq & ! tempory temp. variables 155 , zqice, zqicetot & ! total heat inside the snow/ice system 156 , zdfrl & ! change in ice concentration 157 , zdvsnvol & ! change in snow volume 158 , zdrfrl1, zdrfrl2 & ! tempory scalars 159 , zihsn, zidhb, zihic, zihe, zihq, ziexp, ziqf, zihnf, zibmlt, ziqr, zihgnew, zind 160 !!---------------------------------------------------------------------- 161 162 !----------------------------------------------------------------------- 163 ! 1. Boundaries conditions for snow/ice system internal temperature 164 ! - If tbif_1d(ji,1) > rt0_snow, tbif_1d(ji,1) = rt0_snow 165 ! - If tbif_1d(ji,2/3) > rt0_ice, tbif_1d(ji,2/3) = rt0_ice 166 ! Computation of energies due to surface and bottom melting 167 !----------------------------------------------------------------------- 168 169 DO ji = kideb , kiut 170 zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 171 zihic = MAX( zzero , SIGN( zone , hicdif - h_ice_1d(ji) ) ) 172 !--computation of energy due to surface melting 173 zqcmlt(ji,1) = ( MAX ( zzero , & 174 & rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 175 !--computation of energy due to bottom melting 176 zqcmlt(ji,2) = ( MAX( zzero , & 177 & rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 178 & + MAX( zzero , & 179 & rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 180 & ) * ( 1.0 - zihic ) 181 !--limitation of snow/ice system internal temperature 182 tbif_1d(ji,1) = MIN( rt0_snow, tbif_1d(ji,1) ) 183 tbif_1d(ji,2) = MIN( rt0_ice , tbif_1d(ji,2) ) 184 tbif_1d(ji,3) = MIN( rt0_ice , tbif_1d(ji,3) ) 185 END DO 186 187 !------------------------------------------- 188 ! 2. Calculate some intermediate variables. 189 !------------------------------------------- 190 191 ! initialisation of the thickness of surface layer 192 zhsu = hnzst 193 194 DO ji = kideb , kiut 195 zind = MAX( zzero , SIGN( zone , zhsu - h_snow_1d(ji) ) ) 196 zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 197 zihsn = MAX( zihsn , zind ) 198 zihic = MAX( zzero , sign( zone , hicdif - h_ice_1d(ji) ) ) 199 ! 2.1. Computation of surface melting temperature 200 !---------------------------------------------------- 201 zind = MAX( zzero , SIGN( zone , -h_snow_1d(ji) ) ) 202 ztsmlt(ji) = ( 1.0 - zind ) * rt0_snow + zind * rt0_ice 203 ! 204 ! 2.2. Effective conductivity of snow and ice 205 !----------------------------------------------- 206 207 !---computation of the correction factor on the thermal conductivity 208 !-- (Morales Maqueda, 1995 ; Fichefet and Morales Maqueda, 1997) 209 zhe = ( rcdsn / ( rcdsn + rcdic ) ) * h_ice_1d(ji) & 210 & + ( rcdic / ( rcdsn + rcdic ) ) * h_snow_1d(ji) 211 zihe = MAX( zzero , SIGN( zone , 2.0 * zhe - thth ) ) 212 zheshth = zhe / thth 213 zghe = ( 1.0 - zihe ) * zheshth * ( 2.0 - zheshth ) & 214 & + zihe * 0.5 * ( 1.5 + LOG( 2.0 * zheshth ) ) 215 216 !---effective conductivities 217 zksn(ji) = zghe * rcdsn 218 zkic(ji) = zghe * rcdic 219 220 ! 221 ! 2.3. Computation of the conductive heat flux from the snow/ice 222 ! system interior toward the top surface 223 !------------------------------------------------------------------ 224 225 !---Thermal conductivity at the mid-point of the first snow/ice system layer 226 zksndh(ji) = ( ( 1.0 - zihsn ) * 2.0 * zksn(ji) + zihsn * 4.0 * zkic(ji) ) & 227 & / ( ( 1.0 - zihsn ) * h_snow_1d(ji) & 228 & + zihsn * ( ( 1.0 + 3.0 * zihic ) * h_ice_1d(ji) & 229 & + 4.0 * zkic(ji)/zksn(ji) * h_snow_1d(ji) ) ) 230 231 !---internal temperature at the mid-point of the first snow/ice system layer 232 ztbif(ji) = ( 1.0 - zihsn ) * tbif_1d(ji,1) & 233 & + zihsn * ( ( 1.0 - zihic ) * tbif_1d(ji,2) & 234 & + zihic * tfu_1d(ji) ) 235 !---conductive heat flux 236 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 237 238 END DO 239 240 !-------------------------------------------------------------------- 241 ! 3. Calculate : 242 ! - fstbif_1d, part of solar radiation absorbing inside the ice 243 ! assuming an exponential absorption (Grenfell and Maykut, 1977) 244 ! - zqmax, maximum energy stored in brine pockets 245 ! - qstbif_1d, total energy stored in brine pockets (updating) 246 !------------------------------------------------------------------- 247 248 DO ji = kideb , kiut 249 zihsn = MAX( zzero , SIGN (zone , -h_snow_1d(ji) ) ) 250 zihic = MAX( zzero , 1.0 - ( h_ice_1d(ji) / zhsu ) ) 251 zind = MAX( zzero , SIGN (zone , hicdif - h_ice_1d(ji) ) ) 252 !--Computation of the fraction of the net shortwave radiation which 253 !--penetrates inside the ice cover ( See Forcat) 254 zi0(ji) = zihsn * ( fr1_i0_1d(ji) + zihic * fr2_i0_1d(ji) ) 255 zexp = MIN( zone , EXP( -1.5 * ( h_ice_1d(ji) - zhsu ) ) ) 256 fstbif_1d(ji) = zi0(ji) * qsr_ice_1d(ji) * zexp 257 !--Computation of maximum energy stored in brine pockets zqmax and update 258 !--the total energy stored in brine pockets, if less than zqmax 259 zqmax(ji) = MAX( zzero , 0.5 * xlic * ( h_ice_1d(ji) - hicmin ) ) 260 zfsab = zi0(ji) * qsr_ice_1d(ji) * ( 1.0 - zexp ) 261 zihq = ( 1.0 - zind ) * MAX(zzero, SIGN( zone , qstbif_1d(ji) - zqmax(ji) ) ) & 262 & + zind * zone 263 qstbif_1d(ji) = ( qstbif_1d(ji) + ( 1.0 - zihq ) * zfsab * rdt_ice ) * swiqst 264 !--fraction of shortwave radiation absorbed at surface 265 ziexp = zihq * zexp + ( 1.0 - zihq ) * ( swiqst + ( 1.0 - swiqst ) * zexp ) 266 z1mi0(ji) = 1.0 - zi0(ji) * ziexp 267 END DO 268 269 !-------------------------------------------------------------------------------- 270 ! 4. Computation of the surface temperature : determined by considering the 271 ! budget of a thin layer of thick. zhsu at the top surface (H. Grenier, 1995) 272 ! and based on a surface energy balance : 273 ! hsu * rcp * dT/dt = Fsr + Fnsr(T) + Fcs(T), 274 ! where - Fsr is the net absorbed solar radiation, 275 ! - Fnsr is the total non solar radiation (incoming and outgoing long-wave, 276 ! sensible and latent heat fluxes) 277 ! - Fcs the conductive heat flux at the top of surface 278 !------------------------------------------------------------------------------ 279 280 ! 4.1. Computation of intermediate values 281 !--------------------------------------------- 282 DO ji = kideb, kiut 283 zrcpdt(ji) = ( rcpsn * MIN( h_snow_1d(ji) , zhsu ) & 284 & + rcpic * MAX( zhsu - h_snow_1d(ji) , zzero ) ) / rdt_ice 285 zts_old(ji) = sist_1d(ji) 286 END DO 287 288 ! 4.2. Computation of surface temperature by expanding the eq. of energy balance 289 ! with Ts = Tp + DT. One obtain , F(Tp) + DT * DF(Tp) = 0 290 ! where - F(Tp) = Fsr + Fnsr(Tp) + Fcs(Tp) 291 ! - DF(Tp)= (dFnsr(Tp)/dT) + (dFcs(Tp)/dT) - hsu*rcp/dt 292 !--------------------------------------------------------------------------------- 293 294 DO ji = kideb, kiut 295 !---computation of the derivative of energy balance function 296 zdfts = zksndh(ji) & ! contribution of the conductive heat flux 297 & + zrcpdt(ji) & ! contribution of hsu * rcp / dt 298 & - dqns_ice_1d (ji) ! contribution of the total non solar radiation 299 !---computation of the energy balance function 300 zfts = - z1mi0 (ji) * qsr_ice_1d(ji) & ! net absorbed solar radiation 301 & - qns_ice_1d(ji) & ! total non solar radiation 302 & - zfcsu (ji) ! conductive heat flux from the surface 303 !---computation of surface temperature increment 304 zdts = -zfts / zdfts 305 !---computation of the new surface temperature 306 sist_1d(ji) = sist_1d(ji) + zdts 307 END DO 308 309 !---------------------------------------------------------------------------- 310 ! 5. Boundary condition at the top surface 311 !-- IF Tsb < Tmelt, Fnet = Fcs (the net heat flux equal the conductive heat flux) 312 ! Otherwise Tsb = Tmelt and Qnet(Tmelt) > 0 313 ! Fnet(Tmelt) is therefore the net surface flux needed for melting 314 !---------------------------------------------------------------------------- 308 315 309 316 310 ! 5.1. Limitation of surface temperature and update total non solar fluxes,311 ! latent heat flux and conductive flux at the top surface312 !----------------------------------------------------------------------317 ! 5.1. Limitation of surface temperature and update total non solar fluxes, 318 ! latent heat flux and conductive flux at the top surface 319 !---------------------------------------------------------------------- 313 320 314 IF ( .NOT. lk_cpl ) THEN ! duplicate the loop for performances issues315 DO ji = kideb, kiut316 sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) )317 qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) )318 qla_ice_1d(ji) = qla_ice_1d(ji) + dqla_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) )319 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) )320 END DO321 ELSE322 DO ji = kideb, kiut323 sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) )324 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) )325 END DO326 ENDIF327 328 ! 5.2. Calculate available heat for surface ablation.329 !---------------------------------------------------------------------330 331 DO ji = kideb, kiut332 zfnet(ji) = qns_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji)333 zfnet(ji) = MAX( zzero , zfnet(ji) )334 zfnet(ji) = zfnet(ji) * MAX( zzero , SIGN( zone , sist_1d(ji) - ztsmlt(ji) ) )335 END DO336 337 !-------------------------------------------------------------------------338 ! 6. Calculate changes in internal temperature due to vertical diffusion339 ! processes. The evolution of this temperature is governed by the one-340 ! dimensionnal heat-diffusion equation.341 ! Given the temperature tbif(1/2/3), at time m we solve a set342 ! of finite difference equations to obtain new tempe. Each tempe is coupled343 ! to the temp. immediatly above and below by heat conduction terms. Thus344 ! we have a set of equations of the form A * T = B, where A is a tridiagonal345 ! matrix, T a vector whose components are the unknown new temp.346 !-------------------------------------------------------------------------321 IF ( .NOT. lk_cpl ) THEN ! duplicate the loop for performances issues 322 DO ji = kideb, kiut 323 sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 324 qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 325 qla_ice_1d(ji) = qla_ice_1d(ji) + dqla_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 326 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 327 END DO 328 ELSE 329 DO ji = kideb, kiut 330 sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 331 zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 332 END DO 333 ENDIF 334 335 ! 5.2. Calculate available heat for surface ablation. 336 !--------------------------------------------------------------------- 337 338 DO ji = kideb, kiut 339 zfnet(ji) = qns_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji) 340 zfnet(ji) = MAX( zzero , zfnet(ji) ) 341 zfnet(ji) = zfnet(ji) * MAX( zzero , SIGN( zone , sist_1d(ji) - ztsmlt(ji) ) ) 342 END DO 343 344 !------------------------------------------------------------------------- 345 ! 6. Calculate changes in internal temperature due to vertical diffusion 346 ! processes. The evolution of this temperature is governed by the one- 347 ! dimensionnal heat-diffusion equation. 348 ! Given the temperature tbif(1/2/3), at time m we solve a set 349 ! of finite difference equations to obtain new tempe. Each tempe is coupled 350 ! to the temp. immediatly above and below by heat conduction terms. Thus 351 ! we have a set of equations of the form A * T = B, where A is a tridiagonal 352 ! matrix, T a vector whose components are the unknown new temp. 353 !------------------------------------------------------------------------- 347 354 348 !--parameter for the numerical methode use to solve the heat-diffusion equation349 !- implicit, explicit or Crank-Nicholson350 zumsb = 1.0 - sbeta351 DO ji = kideb, kiut352 zidsn(ji) = MAX ( zzero, SIGN( zone, hsndif - h_snow_1d(ji) ) )353 z1midsn(ji) = 1.0 - zidsn(ji)354 zihic = MAX ( zzero, SIGN( zone, hicdif - h_ice_1d(ji) ) )355 zidsnic(ji) = zidsn(ji) * zihic356 zfcsudt(ji) = zfcsu(ji) * rdt_ice357 END DO355 !--parameter for the numerical methode use to solve the heat-diffusion equation 356 !- implicit, explicit or Crank-Nicholson 357 zumsb = 1.0 - sbeta 358 DO ji = kideb, kiut 359 zidsn(ji) = MAX ( zzero, SIGN( zone, hsndif - h_snow_1d(ji) ) ) 360 z1midsn(ji) = 1.0 - zidsn(ji) 361 zihic = MAX ( zzero, SIGN( zone, hicdif - h_ice_1d(ji) ) ) 362 zidsnic(ji) = zidsn(ji) * zihic 363 zfcsudt(ji) = zfcsu(ji) * rdt_ice 364 END DO 358 365 359 DO ji = kideb, kiut360 361 ! 6.1 Calculate intermediate variables.362 !----------------------------------------363 364 !--conductivity at the snow surface365 zkhsn = 2.0 * zksn(ji) * rdt_ice / rcpsn366 !--conductivity at the ice surface367 zkhic = 4.0 * zkic(ji) * rdt_ice / MAX( h_ice_1d(ji) * h_ice_1d(ji) * rcpic , epsi20 )368 !--conductivity at the snow/ice interface369 zkint = 4.0 * zksn(ji) * zkic(ji) &370 & / ( zksn(ji) * h_ice_1d(ji) + 2.0 * zkic(ji) * h_snow_1d(ji) * z1midsn(ji))371 zkhsnint = zkint * rdt_ice / rcpsn372 zkhicint = zkint * 2.0 * rdt_ice / MAX( h_ice_1d(ji) * rcpic , epsi20 )373 374 ! 6.2. Fulfill the linear system matrix.375 !-----------------------------------------366 DO ji = kideb, kiut 367 368 ! 6.1 Calculate intermediate variables. 369 !---------------------------------------- 370 371 !--conductivity at the snow surface 372 zkhsn = 2.0 * zksn(ji) * rdt_ice / rcpsn 373 !--conductivity at the ice surface 374 zkhic = 4.0 * zkic(ji) * rdt_ice / MAX( h_ice_1d(ji) * h_ice_1d(ji) * rcpic , epsi20 ) 375 !--conductivity at the snow/ice interface 376 zkint = 4.0 * zksn(ji) * zkic(ji) & 377 & / ( zksn(ji) * h_ice_1d(ji) + 2.0 * zkic(ji) * h_snow_1d(ji) * z1midsn(ji)) 378 zkhsnint = zkint * rdt_ice / rcpsn 379 zkhicint = zkint * 2.0 * rdt_ice / MAX( h_ice_1d(ji) * rcpic , epsi20 ) 380 381 ! 6.2. Fulfill the linear system matrix. 382 !----------------------------------------- 376 383 !$$$ zplediag(1) = 1 + sbeta * z1midsn(ji) * ( zkhsn + zkhsnint ) 377 zplediag(1) = zidsn(ji) + z1midsn(ji) * h_snow_1d(ji) &378 & + sbeta * z1midsn(ji) * zkhsnint379 zplediag(2) = 1 + sbeta * ( z1midsn(ji) * zkhicint + zkhic )380 zplediag(3) = 1 + 3.0 * sbeta * zkhic381 382 zsubdiag(1) = 0.e0383 zsubdiag(2) = -1.e0 * z1midsn(ji) * sbeta * zkhicint384 zsubdiag(3) = -1.e0 * sbeta * zkhic385 386 zsupdiag(1) = -1.e0 * z1midsn(ji) * sbeta * zkhsnint387 zsupdiag(2) = zsubdiag(3)388 zsupdiag(3) = 0.e0389 390 ! 6.3. Fulfill the idependent term vector.391 !-------------------------------------------384 zplediag(1) = zidsn(ji) + z1midsn(ji) * h_snow_1d(ji) & 385 & + sbeta * z1midsn(ji) * zkhsnint 386 zplediag(2) = 1 + sbeta * ( z1midsn(ji) * zkhicint + zkhic ) 387 zplediag(3) = 1 + 3.0 * sbeta * zkhic 388 389 zsubdiag(1) = 0.e0 390 zsubdiag(2) = -1.e0 * z1midsn(ji) * sbeta * zkhicint 391 zsubdiag(3) = -1.e0 * sbeta * zkhic 392 393 zsupdiag(1) = -1.e0 * z1midsn(ji) * sbeta * zkhsnint 394 zsupdiag(2) = zsubdiag(3) 395 zsupdiag(3) = 0.e0 396 397 ! 6.3. Fulfill the idependent term vector. 398 !------------------------------------------- 392 399 393 400 !$$$ zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) * & … … 395 402 !$$$ & - zumsb * ( zkhsn * tbif_1d(ji,1) 396 403 !$$$ & + zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) ) ) 397 zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) * &398 & ( h_snow_1d(ji) * tbif_1d(ji,1) - ( zfcsudt(ji) / rcpsn ) &399 & - zumsb * zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) )400 401 zsmbr(2) = tbif_1d(ji,2) &402 & - zidsn(ji) * ( 1.0 - zidsnic(ji) ) &403 & * ( zfcsudt(ji) / MAX( h_ice_1d(ji) * rcpic , epsi20 ) ) &404 & + zumsb * ( zkhicint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) &405 & - zkhic * ( tbif_1d(ji,2) - tbif_1d(ji,3) ) )406 407 zsmbr(3) = tbif_1d(ji,3) &408 & + zkhic * ( 2.0 * tfu_1d(ji) &409 & + zumsb * ( tbif_1d(ji,2) - 3.0 * tbif_1d(ji,3) ) )410 411 ! 6.4. Solve the system (Gauss elimination method).412 !----------------------------------------------------413 414 zpiv1 = zsubdiag(2) / zplediag(1)415 zb2 = zplediag(2) - zpiv1 * zsupdiag(1)416 zd2 = zsmbr(2) - zpiv1 * zsmbr(1)417 418 zpiv2 = zsubdiag(3) / zb2419 zb3 = zplediag(3) - zpiv2 * zsupdiag(2)420 zd3 = zsmbr(3) - zpiv2 * zd2421 422 tbif_1d(ji,3) = zd3 / zb3423 tbif_1d(ji,2) = ( zd2 - zsupdiag(2) * tbif_1d(ji,3) ) / zb2424 tbif_1d(ji,1) = ( zsmbr(1) - zsupdiag(1) * tbif_1d(ji,2) ) / zplediag(1)425 426 !--- taking into account the particular case of zidsnic(ji) = 1427 ztint = ( zkic(ji) * h_snow_1d(ji) * tfu_1d (ji) &428 & + zksn(ji) * h_ice_1d(ji) * sist_1d(ji) ) &429 & / ( zkic(ji) * h_snow_1d(ji) + zksn(ji) * h_ice_1d(ji) )430 431 tbif_1d(ji,1) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,1) &432 & + zidsnic(ji) * ( ztint + sist_1d(ji) ) / 2.0433 tbif_1d(ji,2) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,2) &434 & + zidsnic(ji) * ( 3.0 * ztint + tfu_1d(ji) ) / 4.0435 tbif_1d(ji,3) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,3) &436 & + zidsnic(ji) * ( ztint + 3.0 * tfu_1d(ji) ) / 4.0437 END DO404 zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) * & 405 & ( h_snow_1d(ji) * tbif_1d(ji,1) - ( zfcsudt(ji) / rcpsn ) & 406 & - zumsb * zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) ) 407 408 zsmbr(2) = tbif_1d(ji,2) & 409 & - zidsn(ji) * ( 1.0 - zidsnic(ji) ) & 410 & * ( zfcsudt(ji) / MAX( h_ice_1d(ji) * rcpic , epsi20 ) ) & 411 & + zumsb * ( zkhicint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) & 412 & - zkhic * ( tbif_1d(ji,2) - tbif_1d(ji,3) ) ) 413 414 zsmbr(3) = tbif_1d(ji,3) & 415 & + zkhic * ( 2.0 * tfu_1d(ji) & 416 & + zumsb * ( tbif_1d(ji,2) - 3.0 * tbif_1d(ji,3) ) ) 417 418 ! 6.4. Solve the system (Gauss elimination method). 419 !---------------------------------------------------- 420 421 zpiv1 = zsubdiag(2) / zplediag(1) 422 zb2 = zplediag(2) - zpiv1 * zsupdiag(1) 423 zd2 = zsmbr(2) - zpiv1 * zsmbr(1) 424 425 zpiv2 = zsubdiag(3) / zb2 426 zb3 = zplediag(3) - zpiv2 * zsupdiag(2) 427 zd3 = zsmbr(3) - zpiv2 * zd2 428 429 tbif_1d(ji,3) = zd3 / zb3 430 tbif_1d(ji,2) = ( zd2 - zsupdiag(2) * tbif_1d(ji,3) ) / zb2 431 tbif_1d(ji,1) = ( zsmbr(1) - zsupdiag(1) * tbif_1d(ji,2) ) / zplediag(1) 432 433 !--- taking into account the particular case of zidsnic(ji) = 1 434 ztint = ( zkic(ji) * h_snow_1d(ji) * tfu_1d (ji) & 435 & + zksn(ji) * h_ice_1d(ji) * sist_1d(ji) ) & 436 & / ( zkic(ji) * h_snow_1d(ji) + zksn(ji) * h_ice_1d(ji) ) 437 438 tbif_1d(ji,1) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,1) & 439 & + zidsnic(ji) * ( ztint + sist_1d(ji) ) / 2.0 440 tbif_1d(ji,2) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,2) & 441 & + zidsnic(ji) * ( 3.0 * ztint + tfu_1d(ji) ) / 4.0 442 tbif_1d(ji,3) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,3) & 443 & + zidsnic(ji) * ( ztint + 3.0 * tfu_1d(ji) ) / 4.0 444 END DO 438 445 439 !----------------------------------------------------------------------440 ! 9. Take into account surface ablation and bottom accretion-ablation.|441 !----------------------------------------------------------------------446 !---------------------------------------------------------------------- 447 ! 9. Take into account surface ablation and bottom accretion-ablation.| 448 !---------------------------------------------------------------------- 442 449 443 !---Snow accumulation in one thermodynamic time step444 zsprecip(kideb:kiut) = sprecip_1d(kideb:kiut) * rdt_ice / rhosn445 446 447 DO ji = kideb, kiut448 449 ! 9.1. Surface ablation and update of snow thickness and qstbif_1d450 !--------------------------------------------------------------------450 !---Snow accumulation in one thermodynamic time step 451 zsprecip(kideb:kiut) = sprecip_1d(kideb:kiut) * rdt_ice / rhosn 452 453 454 DO ji = kideb, kiut 455 456 ! 9.1. Surface ablation and update of snow thickness and qstbif_1d 457 !-------------------------------------------------------------------- 451 458 452 459 !-------------------------------------------------------------------------- … … 752 759 qstbif_1d(ji) = zdrfrl2 * qstbif_1d(ji) 753 760 frld_1d(ji) = zfrld_1d(ji) 754 ! 755 END DO 756 ! 761 ! 762 END DO 763 ! 764 END SUBROUTINE lim_thd_zdf_2 765 766 #else 767 !!---------------------------------------------------------------------- 768 !! Default Option NO sea-ice model 769 !!---------------------------------------------------------------------- 770 CONTAINS 771 SUBROUTINE lim_thd_zdf_2 ! Empty routine 757 772 END SUBROUTINE lim_thd_zdf_2 758 759 #else760 !!----------------------------------------------------------------------761 !! Default Option NO sea-ice model762 !!----------------------------------------------------------------------763 773 #endif 764 774 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limtrp_2.F90
r1855 r1857 4 4 !! LIM 2.0 transport ice model : sea-ice advection/diffusion 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (LIM) Original code7 !! 1.0 ! 2001-05 (G. Madec, R. Hordoir) opa norm8 !! 2.0 ! 2004-01 (G. Madec, C. Ethe) F90, mpp9 !!----------------------------------------------------------------------10 6 #if defined key_lim2 11 7 !!---------------------------------------------------------------------- … … 15 11 !! lim_trp_init_2 : initialization and namelist read 16 12 !!---------------------------------------------------------------------- 13 !! * Modules used 17 14 USE phycst 18 15 USE dom_oce … … 29 26 PRIVATE 30 27 31 PUBLIC lim_trp_2 ! called by sbc_ice_lim_2 32 33 REAL(wp), PUBLIC :: bound = 0.e0 !: boundary condit. (0.0 no-slip, 1.0 free-slip) 34 35 REAL(wp) :: epsi06 = 1.e-06 ! constant values 36 REAL(wp) :: epsi03 = 1.e-03 ! 37 REAL(wp) :: epsi16 = 1.e-16 ! 38 REAL(wp) :: rzero = 0.e0 ! 39 REAL(wp) :: rone = 1.e0 ! 28 !! * Routine accessibility 29 PUBLIC lim_trp_2 ! called by sbc_ice_lim_2 30 31 !! * Shared module variables 32 REAL(wp), PUBLIC :: & !: 33 bound = 0.e0 !: boundary condit. (0.0 no-slip, 1.0 free-slip) 34 35 !! * Module variables 36 REAL(wp) :: & ! constant values 37 epsi06 = 1.e-06 , & 38 epsi03 = 1.e-03 , & 39 epsi16 = 1.e-16 , & 40 rzero = 0.e0 , & 41 rone = 1.e0 40 42 41 43 !! * Substitution 42 44 # include "vectopt_loop_substitute.h90" 43 45 !!---------------------------------------------------------------------- 44 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)46 !! LIM 2.0, UCL-LOCEAN-IPSL (2005) 45 47 !! $Id$ 46 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)48 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 47 49 !!---------------------------------------------------------------------- 48 50 … … 60 62 !! 61 63 !! ** action : 64 !! 65 !! History : 66 !! 1.0 ! 00-01 (LIM) Original code 67 !! ! 01-05 (G. Madec, R. Hordoir) opa norm 68 !! 2.0 ! 04-01 (G. Madec, C. Ethe) F90, mpp 62 69 !!--------------------------------------------------------------------- 63 70 INTEGER, INTENT(in) :: kt ! number of iteration 64 !! 65 INTEGER :: ji, jj, jk ! dummy loop indices 66 INTEGER :: initad ! number of sub-timestep for the advection 67 REAL(wp) :: zindb , zacrith, zindsn , zignm , zvbord, zrtt, ztic1, zusnit 68 REAL(wp) :: zindic, zusvosn, zusvoic, zindhe, zcfl , ztsn, ztic2 69 REAL(wp), DIMENSION(jpi,jpj) :: zui_u , zs0sn, zs0c0, zs0a , zsm ! 2D workspace 70 REAL(wp), DIMENSION(jpi,jpj) :: zvi_v , zs0st, zs0c1, zs0c2, zs0ice ! - - 71 72 INTEGER :: ji, jj, jk, & ! dummy loop indices 73 & initad ! number of sub-timestep for the advection 74 75 REAL(wp) :: & 76 zindb , & 77 zacrith, & 78 zindsn , & 79 zindic , & 80 zusvosn, & 81 zusvoic, & 82 zignm , & 83 zindhe , & 84 zvbord , & 85 zcfl , & 86 zusnit , & 87 zrtt, ztsn, ztic1, ztic2 88 89 REAL(wp), DIMENSION(jpi,jpj) :: & ! temporary workspace 90 zui_u , zvi_v , zsm , & 91 zs0ice, zs0sn , zs0a , & 92 zs0c0 , zs0c1 , zs0c2 , & 93 zs0st 71 94 !--------------------------------------------------------------------- 72 95 … … 75 98 zsm(:,:) = area(:,:) 76 99 77 ! !-------------------------------------! 78 IF( ln_limdyn ) THEN ! Advection of sea ice properties ! 79 ! !-------------------------------------! 100 IF( ln_limdyn ) THEN 101 !-------------------------------------! 102 ! Advection of sea ice properties ! 103 !-------------------------------------! 80 104 81 105 ! ice velocities at ocean U- and V-points (zui_u,zvi_v) … … 89 113 END DO 90 114 END DO 91 CALL lbc_lnk( zui_u, 'U', -1. ) ; CALL lbc_lnk( zvi_v, 'V', -1. ) ! Lateral boundary conditions 115 ! Lateral boundary conditions on zui_u, zvi_v 116 CALL lbc_lnk( zui_u, 'U', -1. ) 117 CALL lbc_lnk( zvi_v, 'V', -1. ) 92 118 93 119 ! CFL test for stability … … 97 123 zcfl = MAX( zcfl, MAXVAL( ABS( zvi_v( : ,1:jpjm1) ) * rdt_ice / e2v( : ,1:jpjm1) ) ) 98 124 99 IF ( lk_mpp ) CALL mpp_max( zcfl)100 101 IF ( zcfl > 0.5 .AND. lwp ) WRITE(numout,*) 'lim_trp_2 : CFL violation atthe ',nday,'th day, cfl = ',zcfl125 IF (lk_mpp ) CALL mpp_max(zcfl) 126 127 IF ( zcfl > 0.5 .AND. lwp ) WRITE(numout,*) 'lim_trp_2 : violation of cfl criterion the ',nday,'th day, cfl = ',zcfl 102 128 103 129 ! content of properties … … 118 144 zusnit = 1.0 / REAL( initad ) 119 145 120 IF ( MOD( nday , 2 ) == 0) THEN146 IF ( MOD( nday , 2 ) == 0) THEN 121 147 DO jk = 1,initad 122 148 CALL lim_adv_x_2( zusnit, zui_u, rone , zsm, zs0ice, sxice, sxxice, syice, syyice, sxyice ) … … 202 228 ! Up-dating and limitation of sea ice properties after transport ! 203 229 ! -------------------------------------------------------------------! 230 231 ! Up-dating and limitation of sea ice properties after transport. 204 232 DO jj = 1, jpj 233 !!!iii zindhe = REAL( MAX( 0, isign(1, jj - njeq ) ) ) !ibug mpp !!bugmpp njeq! 205 234 zindhe = MAX( 0.e0, SIGN( 1.e0, fcor(1,jj) ) ) ! = 0 for SH, =1 for NH 206 235 DO ji = 1, jpi 207 ! ! Recover mean values over the grid squares. 236 237 ! Recover mean values over the grid squares. 208 238 zs0sn (ji,jj) = MAX( rzero, zs0sn (ji,jj)/area(ji,jj) ) 209 239 zs0ice(ji,jj) = MAX( rzero, zs0ice(ji,jj)/area(ji,jj) ) … … 213 243 zs0c2 (ji,jj) = MAX( rzero, zs0c2 (ji,jj)/area(ji,jj) ) 214 244 zs0st (ji,jj) = MAX( rzero, zs0st (ji,jj)/area(ji,jj) ) 215 ! ! Recover in situ values. 245 246 ! Recover in situ values. 216 247 zindb = MAX( rzero, SIGN( rone, zs0a(ji,jj) - epsi06 ) ) 217 248 zacrith = 1.0 - ( zindhe * acrit(1) + ( 1.0 - zindhe ) * acrit(2) ) … … 231 262 zrtt = 173.15 * rone 232 263 ztsn = zignm * tbif(ji,jj,1) & 233 &+ ( 1.0 - zignm ) * MIN( MAX( zrtt, rt0_snow * zusvosn * zs0c0(ji,jj)) , tfu(ji,jj) )264 + ( 1.0 - zignm ) * MIN( MAX( zrtt, rt0_snow * zusvosn * zs0c0(ji,jj)) , tfu(ji,jj) ) 234 265 ztic1 = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c1(ji,jj) ) , tfu(ji,jj) ) 235 266 ztic2 = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c2(ji,jj) ) , tfu(ji,jj) ) 236 !267 237 268 tbif(ji,jj,1) = zindsn * ztsn + ( 1.0 - zindsn ) * tfu(ji,jj) 238 269 tbif(ji,jj,2) = zindic * ztic1 + ( 1.0 - zindic ) * tfu(ji,jj) … … 241 272 END DO 242 273 END DO 243 !274 244 275 ENDIF 245 !276 246 277 END SUBROUTINE lim_trp_2 247 278 … … 257 288 !! 258 289 !! ** input : Namelist namicetrp 290 !! 291 !! history : 292 !! 2.0 ! 03-08 (C. Ethe) Original code 259 293 !!------------------------------------------------------------------- 260 294 NAMELIST/namicetrp/ bound 261 295 !!------------------------------------------------------------------- 262 ! 263 REWIND ( numnam_ice ) ! Read Namelist namicetrp 296 297 ! Read Namelist namicetrp 298 REWIND ( numnam_ice ) 264 299 READ ( numnam_ice , namicetrp ) 265 IF(lwp) THEN ! control print300 IF(lwp) THEN 266 301 WRITE(numout,*) 267 302 WRITE(numout,*) 'lim_trp_init_2 : Ice parameters for advection ' … … 269 304 WRITE(numout,*) ' boundary conditions (0. no-slip, 1. free-slip) bound = ', bound 270 305 ENDIF 271 !306 272 307 END SUBROUTINE lim_trp_init_2 273 308 … … 276 311 !! Default option Empty Module No sea-ice model 277 312 !!---------------------------------------------------------------------- 313 CONTAINS 314 SUBROUTINE lim_trp_2 ! Empty routine 315 END SUBROUTINE lim_trp_2 278 316 #endif 279 317 -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limwri_2.F90
r1855 r1857 12 12 !! 'key_lim2' LIM 2.0 sea-ice model 13 13 !!---------------------------------------------------------------------- 14 !! lim_wri_2 : write of the diagnostics variables in ouput file 15 !! lim_wri_init_2 : initialization and namelist read 16 !! lim_wri_state_2 : write for initial state (output.init.nc if ninist=1) or/and in the abort file 14 !!---------------------------------------------------------------------- 15 !! lim_wri_2 : write of the diagnostics variables in ouput file 16 !! lim_wri_init_2 : initialization and namelist read 17 !! lim_wri_state_2 : write for initial state or/and abandon: 18 !! > output.init.nc (if ninist = 1 in namelist) 19 !! > output.abort.nc 17 20 !!---------------------------------------------------------------------- 18 21 USE phycst … … 39 42 INTEGER, PARAMETER :: jpnoumax = 40 ! maximum number of variable for ice output 40 43 INTEGER :: noumef ! number of fields 41 REAL(wp) , DIMENSION(jpnoumax) :: cmulti 42 REAL(wp) , DIMENSION(jpnoumax) ::cadd ! additive constant44 REAL(wp) , DIMENSION(jpnoumax) :: cmulti , & ! multiplicative constant 45 & cadd ! additive constant 43 46 CHARACTER(len = 35), DIMENSION(jpnoumax) :: titn ! title of the field 44 47 CHARACTER(len = 8 ), DIMENSION(jpnoumax) :: nam ! name of the field … … 49 52 INTEGER , DIMENSION( jpij ) :: ndex51 ! ???? 50 53 51 REAL(wp) :: epsi16 = 1.e-16 ! constant values 52 REAL(wp) :: rzero = 0.e0 ! 53 REAL(wp) :: rone = 1.e0 ! 54 REAL(wp) :: & ! constant values 55 epsi16 = 1.e-16 , & 56 zzero = 0.e0 , & 57 zone = 1.e0 54 58 55 59 !! * Substitutions 56 60 # include "vectopt_loop_substitute.h90" 57 61 !!---------------------------------------------------------------------- 58 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)62 !! LIM 2.0, UCL-LOCEAN-IPSL (2006) 59 63 !! $Id$ 60 64 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 77 81 !! 78 82 !! ** Method : computes the average of some variables and write 79 !! it in the NetCDF ouput files 80 !! CAUTION: the sea-ice time-step must be an integer fraction of a day 83 !! it in the NetCDF ouput files 84 !! CAUTION: the sea-ice time-step must be an integer fraction 85 !! of a day 81 86 !!------------------------------------------------------------------- 82 87 INTEGER, INTENT(in) :: kt ! number of iteration 83 88 !! 84 INTEGER :: ji, jj, jf ! dummy loop indices89 INTEGER :: ji, jj, jf ! dummy loop indices 85 90 CHARACTER(len = 40) :: clhstnam, clop 86 REAL(wp) :: zsto, zjulian, zout 87 REAL(wp) ::zindh, zinda, zindb, ztmu91 REAL(wp) :: zsto, zjulian, zout, & ! temporary scalars 92 & zindh, zinda, zindb, ztmu 88 93 REAL(wp), DIMENSION(1) :: zdept 89 94 REAL(wp), DIMENSION(jpi,jpj) :: zfield … … 127 132 DO jj = 2 , jpjm1 128 133 DO ji = 1 , jpim1 ! NO vector opt. 129 zindh = MAX( rzero , SIGN( rone , hicif(ji,jj) * (1.0 - frld(ji,jj) ) - 0.10 ) )130 zinda = MAX( rzero , SIGN( rone , ( 1.0 - frld(ji,jj) ) - 0.10 ) )134 zindh = MAX( zzero , SIGN( zone , hicif(ji,jj) * (1.0 - frld(ji,jj) ) - 0.10 ) ) 135 zinda = MAX( zzero , SIGN( zone , ( 1.0 - frld(ji,jj) ) - 0.10 ) ) 131 136 zindb = zindh * zinda 132 ztmu = MAX( 0.5 * rone , ( tmu(ji,jj) + tmu(ji+1,jj) + tmu(ji,jj+1) + tmu(ji+1,jj+1) ) )137 ztmu = MAX( 0.5 * zone , ( tmu(ji,jj) + tmu(ji+1,jj) + tmu(ji,jj+1) + tmu(ji+1,jj+1) ) ) 133 138 zcmo(ji,jj,1) = hsnif (ji,jj) 134 139 zcmo(ji,jj,2) = hicif (ji,jj) … … 138 143 zcmo(ji,jj,6) = fbif (ji,jj) 139 144 zcmo(ji,jj,7) = zindb * ( u_ice(ji,jj ) * tmu(ji,jj ) + u_ice(ji+1,jj ) * tmu(ji+1,jj ) & 140 &+ u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) &141 &/ ztmu145 + u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 146 / ztmu 142 147 143 148 zcmo(ji,jj,8) = zindb * ( v_ice(ji,jj ) * tmu(ji,jj ) + v_ice(ji+1,jj ) * tmu(ji+1,jj ) & 144 &+ v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) &145 &/ ztmu149 + v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 150 / ztmu 146 151 zcmo(ji,jj,9) = sst_m(ji,jj) 147 152 zcmo(ji,jj,10) = sss_m(ji,jj) … … 195 200 !! ** input : Namelist namicewri 196 201 !!------------------------------------------------------------------- 197 INTEGER :: jf ! dummy loop indices202 INTEGER :: nf ! ??? 198 203 TYPE FIELD 199 204 CHARACTER(len = 35) :: ztitle … … 204 209 REAL :: zcadd 205 210 END TYPE FIELD 206 TYPE(FIELD) :: field_1 , field_2 , field_3 , field_4 , field_5 , field_6 207 TYPE(FIELD) :: field_7 , field_8 , field_9 , field_10, field_11, field_12 208 TYPE(FIELD) :: field_13, field_14, field_15, field_16, field_17, field_18, field_19 211 TYPE(FIELD) :: & 212 field_1 , field_2 , field_3 , field_4 , field_5 , field_6 , & 213 field_7 , field_8 , field_9 , field_10, field_11, field_12, & 214 field_13, field_14, field_15, field_16, field_17, field_18, & 215 field_19 209 216 TYPE(FIELD) , DIMENSION(jpnoumax) :: zfield 210 !! 211 NAMELIST/namiceout/ noumef, & 212 field_1 , field_2 , field_3 , field_4 , field_5 , field_6 , & 213 field_7 , field_8 , field_9 , field_10, field_11, field_12, & 214 field_13, field_14, field_15, field_16, field_17, field_18, field_19 217 218 NAMELIST/namiceout/ noumef, & 219 field_1 , field_2 , field_3 , field_4 , field_5 , field_6 , & 220 field_7 , field_8 , field_9 , field_10, field_11, field_12, & 221 field_13, field_14, field_15, field_16, field_17, field_18, & 222 field_19 215 223 !!------------------------------------------------------------------- 216 224 … … 238 246 zfield(19) = field_19 239 247 240 DO jf = 1, noumef241 titn ( jf) = zfield(jf)%ztitle242 nam ( jf) = zfield(jf)%zname243 uni ( jf) = zfield(jf)%zunit244 nc ( jf) = zfield(jf)%znc245 cmulti( jf) = zfield(jf)%zcmulti246 cadd ( jf) = zfield(jf)%zcadd248 DO nf = 1, noumef 249 titn (nf) = zfield(nf)%ztitle 250 nam (nf) = zfield(nf)%zname 251 uni (nf) = zfield(nf)%zunit 252 nc (nf) = zfield(nf)%znc 253 cmulti(nf) = zfield(nf)%zcmulti 254 cadd (nf) = zfield(nf)%zcadd 247 255 END DO 248 256 … … 254 262 WRITE(numout,*) ' title name unit Saving (1/0) ', & 255 263 & ' multiplicative constant additive constant ' 256 DO jf = 1 , noumef257 WRITE(numout,*) ' ', titn( jf), ' ', nam(jf),' ', uni(jf),' ', nc(jf),' ', cmulti(jf), &258 & ' ', cadd( jf)264 DO nf = 1 , noumef 265 WRITE(numout,*) ' ', titn(nf), ' ', nam(nf),' ', uni(nf),' ', nc(nf),' ', cmulti(nf), & 266 & ' ', cadd(nf) 259 267 END DO 260 268 ENDIF … … 273 281 !! Used to find errors in the initial state or save the last 274 282 !! ocean state in case of abnormal end of a simulation 283 !! 284 !! History : 285 !! 2.0 ! 2009-06 (B. Lemaire) 275 286 !!---------------------------------------------------------------------- 276 287 INTEGER, INTENT( in ) :: kt ! ocean time-step index) -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_3/limtab.F90
r1855 r1857 2 2 !!====================================================================== 3 3 !! *** MODULE limtab *** 4 !! LIM-3 ice model : transform 1D (2D) array to a 2D (1D) array4 !! transform 1D (2D) array to a 2D (1D) table 5 5 !!====================================================================== 6 6 #if defined key_lim3 7 7 !!---------------------------------------------------------------------- 8 !! 'key_lim3' : LIM 2.0sea-ice model8 !! 'key_lim3' LIM3 sea-ice model 9 9 !!---------------------------------------------------------------------- 10 10 !! tab_2d_1d : 2-D to 1-D 11 11 !! tab_1d_2d : 1-D to 2-D 12 12 !!---------------------------------------------------------------------- 13 !! * Modules used 13 14 USE par_kind 14 15 … … 16 17 PRIVATE 17 18 18 PUBLIC tab_2d_1d ! called by lim_thd 19 PUBLIC tab_1d_2d ! called by lim_thd 19 !! * Routine accessibility 20 PUBLIC tab_2d_1d ! called by lim_ther 21 PUBLIC tab_1d_2d ! called by lim_ther 20 22 21 23 !!---------------------------------------------------------------------- 22 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010)24 !! LIM 3.0, UCL-ASTR-LOCEAN-IPSL (2008) 23 25 !! $Id$ 24 26 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 26 28 CONTAINS 27 29 28 SUBROUTINE tab_2d_1d( kdim1d, ptab1d, ptab2d, kdim2d_i, kdim2d_j, ptab_ind ) 29 !!------------------------------------------------------------------- 30 INTEGER , INTENT(in ) :: kdim1d 31 INTEGER , INTENT(in ) :: kdim2d_i, kdim2d_j 32 REAL(wp), INTENT(in ), DIMENSION(kdim2d_i, kdim2d_j) :: ptab2d 33 INTEGER , INTENT(in ), DIMENSION(kdim1d) :: ptab_ind 34 REAL(wp), INTENT( out), DIMENSION(kdim1d) :: ptab1d 35 !! 36 INTEGER :: jn , jid, jjd ! dummy loop indices 37 !!------------------------------------------------------------------- 38 ! 39 DO jn = 1, kdim1d 40 jid = MOD( ptab_ind(jn) - 1, kdim2d_i ) + 1 41 jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 42 ptab1d(jn) = ptab2d(jid,jjd) 43 END DO 44 ! 30 SUBROUTINE tab_2d_1d ( ndim1d, tab1d, tab2d, ndim2d_x, ndim2d_y, tab_ind ) 31 32 INTEGER, INTENT(in) :: & 33 ndim1d, ndim2d_x, ndim2d_y 34 35 REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT(in) :: & 36 tab2d 37 38 INTEGER, DIMENSION ( ndim1d), INTENT ( in) :: & 39 tab_ind 40 41 REAL(wp), DIMENSION(ndim1d), INTENT ( out) :: & 42 tab1d 43 44 INTEGER :: & 45 jn , jid, jjd 46 47 DO jn = 1, ndim1d 48 jid = MOD( tab_ind(jn) - 1, ndim2d_x ) + 1 49 jjd = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 50 tab1d( jn) = tab2d( jid, jjd) 51 END DO 52 45 53 END SUBROUTINE tab_2d_1d 46 54 47 55 48 SUBROUTINE tab_1d_2d( kdim1d, ptab2d, ptab_ind, ptab1d, kdim2d_i, kdim2d_j ) 49 !!------------------------------------------------------------------- 50 INTEGER , INTENT(in ) :: kdim1d 51 INTEGER , INTENT(in ) :: kdim2d_i, kdim2d_j 52 INTEGER , INTENT(in ), DIMENSION(kdim1d) :: ptab_ind 53 REAL(wp), INTENT(in ), DIMENSION(kdim1d) :: ptab1d 54 REAL(wp), INTENT( out), DIMENSION(kdim2d_i, kdim2d_j) :: ptab2d 55 !! 56 INTEGER :: jn, jid, jjd ! dummy loop indices 57 !!------------------------------------------------------------------- 58 ! 59 DO jn = 1, kdim1d 60 jid = MOD( ptab_ind(jn) - 1, kdim2d_i) + 1 61 jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 62 ptab2d(jid,jjd) = ptab1d(jn) 56 SUBROUTINE tab_1d_2d ( ndim1d, tab2d, tab_ind, tab1d, ndim2d_x, ndim2d_y ) 57 58 INTEGER, INTENT ( in) :: & 59 ndim1d, ndim2d_x, ndim2d_y 60 61 INTEGER, DIMENSION (ndim1d) , INTENT (in) :: & 62 tab_ind 63 64 REAL(wp), DIMENSION(ndim1d), INTENT (in) :: & 65 tab1d 66 67 REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT ( out) :: & 68 tab2d 69 70 INTEGER :: & 71 jn, jid, jjd 72 73 DO jn = 1, ndim1d 74 jid = MOD( tab_ind(jn) - 1, ndim2d_x) + 1 75 jjd = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 76 tab2d(jid, jjd) = tab1d( jn) 63 77 END DO 64 ! 78 65 79 END SUBROUTINE tab_1d_2d 66 80 67 #else68 !!----------------------------------------------------------------------69 !! Default option Dummy module NO LIM 2.0 sea-ice model70 !!----------------------------------------------------------------------71 81 #endif 72 73 !!======================================================================74 82 END MODULE limtab -
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/OPA_SRC/SBC/sbc_ice.F90
r1855 r1857 41 41 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl) :: alb_ice !: albedo of ice 42 42 43 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau_ice !: u-stress over ice (I-point for LIM2 or U,V-point for LIM3) [N/m2] 44 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: vtau_ice !: v-stress over ice (I-point for LIM2 or U,V-point for LIM3) [N/m2] 43 45 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr1_i0 !: 1st fraction of sol. rad. which penetrate inside the ice cover 44 46 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr2_i0 !: 2nd fraction of sol. rad. which penetrate inside the ice cover 45 47 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: emp_ice !: solid freshwater budget over ice: sublivation - snow 46 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau_ice !: u-stress over ice (I-point for LIM2 or U,V-point for LIM3) [N/m2]47 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: vtau_ice !: v-stress over ice (I-point for LIM2 or U,V-point for LIM3) [N/m2]48 48 49 49 # if defined key_lim3
Note: See TracChangeset
for help on using the changeset viewer.