source: trunk/LMDZ.VENUS/libf/phyvenus/profile.F @ 238

Last change on this file since 238 was 3, checked in by slebonnois, 15 years ago

Creation de repertoires:

  • chantiers : pour communiquer sur nos projets de modifs
  • documentation : pour stocker les docs

Ajout de:

  • libf/phytitan : physique de Titan
  • libf/chimtitan: chimie de Titan
  • libf/phyvenus : physique de Venus
File size: 5.9 KB
Line 
1      SUBROUTINE profile(unit,nlev,dzst,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 (nlev=llm+1, surf + couches 1 a llm)
12c     dzst    dz/T (avec dz = epaisseur de la couche en m)
13c     ichoice choix de l'atmosphere:
14c             1 Temperature constante
15c             2 profil VIRA lisse
16c             3
17c             4
18c             5
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 dzst(nlev),temp(nlev)
41
42c   local:
43c   ------
44
45      INTEGER il,ichoice,isin,iter
46      REAL pi
47      REAL tref,t1,t2,t3,ww
48      REAL pic,largeur
49      REAL hauteur,tmp
50      REAL zkm(nlev)    ! altitude en km
51
52      isin = 0
53
54c-----------------------------------------------------------------------
55c   choix du profil:
56c   ----------------
57
58c la lecture se fait dans le testphys1d.def, ouvert par testphys1d.F
59      READ(unit,*)
60      READ(unit,*)
61      READ(unit,*)
62      READ(unit,*) ichoice
63      READ(unit,*) tref
64      READ(unit,*) isin
65      READ(unit,*) pic
66      READ(unit,*) largeur
67      READ(unit,*) hauteur
68      CLOSE(unit)
69
70c-----------------------------------------------------------------------
71c   ichoice=1 temperature constante:
72c   --------------------------------
73
74      IF(ichoice.EQ.1) THEN
75         temp(1) = tref
76         zkm(1)  = 0.0
77         DO il=2,nlev
78            temp(il)= tref
79            zkm(il) = zkm(il-1)+temp(il)*dzst(il)/1000.
80         ENDDO
81
82c-----------------------------------------------------------------------
83c   ichoice=2 VIRA lisse:
84c   ---------------------
85
86      ELSE IF(ichoice.EQ.2) THEN
87         temp(1) = 735.
88         zkm(1)  = 0.0
89         DO il=2,nlev
90            zkm(il) = zkm(il-1)+temp(il-1)*dzst(il)/1000. ! approx avec T(l-1)
91            if(zkm(il).lt.60.) then
92              temp(il)=735.-7.95*zkm(il)
93            else
94              temp(il)=AMAX1(258.-3.*(zkm(il)-60.),168.)
95            endif
96            zkm(il) = zkm(il-1)+(temp(il-1)+temp(il))/2.*dzst(il)/1000.
97         ENDDO
98
99c-----------------------------------------------------------------------
100c   ichoice=3 VIRA lisse - 135K:
101c   ----------------------------
102
103      ELSE IF(ichoice.EQ.3) THEN
104         temp(1) = 600.
105         zkm(1)  = 0.0
106         DO il=2,nlev
107            zkm(il) = zkm(il-1)+temp(il-1)*dzst(il)/1000. ! approx avec T(l-1)
108            temp(il)=AMAX1(600.-7.95*zkm(il),168.)
109            zkm(il) = zkm(il-1)+(temp(il-1)+temp(il))/2.*dzst(il)/1000.
110         ENDDO
111
112c-----------------------------------------------------------------------
113c   ichoice=4 :
114c   ------------------
115
116      ELSE IF(ichoice.EQ.4) THEN
117         print*,"Cas non defini..."
118         print*,"Stop dans profile.F"
119         STOP
120
121c-----------------------------------------------------------------------
122c   ichoice=5 :
123c   ----------------
124
125      ELSE IF(ichoice.EQ.5) THEN
126         print*,"Cas non defini..."
127         print*,"Stop dans profile.F"
128         STOP
129
130c-----------------------------------------------------------------------
131c   ichoice=6
132c   ---------
133
134      ELSE IF(ichoice.EQ.6) THEN
135      temp(1) = tref
136      zkm(1)  = 0.0
137      DO il=2,nlev
138        tmp=il-pic
139        temp(il)= tref + hauteur*exp(-tmp*tmp/largeur/largeur)
140        zkm(il) = zkm(il-1)+temp(il)*dzst(il)/1000.
141      ENDDO
142
143
144c-----------------------------------------------------------------------
145c   ichoice=7
146c   ---------
147
148      ELSE IF(ichoice.EQ.7) THEN
149      temp(1) = tref
150      zkm(1)  = 0.0
151      DO il=2,nlev
152        zkm(il) = zkm(il-1)+tref*dzst(il)/1000. ! approx
153        tmp=zkm(il)-pic
154        temp(il)= tref + hauteur*exp(-tmp*tmp*4/largeur/largeur)
155        zkm(il) = zkm(il-1)+(temp(il-1)+temp(il))/2.*dzst(il)/1000.
156      ENDDO
157
158c-----------------------------------------------------------------------
159c   ichoice=8
160c   ---------
161
162      ELSE IF(ichoice.GE.8) THEN
163      OPEN(11,file='profile',status='old',form='formatted',err=101)
164      DO il=1,nlev
165        READ (11,*) temp(il)
166      ENDDO
167      zkm(1) = 0.0
168      DO il=2,nlev
169        zkm(il) = zkm(il-1)+(temp(il-1)+temp(il))/2.*dzst(il)/1000.
170      ENDDO
171
172      GOTO 201
173101   STOP'fichier profile inexistant'
174201   CONTINUE
175      CLOSE(10)
176
177c-----------------------------------------------------------------------
178
179      ENDIF
180
181c-----------------------------------------------------------------------
182c   rajout eventuel d'une perturbation:
183c   -----------------------------------
184
185      IF(isin.EQ.1) THEN
186         pi=2.*ASIN(1.)
187         DO il=1,nlev
188        temp(il)=temp(il)+(1.-1000./(1000+zkm(il)*zkm(il)))*(
189     s      6.*SIN(zkm(il)*pi/6.)+9.*SIN(zkm(il)*pi/10.3) )
190         ENDDO
191      ENDIF
192
193
194c-----------------------------------------------------------------------
195c   Ecriture du profil de temperature dans un fichier profile.out
196c   -------------------------------------------------------------
197
198
199c     OPEN(12,file='profile.out')
200c         DO il=1,nlev
201c            write(12,*) temp(il)
202c           write(12,*) temp(il),zkm(il)
203c         ENDDO
204c     CLOSE(12)
205
206      RETURN
207      END
Note: See TracBrowser for help on using the repository browser.