Context Navigation

defrun_new.F @ 999

Last change on this file since 999 was 999, checked in by tnavarro, 11 years ago
Possibility to store multiple initial states in one start/startfi. This is RETROCOMPATIBLE. New option ecrithist in run.def to write data in start/startfi every ecrithist dynamical timestep. New option timestart in run.def to initialize the GCM with the time timestart stored in start
File size: 19.2 KB

Rev	Line
[38]	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:"
[791]	97	nday_r=1 ! default value
	98	call getin("nday",nday_r)
	99	WRITE(tapeout,*)" nday = ",nday_r
[828]	100
	101	ndynstep=-9999 ! default value
	102	call getin("ndynstep",ndynstep)
	103	if (ndynstep .gt. 0) then
	104	WRITE(tapeout,*) ""
	105	WRITE(tapeout,*) "Number of dynamical steps to run:"
	106	WRITE(tapeout,*) " ndynstep = ",ndynstep
	107	WRITE(tapeout,*) " nday value is now discarded "
	108	endif
[38]	109
	110	WRITE(tapeout,*) ""
	111	WRITE(tapeout,*) "Number of dynamical steps per day:",
	112	& "(should be a multiple of iperiod)"
	113	day_step=960 ! default value
	114	call getin("day_step",day_step)
	115	WRITE(tapeout,*)" day_step = ",day_step
	116
	117	WRITE(tapeout,*) ""
	118	WRITE(tapeout,*) "periode pour le pas Matsuno (en pas)"
	119	iperiod=5 ! default value
	120	call getin("iperiod",iperiod)
[791]	121	! verify that day_step is a multiple of iperiod
	122	if (((1.*day_step)/iperiod).ne.(day_step/iperiod)) then
	123	write(tapeout,*)" Error! iperiod must be such that",
	124	& " day_step is a multiple of iperiod, but iperiod=",
	125	& iperiod," and day_step=",day_step
	126	else
	127	WRITE(tapeout,*)" iperiod = ",iperiod
	128	endif
	129
[38]	130	WRITE(tapeout,*) ""
	131	WRITE(tapeout,*) "periode de sortie des variables de ",
	132	& "controle (en pas)"
	133	iconser=120 ! default value
	134	call getin("iconser",iconser)
	135	WRITE(tapeout,*)" iconser = ",iconser
	136
	137	WRITE(tapeout,*) ""
	138	WRITE(tapeout,*) "periode de la dissipation (en pas)"
	139	idissip=5 ! default value
	140	call getin("idissip",idissip)
	141	WRITE(tapeout,*)" idissip = ",idissip
	142
	143	ccc .... P. Le Van , modif le 29/04/97 .pour la dissipation ...
	144	ccc
	145	WRITE(tapeout,*) ""
	146	WRITE(tapeout,*) "choix de l'operateur de dissipation ",
	147	& "(star ou non star )"
	148	lstardis=.true. ! default value
	149	call getin("lstardis",lstardis)
	150	WRITE(tapeout,*)" lstardis = ",lstardis
	151
	152	WRITE(tapeout,*) ""
	153	WRITE(tapeout,*) "avec ou sans coordonnee hybrides"
	154	hybrid=.true. ! default value
	155	call getin("hybrid",hybrid)
	156	WRITE(tapeout,*)" hybrid = ",hybrid
	157
	158	WRITE(tapeout,*) ""
	159	WRITE(tapeout,*) "nombre d'iterations de l'operateur de ",
	160	& "dissipation gradiv "
	161	nitergdiv=1 ! default value
	162	call getin("nitergdiv",nitergdiv)
	163	WRITE(tapeout,*)" nitergdiv = ",nitergdiv
	164
	165	WRITE(tapeout,*) ""
	166	WRITE(tapeout,*) "nombre d'iterations de l'operateur de ",
	167	& "dissipation nxgradrot"
	168	nitergrot=2 ! default value
	169	call getin("nitergrot",nitergrot)
	170	WRITE(tapeout,*)" nitergrot = ",nitergrot
	171
	172	WRITE(tapeout,*) ""
	173	WRITE(tapeout,*) "nombre d'iterations de l'operateur de ",
	174	& "dissipation divgrad"
	175	niterh=2 ! default value
	176	call getin("niterh",niterh)
	177	WRITE(tapeout,*)" niterh = ",niterh
	178
	179	WRITE(tapeout,*) ""
	180	WRITE(tapeout,*) "temps de dissipation des plus petites ",
	181	& "long.d ondes pour u,v (gradiv)"
	182	tetagdiv=4000. ! default value
	183	call getin("tetagdiv",tetagdiv)
	184	WRITE(tapeout,*)" tetagdiv = ",tetagdiv
	185
	186	WRITE(tapeout,*) ""
	187	WRITE(tapeout,*) "temps de dissipation des plus petites ",
	188	& "long.d ondes pour u,v(nxgradrot)"
	189	tetagrot=5000. ! default value
	190	call getin("tetagrot",tetagrot)
	191	WRITE(tapeout,*)" tetagrot = ",tetagrot
	192
	193	WRITE(tapeout,*) ""
	194	WRITE(tapeout,*) "temps de dissipation des plus petites ",
	195	& "long.d ondes pour h ( divgrad)"
	196	tetatemp=5000. ! default value
	197	call getin("tetatemp",tetatemp)
	198	WRITE(tapeout,*)" tetatemp = ",tetatemp
	199
	200	WRITE(tapeout,*) ""
	201	WRITE(tapeout,*) "coefficient pour gamdissip"
	202	coefdis=0. ! default value
	203	call getin("coefdis",coefdis)
	204	WRITE(tapeout,*)" coefdis = ",coefdis
	205	c
	206	c ...............................................................
	207
	208	WRITE(tapeout,*) ""
	209	WRITE(tapeout,*) "choix du shema d'integration temporelle ",
	210	& "(true==Matsuno ou false==Matsuno-leapfrog)"
	211	purmats=.false. ! default value
	212	call getin("purmats",purmats)
	213	WRITE(tapeout,*)" purmats = ",purmats
	214
	215	WRITE(tapeout,*) ""
	216	WRITE(tapeout,*) "avec ou sans physique"
	217	physic=.true. ! default value
	218	call getin("physic",physic)
	219	WRITE(tapeout,*)" physic = ",physic
	220
	221	WRITE(tapeout,*) ""
	222	WRITE(tapeout,*) "periode de la physique (en pas)"
	223	iphysiq=20 ! default value
	224	call getin("iphysiq",iphysiq)
[791]	225	! verify that day_step is a multiple of iphysiq
	226	if (((1.*day_step)/iphysiq).ne.(day_step/iphysiq)) then
	227	write(tapeout,*)" Error! iphysiq must be such that",
	228	& " day_step is a multiple of iphysiq, but iphysiq=",
	229	& iphysiq," and day_step=",day_step
	230	else
	231	WRITE(tapeout,*)" iphysiq = ",iphysiq
	232	endif
[38]	233
	234	WRITE(tapeout,*) ""
	235	WRITE(tapeout,*) "choix d'une grille reguliere"
	236	grireg=.true.
	237	call getin("grireg",grireg)
	238	WRITE(tapeout,*)" grireg = ",grireg
	239
	240	ccc .... P.Le Van, ajout le 03/01/96 pour l'ecriture phys ...
	241	c
	242	WRITE(tapeout,*) ""
	243	WRITE(tapeout,*) "frequence (en pas) de l'ecriture ",
	244	& "du fichier diagfi.nc"
	245	ecritphy=240
	246	call getin("ecritphy",ecritphy)
[791]	247	! verify that ecriphy is indeed a multiple of iphysiq
	248	if (((1.*ecritphy)/iphysiq).ne.(ecritphy/iphysiq)) then
	249	write(tapeout,*)" Error! ecritphy must be a multiple",
	250	& " of iphysiq, but ecritphy=",ecritphy," and iphysiq=",
	251	& iphysiq
	252	else
	253	WRITE(tapeout,*)" ecritphy = ",ecritphy
	254	endif
[999]	255
	256	ccc .... T.Navarro, read start time (06/2013) ...
	257	c
	258	WRITE(tapeout,*) ""
	259	WRITE(tapeout,*) "date de debut dans le fichier start.nc"
	260	timestart=-9999
	261	call getin("timestart",timestart)
	262	WRITE(tapeout,*)" timestart = ",timestart
[38]	263
[999]	264	ccc .... T.Navarro, start output (01/2013) ...
	265	c
	266	WRITE(tapeout,*) ""
	267	WRITE(tapeout,*) "frequence (en pas) de l'ecriture ",
	268	& "du fichier start.nc"
	269	ecritstart=0
	270	call getin("ecritstart",ecritstart)
	271	! verify that ecritstart is indeed a multiple of iphysiq
	272	if ( ((real(ecritstart))/iphysiq)
	273	& .ne.(real(ecritstart)/iphysiq) ) then
	274	write(tapeout,*)" Error! ecritstart must be a multiple",
	275	& " of iphysiq, but ecritstart=",ecritstart," and iphysiq=",
	276	& iphysiq
	277	else
	278	WRITE(tapeout,*)" ecritstart = ",ecritstart
	279	endif
	280
[38]	281	ccc .... P. Le Van , ajout le 7/03/95 .pour le zoom ...
	282	c ......... ( modif le 17/04/96 ) .........
	283	c
	284	if (.not.etatinit ) then
	285
	286	clonn=63.
	287	call getin("clon",clonn)
	288
	289	IF( ABS(clon - clonn).GE. 0.001 ) THEN
	290	PRINT *,' La valeur de clon passee par run.def est '
	291	* ,'differente de celle lue sur le fichier start '
	292	STOP
	293	ENDIF
	294	c
	295	c
	296	clatt=0.
	297	call getin("clat",clatt)
	298
	299	IF( ABS(clat - clatt).GE. 0.001 ) THEN
	300	PRINT *,' La valeur de clat passee par run.def est '
	301	* ,'differente de celle lue sur le fichier start '
	302	STOP
	303	ENDIF
	304
	305	grossismxx=1.0
	306	call getin("grossismx",grossismxx)
	307
	308	if(grossismxx.eq.0) then
	309	write(,)
	310	write(,)'ERREUR : dans run.def, grossismx =0'
	311	write(,)'Attention a ne pas utiliser une version de'
	312	write(,)'run.def avant le nouveau zoom LMDZ2.3 (06/2000)'
	313	write(,)'(Il faut ajouter grossismx,dzoomx,etc... a la'
	314	write(,)'place de alphax, alphay. cf. dyn3d). '
	315	write(,)
	316	stop
	317	end if
	318
	319	IF( ABS(grossismx - grossismxx).GE. 0.001 ) THEN
	320	PRINT *,' La valeur de grossismx passee par run.def est '
	321	* ,'differente de celle lue sur le fichier start =',
	322	* grossismx
	323	if (grossismx.eq.0) then
	324	write(,) 'OK,Normal : c est un vieux start'
	325	* , 'd avant le nouveau zoom LMDZ2.3 (06/2000)'
	326	grossismx=grossismxx
	327	else
	328	STOP
	329	endif
	330	ENDIF
	331
	332	grossismyy=1.0
	333	call getin("grossismy",grossismyy)
	334
	335	IF( ABS(grossismy - grossismyy).GE. 0.001 ) THEN
	336	PRINT *,' La valeur de grossismy passee par run.def est '
	337	* ,'differente de celle lue sur le fichier start =',
	338	* grossismy
	339	if (grossismy.eq.0) then
	340	write(,) 'OK,Normal : c est un vieux start'
	341	* , 'd avant le nouveau zoom LMDZ2.3 (06/2000)'
	342	grossismy=grossismyy
	343	else
	344	STOP
	345	endif
	346	ENDIF
	347
	348
	349	IF( grossismx.LT.1. ) THEN
	350	PRINT ,' ATTENTION !! grossismx < 1 . * '
	351	STOP
	352	ELSE
	353	alphax = 1. - 1./ grossismx
	354	ENDIF
	355
	356	IF( grossismy.LT.1. ) THEN
	357	PRINT ,' ATTENTION !! grossismy < 1 . * '
	358	STOP
	359	ELSE
	360	alphay = 1. - 1./ grossismy
	361	ENDIF
	362
	363	PRINT *,' '
	364	PRINT *,' --> In defrun: alphax alphay ',alphax,alphay
	365	PRINT *,' '
	366	c
	367	fxyhypbb=.false.
	368	call getin("fxyhypbb",fxyhypbb)
	369
	370	IF( .NOT.fxyhypb ) THEN
	371	IF( fxyhypbb ) THEN
	372	PRINT ,' ***** PBS DANS DEFRUN ****** '
	373	PRINT ,' ** fxyhypb lu sur le fichier start est F ',
	374	* 'alors qu il est T sur run.def ***'
	375	STOP
	376	ENDIF
	377	ELSE
	378	IF( .NOT.fxyhypbb ) THEN
	379	PRINT ,' ***** PBS DANS DEFRUN ****** '
	380	PRINT ,' ** fxyhypb lu sur le fichier start est T ',
	381	* 'alors qu il est F sur run.def **** '
	382	STOP
	383	ENDIF
	384	ENDIF
	385	dzoomxx=0.0
	386	call getin("dzoomx",dzoomxx)
	387
	388	IF( fxyhypb ) THEN
	389	IF( ABS(dzoomx - dzoomxx).GE. 0.001 ) THEN
	390	PRINT *,' La valeur de dzoomx passee par run.def est '
	391	* ,'differente de celle lue sur le fichier start '
	392	STOP
	393	ENDIF
	394	ENDIF
	395
	396	dzoomyy=0.0
	397	call getin("dzoomy",dzoomyy)
	398
	399	IF( fxyhypb ) THEN
	400	IF( ABS(dzoomy - dzoomyy).GE. 0.001 ) THEN
	401	PRINT *,' La valeur de dzoomy passee par run.def est '
	402	* ,'differente de celle lue sur le fichier start '
	403	STOP
	404	ENDIF
	405	ENDIF
	406
	407	tauxx=2.0
	408	call getin("taux",tauxx)
	409
	410	tauyy=2.0
	411	call getin("tauy",tauyy)
	412
	413	IF( fxyhypb ) THEN
	414	IF( ABS(taux - tauxx).GE. 0.001 ) THEN
	415	WRITE(6,*)' La valeur de taux passee par run.def est',
	416	* 'differente de celle lue sur le fichier start '
	417	CALL ABORT
	418	ENDIF
	419
	420	IF( ABS(tauy - tauyy).GE. 0.001 ) THEN
	421	WRITE(6,*)' La valeur de tauy passee par run.def est',
	422	* 'differente de celle lue sur le fichier start '
	423	CALL ABORT
	424	ENDIF
	425	ENDIF
	426
	427	ELSE ! Below, case when etainit=.true.
	428
	429	WRITE(tapeout,*) ""
	430	WRITE(tapeout,*) "longitude en degres du centre du zoom"
	431	clon=63. ! default value
	432	call getin("clon",clon)
	433	WRITE(tapeout,*)" clon = ",clon
	434
	435	c
	436	WRITE(tapeout,*) ""
	437	WRITE(tapeout,*) "latitude en degres du centre du zoom "
	438	clat=0. ! default value
	439	call getin("clat",clat)
	440	WRITE(tapeout,*)" clat = ",clat
	441
	442	WRITE(tapeout,*) ""
	443	WRITE(tapeout,*) "facteur de grossissement du zoom,",
	444	& " selon longitude"
	445	grossismx=1.0 ! default value
	446	call getin("grossismx",grossismx)
	447	WRITE(tapeout,*)" grossismx = ",grossismx
	448
	449	WRITE(tapeout,*) ""
	450	WRITE(tapeout,*) "facteur de grossissement du zoom ,",
	451	& " selon latitude"
	452	grossismy=1.0 ! default value
	453	call getin("grossismy",grossismy)
	454	WRITE(tapeout,*)" grossismy = ",grossismy
	455
	456	IF( grossismx.LT.1. ) THEN
	457	PRINT ,' ATTENTION !! grossismx < 1 . * '
	458	STOP
	459	ELSE
	460	alphax = 1. - 1./ grossismx
	461	ENDIF
	462
	463	IF( grossismy.LT.1. ) THEN
	464	PRINT ,' ATTENTION !! grossismy < 1 . * '
	465	STOP
	466	ELSE
	467	alphay = 1. - 1./ grossismy
	468	ENDIF
	469
	470	PRINT *,' Defrun alphax alphay ',alphax,alphay
	471	c
	472	WRITE(tapeout,*) ""
	473	WRITE(tapeout,*) "Fonction f(y) hyperbolique si = .true.",
	474	& ", sinon sinusoidale"
	475	fxyhypb=.false. ! default value
	476	call getin("fxyhypb",fxyhypb)
	477	WRITE(tapeout,*)" fxyhypb = ",fxyhypb
	478
	479	WRITE(tapeout,*) ""
	480	WRITE(tapeout,*) "extension en longitude de la zone du zoom",
	481	& " (fraction de la zone totale)"
	482	dzoomx=0. ! default value
	483	call getin("dzoomx",dzoomx)
	484	WRITE(tapeout,*)" dzoomx = ",dzoomx
	485
	486	WRITE(tapeout,*) ""
	487	WRITE(tapeout,*) "extension en latitude de la zone du zoom",
	488	& " (fraction de la zone totale)"
	489	dzoomy=0. ! default value
	490	call getin("dzoomy",dzoomy)
	491	WRITE(tapeout,*)" dzoomy = ",dzoomy
	492
	493	WRITE(tapeout,*) ""
	494	WRITE(tapeout,*) "raideur du zoom en X"
	495	taux=2. ! default value
	496	call getin("taux",taux)
	497	WRITE(tapeout,*)" taux = ",taux
	498
	499	WRITE(tapeout,*) ""
	500	WRITE(tapeout,*) "raideur du zoom en Y"
	501	tauy=2.0 ! default value
	502	call getin("tauy",tauy)
	503	WRITE(tapeout,*)" tauy = ",tauy
	504
	505	END IF ! of if (.not.etatinit )
	506
	507	WRITE(tapeout,*) ""
	508	WRITE(tapeout,*) "Avec sponge layer"
	509	callsponge=.true. ! default value
	510	call getin("callsponge",callsponge)
	511	WRITE(tapeout,*)" callsponge = ",callsponge
	512
	513	WRITE(tapeout,*) ""
[575]	514	WRITE(tapeout,*) "Sponge: number of layers over which",
	515	& " sponge extends"
	516	nsponge=3 ! default value
	517	call getin("nsponge",nsponge)
	518	WRITE(tapeout,*)" nsponge = ",nsponge
	519
	520	WRITE(tapeout,*)""
	521	WRITE(tapeout,*)"Sponge mode: (forcing is towards ..."
	522	WRITE(tapeout,*)" over upper nsponge layers)"
	523	WRITE(tapeout,*)" 0: (h=hmean,u=v=0)"
	524	WRITE(tapeout,*)" 1: (h=hmean,u=umean,v=0)"
	525	WRITE(tapeout,*)" 2: (h=hmean,u=umean,v=vmean)"
[38]	526	mode_sponge=2 ! default value
	527	call getin("mode_sponge",mode_sponge)
	528	WRITE(tapeout,*)" mode_sponge = ",mode_sponge
	529
	530	WRITE(tapeout,*) ""
[575]	531	WRITE(tapeout,*) "Sponge: characteristice time scale tetasponge"
	532	WRITE(tapeout,*) "(seconds) at topmost layer (time scale then "
	533	WRITE(tapeout,*) " doubles with decreasing layer index)."
	534	tetasponge=30000.0
[38]	535	call getin("tetasponge",tetasponge)
	536	WRITE(tapeout,*)" tetasponge = ",tetasponge
	537
	538
	539	WRITE(tapeout,*) '-----------------------------------------------'
	540	WRITE(tapeout,*) ' '
	541	WRITE(tapeout,*) ' '
	542	c
	543
	544	c Unlike on Earth (cf LMDZ2.2) , always a regular grid on Mars :
	545	ysinus = .false. !Mars Mettre a jour
	546
	547
	548	WRITE(tapeout,*) '-----------------------------------------------'
	549	WRITE(tapeout,*) ' '
	550	WRITE(tapeout,*) ' '
	551	cc
	552	ELSE
	553	write(tapeerr,*) ' WHERE IS run.def ? WE NEED IT !!!!!!!!!!!!!!'
	554	stop
	555	ENDIF ! of IF(ierr.eq.0)
	556
	557	c Test sur le zoom
	558
	559	if((grossismx.eq.1).and.(grossismy.eq.1)) then
	560	c Pas de zoom :
	561	write(tapeout,*) 'No zoom ? -> fxyhypb set to False'
	562	& ,' (It runs better that way)'
	563	fxyhypb = .false.
	564	else
	565	c Avec Zoom
	566	if (.not.fxyhypb) stop 'With zoom, fxyhypb should be set to T
	567	&in run.def for this version... -> STOP ! '
	568	end if
	569
	570	RETURN
	571	END

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

Download in other formats:

Original Format