1 | !> \file old-dragging_heino_mod.f90 |
---|
2 | !! Module qui définit la loi de glissement dans les expériences Heino. |
---|
3 | !< |
---|
4 | |
---|
5 | !> \namespace sliding_heino |
---|
6 | !! Module qui definit la loi de glissement dans les experiences Heino. |
---|
7 | !! \author Catherine |
---|
8 | !! \date ... |
---|
9 | !! @note Used module |
---|
10 | !! @note - use module3D_phy |
---|
11 | !< |
---|
12 | |
---|
13 | |
---|
14 | module sliding_heino |
---|
15 | |
---|
16 | ! Défini la loi de glissement dans les expériences Heino |
---|
17 | |
---|
18 | USE module3D_phy |
---|
19 | |
---|
20 | |
---|
21 | implicit none |
---|
22 | |
---|
23 | |
---|
24 | integer,dimension(nx,ny) :: mksedim !< masque des regions sediment-hard rock |
---|
25 | integer,dimension(nx,ny) :: mkxsedim !< décalle sur les demi mailles |
---|
26 | integer,dimension(nx,ny) :: mkysedim !< décalle sur les demi mailles |
---|
27 | real :: Cr !< coefficient glissement pour hard rock |
---|
28 | real :: Cs !< coefficient glissement pour sediment |
---|
29 | real :: coefmax |
---|
30 | real :: coefslid |
---|
31 | real :: coefdrag |
---|
32 | logical,dimension(nx,ny) :: gzmx_heino !< pour transporter la valeur dans flottab |
---|
33 | logical,dimension(nx,ny) :: gzmy_heino |
---|
34 | real :: longslope |
---|
35 | |
---|
36 | |
---|
37 | |
---|
38 | contains |
---|
39 | !------------------------------------------------------------------------------- |
---|
40 | |
---|
41 | !> SUBROUTINE: init_sliding |
---|
42 | !! Cette routine fait l'initialisation du glissement et est appelée juste apres |
---|
43 | !! la lecture du masque |
---|
44 | !< |
---|
45 | subroutine init_sliding |
---|
46 | |
---|
47 | |
---|
48 | if (itracebug.eq.1) call tracebug(' Heino: entree dans routine init_sliding') |
---|
49 | |
---|
50 | ! Cette routine fait l'initialisation du glissement et est appelée juste apres |
---|
51 | ! la lecture du masque |
---|
52 | |
---|
53 | |
---|
54 | ! ecriture dans le fichier parametres |
---|
55 | write(42,*)'glissement Heino - initialisation masque' |
---|
56 | write(42,*)'loi glissement=',loigliss |
---|
57 | write(42,*)'----------------------------------------------------------' |
---|
58 | |
---|
59 | |
---|
60 | ! Le masque est égal au mk0 lu du début, |
---|
61 | |
---|
62 | mksedim(:,:)=mk0(:,:) |
---|
63 | where (mk0>1) mk0=1 |
---|
64 | |
---|
65 | !do j=ny, 1, -1 |
---|
66 | ! write(6,'(81i1)') (mk0(i,j), i=1,nx) |
---|
67 | !end do |
---|
68 | !write(6,*) mksedim(41,41) |
---|
69 | !do j=ny, 1, -1 |
---|
70 | ! write(6,'(81i1)') (mksedim(i,j), i=1,nx) |
---|
71 | !end do |
---|
72 | |
---|
73 | |
---|
74 | |
---|
75 | ! calcul du masque sur les demi mailles : frontiere -> sedim |
---|
76 | mkxsedim(:,:)=max(mksedim(:,:),eoshift(mksedim(:,:),shift=-1,boundary=0,dim=1)) |
---|
77 | mkysedim(:,:)=max(mksedim(:,:),eoshift(mksedim(:,:),shift=-1,boundary=0,dim=2)) |
---|
78 | |
---|
79 | |
---|
80 | |
---|
81 | Cr=1.e5 ! en a-1 coefficient hard rock (loi en puissance 3) |
---|
82 | Cs=500. ! en a-1 coefficient sediment (lineaire) |
---|
83 | |
---|
84 | write(42,*)'coefficients: Cr=',Cr,' Cs=',Cs,' (a-1)' |
---|
85 | write(42,*)'glissement seulement si les 2 voisins sont au point de fusion' |
---|
86 | |
---|
87 | if ((loigliss.eq.5).and.(Cs.gt.1.)) then |
---|
88 | coefdrag=rog/Cs |
---|
89 | Cs=0. |
---|
90 | write(42,*)'passage en stream si zone sediment + fusion' |
---|
91 | write(42,*)'coefficient de frottement dans la zone sediment=',coefdrag |
---|
92 | endif |
---|
93 | |
---|
94 | return |
---|
95 | end subroutine init_sliding |
---|
96 | |
---|
97 | !________________________________________________________________________________ |
---|
98 | |
---|
99 | !> SUBROUTINE: sliding |
---|
100 | !! Cette routine calcule le terme ddbx et ddby dans le cas Heino |
---|
101 | !< |
---|
102 | |
---|
103 | subroutine sliding |
---|
104 | |
---|
105 | ! cette routine calcule le terme ddbx et ddby dans le cas Heino |
---|
106 | ! ------------------------------------------------- |
---|
107 | ! GLISSEMENT |
---|
108 | ! DDBX et DDBY termes dus au glissement |
---|
109 | ! relation avec la vitesse de glissement UXB et UYB |
---|
110 | ! UXB=-DDBX*SDX et UYB=-DDBY*SDY |
---|
111 | ! ------------------------------------------------- |
---|
112 | |
---|
113 | |
---|
114 | if (itracebug.eq.1) call tracebug(' Heino: entree dans routine sliding') |
---|
115 | |
---|
116 | coefmax=1. |
---|
117 | |
---|
118 | slidy: do j=2,NY |
---|
119 | do I=1,NX |
---|
120 | pose_y: if (.not.flgzmy(i,j)) then |
---|
121 | |
---|
122 | if (HMY(I,J).lt.1.) then ! glace tres fine |
---|
123 | ddby(i,j)=0. |
---|
124 | |
---|
125 | ! else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).and.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then |
---|
126 | ! ddby(i,j)=0. ! pas de glissement seulement |
---|
127 | ! ! si les 2 voisins sont froids |
---|
128 | |
---|
129 | else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then |
---|
130 | ddby(i,j)=0. ! pas de glissement si au moins un seul |
---|
131 | ! des 2 voisins est froid |
---|
132 | |
---|
133 | else ! glissement |
---|
134 | coefslid=(hwater(i,j)+hwater(i,j-1))*0.5 |
---|
135 | |
---|
136 | coefslid=max(coefslid,0.) |
---|
137 | coefslid=min(coefslid/coefmax,1.) |
---|
138 | coefslid=1. |
---|
139 | |
---|
140 | if (mkysedim(i,j).eq.1) then ! loi hard rock |
---|
141 | ddby(i,j)=Cr*slope2my(i,j)*Hmy(i,j)*coefslid |
---|
142 | |
---|
143 | else if (mkysedim(i,j).eq.2) then ! loi sediment |
---|
144 | ddby(i,j)=Cs*Hmy(i,j)*coefslid |
---|
145 | |
---|
146 | ! if (abs(sdy(i,j)).gt.1.e-8) then |
---|
147 | ! j_1=max(j-1,1) |
---|
148 | ! j1=min(j+1,ny) |
---|
149 | ! longslope=(sdy(i,j_1)+sdy(i,j)+sdy(i,j1))/3. |
---|
150 | ! ddby(i,j)=ddby(i,j)/sdy(i,j)*longslope |
---|
151 | ! endif |
---|
152 | |
---|
153 | else |
---|
154 | ddby(i,j)=0. |
---|
155 | endif |
---|
156 | |
---|
157 | endif |
---|
158 | |
---|
159 | |
---|
160 | endif pose_y |
---|
161 | end do |
---|
162 | end do slidy |
---|
163 | |
---|
164 | |
---|
165 | slidx: do j=1,NY |
---|
166 | do I=2,NX |
---|
167 | pose_x: if (.not.flgzmx(i,j)) then |
---|
168 | |
---|
169 | if (HMx(I,J).lt.1.) then ! glace tres fine |
---|
170 | ddbx(i,j)=0. |
---|
171 | |
---|
172 | ! else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).and.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then |
---|
173 | ! ddbx(i,j)=0. ! pas de glissement seulement |
---|
174 | ! ! si les 2 voisins sont froids |
---|
175 | |
---|
176 | else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then |
---|
177 | ddbx(i,j)=0. ! pas de glissement si au moins un seul |
---|
178 | ! des 2 voisins est froid |
---|
179 | |
---|
180 | |
---|
181 | |
---|
182 | else ! glissement |
---|
183 | |
---|
184 | coefslid=(hwater(i,j)+hwater(i-1,j))*0.5 |
---|
185 | |
---|
186 | coefslid=max(coefslid,0.) |
---|
187 | coefslid=min(coefslid/coefmax,1.) |
---|
188 | coefslid=1. |
---|
189 | |
---|
190 | if (mkxsedim(i,j).eq.1) then ! loi hard rock |
---|
191 | ddbx(i,j)=Cr*slope2mx(i,j)*Hmx(i,j)*coefslid |
---|
192 | |
---|
193 | else if (mkxsedim(i,j).eq.2) then ! loi sediment |
---|
194 | ddbx(i,j)=Cs*Hmx(i,j)*coefslid |
---|
195 | |
---|
196 | ! if (abs(sdx(i,j)).gt.1.e-8) then |
---|
197 | ! i_1=max(i-1,1) |
---|
198 | ! i1=min(i+1,nx) |
---|
199 | ! longslope=(sdx(i_1,j)+sdx(i,j)+sdx(i1,j))/3. |
---|
200 | ! ddbx(i,j)=ddbx(i,j)/sdx(i,j)*longslope |
---|
201 | ! endif |
---|
202 | |
---|
203 | |
---|
204 | else |
---|
205 | ddbx(i,j)=0. |
---|
206 | endif |
---|
207 | |
---|
208 | endif |
---|
209 | |
---|
210 | endif pose_x |
---|
211 | end do |
---|
212 | end do slidx |
---|
213 | |
---|
214 | return |
---|
215 | end subroutine sliding |
---|
216 | |
---|
217 | !------------------------------------------------------------------------- |
---|
218 | !>SUBROUTINE: dragging |
---|
219 | !! Definit le basal drag |
---|
220 | !< |
---|
221 | subroutine dragging |
---|
222 | |
---|
223 | ! dans la zone sediment : appele si loigliss=5, |
---|
224 | ! defini le basal drag |
---|
225 | |
---|
226 | flgzmx(:,:)=.false. |
---|
227 | flgzmy(:,:)=.false. |
---|
228 | |
---|
229 | betamx(:,:)=1.e5 |
---|
230 | betamy(:,:)=1.e5 |
---|
231 | |
---|
232 | froty: do j=2,ny |
---|
233 | do i=1,nx |
---|
234 | |
---|
235 | if (mkysedim(i,j).eq.2) then ! loi sediment |
---|
236 | |
---|
237 | ! seulement si au dessus du point de fusion |
---|
238 | if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then |
---|
239 | ! base froide |
---|
240 | else |
---|
241 | betamy(:,:)=coefdrag |
---|
242 | flgzmy(i,j)=.true. |
---|
243 | endif |
---|
244 | endif |
---|
245 | |
---|
246 | end do |
---|
247 | end do froty |
---|
248 | |
---|
249 | |
---|
250 | frotx: do j=1,ny |
---|
251 | do i=2,nx |
---|
252 | |
---|
253 | if (mkxsedim(i,j).eq.2) then ! loi sediment |
---|
254 | |
---|
255 | ! seulement si au dessus du point de fusion |
---|
256 | if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then |
---|
257 | ! base froide |
---|
258 | else |
---|
259 | betamx(:,:)=coefdrag |
---|
260 | flgzmx(i,j)=.true. |
---|
261 | endif |
---|
262 | |
---|
263 | endif |
---|
264 | end do |
---|
265 | end do frotx |
---|
266 | |
---|
267 | gzmx(:,:)=flgzmx(:,:) |
---|
268 | gzmy(:,:)=flgzmy(:,:) |
---|
269 | gzmx_heino(:,:)=gzmx(:,:) |
---|
270 | gzmy_heino(:,:)=gzmy(:,:) |
---|
271 | |
---|
272 | return |
---|
273 | end subroutine dragging |
---|
274 | |
---|
275 | |
---|
276 | |
---|
277 | END MODULE sliding_heino |
---|
278 | |
---|