[41] | 1 | !> \file dragging_neff_contmaj_mod.f90 |
---|
| 2 | !! Defini les zones de stream avec 3 criteres |
---|
| 3 | !< |
---|
| 4 | |
---|
| 5 | !> \namespace dragging_neff_contmaj |
---|
| 6 | !! Defini les zones de stream avec 3 criteres |
---|
| 7 | !! \author ... |
---|
| 8 | !! \date ... |
---|
| 9 | !! @note Les trois criteres sont: |
---|
| 10 | !! @note * un critere sur la pression effective |
---|
| 11 | !! @note * un critere de continuite depuis la cote |
---|
| 12 | !! @note * un critere sur la courbure du socle (si negatif, vallees) |
---|
| 13 | !! @note Used module |
---|
| 14 | !! @note - use module3D_phy |
---|
| 15 | !! @note - use param_phy_mod |
---|
| 16 | !< |
---|
| 17 | |
---|
[59] | 18 | module dragging_neff_slope |
---|
[41] | 19 | |
---|
| 20 | ! Defini les zones de stream avec : |
---|
| 21 | ! * un critere sur la pression effective |
---|
| 22 | ! * un critere de continuite depuis la cote |
---|
| 23 | ! * un critere sur la courbure du socle (si negatif, vallees) |
---|
| 24 | |
---|
| 25 | use module3d_phy |
---|
| 26 | use param_phy_mod |
---|
[59] | 27 | use interface_input |
---|
| 28 | use io_netcdf_grisli |
---|
[168] | 29 | use fake_beta_iter_vitbil_mod |
---|
[41] | 30 | |
---|
| 31 | implicit none |
---|
[66] | 32 | |
---|
[41] | 33 | logical,dimension(nx,ny) :: fleuve |
---|
| 34 | logical,dimension(nx,ny) :: cote |
---|
[66] | 35 | logical,dimension(nx,ny) :: slowssamx ! not actual stream, but ssa used as Ub |
---|
| 36 | logical,dimension(nx,ny) :: slowssamy ! not actual stream, but ssa used as Ub |
---|
| 37 | logical,dimension(nx,ny) :: slowssa ! not actual stream, but ssa used as Ub |
---|
[41] | 38 | |
---|
[59] | 39 | real,dimension(nx,ny) :: hires_slope ! slope comupted on a high resolution topography |
---|
| 40 | character(len=100) :: slope_fich ! fichier grille |
---|
| 41 | character(len=100) :: file_ncdf !< fichier netcdf issue des fichiers .dat |
---|
[41] | 42 | |
---|
| 43 | real :: valmax |
---|
| 44 | integer :: imax,jmax |
---|
| 45 | integer :: i_moins1,i_plus1,j_moins1,j_plus1 |
---|
| 46 | integer :: lmax=20 |
---|
| 47 | integer :: idep,jdep,iloc,jloc |
---|
| 48 | |
---|
| 49 | |
---|
| 50 | real :: tostick ! pour la glace posee |
---|
| 51 | real :: betamin ! betamin : (Pa) frottement mini sous les streams |
---|
| 52 | real :: tob_ile ! pour les iles |
---|
| 53 | real :: cry_lim=50. ! courbure limite pour le suivi des fleuves |
---|
| 54 | |
---|
| 55 | real,dimension(nx,ny) :: neff ! pression effective noeuds majeurs |
---|
| 56 | real :: seuil_neff ! seuil sur la pression effective pour avoir glissement en zone sediment |
---|
| 57 | real :: coef_gz ! coef frottement zones stream std (10) |
---|
| 58 | real :: coef_ile ! coef frottement zones iles (0.1) |
---|
| 59 | |
---|
| 60 | real, dimension(nx,ny) :: Vcol_x !< uniquement pour compatibilite avec spinup cat |
---|
| 61 | real, dimension(nx,ny) :: Vcol_y !< uniquement pour compatibilite avec spinup cat |
---|
| 62 | real, dimension(nx,ny) :: Vsl_x !< uniquement pour compatibilite avec spinup cat |
---|
| 63 | real, dimension(nx,ny) :: Vsl_y !< uniquement pour compatibilite avec spinup cat |
---|
| 64 | logical :: corr_def = .false. !< for deformation correction, pour compatibilite beta |
---|
| 65 | |
---|
| 66 | contains |
---|
| 67 | !------------------------------------------------------------------------------- |
---|
| 68 | |
---|
| 69 | !> SUBROUTINE: init_dragging |
---|
| 70 | !! Cette routine fait l'initialisation du dragging. |
---|
| 71 | !> |
---|
| 72 | subroutine init_dragging ! Cette routine fait l'initialisation du dragging. |
---|
| 73 | |
---|
| 74 | implicit none |
---|
| 75 | |
---|
[66] | 76 | ! local variables, defining the rugosity-enhanced flow |
---|
| 77 | real :: expo_slope |
---|
| 78 | real :: pente_min, pente_max |
---|
[41] | 79 | |
---|
[66] | 80 | namelist/drag_neff_slope/cf,betamax,betamin,toblim,tostick,seuil_neff,coef_gz,coef_ile,slope_fich,expo_slope,pente_min,pente_max |
---|
[41] | 81 | |
---|
[66] | 82 | if (itracebug.eq.1) call tracebug(' dragging avec neff et slope') |
---|
| 83 | |
---|
| 84 | |
---|
[41] | 85 | ! formats pour les ecritures dans 42 |
---|
| 86 | 428 format(A) |
---|
| 87 | |
---|
| 88 | ! lecture des parametres du run block drag neff |
---|
| 89 | !-------------------------------------------------------------------- |
---|
| 90 | |
---|
| 91 | rewind(num_param) ! pour revenir au debut du fichier param_list.dat |
---|
[59] | 92 | read(num_param,drag_neff_slope) |
---|
[41] | 93 | |
---|
| 94 | write(num_rep_42,428)'!___________________________________________________________' |
---|
[66] | 95 | write(num_rep_42,428) '&drag_neff_slope ! nom du bloc dragging neff slope' |
---|
[41] | 96 | write(num_rep_42,*) |
---|
[66] | 97 | write(num_rep_42,*) 'cf = ',cf |
---|
| 98 | write(num_rep_42,*) 'betamax = ', betamax |
---|
| 99 | write(num_rep_42,*) 'betamin = ', betamin |
---|
| 100 | write(num_rep_42,*) 'toblim = ', toblim |
---|
[59] | 101 | write(num_rep_42,*) 'tostick = ', tostick |
---|
[66] | 102 | write(num_rep_42,*) 'seuil_neff = ', seuil_neff |
---|
| 103 | write(num_rep_42,*) 'coef_gz = ', coef_gz |
---|
| 104 | write(num_rep_42,*) 'coef_ile = ', coef_ile |
---|
| 105 | write(num_rep_42,'(A,A)') 'slope_fich = ', slope_fich |
---|
| 106 | write(num_rep_42,*) 'expo_slope = ', expo_slope |
---|
| 107 | write(num_rep_42,*) 'pente_min = ', pente_min |
---|
| 108 | write(num_rep_42,*) 'pente_max = ', pente_max |
---|
[41] | 109 | write(num_rep_42,*)'/' |
---|
| 110 | write(num_rep_42,428) '! cf coefficient de la loi de frottement fonction Neff' |
---|
| 111 | write(num_rep_42,428) '! betamax : (Pa) frottement maxi sous les streams ' |
---|
| 112 | write(num_rep_42,428) '! betamin : (Pa) frottement mini sous les streams ' |
---|
| 113 | write(num_rep_42,428) '! toblim : (Pa) pour les iles ' |
---|
[59] | 114 | write(num_rep_42,428) '! tostick : (Pa) pour les points non flgzmx ' |
---|
[41] | 115 | write(num_rep_42,428) '! seuil_neff (Pa) seuil sur la pression effective pour avoir glissement' |
---|
| 116 | write(num_rep_42,428) '! coef_gz : coef frottement zones stream std' |
---|
| 117 | write(num_rep_42,428) '! coef_ile : coef frottement zones iles' |
---|
[66] | 118 | write(num_rep_42,428) '! slope_fich : fichier pente bedrock' |
---|
| 119 | write(num_rep_42,428) '! expo_slope : exposant fonction pente' |
---|
| 120 | write(num_rep_42,428) '! pente_min : no flow enhancement below this' |
---|
| 121 | write(num_rep_42,428) '! pente_max : maximum flow enhancement above this' |
---|
[41] | 122 | write(num_rep_42,*) |
---|
| 123 | |
---|
[168] | 124 | inv_beta=0 |
---|
[41] | 125 | tob_ile=betamax/2. |
---|
| 126 | moteurmax=toblim |
---|
| 127 | |
---|
[59] | 128 | ! betamax_2d depends on sub-grid slopes. |
---|
| 129 | slope_fich=trim(dirnameinp)//trim(slope_fich) |
---|
| 130 | call lect_input(1,'slope',1,hires_slope(:,:),slope_fich,file_ncdf) |
---|
| 131 | |
---|
| 132 | ! from slopes, we create an index between 0 & 1 |
---|
| 133 | ! 1 is very mountainous, 0 is flat |
---|
[66] | 134 | hires_slope(:,:)=(max(min(hires_slope(:,:),pente_max),pente_min)-pente_min)/(pente_max-pente_min) |
---|
[59] | 135 | ! now we compute the actual betamax used by the remplimat routines |
---|
| 136 | ! /!\ the slope is used to modify the beta where we have a temperate base, |
---|
| 137 | ! but NO ice stream... -> Slow SSA zone (SSA used to compute Ub) |
---|
[66] | 138 | betamax_2d(:,:)=max ( tostick * (1. - (1 - betamax / tostick ) * hires_slope(:,:)**(1./expo_slope)) , betamax ) |
---|
| 139 | !do j=1,ny |
---|
| 140 | ! do i=1,nx |
---|
| 141 | ! write(18745,*) betamax_2d (i,j) |
---|
| 142 | ! enddo |
---|
| 143 | !enddo |
---|
[59] | 144 | |
---|
[41] | 145 | !------------------------------------------------------------------- |
---|
| 146 | ! masque stream |
---|
| 147 | mstream_mx(:,:)=1 |
---|
| 148 | mstream_my(:,:)=1 |
---|
| 149 | |
---|
| 150 | |
---|
| 151 | ! coefficient permettant de modifier le basal drag. |
---|
| 152 | drag_mx(:,:)=1. |
---|
| 153 | drag_my(:,:)=1. |
---|
| 154 | |
---|
| 155 | where (socle_cry(:,:).lt.cry_lim) |
---|
| 156 | debug_3D(:,:,25)=1 |
---|
| 157 | elsewhere |
---|
| 158 | debug_3D(:,:,25)=0 |
---|
| 159 | endwhere |
---|
| 160 | |
---|
| 161 | |
---|
| 162 | return |
---|
| 163 | end subroutine init_dragging |
---|
| 164 | !________________________________________________________________________________ |
---|
| 165 | |
---|
| 166 | !> SUBROUTINE: dragging |
---|
| 167 | !! Defini la localisation des streams et le frottement basal |
---|
| 168 | !> |
---|
| 169 | !------------------------------------------------------------------------- |
---|
| 170 | subroutine dragging ! defini la localisation des streams et le frottement basal |
---|
[76] | 171 | !$ USE OMP_LIB |
---|
[41] | 172 | |
---|
| 173 | ! les iles n'ont pas de condition neff mais ont la condition toblim |
---|
| 174 | ! (ce bloc etait avant dans flottab) |
---|
| 175 | |
---|
[76] | 176 | !$OMP PARALLEL |
---|
| 177 | !$OMP DO |
---|
[41] | 178 | do j=2,ny |
---|
| 179 | do i=2,nx |
---|
| 180 | ilemx(i,j)=ilemx(i,j).and.(abs(rog*Hmx(i,j)*sdx(i,j)).lt.toblim) |
---|
| 181 | ilemy(i,j)=ilemy(i,j).and.(abs(rog*Hmy(i,j)*sdy(i,j)).lt.toblim) |
---|
| 182 | end do |
---|
| 183 | end do |
---|
[76] | 184 | !$OMP END DO |
---|
[41] | 185 | |
---|
| 186 | ! le calcul des fleuves se fait sur les mailles majeures en partant de la cote |
---|
| 187 | |
---|
| 188 | |
---|
| 189 | |
---|
| 190 | |
---|
| 191 | ! le gzmx initial a ete calculé dans flottab pour les points a cheval sur la grounding line |
---|
| 192 | |
---|
| 193 | ! Partant de la cote, le fleuve est determine en remontant vers le point de hauteur d'eau max |
---|
| 194 | ! tant que ce point est superieur a hwatermax. |
---|
| 195 | |
---|
| 196 | ! Attention : ce systeme ne permet pas d'avoir des fleuves convergents |
---|
| 197 | ! et les fleuves n'ont une épaisseur que de 1 noeud. |
---|
[76] | 198 | !$OMP WORKSHARE |
---|
[41] | 199 | fleuvemx(:,:)=.false. |
---|
| 200 | fleuvemy(:,:)=.false. |
---|
| 201 | fleuve(:,:)=.false. |
---|
| 202 | cote(:,:)=.false. |
---|
[76] | 203 | !$OMP END WORKSHARE |
---|
[41] | 204 | |
---|
| 205 | ! detection des cotes |
---|
[76] | 206 | !$OMP DO |
---|
[41] | 207 | do j=2,ny-1 |
---|
| 208 | do i=2,nx-1 |
---|
| 209 | if ((.not.flot(i,j)).and. & |
---|
| 210 | ((flot(i+1,j)).or.(flot(i,j+1)).or.(flot(i-1,j)).or.(flot(i,j-1)))) then |
---|
| 211 | cote(i,j)=.true. |
---|
| 212 | endif |
---|
| 213 | end do |
---|
| 214 | end do |
---|
[76] | 215 | !$OMP END DO |
---|
[41] | 216 | |
---|
| 217 | ! calcul de neff (pression effective sur noeuds majeurs) |
---|
[76] | 218 | !$OMP DO |
---|
[41] | 219 | do j=1,ny-1 |
---|
| 220 | do i=1,nx-1 |
---|
| 221 | neff(i,j)=(neffmx(i,j)+neffmx(i+1,j)+neffmy(i,j)+neffmy(i,j+1))/4 |
---|
| 222 | enddo |
---|
| 223 | enddo |
---|
[76] | 224 | !$OMP END DO |
---|
[41] | 225 | !aurel, for the borders: |
---|
[76] | 226 | !$OMP WORKSHARE |
---|
[41] | 227 | neff(:,ny)=1.e5 |
---|
| 228 | neff(nx,:)=1.e5 |
---|
| 229 | |
---|
[76] | 230 | |
---|
[41] | 231 | ! aurel, we add the neff threshold: |
---|
[56] | 232 | where ((neff(:,:).le.seuil_neff).and.(.not.flot(:,:)).and.(H(:,:).gt.1.)) fleuve(:,:)=.true. |
---|
[76] | 233 | !$OMP END WORKSHARE |
---|
[41] | 234 | |
---|
[76] | 235 | !$OMP DO |
---|
[41] | 236 | do j=1,ny-1 |
---|
| 237 | do i=1,nx-1 |
---|
| 238 | if (fleuve(i,j)) then |
---|
| 239 | fleuvemx(i,j)=.true. |
---|
| 240 | fleuvemx(i+1,j)=.true. |
---|
| 241 | fleuvemy(i,j)=.true. |
---|
| 242 | fleuvemy(i,j+1)=.true. |
---|
| 243 | end if |
---|
| 244 | end do |
---|
| 245 | end do |
---|
[76] | 246 | !$OMP END DO |
---|
[41] | 247 | |
---|
[59] | 248 | ! we look for the warm based points that will not be treated as stream (ub from SSA): |
---|
[76] | 249 | !$OMP WORKSHARE |
---|
[59] | 250 | slowssa(:,:)=.false. |
---|
| 251 | slowssamx(:,:)=.false. |
---|
| 252 | slowssamy(:,:)=.false. |
---|
[76] | 253 | !$OMP END WORKSHARE |
---|
| 254 | !$OMP DO |
---|
[59] | 255 | do j=1,ny |
---|
| 256 | do i=1,nx |
---|
| 257 | !if ((not(flot(i,j))).and.(hwater(i,j).gt.0.1)) slowssa(i,j)=.true. |
---|
[84] | 258 | if ((.not.(flot(i,j))).and.(ibase(i,j).ne.1).and.(H(i,j).gt.1.)) slowssa(i,j)=.true. |
---|
[59] | 259 | end do |
---|
| 260 | end do |
---|
[76] | 261 | !$OMP END DO |
---|
| 262 | !$OMP DO |
---|
[59] | 263 | do j=1,ny-1 |
---|
| 264 | do i=1,nx-1 |
---|
| 265 | if (slowssa(i,j)) then |
---|
| 266 | slowssamx(i,j)=.true. |
---|
| 267 | slowssamx(i+1,j)=.true. |
---|
| 268 | slowssamy(i,j)=.true. |
---|
| 269 | slowssamy(i,j+1)=.true. |
---|
| 270 | end if |
---|
| 271 | end do |
---|
| 272 | end do |
---|
[76] | 273 | !$OMP END DO |
---|
| 274 | !$OMP DO |
---|
[41] | 275 | do j=1,ny |
---|
| 276 | do i=1,nx |
---|
| 277 | |
---|
[59] | 278 | ! the actual streams and the warm based points are gzmx now: |
---|
| 279 | if ( ((.not.ilemx(i,j)).and.(fleuvemx(i,j))) .or. ((.not.ilemx(i,j)).and.(slowssamx(i,j))) ) gzmx(i,j)=.true. |
---|
| 280 | |
---|
| 281 | |
---|
[41] | 282 | ! calcul du frottement basal (ce bloc etait avant dans neffect) |
---|
| 283 | |
---|
| 284 | if (cotemx(i,j)) then ! point cotier |
---|
| 285 | betamx(i,j)=cf*neffmx(i,j) |
---|
| 286 | betamx(i,j)=min(betamx(i,j),betamax) |
---|
| 287 | |
---|
[59] | 288 | else if ((gzmx(i,j)).and.(.not.cotemx(i,j))) then ! tous les points SSA |
---|
| 289 | |
---|
| 290 | if (fleuvemx(i,j)) then ! the actual streams |
---|
| 291 | betamx(i,j)=cf*neffmx(i,j)*coef_gz |
---|
| 292 | betamx(i,j)=min(betamx(i,j),betamax) |
---|
| 293 | betamx(i,j)=max(betamx(i,j),betamin) |
---|
| 294 | |
---|
| 295 | if (cf*neffmx(i,j).gt.betamax*2.) then ! a stream that's becoming grounded... |
---|
| 296 | if (slowssamx(i,j)) then |
---|
[63] | 297 | betamx(i,j)=betamax_2d(i,j) |
---|
[59] | 298 | else |
---|
| 299 | gzmx(i,j)=.false. |
---|
[63] | 300 | betamx(i,j)=tostick |
---|
[59] | 301 | endif |
---|
| 302 | endif |
---|
| 303 | |
---|
| 304 | else ! not an actual stream, SSA is used to compute Ub |
---|
| 305 | betamx(i,j)=betamax_2d(i,j) |
---|
[41] | 306 | endif |
---|
[59] | 307 | |
---|
[41] | 308 | else if (ilemx(i,j)) then |
---|
| 309 | betamx(i,j)=cf*neffmx(i,j)*coef_ile |
---|
| 310 | betamx(i,j)=min(betamx(i,j),tob_ile) |
---|
| 311 | else if (flotmx(i,j)) then ! flottant ou ile |
---|
| 312 | betamx(i,j)=0. |
---|
| 313 | else ! grounded, SIA |
---|
[63] | 314 | betamx(i,j)=tostick ! frottement glace posee (1 bar) |
---|
[41] | 315 | endif |
---|
| 316 | |
---|
| 317 | end do |
---|
| 318 | end do |
---|
[76] | 319 | !$OMP END DO |
---|
[41] | 320 | |
---|
[76] | 321 | !$OMP DO |
---|
[59] | 322 | do j=1,ny |
---|
| 323 | do i=1,nx |
---|
[41] | 324 | |
---|
[59] | 325 | ! the actual streams and the warm based points are gzmx now: |
---|
| 326 | if ( ((.not.ilemy(i,j)).and.(fleuvemy(i,j))) .or. ((.not.ilemy(i,j)).and.(slowssamy(i,j))) ) gzmy(i,j)=.true. |
---|
[41] | 327 | |
---|
[59] | 328 | |
---|
| 329 | ! calcul du frottement basal (ce bloc etait avant dans neffect) |
---|
[41] | 330 | |
---|
| 331 | if (cotemy(i,j)) then ! point cotier |
---|
| 332 | betamy(i,j)=cf*neffmy(i,j) |
---|
| 333 | betamy(i,j)=min(betamy(i,j),betamax) |
---|
| 334 | |
---|
[59] | 335 | else if ((gzmy(i,j)).and.(.not.cotemy(i,j))) then ! tous les points SSA |
---|
| 336 | |
---|
| 337 | if (fleuvemy(i,j)) then ! the actual streams |
---|
| 338 | betamy(i,j)=cf*neffmy(i,j)*coef_gz |
---|
| 339 | betamy(i,j)=min(betamy(i,j),betamax) |
---|
| 340 | betamy(i,j)=max(betamy(i,j),betamin) |
---|
| 341 | |
---|
| 342 | if (cf*neffmy(i,j).gt.betamax*2.) then ! a stream that's becoming grounded... |
---|
| 343 | if (slowssamy(i,j)) then |
---|
[63] | 344 | betamy(i,j)=betamax_2d(i,j) |
---|
[59] | 345 | else |
---|
| 346 | gzmy(i,j)=.false. |
---|
[63] | 347 | betamy(i,j)=tostick |
---|
[59] | 348 | endif |
---|
| 349 | endif |
---|
| 350 | |
---|
| 351 | else ! not an actual stream, SSA is used to compute Ub |
---|
| 352 | betamy(i,j)=betamax_2d(i,j) |
---|
[41] | 353 | endif |
---|
[59] | 354 | |
---|
[41] | 355 | else if (ilemy(i,j)) then |
---|
[59] | 356 | betamy(i,j)=cf*neffmy(i,j)*coef_ile |
---|
| 357 | betamy(i,j)=min(betamy(i,j),tob_ile) |
---|
[41] | 358 | else if (flotmy(i,j)) then ! flottant ou ile |
---|
| 359 | betamy(i,j)=0. |
---|
[59] | 360 | else ! grounded, SIA |
---|
[63] | 361 | betamy(i,j)=tostick ! frottement glace posee (1 bar) |
---|
[41] | 362 | endif |
---|
[59] | 363 | |
---|
[41] | 364 | end do |
---|
| 365 | end do |
---|
[76] | 366 | !$OMP END DO |
---|
[41] | 367 | |
---|
[76] | 368 | !$OMP DO |
---|
[59] | 369 | do j=2,ny-1 |
---|
| 370 | do i=2,nx-1 |
---|
| 371 | beta_centre(i,j) = ((betamx(i,j)+betamx(i+1,j)) & |
---|
| 372 | + (betamy(i,j)+betamy(i,j+1)))*0.25 |
---|
| 373 | end do |
---|
| 374 | end do |
---|
[76] | 375 | !$OMP END DO |
---|
| 376 | !$OMP END PARALLEL |
---|
[59] | 377 | |
---|
[76] | 378 | !~ where (fleuve(:,:)) |
---|
| 379 | !~ debug_3D(:,:,1)=1 |
---|
| 380 | !~ elsewhere (flot(:,:)) |
---|
| 381 | !~ debug_3D(:,:,1)=2 |
---|
| 382 | !~ elsewhere |
---|
| 383 | !~ debug_3D(:,:,1)=0 |
---|
| 384 | !~ endwhere |
---|
| 385 | !~ |
---|
| 386 | !~ where (cote(:,:)) |
---|
| 387 | !~ debug_3D(:,:,2)=1 |
---|
| 388 | !~ elsewhere |
---|
| 389 | !~ debug_3D(:,:,2)=0 |
---|
| 390 | !~ endwhere |
---|
| 391 | !~ |
---|
| 392 | !~ where (fleuvemx(:,:)) |
---|
| 393 | !~ debug_3D(:,:,3)=1 |
---|
| 394 | !~ elsewhere |
---|
| 395 | !~ debug_3D(:,:,3)=0 |
---|
| 396 | !~ endwhere |
---|
| 397 | !~ |
---|
| 398 | !~ where (flotmx(:,:)) |
---|
| 399 | !~ debug_3D(:,:,3)=-1 |
---|
| 400 | !~ endwhere |
---|
[59] | 401 | |
---|
[41] | 402 | !_____________________ |
---|
| 403 | |
---|
| 404 | |
---|
[76] | 405 | !~ where (fleuvemy(:,:)) |
---|
| 406 | !~ debug_3D(:,:,4)=1 |
---|
| 407 | !~ elsewhere |
---|
| 408 | !~ debug_3D(:,:,4)=0 |
---|
| 409 | !~ endwhere |
---|
| 410 | !~ |
---|
| 411 | !~ where (flotmy(:,:)) |
---|
| 412 | !~ debug_3D(:,:,4)=-1 |
---|
| 413 | !~ end where |
---|
| 414 | !~ |
---|
| 415 | !~ debug_3D(:,:,23)= abs(RO*G*HMX(:,:)*sdx(:,:)*1.e-5) |
---|
| 416 | !~ debug_3D(:,:,24)= abs(RO*G*HMY(:,:)*sdy(:,:)*1.e-5) |
---|
[41] | 417 | |
---|
| 418 | return |
---|
| 419 | end subroutine dragging |
---|
| 420 | |
---|
[59] | 421 | end module dragging_neff_slope |
---|
[41] | 422 | |
---|