source: trunk/SOURCES/Temperature-routines/temp_col.f90

Last change on this file was 102, checked in by dumas, 7 years ago

Update closed water cycle | H = 1m on sea suppressed
File size: 6.5 KB
Line 
1!> \file temp_col.f90
2!! Subroutines for vertical temperature calculation
3!<
4
5!> \namespace temp_col
6!!  Subroutines for vertical temperature calculation
7!<
8
9!> SUBROUTINE: Temp_column
10!! Subroutine for vertical temperature calculation
11!! Used modules:
12!! - use Icetemp_declar
13!! - use Tridiagmod
14!>
15
16Subroutine  Temp_col(Nz,Nzm,Nn,Newtempcol,Ctij,Flotij,Hij,Tbmerij,Ibase,Tij,Tpmpij,Aa,Bb,Cc,Rr,Hh,Ncond,Dee,Cm, &
17                                Dzm,Phidij,Ghfij,Ifail1,Tbdotij,DTT,Dou,Tsij)
18
19!  Use Icetemp_declar
20  Use Tridiagmod
21 
22 
23  Implicit None
24  integer,intent(in) :: Nz,Nzm,Nn
25  Real,Dimension(Nz+Nzm),Intent(Out)  :: Newtempcol        !< Tableau 1d Vert. Des Temperatures En Sortie
26  Real,Dimension(Nz),Intent(In)       :: Ctij                !< Tableau 1d Vert. Conductivite Thermique
27  Logical,Intent(In)                  :: Flotij              !< Vrai Si Flotijtant (Test D'Archimede) 'O'
28  Real,Intent(In)                     :: Hij                 !< Ice Thickness  'O'
29  Real,Intent(In)                     :: Tbmerij             !< Temperature De La Mer A La Base De L'Ice Shelf
30  integer,intent(inout)               :: Ibase
31  real,dimension(Nz+Nzm),intent(inout) :: Tij             !< Temperature sur la colonne
32  real,dimension(Nz),intent(in)       :: Tpmpij               !< pressure melting point temperature in ice sheet 'o'
33  Real,Dimension(Nn),intent(inout) :: Aa    !< Work Arrays For Tridiag  !Dim Nn
34  Real,Dimension(Nn),intent(inout) :: Bb    !< Work Arrays For Tridiag  !Dim Nn
35  Real,Dimension(Nn),intent(inout) :: Cc    !< Work Arrays For Tridiag  !Dim Nn
36  Real,Dimension(Nn),intent(inout) :: Rr    !< Work Arrays For Tridiag  !Dim Nn
37  Real,Dimension(Nn),intent(inout) :: Hh    !< Work Arrays For Tridiag  !Dim Nn
38  Integer,intent(in)               :: Ncond !< Switch : Ncond=1 Thermal Conductivity In Mantle
39  real, intent(in) :: Dee
40  real,intent(in) :: Cm
41  real,intent(in) :: Dzm
42  real,intent(inout) :: Phidij     !< flux de chaleur lie a la deformation et glissement basal
43  real,intent(in) :: Ghfij      !< flux geothermique
44  integer,intent(out) :: Ifail1 !< Diagnostique erreur dans Tridiag
45  real, intent(out) :: Tbdotij  !< variation in time of basal temperature
46  real, intent(in) :: DTT
47  Real, intent(inout) :: Dou
48  real, intent(in) :: Tsij        !< surface ice temperature  'o'
49 
50 
51 
52  ! variables locales
53  integer :: K
54  Real :: Dzi,Tss
55  Real,Dimension(Nzm+1) :: Abis,Bbis,Cbis,Rbis,Hbis   
56
57  If (Hij.Gt.10.) Then
58
59     ! Conditions Aux Limites a La Base De La Glace
60     !________________________________________________
61
62     ! Base Froide
63     If ( ((Ibase.Eq.1).Or.(Ibase.Eq.4) &
64          .Or. ((Ibase.Eq.5).And.(Tij(Nz).Lt.Tpmpij(Nz)))) &
65          .And..Not.Flotij) Then
66
67        Ibase=1
68        Bb(Nz)=1.
69
70        If (Ncond.Eq.1) Then ! Avec Socle
71           Dzi=Hij*Dee*Cm/Ctij(Nz)
72           Aa(Nz)=-Dzm/(Dzm+Dzi)
73           Cc(Nz)=-Dzi/(Dzm+Dzi)
74           Rr(Nz)=Hij*Dee*Phidij*Dzm/Ctij(Nz)/(Dzm+Dzi)
75        Else                 ! Sans Socle
76           Aa(Nz)=-1.
77           Cc(Nz)=0.
78           Rr(Nz)=-(Ghfij-Phidij)/Ctij(Nz)*Hij*Dee
79        Endif
80
81     Else
82
83        ! Base Temperee Ou Shelf
84        ! ----------------------               
85        If (Ibase.Eq.5)  Ibase=2
86        Ibase=Max(Ibase,2)
87        Aa(Nz)=0.
88        Bb(Nz)=1.
89        Cc(Nz)=0.
90
91        If (.Not.Flotij) Then
92           Rr(Nz)=Tpmpij(Nz)
93        Else
94           Rr(Nz)=Tbmerij
95        End If
96
97     Endif
98
99
100     If (Ncond.Eq.1) Then ! Avec Socle
101        Do K=Nz+1,nz+nzm-1
102           Rr(K)=Tij(K)
103        End Do
104        Call Tridiag (Aa,Bb,Cc,Rr,Hh,nz+nzm,Ifail1)
105
106     Else                 ! Sans Socle
107
108        Call Tridiag (Aa,Bb,Cc,Rr,Hh,Nz,Ifail1)
109     Endif
110
111
112     ! Temperature Dans Le Socle Si Ncond=0, Lineaire Avec Le Gradient Ghf
113     If (Ncond.Eq.0) Then
114        Do K=Nz+1,nz+nzm
115           Hh(K)=Hh(Nz)-Dzm*(K-Nz)*Ghfij/Cm
116        End Do
117     Endif
118
119     Do K=1,nz+nzm
120        Newtempcol(K)=Hh(K)
121     End Do
122     Tbdotij=(Newtempcol(Nz)-Tij(Nz))/Dtt
123
124     !      ------------------------------------------------- Glace Tres Fine (H<10m Ou H=0)
125  Else If (Hij.Le.10.) Then
126
127     !       Pour Eviter Des Problemes Lorsque H Est Inferieur A 10 M
128     !       Profil Lineaire, Pente=Flux Geothermique
129     !       Le Cas Sans Glace Est Traite De La Meme Facon  (Dou=0)
130
131     ! ........................................ Avec Calcul Dans Le Socle
132     If (Ncond.Eq.1.And..Not.Flotij) Then   
133        If (Hij.Gt.0.) Then
134           ! Gradient Dans Le Socle
135           Dou=(Tij(Nz+1)-Tij(Nz))/Dzm*Cm     
136           Dou=Dou/Ctij(Nz)*Dee*Hij
137
138           Tss=Min(0.,Tsij)
139           Do K=1,Nz
140              Newtempcol(K)=Tss+Dou*(K-1.)
141           End Do
142
143        Else
144           ! Pas De Glace
145           Tss=Tsij
146           Do K=1,Nz
147              Newtempcol(K)=Tss
148           End Do
149        End If
150
151        ! Calcul Dans Le Socle Meme S'Il N'Y A Pas De Glace
152        Rr(Nz)=Newtempcol(Nz)
153        Aa(Nz)=0
154        Bb(Nz)=1
155        Cc(Nz)=0
156
157        Do K=Nz+1,nz+nzm-1
158           Rr(K)=Tij(K)
159        End Do
160
161        ! Creation De Nouveaux Tableaux Juste Pour Tridiag     
162        Do K=1,Nzm+1
163           Abis(K)=Aa(Nz-1+K)
164           Bbis(K)=Bb(Nz-1+K)
165           Cbis(K)=Cc(Nz-1+K)
166           Rbis(K)=Rr(Nz-1+K)
167        End Do
168
169        Call Tridiag (Abis,Bbis,Cbis,Rbis,Hbis,Nzm+1,Ifail1)
170
171        Do K=1,Nzm+1
172           Hh(Nz-1+K)=Hbis(K)
173        End Do
174
175        ! ........................................ Sans Calcul Dans Le Socle
176     Else
177        ! Calotte Posee                   
178        If ((Hij.Gt.0.).And..Not.Flotij) Then
179           Dou=-Ghfij/Ctij(Nz)*Dee*Hij
180           Tss=Min(0.,Tsij)
181           Do K=1,Nz
182              Newtempcol(K)=Tss+Dou*(K-1.)
183           End Do
184
185           !  Shelf
186        Else If ((Hij.Gt.0.).And.Flotij) Then
187           Tss=Min(0.,Tsij)
188           Dou=(Tbmerij-Tss)*Dee
189           Do K=1,Nz
190              Newtempcol(K)=Tss+Dou*(K-1.)
191           End Do
192
193           ! Pas De Glace
194        Else
195           Tss=Tsij     
196           Do K=1,Nz
197              Newtempcol(K)=Tss
198           End Do
199        End If
200
201        ! Temperature Dans Le Socle, Lineaire Avec Le Gradient Ghf
202        If (.Not.Flotij) Then
203           Do K=Nz+1,nz+nzm
204              Hh(K)=Newtempcol(Nz)-Dzm*(K-Nz)*Ghfij/Cm
205           End Do
206        Else
207           Do K=Nz+1,nz+nzm
208              Hh(K)=Tbmerij-Dzm*(K-Nz)*Ghfij/Cm
209           End Do
210        End If
211
212     Endif
213     ! Fin Du Test Socle
214
215
216     Do K=Nz+1,nz+nzm
217        Newtempcol(K)=Hh(K)
218     End Do
219
220     Tbdotij=(Newtempcol(Nz)-Tij(Nz))/Dtt
221     Ibase=5
222     Phidij=0.
223
224
225  End If
226
227End Subroutine Temp_col
Note: See TracBrowser for help on using the repository browser.