Changeset 1608 for trunk/NEMO/LIM_SRC_3/limmsh.F90
- Timestamp:
- 2009-08-12T10:05:30+02:00 (15 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/LIM_SRC_3/limmsh.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/LIM_SRC_3/limmsh.F90
r1156 r1608 4 4 !! LIM ice model : definition of the ice mesh parameters 5 5 !!====================================================================== 6 !! History : 3.2 ! 2008-01 (NEMO team) LIM-3: adaptation from LIM-2 7 !!---------------------------------------------------------------------- 6 8 #if defined key_lim3 7 9 !!---------------------------------------------------------------------- … … 10 12 !! lim_msh : definition of the ice mesh 11 13 !!---------------------------------------------------------------------- 12 !! * Modules used 13 USE phycst 14 USE dom_oce 15 USE dom_ice 16 USE lbclnk 17 USE in_out_manager 14 USE phycst ! physical constants 15 USE dom_oce ! ocean domain 16 USE dom_ice ! sea-ice domain 17 USE in_out_manager ! I/O manager 18 USE lbclnk ! 18 19 19 20 IMPLICIT NONE 20 21 PRIVATE 21 22 22 !! * Accessibility 23 PUBLIC lim_msh ! routine called by ice_ini.F90 23 PUBLIC lim_msh ! routine called by ice_ini.F90 24 24 25 25 !!---------------------------------------------------------------------- 26 !! LIM 3.0, UCL-ASTR-LOCEAN-IPSL (2008)26 !! NEMO/LIM 3.2, UCL-ASTR-LOCEAN-IPSL (2009) 27 27 !! $Id$ 28 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt28 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 29 29 !!---------------------------------------------------------------------- 30 30 … … 42 42 !! - Initialization of the ice masks (tmsk, umsk) 43 43 !! 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) 44 !! Reference : Deleersnijder et al. Ocean Modelling 100, 7-10 49 45 !!--------------------------------------------------------------------- 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 46 INTEGER :: ji, jj ! dummy loop indices 47 REAL(wp) :: zusden ! temporary scalar 60 48 !!--------------------------------------------------------------------- 61 49 62 50 IF(lwp) THEN 63 51 WRITE(numout,*) 64 WRITE(numout,*) 'lim_msh : LIM sea-ice model, mesh initialization'52 WRITE(numout,*) 'lim_msh : LIM-3 sea-ice model, mesh initialization' 65 53 WRITE(numout,*) '~~~~~~~' 66 54 ENDIF 67 55 68 !---------------------------------------------------------- 69 ! Initialization of local and some global (common) variables 70 !------------------------------------------------------------------ 71 72 njeq = INT( jpj / 2 ) !i bug mpp potentiel 56 ! !== coriolis factor & Equator position ==! 57 njeq = INT( jpj / 2 ) 73 58 njeqm1 = njeq - 1 74 75 fcor(:,:) = 2. * omega * SIN( gphit(:,:) * rad ) ! 76 59 ! 60 fcor(:,:) = 2. * omega * SIN( gphit(:,:) * rad ) ! coriolis factor 61 ! 77 62 IF( fcor(1,1) * fcor(1,nlcj) < 0.e0 ) THEN ! local domain include both hemisphere 78 63 l_jeq = .TRUE. … … 91 76 IF(lwp ) WRITE(numout,*) ' the model domain is entirely in the northern hemisphere: njeq = ', njeq 92 77 ENDIF 93 78 ! 94 79 njeqm1 = njeq - 1 95 80 96 81 97 ! For each grid, definition of geometric tables 98 !------------------------------------------------------------------ 99 100 !------------------- 101 ! Conventions : ! 102 !------------------- 103 ! indices 1 \ 2 <-> localisation in the 2 direction x \ y 104 ! 3rd indice <-> localisation on the mesh : 105 ! 0 = Centre ; 1 = corner W x(i-1/2) ; 2 = corner S y(j-1/2) ; 106 ! 3 = corner SW x(i-1/2),y(j-1/2) 107 !------------------- 108 !!ibug ??? 109 akappa(:,:,:,:) = 0.e0 82 ! !== metric coefficients for sea ice dynamic ==! 110 83 wght(:,:,:,:) = 0.e0 111 alambd(:,:,:,:,:,:) = 0.e0 112 tmu(:,:) = 0.e0 113 tmv(:,:) = 0.0e0 ! CGrid EVP 114 !!i 115 ! metric coefficients for sea ice dynamic 116 !---------------------------------------- 117 ! ! akappa 118 DO jj = 2, jpj 119 zd1d2(:,jj) = e1v(:,jj) - e1v(:,jj-1) 120 END DO 121 CALL lbc_lnk( zd1d2, 'T', -1. ) 122 123 DO ji = 2, jpi 124 zd2d1(ji,:) = e2u(ji,:) - e2u(ji-1,:) 125 END DO 126 CALL lbc_lnk( zd2d1, 'T', -1. ) 127 128 akappa(:,:,1,1) = 1.0 / ( 2.0 * e1t(:,:) ) 129 akappa(:,:,1,2) = zd1d2(:,:) / ( 4.0 * e1t(:,:) * e2t(:,:) ) 130 akappa(:,:,2,1) = zd2d1(:,:) / ( 4.0 * e1t(:,:) * e2t(:,:) ) 131 akappa(:,:,2,2) = 1.0 / ( 2.0 * e2t(:,:) ) 132 133 ! ! weights (wght) 84 !!gm Optimisation : wght to be defined at F-point, not I-point and change in limrhg 134 85 DO jj = 2, jpj 135 86 DO ji = 2, jpi 136 zusden = 1. / ( ( e1t(ji,jj) + e1t(ji-1,jj ) ) &137 & * ( e2t(ji,jj) + e2t(ji ,jj-1) ) )87 zusden = 1.e0 / ( ( e1t(ji,jj) + e1t(ji-1,jj ) ) & 88 & * ( e2t(ji,jj) + e2t(ji ,jj-1) ) ) 138 89 wght(ji,jj,1,1) = zusden * e1t(ji ,jj) * e2t(ji,jj ) 139 90 wght(ji,jj,1,2) = zusden * e1t(ji ,jj) * e2t(ji,jj-1) … … 146 97 CALL lbc_lnk( wght(:,:,2,1), 'I', 1. ) ! but it is never used 147 98 CALL lbc_lnk( wght(:,:,2,2), 'I', 1. ) 99 !!gm end 148 100 149 ! Coefficients for divergence of the stress tensor 150 !------------------------------------------------- 151 152 DO jj = 2, jpj 153 DO ji = 2, jpi 154 zh1p = e1t(ji ,jj ) * wght(ji,jj,2,2) & 155 & + e1t(ji-1,jj ) * wght(ji,jj,1,2) & 156 & + e1t(ji ,jj-1) * wght(ji,jj,2,1) & 157 & + e1t(ji-1,jj-1) * wght(ji,jj,1,1) 158 159 zh2p = e2t(ji ,jj ) * wght(ji,jj,2,2) & 160 & + e2t(ji-1,jj ) * wght(ji,jj,1,2) & 161 & + e2t(ji ,jj-1) * wght(ji,jj,2,1) & 162 & + e2t(ji-1,jj-1) * wght(ji,jj,1,1) 163 164 zusden = 1.0 / MAX( zh1p * zh2p * 4.e0 , 1.e-20 ) 165 zusden2 = zusden * 2.0 166 167 zd1d2p = zusden * 0.5 * ( -e1t(ji-1,jj-1) + e1t(ji-1,jj ) - e1t(ji,jj-1) + e1t(ji ,jj) ) 168 zd2d1p = zusden * 0.5 * ( e2t(ji ,jj-1) - e2t(ji-1,jj-1) + e2t(ji,jj ) - e2t(ji-1,jj) ) 169 170 alambd(ji,jj,2,2,2,1) = zusden2 * e2t(ji ,jj-1) 171 alambd(ji,jj,2,2,2,2) = zusden2 * e2t(ji ,jj ) 172 alambd(ji,jj,2,2,1,1) = zusden2 * e2t(ji-1,jj-1) 173 alambd(ji,jj,2,2,1,2) = zusden2 * e2t(ji-1,jj ) 174 175 alambd(ji,jj,1,1,2,1) = zusden2 * e1t(ji ,jj-1) 176 alambd(ji,jj,1,1,2,2) = zusden2 * e1t(ji ,jj ) 177 alambd(ji,jj,1,1,1,1) = zusden2 * e1t(ji-1,jj-1) 178 alambd(ji,jj,1,1,1,2) = zusden2 * e1t(ji-1,jj ) 179 180 alambd(ji,jj,1,2,2,1) = zd1d2p 181 alambd(ji,jj,1,2,2,2) = zd1d2p 182 alambd(ji,jj,1,2,1,1) = zd1d2p 183 alambd(ji,jj,1,2,1,2) = zd1d2p 184 185 alambd(ji,jj,2,1,2,1) = zd2d1p 186 alambd(ji,jj,2,1,2,2) = zd2d1p 187 alambd(ji,jj,2,1,1,1) = zd2d1p 188 alambd(ji,jj,2,1,1,2) = zd2d1p 101 ! !== ice masks ==! 102 tms(:,:) = tmask(:,:,1) ! ice T-point : use surface tmask 103 tmu(:,:) = umask(:,:,1) ! ice U-point : use surface umask (C-grid EVP) 104 tmv(:,:) = vmask(:,:,1) ! ice V-point : use surface vmask (C-grid EVP) 105 DO jj = 1, jpjm1 ! ice F-point : recompute fmask (due to nn_shlat) 106 DO ji = 1 , jpim1 107 tmf(ji,jj) = tms(ji,jj) * tms(ji+1,jj) * tms(ji,jj+1) * tms(ji+1,jj+1) 189 108 END DO 190 109 END DO 110 CALL lbc_lnk( tmf(:,:), 'F', 1. ) ! lateral boundary conditions 191 111 192 CALL lbc_lnk( alambd(:,:,2,2,2,1), 'I', 1. ) ! CAUTION: even with the lbc_lnk at ice U-V point 193 CALL lbc_lnk( alambd(:,:,2,2,2,2), 'I', 1. ) ! the value of wght at jpj is wrong 194 CALL lbc_lnk( alambd(:,:,2,2,1,1), 'I', 1. ) ! but it is never used 195 CALL lbc_lnk( alambd(:,:,2,2,1,2), 'I', 1. ) ! 196 197 CALL lbc_lnk( alambd(:,:,1,1,2,1), 'I', 1. ) ! CAUTION: idem 198 CALL lbc_lnk( alambd(:,:,1,1,2,2), 'I', 1. ) ! 199 CALL lbc_lnk( alambd(:,:,1,1,1,1), 'I', 1. ) ! 200 CALL lbc_lnk( alambd(:,:,1,1,1,2), 'I', 1. ) ! 201 202 CALL lbc_lnk( alambd(:,:,1,2,2,1), 'I', 1. ) ! CAUTION: idem 203 CALL lbc_lnk( alambd(:,:,1,2,2,2), 'I', 1. ) ! 204 CALL lbc_lnk( alambd(:,:,1,2,1,1), 'I', 1. ) ! 205 CALL lbc_lnk( alambd(:,:,1,2,1,2), 'I', 1. ) ! 206 207 CALL lbc_lnk( alambd(:,:,2,1,2,1), 'I', 1. ) ! CAUTION: idem 208 CALL lbc_lnk( alambd(:,:,2,1,2,2), 'I', 1. ) ! 209 CALL lbc_lnk( alambd(:,:,2,1,1,1), 'I', 1. ) ! 210 CALL lbc_lnk( alambd(:,:,2,1,1,2), 'I', 1. ) ! 211 212 213 ! Initialization of ice masks 214 !---------------------------- 215 216 tms(:,:) = tmask(:,:,1) ! ice T-point : use surface tmask 217 218 ! tmu(:,1) = 0.e0 219 ! tmu(1,:) = 0.e0 220 ! tmv(:,1) = 0.e0 221 ! tmv(1,:) = 0.e0 222 223 DO jj = 1, jpj - 1 224 DO ji = 1 , jpi - 1 225 tmu(ji,jj) = tms(ji,jj) * tms(ji+1,jj) 226 tmv(ji,jj) = tms(ji,jj) * tms(ji,jj+1) 227 tmf(ji,jj) = tms(ji,jj) * tms(ji+1,jj) * tms(ji,jj+1) * & 228 tms(ji+1,jj+1) 229 END DO 230 END DO 231 232 !--lateral boundary conditions 233 CALL lbc_lnk( tmu(:,:), 'U', 1. ) 234 CALL lbc_lnk( tmv(:,:), 'V', 1. ) 235 CALL lbc_lnk( tmf(:,:), 'F', 1. ) 236 237 ! unmasked and masked area of T-grid cell 112 ! !== unmasked and masked area of T-grid cell 238 113 area(:,:) = e1t(:,:) * e2t(:,:) 239 114 ! 240 115 END SUBROUTINE lim_msh 241 116
Note: See TracChangeset
for help on using the changeset viewer.
