source: trunk/LMDZ.MARS/libf/phymars/dyn1d/profile_temp_mod.F90 @ 4063

Last change on this file since 4063 was 3918, checked in by jbclement, 5 months ago

Mars PCM 1D:

  • "profile.F" is transformed into a module "profile_temp_mod.F90" with F90 format. The file "profile" where the temperature is read is renamed into "profile_temp" to avoid any confusion.
  • Small corrections and improvements throughout the code and scripts in the deftank.

JBC
File size: 6.4 KB
Line 
1MODULE profile_temp_mod
2
3implicit none
4
5CONTAINS
6
7SUBROUTINE profile_temp(nlev,zkm,temp)
8
9use ioipsl_getincom, only: getin ! To use 'getin'
10use comcstfi_h,      only: pi
11     
12IMPLICIT NONE
13
14!=======================================================================
15!     Subroutine utilisee dans le modele 1-D "testphys1d"
16!     pour l'initialisation du profil atmospherique
17!=======================================================================
18!
19!   differents profils d'atmospheres. T=f(z)
20!   entree:
21!     nlev    nombre de niveaux
22!     zkm     alititudes en km
23!     ichoice choix de l'atmosphere:
24!             1 Temperature constante
25!             2 profil Savidjari
26!             3 Lindner (profil polaire)
27!             4 Inversion pour Francois
28!             5 Seiff (moyen)
29!             6 T constante + perturbation gauss (level) (christophe 10/98)
30!             7 T constante + perturbation gauss    (km) (christophe 10/98)
31!             8 Lecture du profile dans un fichier ASCII (profile_temp)
32!     tref    temperature de reference
33!     isin    ajout d'une perturbation (isin=1)
34!     pic     pic perturbation gauss pour ichoice = 6 ou 7
35!     largeur largeur de la perturbation gauss pour ichoice = 6 ou 7
36!     hauteur hauteur de la perturbation gauss pour ichoice = 6 ou 7
37!
38!   sortie:
39!     temp    temperatures en K
40!     
41!=======================================================================
42!-----------------------------------------------------------------------
43!   declarations:
44!   -------------
45
46!   arguments:
47!   ----------
48
49INTEGER,               intent(in) :: nlev
50REAL, dimension(nlev), intent(in) :: zkm
51REAL, dimension(nlev), intent(out) :: temp
52
53!   local:
54!   ------
55
56INTEGER :: il,ichoice,iseiff,isin,iter,ierr
57INTEGER, PARAMETER :: nseiff = 37
58REAL :: tref,t1,t2,t3,ww
59REAL, dimension(nseiff) :: tseiff
60DATA tseiff/214. ,213.8,213.4,212.4,209.3,205. ,201.4,197.8, &
61            194.6,191.4,188.2,185.2,182.5,180. ,177.5,175  , &
62            172.5,170. ,167.5,164.8,162.4,160. ,158. ,156. , &
63            154.1,152.2,150.3,148.7,147.2,145.7,144.2,143. , &
64            142. ,141. ,140  ,139.5,139./
65REAL :: pic,largeur
66REAL :: hauteur,tmp
67
68!-----------------------------------------------------------------------
69!   read input profile type:
70!--------------------------
71
72ichoice=1 ! default value for ichoice
73call getin("ichoice",ichoice)
74tref=200 ! default value for tref
75call getin("tref",tref)
76isin=0 ! default value for isin (=0 means no perturbation)
77call getin("isin",isin)
78pic=26.522 ! default value for pic
79call getin("pic",pic)
80largeur=10 ! default value for largeur
81call getin("largeur",largeur)
82hauteur=30 ! default value for hauteur
83call getin("hauteur",hauteur)
84
85!-----------------------------------------------------------------------
86!   ichoice=1 temperature constante:
87!   --------------------------------
88
89IF(ichoice.EQ.1) THEN
90    DO il=1,nlev
91        temp(il)=tref
92    ENDDO
93
94!-----------------------------------------------------------------------
95!   ichoice=2 savidjari:
96!   --------------------
97
98ELSE IF(ichoice.EQ.2) THEN
99    DO il=1,nlev
100        temp(il)=AMAX1(219.-2.5*zkm(il),140.)
101    ENDDO
102
103!-----------------------------------------------------------------------
104!   ichoice=3 Lindner:
105!   ------------------
106
107ELSE IF(ichoice.EQ.3) THEN
108    DO il=1,nlev
109        IF(zkm(il).LT.2.5) THEN
110            temp(il)=150.+30.*zkm(il)/2.5
111        ELSE IF(zkm(il).LT.5.) THEN
112            temp(il)=180.
113        ELSE
114            temp(il)=AMAX1(180.-2.*(zkm(il)-5.),130.)
115        ENDIF
116    ENDDO
117
118!-----------------------------------------------------------------------
119!   ichoice=4 Inversion pour Francois:
120!   ------------------
121
122ELSE IF(ichoice.EQ.4) THEN
123    DO il=1,nlev
124        IF(zkm(il).LT.20.) THEN
125            temp(il)=135.
126        ELSE
127            temp(il)=AMIN1(135.+5.*(zkm(il)-20.),200.)
128        ENDIF
129    ENDDO
130
131!-----------------------------------------------------------------------
132!   ichoice=5 Seiff:
133!   ----------------
134
135ELSE IF(ichoice.EQ.5) THEN
136    DO il=1,nlev
137        iseiff=INT(zkm(il)/2.)+1
138        IF(iseiff.LT.nseiff-1) THEN
139            temp(il)=tseiff(iseiff)+(zkm(il)-2.*(iseiff-1))*(tseiff(iseiff+1)-tseiff(iseiff))/2.
140        ELSE
141            temp(il)=tseiff(nseiff)
142        ENDIF
143    ENDDO
144! IF(ichoice.EQ.6) THEN
145!    DO iter=1,6
146!        t2=temp(1)
147!        t3=temp(2)
148!        DO il=2,nlev-1
149!            t1=t2
150!            t2=t3
151!            t3=temp(il+1)
152!            ww=tanh(zkm(il)/20.)
153!            ww=ww*ww*ww
154!            temp(il)=t2+ww*.5*(t1-2.*t2+t3)
155!        ENDDO
156!    ENDDO
157!ENDIF
158
159!-----------------------------------------------------------------------
160!   ichoice=6
161!   ---------
162
163ELSE IF(ichoice.EQ.6) THEN
164    DO il=1,nlev
165        tmp=il-pic
166        temp(il)=tref + hauteur*exp(-tmp*tmp/largeur/largeur)
167    ENDDO
168
169!-----------------------------------------------------------------------
170!   ichoice=7
171!   ---------
172
173ELSE IF(ichoice.EQ.7) THEN
174    DO il=1,nlev
175        tmp=zkm(il)-pic
176        temp(il)=tref + hauteur*exp(-tmp*tmp*4/largeur/largeur)
177    ENDDO
178
179!-----------------------------------------------------------------------
180!   ichoice=8
181!   ---------
182
183ELSE IF(ichoice.GE.8) THEN
184    OPEN(11,file='profile_temp',status='old',form='formatted',iostat=ierr)
185    if (ierr == 0) then
186        DO il=1,nlev
187            READ (11,*) temp(il)
188        ENDDO
189    else
190        error stop 'File "profile_temp" not found!'
191    endif
192    CLOSE(11)
193
194!-----------------------------------------------------------------------
195ENDIF
196
197!-----------------------------------------------------------------------
198!   rajout eventuel d'une perturbation:
199!   -----------------------------------
200
201IF (isin.EQ.1) THEN
202    DO il=1,nlev
203!        if (nlev.EQ.501) then
204!            if (zkm(il).LE.70.5) then
205!                temp(il)=temp(il)+(1.-1000./(1000+zkm(il)*zkm(il)))*( 6.*SIN(zkm(il)*pi/6.)+9.*SIN(zkm(il)*pi/10.3) )
206!            endif
207!        else
208        temp(il)=temp(il)+(1.-1000./(1000+zkm(il)*zkm(il)))*( 6.*SIN(zkm(il)*pi/6.)+9.*SIN(zkm(il)*pi/10.3) )
209!        endif
210    ENDDO
211ENDIF
212
213!-----------------------------------------------------------------------
214!   Ecriture du profil de temperature dans un fichier profile_temp.out
215!   ------------------------------------------------------------------
216
217OPEN(12,file='profile_temp.out',status='unknown',form='formatted')
218    DO il=1,nlev
219        write(12,*) temp(il)
220    ENDDO
221CLOSE(12)
222
223END SUBROUTINE profile_temp
224
225END MODULE profile_temp_mod
Note: See TracBrowser for help on using the repository browser.