source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/profile.F @ 313

Last change on this file since 313 was 313, checked in by ymipsl, 10 years ago
  • implement splitting of XIOS file for lmdz physics
  • Termination is done properly in parallel by calling MPI_ABORT instead of abort or stop

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