source: trunk/SOURCES/deformation_mod-0.3.f90 @ 4

Last change on this file since 4 was 4, checked in by dumas, 10 years ago

initial import GRISLI trunk
File size: 7.5 KB
Line 
1!> \file deformation_mod-0.3.f90
2!! C'est ce module qui doit etre selectionne pour faire le calcul de la
3!! deformation de la glace en utilisant une loi de deformation de corps
4!! visqueux non newtonnien.
5!<
6
7!> \namespace deformation_mod
8!!C'est ce module qui doit etre selectionne pour faire le calcul de la
9!!deformation de la glace en utilisant une loi de deformation de corps
10!!visqueux non newtonnien.
11!! \author CatRitz
12!! \date decmebre 1999
13!!@note Ce module contient la routine d'initialisation des variables utilisees
14!! pour la deformation.
15!! @note Used modules
16!! @note   - use deform_declar
17!<
18
19
20module DEFORMATION_MOD
21
22  use DEFORM_DECLAR
23
24
25CONTAINS
26
27
28  !> SUBROUTINE: init_deformation
29  !!Routine qui initialise les variables et alloue les tableaux qui
30  !!@note interviennent dans la loi de deformation
31  !!@note Used modules
32  !!@note use module3D_phy
33  !>
34  SUBROUTINE INIT_DEFORMATION 
35
36    USE module3D_phy
37
38    implicit none
39
40    REAL,dimension(NZ) :: EDECAL ! tableau de travail (decalage de E de 1 indice)
41
42    ! tous les tableaux qui dependent de la loi de deformation sont allocatable.
43    ! L'allocation se fait juste apres avoir determiner le nombre d'exposants
44    ! differents utilises.
45
46
47    !************  allocation des tableaux *********************
48
49    if (.not.allocated(BTT)) THEN
50       allocate(BTT(NX,NY,NZ,n1poly:n2poly),stat=err)
51       if (err/=0) then
52          print *,"Erreur a l'allocation du tableau BTT",err
53          stop 4
54       end if
55    end if
56
57    if (.not.allocated(SA)) THEN
58       allocate(SA(NX,NY,NZ,n1poly:n2poly),stat=err)
59       if (err/=0) then
60          print *,"Erreur a l'allocation du tableau SA",err
61          stop 4
62       end if
63    end if
64
65    if (.not.allocated(SA_mx)) THEN
66       allocate(SA_mx(NX,NY,NZ,n1poly:n2poly),stat=err)
67       if (err/=0) then
68          print *,"Erreur a l'allocation du tableau SA_mx",err
69          stop 4
70       end if
71    end if
72
73    if (.not.allocated(SA_my)) THEN
74       allocate(SA_my(NX,NY,NZ,n1poly:n2poly),stat=err)
75       if (err/=0) then
76          print *,"Erreur a l'allocation du tableau SA_my",err
77          stop 4
78       end if
79    end if
80
81    if (.not.allocated(S2A)) THEN
82       allocate(S2A(NX,NY,NZ,n1poly:n2poly),stat=err)
83       if (err/=0) then
84          print *,"Erreur a l'allocation du tableau S2A",err
85          stop 4
86       end if
87    end if
88
89    if (.not.allocated(S2A_mx)) THEN
90       allocate(S2A_mx(NX,NY,NZ,n1poly:n2poly),stat=err)
91       if (err/=0) then
92          print *,"Erreur a l'allocation du tableau S2A_mx",err
93          stop 4
94       end if
95    end if
96
97    if (.not.allocated(S2A_my)) THEN
98       allocate(S2A_my(NX,NY,NZ,n1poly:n2poly),stat=err)
99       if (err/=0) then
100          print *,"Erreur a l'allocation du tableau S2A_my",err
101          stop 4
102       end if
103    end if
104
105
106    if (.not.allocated(SF)) THEN
107       allocate(SF(n1poly:n2poly),stat=err)
108       if (err/=0) then
109          print *,"Erreur a l'allocation du tableau SF",err
110          stop 4
111       end if
112    end if
113
114    if (.not.allocated(Q1)) THEN
115       allocate(Q1(n1poly:n2poly),stat=err)
116       if (err/=0) then
117          print *,"Erreur a l'allocation du tableau Q1",err
118          stop 4
119       end if
120    end if
121
122    if (.not.allocated(Q2)) THEN
123       allocate(Q2(n1poly:n2poly),stat=err)
124       if (err/=0) then
125          print *,"Erreur a l'allocation du tableau Q2",err
126          stop 4
127       end if
128    end if
129
130    if (.not.allocated(BAT1)) THEN
131       allocate(BAT1(n1poly:n2poly),stat=err)
132       if (err/=0) then
133          print *,"Erreur a l'allocation du tableau BAT1",err
134          stop 4
135       end if
136    end if
137
138    if (.not.allocated(BAT2)) THEN
139       allocate(BAT2(n1poly:n2poly),stat=err)
140       if (err/=0) then
141          print *,"Erreur a l'allocation du tableau BAT2",err
142          stop 4
143       end if
144    end if
145
146    if (.not.allocated(glen)) THEN
147       allocate(glen(n1poly:n2poly),stat=err)
148       if (err/=0) then
149          print *,"Erreur a l'allocation du tableau glen",err
150          stop 4
151       end if
152    end if
153
154    if (.not.allocated(TTRANS)) THEN
155       allocate(TTRANS(n1poly:n2poly),stat=err)
156       if (err/=0) then
157          print *,"Erreur a l'allocation du tableau TTRANS",err
158          stop 4
159       end if
160    end if
161    if (.not.allocated(DDX)) THEN
162       allocate(DDX(NX,NY,n1poly:n2poly),stat=err)
163       if (err/=0) then
164          print *,"Erreur a l'allocation du tableau DDX",err
165          stop 4
166       end if
167    end if
168
169    if (.not.allocated(DDY)) THEN
170       allocate(DDY(NX,NY,n1poly:n2poly),stat=err)
171       if (err/=0) then
172          print *,"Erreur a l'allocation du tableau DDY",err
173          stop 4
174       end if
175    end if
176
177    !************* FIN DE L'ALLOCATION DES TABLEAUX *************
178
179
180    !************** COEFFICIENTS DES LOIS DE COMPORTEMENT ********
181    glen(1) = 3. ! LOI DE GLEN
182    glen(2) = 1. ! LOI LINEAIRE
183
184    !********* INITIALISATION DES LOIS DE DEFORMATION **************
185
186    ! ************** FLOW LAW COEFFICIENT LOI DE GLEN **************
187    !     softening factor for flow law
188    if ((GEOPLACE.eq.'eismint').and.(ICOUPLE.eq.2)) then
189       if (IMARGIN.eq.0) then
190          SF(1)=3.52 ! coefficient multiplicateur
191       else
192          SF(1)=1.745
193       endif
194    else if ((GEOPLACE(1:6).eq.'marine') &
195         .and.(GEOPLACE.ne.'marine0')) then
196
197       !        SF(1)=6.
198       SF(1)=1.
199       !        SF(1)=0.4
200    else
201       SF(1)=3.
202    endif
203    !     arrhenius constants /Pa3/a, J/mol
204    !      A0=1.E-16  initialise dans flowlawcoefbis
205
206    if (icouple.eq.4) then
207       !       idem these catritz
208       SF(1)=3.
209       !         SF(1)=1.
210       !        SF(1)=0.4
211       !        BAT1(1)=0.
212       BAT1(1)=0.166e-15 ! fluidite Newtonienne T < TTRANS(1)
213       Q1(1)=78.2e+3        ! energie d'activation T < TTRANS(1)
214       BAT2(1)=0.2e-15   ! fluidite Newtonienne T > TTRANS(1)
215       Q2(1)=95.45e+3       ! energie d'activation T > TTRANS(1) 
216       !       temperature de transition entre les deux lois en deg Celsius
217       TTRANS(1)=-6.5
218    else
219       !       idem P. Huybrechts
220       SF(1)=5.   ! pour run Antarctique Grounded ou level 3a
221       BAT1(1)=7.65e-17
222       Q1(1)=60.e+3
223       BAT2(1)=0.28642e-15
224       Q2(1)=139.e+3
225       !     temperature de transition entre les deux lois en deg Celsius
226       TTRANS(1)=-10.
227    endif
228
229
230    !************* LINEAR FLOW LAW COEFFICIENT ***************
231
232    !     temperature de transition entre les deux lois en deg Celsius
233    TTRANS(2)=-10.
234    !  SF(2)=0.
235    SF(2)=1.
236    !        SF(2)=3.
237    !        SF(2)=6.    ! coefficient multiplicateur
238
239    !     Pour les temperatures  inferieures a TTRPHI
240
241    BAT1(2)=8.313e-8   ! fluidite Newtonienne
242    Q1(2)=40.e+3            ! energie d'activation
243
244    !     Pour les temperatures  superieures a TTRPHI
245    BAT2(2)=0.      ! fluidite Newtonienne
246    BAT2(2)=BAT1(2)  ! attention rajouté pour l'ice shelf
247
248    !  Q2(2)=SPHI*8.313e-8
249    Q2(2)=60.e+3    ! energie d'activation
250
251    !********* gas constant (J/mol/K) ************
252    RGAS=8.314
253
254
255    !************* INITIALISATION DE DDX ET DDY *****************
256    DDX(:,:,:)=0.
257    DDY(:,:,:)=0.
258
259
260    write(num_rep_42,*) 'Loi de deformation   n, sf, ttrans, bat1, Q1, bat2, Q2'
261    do iglen=1,2
262       write(num_rep_42,fmt=123) int(glen(iglen)), sf(iglen), ttrans(iglen),bat1(iglen), & 
263            Q1(iglen), bat2(iglen), Q2(iglen)
264123    format(i2,1x,f0.2,1x,f0.3,1x,4(es10.3,1x))
265       ! application des sf
266       bat1(iglen)=bat1(iglen)*sf(iglen)
267       bat2(iglen)=bat2(iglen)*sf(iglen)
268    end do
269
270
271
272
273  END SUBROUTINE INIT_DEFORMATION
274
275
276                                                END MODULE DEFORMATION_MOD
Note: See TracBrowser for help on using the repository browser.