Context Navigation

defrun_new.F @ 799

Last change on this file since 799 was 791, checked in by emillour, 12 years ago

Mars GCM:

Adapted code so that it can run fractions of days: e.g. if "nday=1.5" in

run.def, then run a sol and a half, "nday=0.75" to run three quarters of
a sol... Of course the fraction should correspond to a number of complete
dynamics/physics cycles. The fraction of the sol that a (re)start.nc
file corresponds to is (read)written as 'Time' in the file.

File size: 18.0 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format