source: LMDZ.3.3/branches/rel-LF/libf/dyn3d/nudge.F @ 242

Last change on this file since 242 was 232, checked in by lmdzadmin, 23 years ago

Merge par rapport a la branche principale
LF

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 7.1 KB
Line 
1c $Header$
2      subroutine nudge(itau,ucov,vcov,teta,masse,ps)
3
4      IMPLICIT NONE
5
6c      ......   Version  du 10/01/98    ..........
7
8c             avec  coordonnees  verticales hybrides
9c   avec nouveaux operat. dissipation * ( gradiv2,divgrad2,nxgraro2 )
10
11c=======================================================================
12c
13c   Auteur:  P. Le Van /L. Fairhead/F.Hourdin
14c   -------
15c
16c   Objet:
17c   ------
18c
19c   GCM LMD nouvelle grille
20c
21c=======================================================================
22
23c   ...  Dans inigeom , nouveaux calculs pour les elongations  cu , cv
24c        et possibilite d'appeler une fonction f(y)  a derivee tangente
25c        hyperbolique a la  place de la fonction a derivee sinusoidale.         
26
27c   ...  Possibilite de choisir le shema de Van-leer pour l'advection de
28c         q  , en faisant iadv = 3  dans   traceur  (29/04/97) .
29c
30c-----------------------------------------------------------------------
31c   Declarations:
32c   -------------
33
34#include "dimensions.h"
35#include "paramet.h"
36#include "comconst.h"
37#include "comdissnew.h"
38#include "comvert.h"
39#include "comgeom.h"
40#include "logic.h"
41#include "temps.h"
42#include "control.h"
43#include "ener.h"
44#include "netcdf.inc"
45#include "description.h"
46#include "serre.h"
47#include "tracstoke.h"
48
49
50c   variables dynamiques
51      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
52      REAL teta(ip1jmp1,llm)                 ! temperature potentielle
53      REAL ps(ip1jmp1)                       ! pression  au sol
54      REAL masse(ip1jmp1,llm)                ! masse d'air
55
56c   variables dynamiques pour les reanalyses.
57      REAL ucovrea1(ip1jmp1,llm),vcovrea1(ip1jm,llm) !vts cov reas
58      REAL tetarea1(ip1jmp1,llm)             ! temp pot  reales
59      REAL masserea1(ip1jmp1,llm)             ! masse
60      REAL psrea1(ip1jmp1)             ! ps
61      REAL ucovrea2(ip1jmp1,llm),vcovrea2(ip1jm,llm) !vts cov reas
62      REAL tetarea2(ip1jmp1,llm)             ! temp pot  reales
63      REAL masserea2(ip1jmp1,llm)             ! masse
64      REAL psrea2(ip1jmp1)             ! ps
65      real alphaT,alphaP,alphau,alphav
66      real dday_step,toto,reste,itau_test
67      INTEGER step_rea,count_no_rea
68
69
70c
71      INTEGER itau,ij,l
72      real ditau,tau,a
73
74      LOGICAL first
75      integer online
76      save first,online
77      data first,online/.true.,1/
78
79      save ucovrea1,vcovrea1,tetarea1,masserea1,psrea1
80      save ucovrea2,vcovrea2,tetarea2,masserea2,psrea2
81
82      save alphaT,alphau,alphav,alphaP,itau_test
83      save step_rea,count_no_rea
84
85c-----------------------------------------------------------------------
86c   initialisations pour la lecture des reanalyses.
87c    alpha determine la part des injections de donnees a chaque etape
88c    alpha=1 signifie pas d'injection
89c    alpha=0 signifie injection totale
90c-----------------------------------------------------------------------
91
92      print*,'ONLINE=',online
93      if(online.eq.-1) then
94          return
95      endif
96
97      if (first) then
98         print*
99     s   ,'1: en-ligne, 0: hors-ligne (x=x_rea), -1: climat (x=x_gcm)'
100cnec         read(*,*) online
101cnec      online=-1
102cnec     print*,'Entrer les constantes de temps de rappel en jours'
103cnec         read(*,*)alphaT
104cnec         read(*,*)alphau
105cnec         read(*,*)alphav
106cnec         read(*,*)alphaP
107             alphaT=0.1
108             alphau=0.1
109             alphav=0.1
110             alphaP=1.e10
111         print*,'alphaT,alphau,alphav,alphaP'
112     s          ,alphaT,alphau,alphav,alphaP
113         if(online.eq.-1) return
114         print*,'alpha rappel pour T, u, v, P ',
115     s    alphaT,alphau,alphav,alphaP
116c   En fait, le coef alpha qu'on utilise est
117c   x_gcm = a * x_gcm + ( 1 - a ) * x_reanalys
118c   on a alors a=1-dt/tau
119c   ou dt est la frequence a laquelle on corrige (ici iperiod*dtvr)
120c   et tau la constante de temps lue. D'ou:
121         if (online.eq.1) then
122            alphaT=1.-iperiod*dtvr/(daysec*alphaT)
123            alphau=1.-iperiod*dtvr/(daysec*alphau)
124            alphav=1.-iperiod*dtvr/(daysec*alphav)
125            alphaP=1.-iperiod*dtvr/(daysec*alphaP)
126         else
127            alphaT=0.
128            alphau=0.
129            alphav=0.
130            alphaP=0.
131c           physic=.false.
132         endif
133         print*,'alpha rappel pour T, u, v, P ',
134     s    alphaT,alphau,alphav,alphaP
135
136         itau_test=1001
137         step_rea=1
138         count_no_rea=0
139
140c    itau_test    montre si l'importation a deja ete faite au rang itau
141
142
143c   Lecture du premier etat des reanalyses.
144         call read_reanalyse(1
145     s   ,ucovrea2,vcovrea2,tetarea2,masserea2,ps,1)
146
147c-----------------------------------------------------------------------
148c   Debut de l'integration temporelle:
149c   ----------------------------------
150
151         first=.false.
152      endif ! first
153c
154C-----------------------------------------------------------------------
155C----- IMPORTATION DES VENTS,PRESSION ET TEMPERATURE REELS:
156C-----------------------------------------------------------------------
157
158      ditau=real(itau)
159      dday_step=real(day_step)
160      write(*,*)'ditau,dday_step'
161      write(*,*)ditau,dday_step
162      toto=4*ditau/dday_step
163      reste=toto-aint(toto)
164c     write(*,*)'toto,reste',toto,reste
165
166      if (reste.eq.0.) then
167        if (itau_test.eq.itau) then
168          write(*,*)'deuxieme passage de advreel a itau=',itau
169          stop
170        else
171           CALL SCOPY( ijmllm ,vcovrea2, 1, vcovrea1 , 1 )
172           CALL SCOPY( ijp1llm,ucovrea2, 1, ucovrea1 , 1 )
173           CALL SCOPY( ijp1llm,tetarea2,1,tetarea1   , 1 )
174c          CALL SCOPY( ijp1llm,masserea2,1,masserea1   , 1 )
175c          CALL SCOPY( ip1jmp1,psrea2, 1,  psrea1, 1 )
176
177          print*,'LECTURE REANALYSES, pas ',step_rea
178     s         ,'apres ',count_no_rea,' non lectures'
179           step_rea=step_rea+1
180           itau_test=itau
181           call read_reanalyse(step_rea
182     s     ,ucovrea2,vcovrea2,tetarea2,masserea2,ps,1)
183        endif
184      else
185        count_no_rea=count_no_rea+1
186      endif
187 
188C-----------------------------------------------------------------------
189c   Corrrections
190c
191c    x_gcm = a * x_gcm + (1-a) * x_reanalyses
192C-----------------------------------------------------------------------
193
194      ditau=real(itau)
195      dday_step=real(day_step)
196
197c     print*,'alpha1',alphaT,alphaP,alphau,alphav
198
199      tau=4*ditau/dday_step
200      tau=tau-aint(tau)
201
202      do l=1,llm
203         do ij=1,ip1jmp1
204            a=(1.-tau)*ucovrea1(ij,l)+tau*ucovrea2(ij,l)
205            ucov(ij,l)=alphau*ucov(ij,l)+(1-alphau)*a
206            a=(1.-tau)*tetarea1(ij,l)+tau*tetarea2(ij,l)
207            teta(ij,l)=alphaT*teta(ij,l)+(1-alphaT)*a
208c           a=(1.-tau)*masserea1(ij,l)+tau*masserea2(ij,l)
209c           masse(ij,l)=alphaP*masse(ij,l)+(1-alphaP)*a
210         enddo
211         do ij=1,ip1jm
212            a=(1.-tau)*vcovrea1(ij,l)+tau*vcovrea2(ij,l)
213            vcov(ij,l)=alphav*vcov(ij,l)+(1-alphav)*a
214         enddo
215      enddo
216
217c     call dump2d(iip1,jjp1,tetarea1,'TETA REA 1     ')
218c     call dump2d(iip1,jjp1,tetarea2,'TETA REA 2     ')
219c     call dump2d(iip1,jjp1,teta,'TETA           ')
220Cmaf  on ne nudge pas sur la pression de surface
221c      do ij=1,ip1jmp1
222c         a=(1.-tau)*psrea1(ij)+tau*psrea2(ij)
223c         ps(ij)=alphaP*ps(ij)+(1-alphaP)*a
224c      enddo
225
226
227      return
228      end
Note: See TracBrowser for help on using the repository browser.