source: trunk/LMDZ.GENERIC/libf/dynphy_lonlat/phystd/defrun_new.F @ 3937

Last change on this file since 3937 was 3893, checked in by gmilcareck, 3 months ago

Remove all "call abort" and "stop" statement in LMDZ.GENERIC and replacing them by call abort_physic().

File size: 18.5 KB
Line 
1      SUBROUTINE defrun_new( tapedef, etatinit )
2#ifndef CPP_PARA
3c
4c-----------------------------------------------------------------------
5c     Auteurs :   L. Fairhead , P. Le Van  .
6c      Modif C. Hourdin F. Forget VERSION MARTIENNE
7c
8c
9c  -------------------------------------------------------------------
10c
11c                    MODIF JUIN 2000 (zoom)
12c       .........     Version  du 29/04/97       ..........
13c
14c   Nouveaux parametres nitergdiv,nitergrot,niterh,tetagdiv,tetagrot,
15c   tetatemp   ajoutes  pour la dissipation   .
16c
17c   Autre parametre ajoute en fin de liste : ** fxyhypb **
18c
19c   Si fxyhypb = .TRUE. , choix de la fonction a derivee tangente hyperb.
20c   Sinon , choix de fxynew  , a derivee sinusoidale  ..
21c
22c   ......  etatinit = . TRUE. si defrun_new  est appele dans NEWSTART
23c   ETAT0_LMD  ou  LIMIT_LMD  pour l'initialisation de start.dat (dic) et
24c   de limit.dat (dic)  ...........
25c   Sinon  etatinit = . FALSE .
26c
27c   Donc etatinit = .F.  si on veut comparer les valeurs de  alphax ,
28c   alphay,clon,clat, fxyhypb  lues sur  le fichier  start  avec
29c   celles passees  par run.def ,  au debut du gcm, apres l'appel a
30c   lectba . 
31c   Ces parametres definissant entre autres la grille et doivent etre
32c   pareils et coherents , sinon il y aura  divergence du gcm .
33c
34c
35c-----------------------------------------------------------------------
36c   Declarations :
37c   --------------
38! to use  'getin'
39      USE ioipsl_getincom
40      use sponge_mod,only: callsponge,nsponge,mode_sponge,tetasponge
41      use control_mod,only: nday, day_step, iperiod, anneeref,
42     &                      iconser, dissip_period, iphysiq
43      USE logic_mod, ONLY: hybrid,purmats,fxyhypb,ysinus,iflag_phys
44      USE serre_mod, ONLY: clon,clat,grossismx,grossismy,dzoomx,dzoomy,
45     .                  alphax,alphay,taux,tauy
46      IMPLICIT NONE
47
48#include "dimensions.h"
49#include "paramet.h"
50!#include "control.h"
51#include "comdissnew.h"
52#include "iniprint.h"
53c
54c   arguments:
55c   ---------
56      LOGICAL  etatinit ! should be .false. for a call from gcm.F
57                        ! and .true. for a call from newstart.F
58      INTEGER  tapedef  ! unit number to assign to 'run.def' file
59c
60c   local variables:
61c   ---------------
62
63      CHARACTER ch1*72,ch2*72,ch3*72,ch4*8 ! to store various strings
64      integer tapeerr ! unit number for error message
65      parameter (tapeerr=0)
66
67c     REAL clonn,clatt,alphaxx,alphayy
68c     LOGICAL  fxyhypbb
69      INTEGER ierr
70      REAL clonn,clatt,grossismxx,grossismyy
71      REAL dzoomxx,dzoomyy,tauxx,tauyy,temp
72      LOGICAL  fxyhypbb, ysinuss
73
74
75c   initialisations:
76c   ----------------
77 
78!      lunout=6
79
80c-----------------------------------------------------------------------
81c  Parametres de controle du run:
82c-----------------------------------------------------------------------
83
84
85!Initialisation des parametres "terrestres", qui ne concernent pas
86!le modele martien et ne sont donc plus lues dans "run.def"
87
88        anneeref=0
89        ! Note: anneref is a common in 'control.h'
90
91      OPEN(tapedef,file='run.def',status='old',form='formatted'
92     .     ,iostat=ierr)
93      CLOSE(tapedef) ! first call to getin will open the file
94
95      !lunout: default unit for the text outputs
96      lunout=6
97      CALL getin('lunout', lunout)
98      IF (lunout /= 5 .and. lunout /= 6) THEN
99        OPEN(UNIT=lunout,FILE='lmdz.out',ACTION='write',
100     &       STATUS='unknown',FORM='formatted')
101      ENDIF
102
103      IF(ierr.EQ.0) THEN ! if file run.def is found
104        WRITE(lunout,*) "DEFRUN_NEW: reading file run.def"
105       
106        WRITE(lunout,*) ""
107        WRITE(lunout,*) "Number of days to run:"
108        nday=1 ! default value
109        call getin("nday",nday)
110        WRITE(lunout,*)" nday = ",nday
111
112        WRITE(lunout,*) ""
113        WRITE(lunout,*) "Number of dynamical steps per day:",
114     & "(should be a multiple of iperiod)"
115        day_step=960 ! default value
116        call getin("day_step",day_step)
117        WRITE(lunout,*)" day_step = ",day_step
118
119        WRITE(lunout,*) ""
120        WRITE(lunout,*) "periode pour le pas Matsuno (en pas)"
121        iperiod=5 ! default value
122        call getin("iperiod",iperiod)
123        WRITE(lunout,*)" iperiod = ",iperiod
124
125        WRITE(lunout,*) ""
126        WRITE(lunout,*) "periode de sortie des variables de ",
127     &  "controle (en pas)"
128        iconser=120 ! default value
129        call getin("iconser",iconser)
130        WRITE(lunout,*)" iconser = ",iconser
131
132        WRITE(lunout,*) ""
133        WRITE(lunout,*) "periode de la dissipation (en pas)"
134        dissip_period=5 ! default value
135        call getin("idissip",dissip_period)
136        ! if there is a "dissip_period" in run.def, it overrides "idissip"
137        call getin("dissip_period",dissip_period)
138        WRITE(lunout,*)" dissip_period = ",dissip_period
139
140ccc  ....   P. Le Van , modif le 29/04/97 .pour la dissipation  ...
141ccc
142        WRITE(lunout,*) ""
143        WRITE(lunout,*) "choix de l'operateur de dissipation ",
144     & "(star ou  non star )"
145        lstardis=.true. ! default value
146        call getin("lstardis",lstardis)
147        WRITE(lunout,*)" lstardis = ",lstardis
148
149        WRITE(lunout,*) ""
150        WRITE(lunout,*) "avec ou sans coordonnee hybrides"
151        hybrid=.true. ! default value
152        call getin("hybrid",hybrid)
153        WRITE(lunout,*)" hybrid = ",hybrid
154
155        WRITE(lunout,*) ""
156        WRITE(lunout,*) "nombre d'iterations de l'operateur de ",
157     & "dissipation   gradiv "
158        nitergdiv=1 ! default value
159        call getin("nitergdiv",nitergdiv)
160        WRITE(lunout,*)" nitergdiv = ",nitergdiv
161
162        WRITE(lunout,*) ""
163        WRITE(lunout,*) "nombre d'iterations de l'operateur de ",
164     & "dissipation  nxgradrot"
165        nitergrot=2 ! default value
166        call getin("nitergrot",nitergrot)
167        WRITE(lunout,*)" nitergrot = ",nitergrot
168
169        WRITE(lunout,*) ""
170        WRITE(lunout,*) "nombre d'iterations de l'operateur de ",
171     & "dissipation  divgrad"
172        niterh=2 ! default value
173        call getin("niterh",niterh)
174        WRITE(lunout,*)" niterh = ",niterh
175
176        WRITE(lunout,*) ""
177        WRITE(lunout,*) "temps de dissipation des plus petites ",
178     & "long.d ondes pour u,v (gradiv)"
179        tetagdiv=4000. ! default value
180        call getin("tetagdiv",tetagdiv)
181        WRITE(lunout,*)" tetagdiv = ",tetagdiv
182
183        WRITE(lunout,*) ""
184        WRITE(lunout,*) "temps de dissipation des plus petites ",
185     & "long.d ondes pour u,v(nxgradrot)"
186        tetagrot=5000. ! default value
187        call getin("tetagrot",tetagrot)
188        WRITE(lunout,*)" tetagrot = ",tetagrot
189
190        WRITE(lunout,*) ""
191        WRITE(lunout,*) "temps de dissipation des plus petites ",
192     & "long.d ondes pour  h ( divgrad)"
193        tetatemp=5000. ! default value
194        call getin("tetatemp",tetatemp)
195        WRITE(lunout,*)" tetatemp = ",tetatemp
196
197        WRITE(lunout,*) ""
198        WRITE(lunout,*) "coefficient pour gamdissip"
199        coefdis=0. ! default value
200        call getin("coefdis",coefdis)
201        WRITE(lunout,*)" coefdis = ",coefdis
202c
203c    ...............................................................
204
205        WRITE(lunout,*) ""
206        WRITE(lunout,*) "choix du shema d'integration temporelle ",
207     & "(true==Matsuno ou false==Matsuno-leapfrog)"
208        purmats=.false. ! default value
209        call getin("purmats",purmats)
210        WRITE(lunout,*)" purmats = ",purmats
211
212        WRITE(lunout,*) ""
213        WRITE(lunout,*) "avec ou sans physique"
214!        physic=.true. ! default value
215!        call getin("physic",physic)
216!        WRITE(lunout,*)" physic = ",physic
217        iflag_phys=1 ! default value
218        call getin("iflag_phys",iflag_phys)
219        WRITE(lunout,*)" iflag_phys = ",iflag_phys
220
221        WRITE(lunout,*) ""
222        WRITE(lunout,*) "periode de la physique (en pas)"
223        iphysiq=20 ! default value
224        call getin("iphysiq",iphysiq)
225        WRITE(lunout,*)" iphysiq = ",iphysiq
226
227!        WRITE(lunout,*) ""
228!        WRITE(lunout,*) "choix d'une grille reguliere"
229!        grireg=.true.
230!        call getin("grireg",grireg)
231!        WRITE(lunout,*)" grireg = ",grireg
232
233ccc   .... P.Le Van, ajout le 03/01/96 pour l'ecriture phys ...
234c
235!        WRITE(lunout,*) ""
236!        WRITE(lunout,*) "frequence (en pas) de l'ecriture ",
237!     & "du fichier diagfi.nc"
238!        ecritphy=240
239!        call getin("ecritphy",ecritphy)
240!        WRITE(lunout,*)" ecritphy = ",ecritphy
241
242ccc  ....   P. Le Van , ajout  le 7/03/95 .pour le zoom ...
243c     .........   (  modif  le 17/04/96 )   .........
244c
245        if (.not.etatinit ) then
246
247           clonn=63.
248           call getin("clon",clonn)
249           
250           IF( ABS(clon - clonn).GE. 0.001 )  THEN
251             call abort_physic("defrun_new",
252     &        "Value mismatch of clon in run.def and start file",1)
253           ENDIF
254c
255c
256           clatt=0.
257           call getin("clat",clatt)
258 
259           IF( ABS(clat - clatt).GE. 0.001 )  THEN
260             call abort_physic("defrun_new",
261     &        "Value mismatch of clat in run.def and start file",1)
262           ENDIF
263
264           grossismxx=1.0
265           call getin("grossismx",grossismxx)
266
267           if(grossismxx.eq.0) then 
268             write(*,*)
269             write(*,*)'ERREUR : dans run.def, grossismx =0'
270             write(*,*)'Attention a ne pas utiliser une version de'
271             write(*,*)'run.def avant le nouveau zoom LMDZ2.3 (06/2000)'
272             write(*,*)'(Il faut ajouter grossismx,dzoomx,etc... a la'
273             write(*,*)'place de alphax, alphay. cf. dyn3d). '
274             write(*,*)
275             call abort_physic("defrun_new",
276     &            "ERREUR : dans run.def, grossismx =0",1)
277           end if
278
279           IF( ABS(grossismx - grossismxx).GE. 0.001 )  THEN
280             PRINT *,' La valeur de grossismx passee par run.def est '
281     *       ,'differente de celle lue sur le fichier start =',
282     *        grossismx
283             if (grossismx.eq.0) then
284                  write(*,*) 'OK,Normal : c est un vieux start'
285     *             , 'd avant le nouveau zoom LMDZ2.3 (06/2000)'
286                 grossismx=grossismxx
287             else
288              call abort_physic("defrun_new",
289     &        "Value mismatch of grossismx in run.def and start file",1)
290             endif
291           ENDIF
292
293           grossismyy=1.0
294           call getin("grossismy",grossismyy)
295
296           IF( ABS(grossismy - grossismyy).GE. 0.001 )  THEN
297             PRINT *,' La valeur de grossismy passee par run.def est '
298     *       ,'differente de celle lue sur le fichier  start =',
299     *        grossismy
300             if (grossismy.eq.0) then
301                  write(*,*) 'OK,Normal : c est un vieux start'
302     *             , 'd avant le nouveau zoom LMDZ2.3 (06/2000)'
303                 grossismy=grossismyy
304             else
305              call abort_physic("defrun_new",
306     &        "Value mismatch of grossismy in run.def and start file",1)
307             endif
308           ENDIF
309
310
311           IF( grossismx.LT.1. )  THEN
312             call abort_physic("defrun_new",
313     &            "ATTENTION !! grossismx < 1",1)
314           ELSE
315             alphax = 1. - 1./ grossismx
316           ENDIF
317
318           IF( grossismy.LT.1. )  THEN
319             call abort_physic("defrun_new",
320     &            "ATTENTION !! grossismy < 1",1)
321           ELSE
322             alphay = 1. - 1./ grossismy
323           ENDIF
324
325           PRINT *,' '
326           PRINT *,' --> In defrun: alphax alphay  ',alphax,alphay
327           PRINT *,' '
328c
329           fxyhypbb=.false.
330           call getin("fxyhypbb",fxyhypbb)
331 
332           IF( .NOT.fxyhypb )  THEN
333             IF( fxyhypbb )     THEN
334                PRINT *,' ********  PBS DANS  DEFRUN  ******** '
335                PRINT *,' *** fxyhypb lu sur le fichier start est F ',
336     *          'alors  qu il est  T  sur  run.def  ***'
337                call abort_physic("defrun_new",
338     &               "fxyhypb=F dans start alors que =T dans run.def",1)
339             ENDIF
340           ELSE
341             IF( .NOT.fxyhypbb )   THEN
342                PRINT *,' ********  PBS DANS  DEFRUN  ******** '
343                PRINT *,' ***  fxyhypb lu sur le fichier start est T ',
344     *         'alors  qu il est  F  sur  run.def  ****  '
345                call abort_physic("defrun_new",
346     &               "fxyhypb=T dans start alors que =F dans run.def",1)
347             ENDIF
348           ENDIF
349           dzoomxx=0.0
350           call getin("dzoomx",dzoomxx)
351
352           IF( fxyhypb )  THEN
353            IF( ABS(dzoomx - dzoomxx).GE. 0.001 )  THEN
354             call abort_physic("defrun_new",
355     &       "Value mismatch of dzoomx between run.def and startfile",1)
356            ENDIF
357           ENDIF
358
359           dzoomyy=0.0
360           call getin("dzoomy",dzoomyy)
361
362           IF( fxyhypb )  THEN
363            IF( ABS(dzoomy - dzoomyy).GE. 0.001 )  THEN
364             call abort_physic("defrun_new",
365     &       "Value mismatch of dzoomy between run.def and startfile",1)
366            ENDIF
367           ENDIF
368
369           tauxx=2.0
370           call getin("taux",tauxx)
371
372           tauyy=2.0
373           call getin("tauy",tauyy)
374
375           IF( fxyhypb )  THEN
376              IF( ABS(taux - tauxx).GE. 0.001 )  THEN
377               WRITE(6,*)' La valeur de taux passee par run.def est',
378     *             'differente de celle lue sur le fichier  start '
379               call abort_physic("defrun_new",
380     &         "Value mismatch of taux between run.def and startfile",1)
381              ENDIF
382
383              IF( ABS(tauy - tauyy).GE. 0.001 )  THEN
384               WRITE(6,*)' La valeur de tauy passee par run.def est',
385     *          'differente de celle lue sur le fichier  start '
386               call abort_physic("defrun_new",
387     &         "Value mismatch of tauy between run.def and startfile",1)
388              ENDIF
389           ENDIF
390 
391        ELSE    ! Below, case when etainit=.true.
392
393           WRITE(lunout,*) ""
394           WRITE(lunout,*) "longitude en degres du centre du zoom"
395           clon=63. ! default value
396           call getin("clon",clon)
397           WRITE(lunout,*)" clon = ",clon
398           
399c
400           WRITE(lunout,*) ""
401           WRITE(lunout,*) "latitude en degres du centre du zoom "
402           clat=0. ! default value
403           call getin("clat",clat)
404           WRITE(lunout,*)" clat = ",clat
405
406           WRITE(lunout,*) ""
407           WRITE(lunout,*) "facteur de grossissement du zoom,",
408     & " selon longitude"
409           grossismx=1.0 ! default value
410           call getin("grossismx",grossismx)
411           WRITE(lunout,*)" grossismx = ",grossismx
412
413           WRITE(lunout,*) ""
414           WRITE(lunout,*) "facteur de grossissement du zoom ,",
415     & " selon latitude"
416           grossismy=1.0 ! default value
417           call getin("grossismy",grossismy)
418           WRITE(lunout,*)" grossismy = ",grossismy
419
420           IF( grossismx.LT.1. )  THEN
421            call abort_physic("defrun_new",
422     &           "ATTENTION !! grossismx < 1",1)
423           ELSE
424             alphax = 1. - 1./ grossismx
425           ENDIF
426
427           IF( grossismy.LT.1. )  THEN
428             call abort_physic("defrun_new",
429     &            "ATTENTION !! grossismy < 1",1)
430           ELSE
431             alphay = 1. - 1./ grossismy
432           ENDIF
433
434           PRINT *,' Defrun  alphax alphay  ',alphax,alphay
435c
436           WRITE(lunout,*) ""
437           WRITE(lunout,*) "Fonction  f(y)  hyperbolique  si = .true.",
438     &  ", sinon  sinusoidale"
439           fxyhypb=.false. ! default value
440           call getin("fxyhypb",fxyhypb)
441           WRITE(lunout,*)" fxyhypb = ",fxyhypb
442
443           WRITE(lunout,*) ""
444           WRITE(lunout,*) "extension en longitude de la zone du zoom",
445     & " (fraction de la zone totale)"
446           dzoomx=0. ! default value
447           call getin("dzoomx",dzoomx)
448           WRITE(lunout,*)" dzoomx = ",dzoomx
449
450           WRITE(lunout,*) ""
451           WRITE(lunout,*) "extension en latitude de la zone du zoom",
452     & " (fraction de la zone totale)"
453           dzoomy=0. ! default value
454           call getin("dzoomy",dzoomy)
455           WRITE(lunout,*)" dzoomy = ",dzoomy
456
457           WRITE(lunout,*) ""
458           WRITE(lunout,*) "raideur du zoom en  X"
459           taux=2. ! default value
460           call getin("taux",taux)
461           WRITE(lunout,*)" taux = ",taux
462
463           WRITE(lunout,*) ""
464           WRITE(lunout,*) "raideur du zoom en  Y"
465           tauy=2.0 ! default value
466           call getin("tauy",tauy)
467           WRITE(lunout,*)" tauy = ",tauy
468
469        END IF ! of if (.not.etatinit )
470
471        WRITE(lunout,*) ""
472        WRITE(lunout,*) "Use a sponge layer?"
473        callsponge=.true. ! default value
474        call getin("callsponge",callsponge)
475        WRITE(lunout,*)" callsponge = ",callsponge
476
477        WRITE(lunout,*) ""
478        WRITE(lunout,*) "Sponge: number of layers over which",
479     &                    " sponge extends"
480        nsponge=3 ! default value
481        call getin("nsponge",nsponge)
482        WRITE(lunout,*)" nsponge = ",nsponge
483
484        WRITE(lunout,*)""
485        WRITE(lunout,*)"Sponge mode: (forcing is towards ..."
486        WRITE(lunout,*)"  over upper nsponge layers)"
487        WRITE(lunout,*)"  0: (h=hmean,u=v=0)"
488        WRITE(lunout,*)"  1: (h=hmean,u=umean,v=0)"
489        WRITE(lunout,*)"  2: (h=hmean,u=umean,v=vmean)"
490        mode_sponge=2 ! default value
491        call getin("mode_sponge",mode_sponge)
492        WRITE(lunout,*)" mode_sponge = ",mode_sponge
493
494        WRITE(lunout,*) ""
495        WRITE(lunout,*) "Sponge: characteristic time scale tetasponge"
496        WRITE(lunout,*) "(seconds) at topmost layer (time scale then "
497        WRITE(lunout,*) " doubles with decreasing layer index)."
498        tetasponge=50000.0
499        call getin("tetasponge",tetasponge)
500        WRITE(lunout,*)" tetasponge = ",tetasponge
501
502
503      WRITE(lunout,*) '-----------------------------------------------'
504      WRITE(lunout,*) ' '
505      WRITE(lunout,*) ' '
506c
507
508c       Unlike on Earth (cf LMDZ2.2) , always a regular grid on Mars :
509        ysinus = .false. !Mars Mettre a jour
510
511
512      WRITE(lunout,*) '-----------------------------------------------'
513      WRITE(lunout,*) ' '
514      WRITE(lunout,*) ' '
515cc
516      ELSE
517        write(tapeerr,*) ' WHERE IS run.def ? WE NEED IT !!!!!!!!!!!!!!'
518        call abort_physic("defrun_new","Unable to find run.def",1)
519      ENDIF ! of IF(ierr.eq.0)
520
521c     Test sur le zoom
522
523      if((grossismx.eq.1).and.(grossismy.eq.1)) then 
524c        Pas de zoom :
525         write(lunout,*) 'No zoom ? -> fxyhypb set to False'
526     &   ,'           (It runs better that way)'
527         fxyhypb = .false.
528      else     
529c        Avec Zoom
530         if (.not.fxyhypb) call abort_physic("defrun_new",
531     &   "With zoom, fxyhypb should be T in run.def for this version",1)
532      end if
533#else
534      call abort_physic("defrun_new",
535     &    "defrun_new should not be used in parallel mode",1)   
536#endif
537! of #ifndef CPP_PARA
538      END
Note: See TracBrowser for help on using the repository browser.