Context Navigation

source: tags/ORCHIDEE_1_9_6/ORCHIDEE_OL/teststomate.f90 @ 881

Last change on this file since 881 was 844, checked in by didier.solyga, 12 years ago
Improved parameters documentation and presentationfor ORCHIDEE_OL. The presentation is the same for all paramters so a script can automatically build a orchidee.default file.
File size: 35.2 KB

Line
1	PROGRAM teststomate
2	!---------------------------------------------------------------------
3	! This program reads a forcing file for STOMATE
4	! which was created with STOMATE coupled to SECHIBA.
5	! It then integrates STOMATE for a
6	! given number of years using this forcing.
7	! The GPP is scaled to the new calculated LAI...
8	! Nothing is done for soil moisture which should
9	! actually evolve with the vegetation.
10	!---------------------------------------------------------------------
11	!- IPSL (2006)
12	!- This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
13	USE netcdf
14	!-
15	USE defprec
16	USE constantes
17	USE pft_parameters
18	USE stomate_data
19	USE ioipsl
20	USE grid
21	USE slowproc
22	USE stomate
23	USE intersurf, ONLY: stom_define_history, stom_ipcc_define_history, intsurf_time, l_first_intersurf, check_time,impose_param
24	USE parallel
25	!-
26	IMPLICIT NONE
27	!-
28	! Declarations
29	!-
30	INTEGER(i_std) :: kjpij,kjpindex
31	REAL(r_std) :: dtradia,dt_force
32
33	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indices
34	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indexveg
35	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: soiltype, veget_x
36	REAL(r_std),DIMENSION(:),ALLOCATABLE :: totfrac_nobio
37	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: frac_nobio
38	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_max_x
39	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lai_x
40	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_force_x
41	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_max_force_x
42	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lai_force_x
43	REAL(r_std),DIMENSION(:),ALLOCATABLE :: t2m,t2m_min,temp_sol
44	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: soiltemp,soilhum
45	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: humrel_x
46	REAL(r_std),DIMENSION(:),ALLOCATABLE :: litterhum
47	REAL(r_std),DIMENSION(:),ALLOCATABLE :: precip_rain,precip_snow
48	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: gpp_x
49	REAL(r_std),DIMENSION(:),ALLOCATABLE :: deadleaf_cover
50	REAL(r_std),DIMENSION(:,:,:),ALLOCATABLE :: assim_param_x
51	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: height_x
52	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: qsintmax_x
53	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: co2_flux
54	REAL(r_std),DIMENSION(:),ALLOCATABLE :: fco2_lu
55
56	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indices_g
57	REAL(r_std),DIMENSION(:),ALLOCATABLE :: x_indices_g
58	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: x_neighbours_g
59
60	INTEGER :: ier,iret
61	INTEGER :: ncfid
62	CHARACTER(LEN=20),SAVE :: thecalendar='noleap'
63
64	LOGICAL :: a_er
65	CHARACTER(LEN=80) :: &
66	& dri_restname_in,dri_restname_out, &
67	& sec_restname_in,sec_restname_out, &
68	& sto_restname_in,sto_restname_out, &
69	& stom_histname, stom_ipcc_histname
70	INTEGER(i_std) :: iim,jjm,llm
71	REAL, ALLOCATABLE, DIMENSION(:,:) :: lon, lat
72	REAL, ALLOCATABLE, DIMENSION(:) :: lev
73	LOGICAL :: rectilinear
74	REAL, ALLOCATABLE, DIMENSION(:) :: lon_rect, lat_rect
75	REAL(r_std) :: dt
76	INTEGER(i_std) :: itau_dep,itau_fin,itau,itau_len,itau_step
77	REAL(r_std) :: date0
78	INTEGER(i_std) :: rest_id_sec,rest_id_sto
79	INTEGER(i_std) :: hist_id_sec,hist_id_sec2,hist_id_stom,hist_id_stom_IPCC
80	CHARACTER(LEN=30) :: time_str
81	REAL(r_std) :: dt_slow_
82	REAL :: hist_days_stom,hist_days_stom_ipcc,hist_dt_stom,hist_dt_stom_ipcc
83	REAL(r_std),ALLOCATABLE,DIMENSION(:) :: hist_PFTaxis
84	REAL(r_std),DIMENSION(10) :: hist_pool_10axis
85	REAL(r_std),DIMENSION(100) :: hist_pool_100axis
86	REAL(r_std),DIMENSION(11) :: hist_pool_11axis
87	REAL(r_std),DIMENSION(101) :: hist_pool_101axis
88	INTEGER :: hist_PFTaxis_id,hist_IPCC_PFTaxis_id,hori_id
89	INTEGER :: hist_pool_10axis_id
90	INTEGER :: hist_pool_100axis_id
91	INTEGER :: hist_pool_11axis_id
92	INTEGER :: hist_pool_101axis_id
93	INTEGER(i_std) :: i,j,iv
94	LOGICAL :: ldrestart_read,ldrestart_write
95	LOGICAL :: l1d
96	INTEGER(i_std),PARAMETER :: nbvarmax=200
97	INTEGER(i_std) :: nbvar
98	CHARACTER(LEN=50),DIMENSION(nbvarmax) :: varnames
99	INTEGER(i_std) :: varnbdim
100	INTEGER(i_std),PARAMETER :: varnbdim_max=20
101	INTEGER,DIMENSION(varnbdim_max) :: vardims
102	CHARACTER(LEN=200) :: taboo_vars
103	REAL(r_std),DIMENSION(1) :: xtmp
104	INTEGER :: nsfm,nsft
105	INTEGER :: iisf,iiisf
106	INTEGER(i_std) :: max_totsize,totsize_1step,totsize_tmp
107
108	INTEGER :: vid
109	CHARACTER(LEN=100) :: forcing_name
110	REAL :: x
111
112	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: var_3d
113	REAL(r_std) :: var_1d(1)
114
115	!-
116	REAL(r_std) :: time_sec,time_step_sec
117	REAL(r_std) :: time_dri,time_step_dri
118	REAL(r_std),DIMENSION(1) :: r1d
119	REAL(r_std) :: julian,djulian
120
121	INTEGER(i_std) :: ji,jv,l
122
123	LOGICAL :: debug
124
125	!---------------------------------------------------------------------
126
127	CALL init_para(.FALSE.)
128	CALL init_timer
129
130	IF (is_root_prc) THEN
131	!-
132	! open STOMATE's forcing file to read some basic info
133	!-
134	forcing_name = 'stomate_forcing.nc'
135	CALL getin ('STOMATE_FORCING_NAME',forcing_name)
136	iret = NF90_OPEN (TRIM(forcing_name),NF90_NOWRITE,forcing_id)
137	IF (iret /= NF90_NOERR) THEN
138	CALL ipslerr (3,'teststomate', &
139	& 'Could not open file : ', &
140	& forcing_name,'(Do you have forget it ?)')
141	ENDIF
142	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'dtradia',dtradia)
143	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'dt_slow',dt_force)
144	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'nsft',x)
145	nsft = NINT(x)
146	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'kjpij',x)
147	kjpij = NINT(x)
148	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'kjpindex',x)
149	nbp_glo = NINT(x)
150	ENDIF
151	CALL bcast(dtradia)
152	CALL bcast(dt_force)
153	CALL bcast(nsft)
154	CALL bcast(nbp_glo)
155	!-
156	write(numout,*) 'ATTENTION',dtradia,dt_force
157	!-
158	! read info about land points
159	!-
160	IF (is_root_prc) THEN
161	a_er=.FALSE.
162	ALLOCATE (indices_g(nbp_glo),stat=ier)
163	a_er = a_er .OR. (ier.NE.0)
164	IF (a_er) THEN
165	CALL ipslerr (3,'teststomate', &
166	& 'PROBLEM WITH ALLOCATION', &
167	& 'for local variables 1','')
168	ENDIF
169	!
170	ALLOCATE (x_indices_g(nbp_glo),stat=ier)
171	a_er = a_er .OR. (ier.NE.0)
172	IF (a_er) THEN
173	CALL ipslerr (3,'teststomate', &
174	& 'PROBLEM WITH ALLOCATION', &
175	& 'for global variables 1','')
176	ENDIF
177	ier = NF90_INQ_VARID (forcing_id,'index',vid)
178	IF (ier .NE. 0) THEN
179	CALL ipslerr (3,'teststomate', &
180	& 'PROBLEM WITH ALLOCATION', &
181	& 'for global variables 1','')
182	ENDIF
183	ier = NF90_GET_VAR (forcing_id,vid,x_indices_g)
184	IF (iret /= NF90_NOERR) THEN
185	CALL ipslerr (3,'teststomate', &
186	& 'PROBLEM WITH variable "index" in file ', &
187	& forcing_name,'(check this file)')
188	ENDIF
189	indices_g(:) = NINT(x_indices_g(:))
190	DEALLOCATE (x_indices_g)
191	ELSE
192	ALLOCATE (indices_g(0))
193	ENDIF
194	!---------------------------------------------------------------------
195	!-
196	! set debug to have more information
197	!-
198	!Config Key = DEBUG_INFO
199	!Config Desc = Flag for debug information
200	!Config If = [-]
201	!Config Def = n
202	!Config Help = This option allows to switch on the output of debug
203	!Config information without recompiling the code.
204	!Config Units = [FLAG]
205	!-
206	debug = .FALSE.
207	CALL getin_p('DEBUG_INFO',debug)
208	!
209	!Config Key = LONGPRINT
210	!Config Desc = ORCHIDEE will print more messages
211	!Config If = [-]
212	!Config Def = n
213	!Config Help = This flag permits to print more debug messages in the run.
214	!Config Units = [FLAG]
215	!
216	long_print = .FALSE.
217	CALL getin_p('LONGPRINT',long_print)
218	!-
219	! activate CO2, STOMATE, but not sechiba
220	!-
221	control%river_routing = .FALSE.
222	control%hydrol_cwrr = .FALSE.
223	control%ok_sechiba = .FALSE.
224	!
225	control%stomate_watchout = .TRUE.
226	control%ok_co2 = .TRUE.
227	control%ok_stomate = .TRUE.
228	!-
229	! is DGVM activated?
230	!-
231	control%ok_dgvm = .FALSE.
232	CALL getin_p('STOMATE_OK_DGVM',control%ok_dgvm)
233	WRITE(numout,*) 'LPJ is activated: ',control%ok_dgvm
234
235	!-
236	! Configuration
237	!-
238	! 1. Number of PFTs defined by the user
239
240	!Config Key = NVM
241	!Config Desc = number of PFTs
242	!Config If = OK_SECHIBA or OK_STOMATE
243	!Config Def = 13
244	!Config Help = The number of vegetation types define by the user
245	!Config Units = [-]
246	!
247	CALL getin_p('NVM',nvm)
248	WRITE(numout,*)'the number of pfts is : ', nvm
249
250	! 2. Should we read the parameters in the run.def file ?
251
252	!Config Key = IMPOSE_PARAM
253	!Config Desc = Do you impose the values of the parameters?
254	!Config if = OK_SECHIBA or OK_STOMATE
255	!Config Def = y
256	!Config Help = This flag can deactivate the reading of some parameters.
257	!Config Useful if you want to use the standard values without commenting the run.def
258	!Config Units = [FLAG]
259	!
260	CALL getin_p('IMPOSE_PARAM',impose_param)
261
262	! 3. Allocate and intialize the pft parameters
263
264	CALL pft_parameters_main(control)
265
266	! 4. Activation sub-models of ORCHIDEE
267
268	CALL activate_sub_models(control)
269
270	! 5. Vegetation configuration (impose_veg, land_use, lcchange...previously in slowproc)
271
272	CALL veget_config
273
274	! 6. Read the parameters in the run.def file
275
276	IF ( impose_param) THEN
277	CALL config_pft_parameters
278	CALL config_stomate_pft_parameters
279	CALL config_co2_parameters
280	CALL config_stomate_parameters
281	ENDIF
282	!-
283	IF (control%ok_dgvm) THEN
284	IF ( impose_param ) THEN
285	CALL config_dgvm_parameters
286	ENDIF
287	ENDIF
288	!-
289	!-
290	! restart files
291	!-
292	IF (is_root_prc) THEN
293	! Sechiba's restart files
294	sec_restname_in = 'sechiba_start.nc'
295	CALL getin('SECHIBA_restart_in',sec_restname_in)
296	WRITE(numout,*) 'SECHIBA INPUT RESTART_FILE: ',TRIM(sec_restname_in)
297	IF ( TRIM(sec_restname_in) .EQ. 'NONE' ) THEN
298	STOP 'Need a restart file for Sechiba'
299	ENDIF
300	sec_restname_out = 'sechiba_rest_out.nc'
301	CALL getin('SECHIBA_rest_out',sec_restname_out)
302	WRITE(numout,*) 'SECHIBA OUTPUT RESTART_FILE: ',TRIM(sec_restname_out)
303	! Stomate's restart files
304	sto_restname_in = 'stomate_start.nc'
305	CALL getin('STOMATE_RESTART_FILEIN',sto_restname_in)
306	WRITE(numout,*) 'STOMATE INPUT RESTART_FILE: ',TRIM(sto_restname_in)
307	sto_restname_out = 'stomate_rest_out.nc'
308	CALL getin('STOMATE_RESTART_FILEOUT',sto_restname_out)
309	WRITE(numout,*) 'STOMATE OUTPUT RESTART_FILE: ',TRIM(sto_restname_out)
310
311	!-
312	! We need to know iim, jjm.
313	! Get them from the restart files themselves.
314	!-
315	iret = NF90_OPEN (sec_restname_in,NF90_NOWRITE,ncfid)
316	IF (iret /= NF90_NOERR) THEN
317	CALL ipslerr (3,'teststomate', &
318	& 'Could not open file : ', &
319	& sec_restname_in,'(Do you have forget it ?)')
320	ENDIF
321	iret = NF90_INQUIRE_DIMENSION (ncfid,1,len=iim_g)
322	iret = NF90_INQUIRE_DIMENSION (ncfid,2,len=jjm_g)
323	iret = NF90_INQ_VARID (ncfid, "time", iv)
324	iret = NF90_GET_ATT (ncfid, iv, 'calendar',thecalendar)
325	iret = NF90_CLOSE (ncfid)
326	i=INDEX(thecalendar,ACHAR(0))
327	IF ( i > 0 ) THEN
328	thecalendar(i:20)=' '
329	ENDIF
330	ENDIF
331	CALL bcast(iim_g)
332	CALL bcast(jjm_g)
333	CALL bcast(thecalendar)
334	!-
335	! calendar
336	!-
337	CALL ioconf_calendar (thecalendar)
338	CALL ioget_calendar (one_year,one_day)
339	!
340	! Parallelization :
341	!
342	CALL init_data_para(iim_g,jjm_g,nbp_glo,indices_g)
343	kjpindex=nbp_loc
344	jjm=jj_nb
345	iim=iim_g
346	kjpij=iim*jjm
347	IF (debug) WRITE(numout,*) "Local grid : ",kjpindex,iim,jjm
348	!-
349	!-
350	! read info about grids
351	!-
352	!-
353	llm=1
354	ALLOCATE(lev(llm))
355	IF (is_root_prc) THEN
356	!-
357	ier = NF90_INQ_VARID (forcing_id,'lalo',vid)
358	ier = NF90_GET_VAR (forcing_id,vid,lalo_g)
359	!-
360	ALLOCATE (x_neighbours_g(nbp_glo,8),stat=ier)
361	ier = NF90_INQ_VARID (forcing_id,'neighbours',vid)
362	ier = NF90_GET_VAR (forcing_id,vid,x_neighbours_g)
363	neighbours_g(:,:) = NINT(x_neighbours_g(:,:))
364	DEALLOCATE (x_neighbours_g)
365	!-
366	ier = NF90_INQ_VARID (forcing_id,'resolution',vid)
367	ier = NF90_GET_VAR (forcing_id,vid,resolution_g)
368	!-
369	ier = NF90_INQ_VARID (forcing_id,'contfrac',vid)
370	ier = NF90_GET_VAR (forcing_id,vid,contfrac_g)
371
372	lon_g(:,:) = zero
373	lat_g(:,:) = zero
374	lev(1) = zero
375	!-
376	CALL restini &
377	& (sec_restname_in, iim_g, jjm_g, lon_g, lat_g, llm, lev, &
378	& sec_restname_out, itau_dep, date0, dt, rest_id_sec)
379	!-
380	IF ( dt .NE. dtradia ) THEN
381	WRITE(numout,*) 'dt',dt
382	WRITE(numout,*) 'dtradia',dtradia
383	CALL ipslerr (3,'teststomate', &
384	& 'PROBLEM with time steps.', &
385	& sec_restname_in,'(dt .NE. dtradia)')
386	ENDIF
387	!-
388	CALL restini &
389	& (sto_restname_in, iim_g, jjm_g, lon_g, lat_g, llm, lev, &
390	& sto_restname_out, itau_dep, date0, dt, rest_id_sto)
391	!-
392	IF ( dt .NE. dtradia ) THEN
393	WRITE(numout,*) 'dt',dt
394	WRITE(numout,*) 'dtradia',dtradia
395	CALL ipslerr (3,'teststomate', &
396	& 'PROBLEM with time steps.', &
397	& sto_restname_in,'(dt .NE. dtradia)')
398	ENDIF
399	ENDIF
400	CALL bcast(rest_id_sec)
401	CALL bcast(rest_id_sto)
402	CALL bcast(itau_dep)
403	CALL bcast(date0)
404	CALL bcast(dt)
405	CALL bcast(lev)
406	!---
407	!--- Create the index table
408	!---
409	!--- This job return a LOCAL kindex
410	!---
411	ALLOCATE (indices(kjpindex),stat=ier)
412	IF (debug .AND. is_root_prc) WRITE(numout,*) "indices_g = ",indices_g(1:nbp_glo)
413	CALL scatter(indices_g,indices)
414	indices(1:kjpindex)=indices(1:kjpindex)-(jj_begin-1)*iim_g
415	IF (debug) WRITE(numout,*) "indices = ",indices(1:kjpindex)
416
417	!---
418	!--- initialize global grid
419	!---
420	CALL init_grid( kjpindex )
421	CALL grid_stuff (nbp_glo, iim_g, jjm_g, lon_g, lat_g, indices_g)
422
423	!---
424	!--- initialize local grid
425	!---
426	jlandindex = (((indices(1:kjpindex)-1)/iim) + 1)
427	if (debug) WRITE(numout,*) "jlandindex = ",jlandindex(1:kjpindex)
428	ilandindex = (indices(1:kjpindex) - (jlandindex(1:kjpindex)-1)*iim)
429	IF (debug) WRITE(numout,*) "ilandindex = ",ilandindex(1:kjpindex)
430	ALLOCATE(lon(iim,jjm))
431	ALLOCATE(lat(iim,jjm))
432	lon=zero
433	lat=zero
434	CALL scatter2D(lon_g,lon)
435	CALL scatter2D(lat_g,lat)
436
437	DO ji=1,kjpindex
438
439	j = jlandindex(ji)
440	i = ilandindex(ji)
441
442	!- Create the internal coordinate table
443	!-
444	lalo(ji,1) = lat(i,j)
445	lalo(ji,2) = lon(i,j)
446	ENDDO
447	CALL scatter(neighbours_g,neighbours)
448	CALL scatter(resolution_g,resolution)
449	CALL scatter(contfrac_g,contfrac)
450	CALL scatter(area_g,area)
451	!-
452	!- Check if we have by any change a rectilinear grid. This would allow us to
453	!- simplify the output files.
454	!
455	rectilinear = .FALSE.
456	IF (is_root_prc) THEN
457	IF ( ALL(lon_g(:,:) == SPREAD(lon_g(:,1), 2, SIZE(lon_g,2))) .AND. &
458	& ALL(lat_g(:,:) == SPREAD(lat_g(1,:), 1, SIZE(lat_g,1))) ) THEN
459	rectilinear = .TRUE.
460	ENDIF
461	ENDIF
462	CALL bcast(rectilinear)
463	IF (rectilinear) THEN
464	ALLOCATE(lon_rect(iim),stat=ier)
465	IF (ier .NE. 0) THEN
466	WRITE (numout,*) ' error in lon_rect allocation. We stop. We need iim words = ',iim
467	STOP 'intersurf_history'
468	ENDIF
469	ALLOCATE(lat_rect(jjm),stat=ier)
470	IF (ier .NE. 0) THEN
471	WRITE (numout,*) ' error in lat_rect allocation. We stop. We need jjm words = ',jjm
472	STOP 'intersurf_history'
473	ENDIF
474	lon_rect(:) = lon(:,1)
475	lat_rect(:) = lat(1,:)
476	ENDIF
477	!-
478	! allocate arrays
479	!-
480	!
481	a_er = .FALSE.
482	ALLOCATE (hist_PFTaxis(nvm),stat=ier)
483	a_er = a_er .OR. (ier.NE.0)
484	ALLOCATE (indexveg(kjpindex*nvm), stat=ier)
485	a_er = a_er .OR. (ier.NE.0)
486	ALLOCATE (soiltype(kjpindex,nstm),stat=ier)
487	a_er = a_er .OR. (ier.NE.0)
488	ALLOCATE (veget_x(kjpindex,nvm),stat=ier)
489	a_er = a_er .OR. (ier.NE.0)
490	ALLOCATE (totfrac_nobio(kjpindex),stat=ier)
491	a_er = a_er .OR. (ier.NE.0)
492	ALLOCATE (frac_nobio(kjpindex,nnobio),stat=ier)
493	a_er = a_er .OR. (ier.NE.0)
494	ALLOCATE (veget_max_x(kjpindex,nvm),stat=ier)
495	a_er = a_er .OR. (ier.NE.0)
496	ALLOCATE (lai_x(kjpindex,nvm),stat=ier)
497	a_er = a_er .OR. (ier.NE.0)
498	ALLOCATE (veget_force_x(kjpindex,nvm),stat=ier)
499	a_er = a_er .OR. (ier.NE.0)
500	ALLOCATE (veget_max_force_x(kjpindex,nvm),stat=ier)
501	a_er = a_er .OR. (ier.NE.0)
502	ALLOCATE (lai_force_x(kjpindex,nvm),stat=ier)
503	a_er = a_er .OR. (ier.NE.0)
504	ALLOCATE (t2m(kjpindex),stat=ier)
505	a_er = a_er .OR. (ier.NE.0)
506	ALLOCATE (t2m_min(kjpindex),stat=ier)
507	a_er = a_er .OR. (ier.NE.0)
508	ALLOCATE (temp_sol(kjpindex),stat=ier)
509	a_er = a_er .OR. (ier.NE.0)
510	ALLOCATE (soiltemp(kjpindex,nbdl),stat=ier)
511	a_er = a_er .OR. (ier.NE.0)
512	ALLOCATE (soilhum(kjpindex,nbdl),stat=ier)
513	a_er = a_er .OR. (ier.NE.0)
514	ALLOCATE (humrel_x(kjpindex,nvm),stat=ier)
515	a_er = a_er .OR. (ier.NE.0)
516	ALLOCATE (litterhum(kjpindex),stat=ier)
517	a_er = a_er .OR. (ier.NE.0)
518	ALLOCATE (precip_rain(kjpindex),stat=ier)
519	a_er = a_er .OR. (ier.NE.0)
520	ALLOCATE (precip_snow(kjpindex),stat=ier)
521	a_er = a_er .OR. (ier.NE.0)
522	ALLOCATE (gpp_x(kjpindex,nvm),stat=ier)
523	a_er = a_er .OR. (ier.NE.0)
524	ALLOCATE (deadleaf_cover(kjpindex),stat=ier)
525	a_er = a_er .OR. (ier.NE.0)
526	ALLOCATE (assim_param_x(kjpindex,nvm,npco2),stat=ier)
527	a_er = a_er .OR. (ier.NE.0)
528	ALLOCATE (height_x(kjpindex,nvm),stat=ier)
529	a_er = a_er .OR. (ier.NE.0)
530	ALLOCATE (qsintmax_x(kjpindex,nvm),stat=ier)
531	a_er = a_er .OR. (ier.NE.0)
532	ALLOCATE (co2_flux(kjpindex,nvm),stat=ier)
533	a_er = a_er .OR. (ier.NE.0)
534	ALLOCATE (fco2_lu(kjpindex),stat=ier)
535	a_er = a_er .OR. (ier.NE.0)
536	IF (a_er) THEN
537	CALL ipslerr (3,'teststomate', &
538	& 'PROBLEM WITH ALLOCATION', &
539	& 'for local variables 1','')
540	ENDIF
541	!
542	! prepare forcing
543	!
544	max_totsize = 50
545	CALL getin_p ('STOMATE_FORCING_MEMSIZE',max_totsize)
546	max_totsize = max_totsize * 1000000
547
548	totsize_1step = SIZE(soiltype(:,3))*KIND(soiltype(:,3)) + &
549	SIZE(humrel_x)*KIND(humrel_x) + &
550	SIZE(litterhum)*KIND(litterhum) + &
551	SIZE(t2m)*KIND(t2m) + &
552	SIZE(t2m_min)*KIND(t2m_min) + &
553	SIZE(temp_sol)*KIND(temp_sol) + &
554	SIZE(soiltemp)*KIND(soiltemp) + &
555	SIZE(soilhum)*KIND(soilhum) + &
556	SIZE(precip_rain)*KIND(precip_rain) + &
557	SIZE(precip_snow)*KIND(precip_snow) + &
558	SIZE(gpp_x)*KIND(gpp_x) + &
559	SIZE(veget_force_x)*KIND(veget_force_x) + &
560	SIZE(veget_max_force_x)*KIND(veget_max_force_x) + &
561	SIZE(lai_force_x)*KIND(lai_force_x)
562	CALL reduce_sum(totsize_1step,totsize_tmp)
563	CALL bcast(totsize_tmp)
564	totsize_1step=totsize_tmp
565
566	! total number of forcing steps
567	IF ( nsft .NE. INT(one_year/(dt_force/one_day)) ) THEN
568	CALL ipslerr (3,'teststomate', &
569	& 'stomate: error with total number of forcing steps', &
570	& 'nsft','teststomate computation different with forcing file value.')
571	ENDIF
572	! number of forcing steps in memory
573	nsfm = MIN(nsft, &
574	& MAX(1,NINT( REAL(max_totsize,r_std) &
575	& /REAL(totsize_1step,r_std))))
576	!-
577	WRITE(numout,*) 'Offline forcing of Stomate:'
578	WRITE(numout,*) ' Total number of forcing steps:',nsft
579	WRITE(numout,*) ' Number of forcing steps in memory:',nsfm
580	!-
581	CALL init_forcing(kjpindex,nsfm,nsft)
582	!-
583	! ensure that we read all new forcing states
584	iisf = nsfm
585	! initialize that table that contains the indices
586	! of the forcing states that will be in memory
587	isf(:) = (/ (i,i=1,nsfm) /)
588
589	nf_written(:) = .FALSE.
590	nf_written(isf(:)) = .TRUE.
591
592	!-
593	! a time step for STOMATE corresponds to itau_step timesteps in SECHIBA
594	!-
595	itau_step = NINT(dt_force/dtradia)
596	IF (debug) WRITE(numout,*) "dtradia, dt_rest, dt_force, itau_step",dtradia, dt, dt_force, itau_step
597	!
598	CALL ioconf_startdate(date0)
599	CALL intsurf_time( itau_dep, date0, dtradia )
600	!-
601	! Length of integration
602	!-
603	WRITE(time_str,'(a)') '1Y'
604	CALL getin_p ('TIME_LENGTH', time_str)
605	! transform into itau
606	CALL tlen2itau(time_str, dt, date0, itau_len)
607	! itau_len*dtradia must be a multiple of dt_force
608	itau_len = NINT( MAX(1.,FLOAT(NINT(itau_len*dtradia/dt_force))) &
609	& *dt_force/dtradia)
610	!-
611	itau_fin = itau_dep+itau_len
612	!-
613	! set up STOMATE history file
614	!-
615	!Config Key = STOMATE_OUTPUT_FILE
616	!Config Desc = Name of file in which STOMATE's output is going to be written
617	!Config If =
618	!Config Def = stomate_history.nc
619	!Config Help = This file is going to be created by the model
620	!Config and will contain the output from the model.
621	!Config This file is a truly COADS compliant netCDF file.
622	!Config It will be generated by the hist software from
623	!Config the IOIPSL package.
624	!Config Units = FILE
625	!-
626	stom_histname='stomate_history.nc'
627	CALL getin_p ('STOMATE_OUTPUT_FILE', stom_histname)
628	WRITE(numout,*) 'STOMATE_OUTPUT_FILE', TRIM(stom_histname)
629	!-
630	!Config Key = STOMATE_HIST_DT
631	!Config Desc = STOMATE history time step
632	!Config If =
633	!Config Def = 10.
634	!Config Help = Time step of the STOMATE history file
635	!Config Units = days [d]
636	!-
637	hist_days_stom = 10.
638	CALL getin_p ('STOMATE_HIST_DT', hist_days_stom)
639	IF ( hist_days_stom == -1. ) THEN
640	hist_dt_stom = -1.
641	WRITE(numout,*) 'output frequency for STOMATE history file (d): one month.'
642	ELSE
643	hist_dt_stom = NINT( hist_days_stom ) * one_day
644	WRITE(numout,*) 'output frequency for STOMATE history file (d): ', &
645	hist_dt_stom/one_day
646	ENDIF
647	!-
648	!-
649	!- initialize
650	WRITE(numout,*) "before histbeg : ",date0,dt
651	IF ( rectilinear ) THEN
652	#ifdef CPP_PARA
653	CALL histbeg(stom_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
654	& itau_dep, date0, dt, hori_id, hist_id_stom, domain_id=orch_domain_id)
655	#else
656	CALL histbeg(stom_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
657	& itau_dep, date0, dt, hori_id, hist_id_stom)
658	#endif
659	ELSE
660	#ifdef CPP_PARA
661	CALL histbeg(stom_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
662	& itau_dep, date0, dt, hori_id, hist_id_stom, domain_id=orch_domain_id)
663	#else
664	CALL histbeg(stom_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
665	& itau_dep, date0, dt, hori_id, hist_id_stom)
666	#endif
667	ENDIF
668	!- define PFT axis
669	hist_PFTaxis = (/ ( REAL(i,r_std), i=1,nvm ) /)
670	!- declare this axis
671	CALL histvert (hist_id_stom, 'PFT', 'Plant functional type', &
672	& '1', nvm, hist_PFTaxis, hist_PFTaxis_id)
673	!- define Pool_10 axis
674	hist_pool_10axis = (/ ( REAL(i,r_std), i=1,10 ) /)
675	!- declare this axis
676	CALL histvert (hist_id_stom, 'P10', 'Pool 10 years', &
677	& '1', 10, hist_pool_10axis, hist_pool_10axis_id)
678
679	!- define Pool_100 axis
680	hist_pool_100axis = (/ ( REAL(i,r_std), i=1,100 ) /)
681	!- declare this axis
682	CALL histvert (hist_id_stom, 'P100', 'Pool 100 years', &
683	& '1', 100, hist_pool_100axis, hist_pool_100axis_id)
684
685	!- define Pool_11 axis
686	hist_pool_11axis = (/ ( REAL(i,r_std), i=1,11 ) /)
687	!- declare this axis
688	CALL histvert (hist_id_stom, 'P11', 'Pool 10 years + 1', &
689	& '1', 11, hist_pool_11axis, hist_pool_11axis_id)
690
691	!- define Pool_101 axis
692	hist_pool_101axis = (/ ( REAL(i,r_std), i=1,101 ) /)
693	!- declare this axis
694	CALL histvert (hist_id_stom, 'P101', 'Pool 100 years + 1', &
695	& '1', 101, hist_pool_101axis, hist_pool_101axis_id)
696
697	!- define STOMATE history file
698	CALL stom_define_history (hist_id_stom, nvm, iim, jjm, &
699	& dt, hist_dt_stom, hori_id, hist_PFTaxis_id, &
700	& hist_pool_10axis_id, hist_pool_100axis_id, &
701	& hist_pool_11axis_id, hist_pool_101axis_id)
702
703	!- end definition
704	CALL histend(hist_id_stom)
705	!-
706	!-
707	! STOMATE IPCC OUTPUTS IS ACTIVATED
708	!-
709	!Config Key = STOMATE_IPCC_OUTPUT_FILE
710	!Config Desc = Name of file in which STOMATE's output is going to be written
711	!Config If =
712	!Config Def = stomate_ipcc_history.nc
713	!Config Help = This file is going to be created by the model
714	!Config and will contain the output from the model.
715	!Config This file is a truly COADS compliant netCDF file.
716	!Config It will be generated by the hist software from
717	!Config the IOIPSL package.
718	!Config Units = FILE
719	!-
720	stom_ipcc_histname='stomate_ipcc_history.nc'
721	CALL getin_p('STOMATE_IPCC_OUTPUT_FILE', stom_ipcc_histname)
722	WRITE(numout,*) 'STOMATE_IPCC_OUTPUT_FILE', TRIM(stom_ipcc_histname)
723	!-
724	!Config Key = STOMATE_IPCC_HIST_DT
725	!Config Desc = STOMATE IPCC history time step
726	!Config If =
727	!Config Def = 0.
728	!Config Help = Time step of the STOMATE IPCC history file
729	!Config Units = days [d]
730	!-
731	hist_days_stom_ipcc = zero
732	CALL getin_p('STOMATE_IPCC_HIST_DT', hist_days_stom_ipcc)
733	IF ( hist_days_stom_ipcc == moins_un ) THEN
734	hist_dt_stom_ipcc = moins_un
735	WRITE(numout,*) 'output frequency for STOMATE IPCC history file (d): one month.'
736	ELSE
737	hist_dt_stom_ipcc = NINT( hist_days_stom_ipcc ) * one_day
738	WRITE(numout,*) 'output frequency for STOMATE IPCC history file (d): ', &
739	hist_dt_stom_ipcc/one_day
740	ENDIF
741
742	! test consistency between STOMATE_IPCC_HIST_DT and DT_SLOW parameters
743	dt_slow_ = one_day
744	CALL getin_p('DT_SLOW', dt_slow_)
745	IF ( hist_days_stom_ipcc > zero ) THEN
746	IF (dt_slow_ > hist_dt_stom_ipcc) THEN
747	WRITE(numout,*) "DT_SLOW = ",dt_slow_," , STOMATE_IPCC_HIST_DT = ",hist_dt_stom_ipcc
748	CALL ipslerr (3,'intsurf_history', &
749	& 'Problem with DT_SLOW > STOMATE_IPCC_HIST_DT','', &
750	& '(must be less or equal)')
751	ENDIF
752	ENDIF
753
754	IF ( hist_dt_stom_ipcc == 0 ) THEN
755	hist_id_stom_ipcc = -1
756	ELSE
757	!-
758	!- initialize
759	IF ( rectilinear ) THEN
760	#ifdef CPP_PARA
761	CALL histbeg(stom_ipcc_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
762	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc,domain_id=orch_domain_id)
763	#else
764	CALL histbeg(stom_ipcc_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
765	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc)
766	#endif
767	ELSE
768	#ifdef CPP_PARA
769	CALL histbeg(stom_ipcc_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
770	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc,domain_id=orch_domain_id)
771	#else
772	CALL histbeg(stom_ipcc_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
773	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc)
774	#endif
775	ENDIF
776	!- declare this axis
777	CALL histvert (hist_id_stom_IPCC, 'PFT', 'Plant functional type', &
778	& '1', nvm, hist_PFTaxis, hist_IPCC_PFTaxis_id)
779
780	!- define STOMATE history file
781	CALL stom_IPCC_define_history (hist_id_stom_IPCC, nvm, iim, jjm, &
782	& dt, hist_dt_stom_ipcc, hori_id, hist_IPCC_PFTaxis_id)
783
784	!- end definition
785	CALL histend(hist_id_stom_IPCC)
786
787	ENDIF
788	!
789	CALL histwrite(hist_id_stom, 'Areas', itau_dep+itau_step, area, kjpindex, indices)
790	IF ( hist_id_stom_IPCC > 0 ) THEN
791	CALL histwrite(hist_id_stom_IPCC, 'Areas', itau_dep+itau_step, area, kjpindex, indices)
792	ENDIF
793	!
794	hist_id_sec = -1
795	hist_id_sec2 = -1
796	!-
797	! first call of slowproc to initialize variables
798	!-
799	itau = itau_dep
800	!
801	DO ji=1,kjpindex
802	DO jv=1,nvm
803	indexveg((jv-1)kjpindex + ji) = indices(ji) + (jv-1)kjpij
804	ENDDO
805	ENDDO
806	!-
807	!MM Problem here with dpu which depends on soil type
808	DO l = 1, nbdl-1
809	! first 2.0 is dpu
810	! second 2.0 is average
811	diaglev(l) = dpu_cste/(2*(nbdl-1) -1) ( ( 2(l-1) -1) + ( 2(l) -1) ) / 2.0
812	ENDDO
813	diaglev(nbdl) = dpu_cste
814	!
815	!-
816	! Read the parameters in the "*.def" files
817	!-
818	!
819	!Config Key = CLAYFRACTION_DEFAULT
820	!Config Desc =
821	!Config If = OK_SECHIBA
822	!Config Def = 0.2
823	!Config Help =
824	!Config Units = [-]
825	CALL getin_p('CLAYFRACTION_DEFAULT',clayfraction_default)
826	!
827	!Config Key = MIN_VEGFRAC
828	!Config Desc = Minimal fraction of mesh a vegetation type can occupy
829	!Config If = OK_SECHIBA
830	!Config Def = 0.001
831	!Config Help =
832	!Config Units = [-]
833	CALL getin_p('MIN_VEGFRAC',min_vegfrac)
834	!
835	!Config Key = STEMPDIAG_BID
836	!Config Desc = only needed for an initial LAI if there is no restart file
837	!Config If = OK_SECHIBA
838	!Config Def = 280.
839	!Config Help =
840	!Config Units = [K]
841	CALL getin_p('STEMPDIAG_BID',stempdiag_bid)
842	!
843	!Config Key = SOILTYPE_DEFAULT
844	!Config Desc = Default soil texture distribution in the following order : sand, loam and clay
845	!Config If = OK_SECHIBA
846	!Config Def = 0.0, 1.0, 0.0
847	!Config Help =
848	!Config Units = [-]
849	CALL getin_p('SOILTYPE_DEFAULT',soiltype_default)
850	!-
851	CALL ioget_expval(val_exp)
852	ldrestart_read = .TRUE.
853	ldrestart_write = .FALSE.
854	!-
855	! read some variables we need from SECHIBA's restart file
856	!-
857	CALL slowproc_main &
858	& (itau, kjpij, kjpindex, dt_force, date0, &
859	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
860	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
861	& t2m, t2m_min, temp_sol, soiltemp, &
862	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
863	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
864	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
865	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
866	! correct date
867	day_counter = one_day - dt_force
868	date=1
869	!-
870	! time loop
871	!-
872	IF (debug) check_time=.TRUE.
873	CALL intsurf_time( itau_dep+itau_step, date0, dtradia )
874	l_first_intersurf=.FALSE.
875	!
876	DO itau = itau_dep+itau_step,itau_fin,itau_step
877	!
878	CALL intsurf_time( itau, date0, dtradia )
879	!
880	!-- next forcing state
881	iisf = iisf+1
882	IF (debug) WRITE(numout,*) "itau,iisf : ",itau,iisf
883	!---
884	IF (iisf .GT. nsfm) THEN
885	!-----
886	!---- we have to read new forcing states from the file
887	!-----
888	!---- determine blocks of forcing states that are contiguous in memory
889	!-----
890	CALL forcing_read(forcing_id,nsfm)
891
892	!--------------------------
893
894	!-----
895	!---- determine which forcing states must be read next time
896	!-----
897	isf(1) = isf(nsfm)+1
898	IF ( isf(1) .GT. nsft ) isf(1) = 1
899	DO iiisf = 2, nsfm
900	isf(iiisf) = isf(iiisf-1)+1
901	IF ( isf(iiisf) .GT. nsft ) isf(iiisf) = 1
902	ENDDO
903	nf_written(isf(:)) = .TRUE.
904	!---- start again at first forcing state
905	iisf = 1
906	ENDIF
907	! Bug here ! soiltype(:,3) != clay
908	! soiltype(:,3) = clay_fm(:,iisf)
909	humrel_x(:,:) = humrel_daily_fm(:,:,iisf)
910	litterhum(:) = litterhum_daily_fm(:,iisf)
911	t2m(:) = t2m_daily_fm(:,iisf)
912	t2m_min(:) = t2m_min_daily_fm(:,iisf)
913	temp_sol(:) = tsurf_daily_fm(:,iisf)
914	soiltemp(:,:) = tsoil_daily_fm(:,:,iisf)
915	soilhum(:,:) = soilhum_daily_fm(:,:,iisf)
916	precip_rain(:) = precip_fm(:,iisf)
917	gpp_x(:,:) = gpp_daily_fm(:,:,iisf)
918	veget_force_x(:,:) = veget_fm(:,:,iisf)
919	veget_max_force_x(:,:) = veget_max_fm(:,:,iisf)
920	lai_force_x(:,:) = lai_fm(:,:,iisf)
921	WHERE ( t2m(:) .LT. ZeroCelsius )
922	precip_snow(:) = precip_rain(:)
923	precip_rain(:) = zero
924	ELSEWHERE
925	precip_snow(:) = zero
926	ENDWHERE
927	!---
928	!-- scale GPP to new lai and veget_max
929	!---
930	WHERE ( lai_x(:,:) .EQ. zero ) gpp_x(:,:) = zero
931	!-- scale GPP to new LAI
932	WHERE (lai_force_x(:,:) .GT. zero )
933	gpp_x(:,:) = gpp_x(:,:)ATAN(2.lai_x(:,:)) &
934	& /ATAN( 2.*MAX(lai_force_x(:,:),0.01))
935	ENDWHERE
936	!-- scale GPP to new veget_max
937	WHERE (veget_max_force_x(:,:) .GT. zero )
938	gpp_x(:,:) = gpp_x(:,:)*veget_max_x(:,:)/veget_max_force_x(:,:)
939	ENDWHERE
940	!---
941	!-- number crunching
942	!---
943	ldrestart_read = .FALSE.
944	ldrestart_write = .FALSE.
945	CALL slowproc_main &
946	& (itau, kjpij, kjpindex, dt_force, date0, &
947	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
948	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
949	& t2m, t2m_min, temp_sol, soiltemp, &
950	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
951	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
952	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
953	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
954	day_counter = one_day - dt_force
955	ENDDO ! end of the time loop
956	!-
957	! write restart files
958	!-
959	IF (is_root_prc) THEN
960	! first, read and write variables that are not managed otherwise
961	taboo_vars = &
962	& '$lat$ $lon$ $lev$ $veget_year$ '// &
963	& '$height$ $day_counter$ $veget$ $veget_max$ $frac_nobio$ '// &
964	& '$lai$ $soiltype_frac$ $clay_frac$ '// &
965	& '$nav_lon$ $nav_lat$ $nav_lev$ $time$ $time_steps$'
966	!-
967	CALL ioget_vname(rest_id_sec, nbvar, varnames)
968	!-
969	DO iv = 1, nbvar
970	!-- check if the variable is to be written here
971	IF (INDEX( taboo_vars,'$'//TRIM(varnames(iv))//'$') .EQ. 0 ) THEN
972	IF (debug) WRITE(numout,*) "restart var : ",TRIM(varnames(iv)),itau_dep,itau_fin
973
974	!---- get variable dimensions, especially 3rd dimension
975	CALL ioget_vdim &
976	& (rest_id_sec, varnames(iv), varnbdim_max, varnbdim, vardims)
977	l1d = ALL(vardims(1:varnbdim) .EQ. 1)
978	!---- read it
979	IF (l1d) THEN
980	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
981	1,1,1,itau_dep,.TRUE.,var_1d)
982	ELSE
983	ALLOCATE(var_3d(nbp_glo,vardims(3)),stat=ier)
984	IF (ier .NE. 0) STOP 'ALLOCATION PROBLEM'
985	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
986	nbp_glo,vardims(3),1,itau_dep,.TRUE.,var_3d, &
987	"gather",nbp_glo,indices_g)
988	ENDIF
989	!---- write it
990	IF (l1d) THEN
991	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
992	1,1,1,itau_fin,var_1d)
993	ELSE
994	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
995	nbp_glo,vardims(3),1,itau_fin,var_3d, &
996	'scatter',nbp_glo,indices_g)
997	DEALLOCATE (var_3d)
998	ENDIF
999	ENDIF
1000	ENDDO
1001	ENDIF
1002	CALL barrier_para()
1003
1004	! call slowproc to write restart files
1005	ldrestart_read = .FALSE.
1006	ldrestart_write = .TRUE.
1007	!-
1008	IF (debug) WRITE(numout,*) "Call slowproc for restart."
1009	CALL slowproc_main &
1010	& (itau_fin, kjpij, kjpindex, dt_force, date0, &
1011	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
1012	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
1013	& t2m, t2m_min, temp_sol, soiltemp, &
1014	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
1015	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
1016	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
1017	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
1018	!-
1019	! close files
1020	!-
1021	IF (is_root_prc) THEN
1022	CALL restclo
1023	IF ( debug ) WRITE(numout,*) 'REST CLOSED'
1024	ENDIF
1025	CALL histclo
1026
1027	IF (is_root_prc) &
1028	ier = NF90_CLOSE (forcing_id)
1029
1030	IF (is_root_prc) THEN
1031	CALL getin_dump()
1032	ENDIF
1033	#ifdef CPP_PARA
1034	CALL MPI_FINALIZE(ier)
1035	#endif
1036	WRITE(numout,*) "End of teststomate."
1037
1038	!---------------
1039	END PROGRAM teststomate
1040	!
1041	!===
1042	!

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

Download in other formats: