[4] | 1 | !> \file dragging_mismip3d_mod.f90 |
---|
| 2 | !< |
---|
| 3 | |
---|
| 4 | !> \namespace dragging_mismip3 |
---|
| 5 | !! Module qui fait le dragging de MISMIP avec un beta non lineraire |
---|
| 6 | !! tau_b = C | U_b |**(m-1) U_b |
---|
| 7 | |
---|
| 8 | !! \author Cat |
---|
| 9 | !! \date 11 mars 2012 |
---|
| 10 | !! @note Used module |
---|
| 11 | !! @note - use module3D_phy |
---|
| 12 | !! @note - use param_phy_mod |
---|
| 13 | !< |
---|
| 14 | |
---|
| 15 | module dragging_mismip3 |
---|
| 16 | |
---|
| 17 | |
---|
| 18 | use module3d_phy |
---|
| 19 | use param_phy_mod |
---|
| 20 | use interface_input |
---|
| 21 | |
---|
| 22 | implicit none |
---|
| 23 | |
---|
| 24 | real, dimension(nx,ny) :: Vcol_x !< uniquement pour compatibilite dans le spinup |
---|
| 25 | real, dimension(nx,ny) :: Vcol_y !< uniquement pour compatibilite |
---|
| 26 | real, dimension(nx,ny) :: Vsl_x !< |
---|
| 27 | real, dimension(nx,ny) :: Vsl_y !< |
---|
| 28 | logical :: corr_def = .False. !< for deformation correction (compatibility) |
---|
| 29 | |
---|
| 30 | |
---|
| 31 | real :: beta_limgz !< when beta gt beta_limgz -> not gzmx |
---|
| 32 | |
---|
| 33 | ! equation glissement tau_b = C | U_b |**(m-1) U_b |
---|
| 34 | |
---|
| 35 | real :: Coef_C ! coefficient du glissement dans equation (C) |
---|
| 36 | real :: exp_m ! exposant du glissement |
---|
| 37 | real :: m_1 ! exposant dans l'equation |
---|
| 38 | real :: lim_unorm ! vitesse correspondante a beta_limgz |
---|
| 39 | integer :: i_dome ! indice du noeud milieu en x |
---|
| 40 | integer :: j_dome ! indice du noeud milieu en y |
---|
| 41 | |
---|
| 42 | contains |
---|
| 43 | !------------------------------------------------------------------------------- |
---|
| 44 | |
---|
| 45 | !> subroutine: init_dragging |
---|
| 46 | !! Cette routine fait l'initialisation du dragging. |
---|
| 47 | !< |
---|
| 48 | |
---|
| 49 | subroutine init_dragging ! Cette routine fait l'initialisation du dragging. |
---|
| 50 | |
---|
| 51 | implicit none |
---|
| 52 | |
---|
| 53 | if (itracebug.eq.1) call tracebug(' entree dans init_dragging de MISMIP') |
---|
| 54 | ! les coefficients sont codes en dur |
---|
| 55 | |
---|
| 56 | Coef_C = 1.e7 ! Coefficient glissement en Pa m-1/3 s-1/3 |
---|
| 57 | exp_m = 1./3. ! exposant du glissement (m) |
---|
| 58 | Coef_C = secyear**(-exp_m) * Coef_C ! en Pa m-1/3 a-1/3 |
---|
| 59 | m_1 = exp_m - 1 ! exposant dans le calcul de beta |
---|
| 60 | |
---|
| 61 | |
---|
| 62 | ! limites sur beta |
---|
| 63 | beta_limgz = 1.e6 ! valeur maxi de beta <--> U = 0 |
---|
| 64 | betamax = beta_limgz |
---|
| 65 | lim_unorm = (beta_limgz/Coef_C)**(1./m_1) ! vitesse correspondante a beta_limgz (de l'ordre de 5 mm/an) |
---|
| 66 | |
---|
| 67 | |
---|
| 68 | ! coefficient permettant de modifier le basal drag. |
---|
| 69 | ! drag est celui fixe au depart qui reflete les proprietes du socle |
---|
| 70 | drag_mx(:,:)=1. |
---|
| 71 | drag_my(:,:)=1. |
---|
| 72 | |
---|
| 73 | ! mstream est celui courant, eventuellement affecte par les conditions basales |
---|
| 74 | mstream_mx(:,:)=1. |
---|
| 75 | mstream_my(:,:)=1. |
---|
| 76 | |
---|
| 77 | ! il faut quand meme au moins un point immobile |
---|
| 78 | ! Ca fait 4 points pour des raisons de symetrie |
---|
| 79 | |
---|
| 80 | i_dome = (nx-1)/2 +1 |
---|
| 81 | j_dome = (ny-1)/2 +1 |
---|
| 82 | |
---|
| 83 | drag_mx(i_dome:i_dome+1,j_dome ) = 0 |
---|
| 84 | drag_my(i_dome ,j_dome:j_dome+1 ) = 0 |
---|
| 85 | mstream_mx(i_dome:i_dome+1,j_dome ) = 0 |
---|
| 86 | mstream_my(i_dome ,j_dome:j_dome+1 ) = 0 |
---|
| 87 | |
---|
| 88 | |
---|
| 89 | ! initialisation de certains tableaux |
---|
| 90 | ilemx(:,:) = .false. |
---|
| 91 | ilemy(:,:) = .false. |
---|
| 92 | |
---|
| 93 | |
---|
| 94 | ! initialisation des vitesses de glissement : lues dans lect_mismip_mod |
---|
| 95 | |
---|
| 96 | ! initialisation du beta en fonction des vitesses |
---|
| 97 | ! |
---|
| 98 | |
---|
| 99 | call drag_non_lin_stag (nx,ny,Uxflgz,Uyflgz,betamx,betamy) |
---|
| 100 | |
---|
| 101 | return |
---|
| 102 | end subroutine init_dragging |
---|
| 103 | |
---|
| 104 | |
---|
| 105 | !________________________________________________________________________________ |
---|
| 106 | |
---|
| 107 | !> subroutine: dragging |
---|
| 108 | !! Defini la localisation des streams et le frottement basal |
---|
| 109 | !! est appele par flottab |
---|
| 110 | !! en entrant dans cette routine les logiques |
---|
| 111 | !! float et flotmx/flotmy sont connus |
---|
| 112 | !< |
---|
| 113 | !------------------------------------------------------------------------- |
---|
| 114 | subroutine dragging ! defini la localisation des streams et le frottement basal |
---|
| 115 | |
---|
| 116 | |
---|
| 117 | if (itracebug.eq.1) call tracebug(' entree dans dragging de MISMIP') |
---|
| 118 | |
---|
| 119 | call drag_non_lin_stag (nx,ny,Uxflgz,Uyflgz,betamx,betamy) |
---|
| 120 | |
---|
| 121 | |
---|
| 122 | ! pas d'iles |
---|
| 123 | ilemx(:,:) = .false. |
---|
| 124 | ilemy(:,:) = .false. |
---|
| 125 | |
---|
| 126 | |
---|
| 127 | where (mstream_mx(:,:).eq.1) ! SSA betamx garde la valeur |
---|
| 128 | gzmx(:,:) = .true. |
---|
| 129 | |
---|
| 130 | |
---|
| 131 | ! betamx(:,:) = beta_limgz ! attention ligne pour debug qui enleve le gissement |
---|
| 132 | |
---|
| 133 | elsewhere |
---|
| 134 | gzmx(:,:) = .false. |
---|
| 135 | betamx(:,:) = beta_limgz |
---|
| 136 | Uxflgz(:,:) = 0. |
---|
| 137 | end where |
---|
| 138 | |
---|
| 139 | where (flotmx(:,:)) ! pas de frottement sous les parties flottantes |
---|
| 140 | betamx(:,:) = 0. |
---|
| 141 | end where |
---|
| 142 | |
---|
| 143 | where (mstream_my(:,:).eq.1) ! SSA betamy garde la valeur |
---|
| 144 | gzmy(:,:) = .true. |
---|
| 145 | |
---|
| 146 | |
---|
| 147 | ! betamy(:,:) = beta_limgz ! attention ligne pour debug qui enleve le gissement |
---|
| 148 | |
---|
| 149 | |
---|
| 150 | |
---|
| 151 | elsewhere |
---|
| 152 | gzmy(:,:) = .false. |
---|
| 153 | betamy(:,:) = beta_limgz |
---|
| 154 | Uyflgz(:,:) = 0. |
---|
| 155 | end where |
---|
| 156 | |
---|
| 157 | where (flotmy(:,:)) ! pas de frottement sous les parties flottantes |
---|
| 158 | betamy(:,:) = 0. |
---|
| 159 | end where |
---|
| 160 | |
---|
| 161 | if (itracebug.eq.1) call tracebug(' fin dragging de MISMIP') |
---|
| 162 | |
---|
| 163 | end subroutine dragging |
---|
| 164 | |
---|
| 165 | !________________________________________________________________________________________ |
---|
| 166 | !< calcule la valeur de beta dependant de la vitesse basale |
---|
| 167 | |
---|
| 168 | |
---|
| 169 | subroutine drag_non_lin_stag (nxx,nyy,Ux_mx,Uy_my,beta_x,beta_y) |
---|
| 170 | |
---|
| 171 | implicit none |
---|
| 172 | |
---|
| 173 | !declarations variables dummy |
---|
| 174 | integer :: nxx,nyy !< dimension des tableaux |
---|
| 175 | real, dimension(nxx,nyy),intent(in) :: Ux_mx !< vitesse de glissement x sur le noeud mineur x |
---|
| 176 | real, dimension(nxx,nyy),intent(in) :: Uy_my !< vitesse de glissement y sur le noeud mineur y |
---|
| 177 | real, dimension(nxx,nyy),intent(out) :: beta_x !< valeur de beta sur le noeud mineur x |
---|
| 178 | real, dimension(nxx,nyy),intent(out) :: beta_y !< valeur de beta sur le noeud mineur y |
---|
| 179 | |
---|
| 180 | ! variables locales |
---|
| 181 | real, dimension(nxx,nyy) :: Uxmy ! vitesse selon x en my |
---|
| 182 | real, dimension(nxx,nyy) :: Uymx ! vitesse selon y en mx |
---|
| 183 | real, dimension(nxx,nyy) :: Unorm_mx !< norme vitesse de glissement sur le noeud mineur x |
---|
| 184 | real, dimension(nxx,nyy) :: Unorm_my !< norme vitesse de glissement sur le noeud mineur y |
---|
| 185 | |
---|
| 186 | |
---|
| 187 | if (itracebug.eq.1) call tracebug(' entree dans drag_non_lin_stag MISMIP') |
---|
| 188 | |
---|
| 189 | |
---|
| 190 | do j=2,nyy-1 |
---|
| 191 | do i=2,nxx |
---|
| 192 | Uymx(i,j) = ((Uy_my(i,j)+Uy_my(i-1,j+1))+(Uy_my(i-1,j)+Uy_my(i,j+1)))/4. ! vitesse selon y en mx |
---|
| 193 | Unorm_mx(i,j) = (Ux_mx(i,j)**2+Uymx(i,j)**2)**0.5 |
---|
| 194 | end do |
---|
| 195 | end do |
---|
| 196 | |
---|
| 197 | |
---|
| 198 | do j=2,nyy |
---|
| 199 | do i=2,nxx-1 |
---|
| 200 | Uxmy(i,j) = ((Ux_mx(i,j)+Ux_mx(i+1,j-1))+(Ux_mx(i,j-1)+Ux_mx(i+1,j)))/4. ! vitesse selon y en mx |
---|
| 201 | Unorm_my(i,j) = (Uy_my(i,j)**2+Uxmy(i,j)**2)**0.5 |
---|
| 202 | end do |
---|
| 203 | end do |
---|
| 204 | |
---|
| 205 | |
---|
| 206 | do j=2,nyy-1 |
---|
| 207 | do i=2,nxx |
---|
| 208 | |
---|
| 209 | if (Unorm_mx(i,j).gt.lim_unorm) then |
---|
| 210 | beta_x(i,j) = Coef_C * Unorm_mx(i,j)**m_1 |
---|
| 211 | else |
---|
| 212 | beta_x(i,j) = beta_limgz |
---|
| 213 | end if |
---|
| 214 | end do |
---|
| 215 | end do |
---|
| 216 | |
---|
| 217 | |
---|
| 218 | do j=2,nyy |
---|
| 219 | do i=2,nxx-1 |
---|
| 220 | |
---|
| 221 | if (Unorm_my(i,j).gt.lim_unorm) then |
---|
| 222 | beta_y(i,j) = Coef_C * Unorm_my(i,j)**m_1 |
---|
| 223 | else |
---|
| 224 | beta_y(i,j) = beta_limgz |
---|
| 225 | end if |
---|
| 226 | end do |
---|
| 227 | end do |
---|
| 228 | |
---|
| 229 | end subroutine drag_non_lin_stag |
---|
| 230 | |
---|
| 231 | end module dragging_mismip3 |
---|
| 232 | |
---|