source: trunk/mesoscale/LMDZ.MARS/libf_gcm/phymars/profile.F @ 57

Last change on this file since 57 was 57, checked in by aslmd, 14 years ago

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