Changeset 76 for trunk/NEMO/LIM_SRC/limdyn.F90
- Timestamp:
- 2004-04-22T14:45:08+02:00 (20 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/LIM_SRC/limdyn.F90
r12 r76 20 20 USE iceini 21 21 USE limistate 22 USE limrhg ! ice rheology22 USE limrhg ! ice rheology 23 23 USE lbclnk 24 USE lib_mpp 24 25 25 26 IMPLICIT NONE … … 30 31 31 32 !! * Module variables 32 REAL(wp) :: rone = 1. 0 ! constant value33 REAL(wp) :: rone = 1.e0 ! constant value 33 34 34 35 !!---------------------------------------------------------------------- … … 55 56 !!--------------------------------------------------------------------- 56 57 !! * Loal variables 57 INTEGER :: ji, jj , &! dummy loop indices58 jhemis ! jhemis = 1 (NH) ; jhemis = -1 (SH)58 INTEGER :: ji, jj ! dummy loop indices 59 INTEGER :: i_j1, i_jpj ! Starting/ending j-indices for rheology 59 60 REAL(wp) :: & 60 ztairx, ztairy, & ! tempory scalars61 zsang , zmod, &62 ztglx , ztgly , &63 zt11, zt12, zt21, zt22 , &64 zustm, zsfrld, zsfrldm4, &61 ztairx, ztairy, & ! tempory scalars 62 zsang , zmod, & 63 ztglx , ztgly , & 64 zt11, zt12, zt21, zt22 , & 65 zustm, zsfrld, zsfrldm4, & 65 66 zu_ice, zv_ice, ztair2 66 67 REAL(wp),DIMENSION(jpj) :: & 68 zind, & ! i-averaged indicator of sea-ice 69 zmsk ! i-averaged of tmask 67 70 !!--------------------------------------------------------------------- 68 69 71 70 72 IF( numit == nstart ) CALL lim_dyn_init ! Initialization (first time-step only) 71 73 72 IF ( l dyn ) THEN74 IF ( ln_limdyn ) THEN 73 75 74 76 ! Mean ice and snow thicknesses. … … 80 82 81 83 ! ! Rheology (ice dynamics) 82 !-- Northern hemisphere ! ======== 83 jhemis = +1 84 CALL lim_rhg( jhemis ) 85 86 !-- Southern hemisphere. 87 jhemis = -1 88 CALL lim_rhg( jhemis ) 89 90 u_ice(:,1) = 0.0 !ibug est-ce vraiment necessaire? 91 v_ice(:,1) = 0.0 92 93 IF( l_ctl .AND. lwp ) THEN 84 ! ! ======== 85 86 ! Define the j-limits where ice rheology is computed 87 ! --------------------------------------------------- 88 89 IF( lk_mpp ) THEN ! mpp: compute over the whole domain 90 i_j1 = 1 91 i_jpj = jpj 92 IF(l_ctl) WRITE(numout,*) 'lim_rhg : i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 93 CALL lim_rhg( i_j1, i_jpj ) 94 95 ELSE ! optimization of the computational area 96 97 DO jj = 1, jpj 98 zind(jj) = SUM( frld (:,jj ) ) ! = FLOAT(jpj) if ocean everywhere on a j-line 99 zmsk(jj) = SUM( tmask(:,jj,1) ) ! = 0 if land everywhere on a j-line 100 !!i write(numout,*) narea, 'limrhg' , jj, zind(jj), zmsk(jj) 101 END DO 102 103 IF( l_jeq ) THEN ! local domain include both hemisphere 104 ! ! Rheology is computed in each hemisphere 105 ! ! only over the ice cover latitude strip 106 ! Northern hemisphere 107 i_j1 = njeq 108 i_jpj = jpj 109 DO WHILE ( i_j1 <= jpj .AND. zind(i_j1) == FLOAT(jpi) .AND. zmsk(i_j1) /=0 ) 110 i_j1 = i_j1 + 1 111 END DO 112 i_j1 = MAX( 1, i_j1-1 ) 113 IF(l_ctl) WRITE(numout,*) 'lim_rhg : NH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 114 115 CALL lim_rhg( i_j1, i_jpj ) 116 117 ! Southern hemisphere 118 i_j1 = 1 119 i_jpj = njeq 120 DO WHILE ( i_jpj >= 1 .AND. zind(i_jpj) == FLOAT(jpi) .AND. zmsk(i_jpj) /=0 ) 121 i_jpj = i_jpj - 1 122 END DO 123 i_jpj = MIN( jpj, i_jpj+2 ) 124 IF(l_ctl) WRITE(numout,*) 'lim_rhg : SH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 125 126 CALL lim_rhg( i_j1, i_jpj ) 127 128 ELSE ! local domain extends over one hemisphere only 129 ! ! Rheology is computed only over the ice cover 130 ! ! latitude strip 131 i_j1 = 1 132 DO WHILE ( i_j1 <= jpj .AND. zind(i_j1) == FLOAT(jpi) .AND. zmsk(i_j1) /=0 ) 133 i_j1 = i_j1 + 1 134 END DO 135 i_j1 = MAX( 1, i_j1-1 ) 136 137 i_jpj = jpj 138 DO WHILE ( i_jpj >= 1 .AND. zind(i_jpj) == FLOAT(jpi) .AND. zmsk(i_jpj) /=0 ) 139 i_jpj = i_jpj - 1 140 END DO 141 i_jpj = MIN( jpj, i_jpj+2) 142 143 IF(l_ctl) WRITE(numout,*) 'lim_rhg : one hemisphere: i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 144 145 CALL lim_rhg( i_j1, i_jpj ) 146 147 ENDIF 148 149 ENDIF 150 151 u_ice(:,1) = 0.e0 !ibug est-ce vraiment necessaire? 152 v_ice(:,1) = 0.e0 153 154 IF(l_ctl) THEN 94 155 WRITE(numout,*) ' lim_dyn : u_oce ', SUM( u_oce ), ' v_oce ', SUM( v_oce ) 95 156 WRITE(numout,*) ' lim_dyn : u_ice ', SUM( u_ice ), ' v_ice ', SUM( v_ice ) … … 99 160 ! ! ================ 100 161 DO jj = 2, jpjm1 101 jhemis = SIGN(1, jj - jeq ) 102 zsang = REAL(jhemis) * sangvg 162 zsang = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 103 163 DO ji = 2, jpim1 104 164 ! computation of wind stress over ocean in X and Y direction … … 155 215 END DO 156 216 157 ELSE ! If no ice dynamics217 ELSE ! no ice dynamics : transmit directly the atmospheric stress to the ocean 158 218 159 219 DO jj = 2, jpjm1 … … 182 242 CALL lbc_lnk( tio_v, 'I', -1. ) ! I-point (i.e. ice U-V point) 183 243 184 IF( l_ctl .AND. lwp) THEN244 IF(l_ctl) THEN 185 245 WRITE(numout,*) ' lim_dyn : tio_u ', SUM( tio_u ), ' tio_v ', SUM( tio_v ) 186 246 WRITE(numout,*) ' lim_dyn : ust2s ', SUM( ust2s ) … … 189 249 END SUBROUTINE lim_dyn 190 250 191 SUBROUTINE lim_dyn_init 251 252 SUBROUTINE lim_dyn_init 192 253 !!------------------------------------------------------------------- 193 254 !! *** ROUTINE lim_dyn_init *** … … 219 280 WRITE(numout,*) 'lim_dyn_init : ice parameters for ice dynamics ' 220 281 WRITE(numout,*) '~~~~~~~~~~~~' 221 WRITE(numout,*) ' tolerance parameter epsd = ', epsd222 WRITE(numout,*) ' coefficient for semi-implicit coriolis alpha = ', alpha223 WRITE(numout,*) ' diffusion constant for dynamics dm = ', dm224 WRITE(numout,*) ' number of sub-time steps for relaxation nbiter = ', nbiter225 WRITE(numout,*) ' maximum number of iterations for relaxation nbitdr = ', nbitdr226 WRITE(numout,*) ' relaxation constant om = ', om227 WRITE(numout,*) ' maximum value for the residual of relaxation resl = ', resl228 WRITE(numout,*) ' drag coefficient for oceanic stress cw = ', cw229 WRITE(numout,*) ' turning angle for oceanic stress angvg = ', angvg230 WRITE(numout,*) ' first bulk-rheology parameter pstar = ', pstar231 WRITE(numout,*) ' second bulk-rhelogy parameter c_rhg = ', c_rhg232 WRITE(numout,*) ' minimun value for viscosity etamn = ', etamn233 WRITE(numout,*) ' creep limit creepl = ', creepl234 WRITE(numout,*) ' eccentricity of the elliptical yield curve ecc = ', ecc235 WRITE(numout,*) ' horizontal diffusivity coeff. for sea-ice ahi0 = ', ahi0282 WRITE(numout,*) ' tolerance parameter epsd = ', epsd 283 WRITE(numout,*) ' coefficient for semi-implicit coriolis alpha = ', alpha 284 WRITE(numout,*) ' diffusion constant for dynamics dm = ', dm 285 WRITE(numout,*) ' number of sub-time steps for relaxation nbiter = ', nbiter 286 WRITE(numout,*) ' maximum number of iterations for relaxation nbitdr = ', nbitdr 287 WRITE(numout,*) ' relaxation constant om = ', om 288 WRITE(numout,*) ' maximum value for the residual of relaxation resl = ', resl 289 WRITE(numout,*) ' drag coefficient for oceanic stress cw = ', cw 290 WRITE(numout,*) ' turning angle for oceanic stress angvg = ', angvg 291 WRITE(numout,*) ' first bulk-rheology parameter pstar = ', pstar 292 WRITE(numout,*) ' second bulk-rhelogy parameter c_rhg = ', c_rhg 293 WRITE(numout,*) ' minimun value for viscosity etamn = ', etamn 294 WRITE(numout,*) ' creep limit creepl = ', creepl 295 WRITE(numout,*) ' eccentricity of the elliptical yield curve ecc = ', ecc 296 WRITE(numout,*) ' horizontal diffusivity coeff. for sea-ice ahi0 = ', ahi0 236 297 ENDIF 237 298
Note: See TracChangeset
for help on using the changeset viewer.