Context Navigation

source: branches/ORCHIDEE_EXT/ORCHIDEE_OL/teststomate.f90 @ 541

Last change on this file since 541 was 508, checked in by didier.solyga, 13 years ago
Correct teststomate bug : double allocation of the variable fco2_lu
File size: 35.1 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	!!$ >> DS for the externalized version
84	REAL(r_std),ALLOCATABLE,DIMENSION(:) :: hist_PFTaxis
85	!!$ >> DS
86	REAL(r_std),DIMENSION(10) :: hist_pool_10axis
87	REAL(r_std),DIMENSION(100) :: hist_pool_100axis
88	REAL(r_std),DIMENSION(11) :: hist_pool_11axis
89	REAL(r_std),DIMENSION(101) :: hist_pool_101axis
90	INTEGER :: hist_PFTaxis_id,hist_IPCC_PFTaxis_id,hori_id
91	INTEGER :: hist_pool_10axis_id
92	INTEGER :: hist_pool_100axis_id
93	INTEGER :: hist_pool_11axis_id
94	INTEGER :: hist_pool_101axis_id
95	INTEGER(i_std) :: i,j,iv
96	LOGICAL :: ldrestart_read,ldrestart_write
97	LOGICAL :: l1d
98	INTEGER(i_std),PARAMETER :: nbvarmax=200
99	INTEGER(i_std) :: nbvar
100	CHARACTER*50,DIMENSION(nbvarmax) :: varnames
101	INTEGER(i_std) :: varnbdim
102	INTEGER(i_std),PARAMETER :: varnbdim_max=20
103	INTEGER,DIMENSION(varnbdim_max) :: vardims
104	CHARACTER*200 :: taboo_vars
105	REAL(r_std),DIMENSION(1) :: xtmp
106	INTEGER :: nsfm,nsft
107	INTEGER :: iisf,iiisf
108	INTEGER(i_std) :: max_totsize,totsize_1step,totsize_tmp
109
110	INTEGER :: vid
111	CHARACTER(LEN=100) :: forcing_name
112	REAL :: x
113
114	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: var_3d
115	REAL(r_std) :: var_1d(1)
116
117	!-
118	REAL(r_std) :: time_sec,time_step_sec
119	REAL(r_std) :: time_dri,time_step_dri
120	REAL(r_std),DIMENSION(1) :: r1d
121	REAL(r_std) :: julian,djulian
122
123	INTEGER(i_std) :: ji,jv,l
124
125	LOGICAL :: debug
126
127	!---------------------------------------------------------------------
128
129	CALL init_para(.FALSE.)
130	CALL init_timer
131
132	IF (is_root_prc) THEN
133	!-
134	! open STOMATE's forcing file to read some basic info
135	!-
136	forcing_name = 'stomate_forcing.nc'
137	CALL getin ('STOMATE_FORCING_NAME',forcing_name)
138	iret = NF90_OPEN (TRIM(forcing_name),NF90_NOWRITE,forcing_id)
139	IF (iret /= NF90_NOERR) THEN
140	CALL ipslerr (3,'teststomate', &
141	& 'Could not open file : ', &
142	& forcing_name,'(Do you have forget it ?)')
143	ENDIF
144	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'dtradia',dtradia)
145	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'dt_slow',dt_force)
146	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'nsft',x)
147	nsft = NINT(x)
148	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'kjpij',x)
149	kjpij = NINT(x)
150	ier = NF90_GET_ATT (forcing_id,NF90_GLOBAL,'kjpindex',x)
151	nbp_glo = NINT(x)
152	ENDIF
153	CALL bcast(dtradia)
154	CALL bcast(dt_force)
155	CALL bcast(nsft)
156	CALL bcast(nbp_glo)
157	!-
158	write(numout,*) 'ATTENTION',dtradia,dt_force
159	!-
160	! read info about land points
161	!-
162	IF (is_root_prc) THEN
163	a_er=.FALSE.
164	ALLOCATE (indices_g(nbp_glo),stat=ier)
165	a_er = a_er .OR. (ier.NE.0)
166	IF (a_er) THEN
167	CALL ipslerr (3,'teststomate', &
168	& 'PROBLEM WITH ALLOCATION', &
169	& 'for local variables 1','')
170	ENDIF
171	!
172	ALLOCATE (x_indices_g(nbp_glo),stat=ier)
173	a_er = a_er .OR. (ier.NE.0)
174	IF (a_er) THEN
175	CALL ipslerr (3,'teststomate', &
176	& 'PROBLEM WITH ALLOCATION', &
177	& 'for global variables 1','')
178	ENDIF
179	ier = NF90_INQ_VARID (forcing_id,'index',vid)
180	IF (ier .NE. 0) THEN
181	CALL ipslerr (3,'teststomate', &
182	& 'PROBLEM WITH ALLOCATION', &
183	& 'for global variables 1','')
184	ENDIF
185	ier = NF90_GET_VAR (forcing_id,vid,x_indices_g)
186	IF (iret /= NF90_NOERR) THEN
187	CALL ipslerr (3,'teststomate', &
188	& 'PROBLEM WITH variable "index" in file ', &
189	& forcing_name,'(check this file)')
190	ENDIF
191	indices_g(:) = NINT(x_indices_g(:))
192	DEALLOCATE (x_indices_g)
193	ELSE
194	ALLOCATE (indices_g(0))
195	ENDIF
196	!---------------------------------------------------------------------
197	!-
198	! set debug to have more information
199	!-
200	!Config Key = DEBUG_INFO
201	!Config Desc = Flag for debug information
202	!Config Def = n
203	!Config Help = This option allows to switch on the output of debug
204	!Config information without recompiling the code.
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 Def = n
212	!Config Help = This flag permits to print more debug messages in the run.
213	!
214	long_print = .FALSE.
215	CALL getin_p('LONGPRINT',long_print)
216
217	!-
218	! activate CO2, STOMATE, but not sechiba
219	!-
220	control%river_routing = .FALSE.
221	control%hydrol_cwrr = .FALSE.
222	control%ok_sechiba = .FALSE.
223	!
224	control%stomate_watchout = .TRUE.
225	control%ok_co2 = .TRUE.
226	control%ok_stomate = .TRUE.
227	!-
228	! is DGVM activated?
229	!-
230	control%ok_dgvm = .FALSE.
231	CALL getin_p('STOMATE_OK_DGVM',control%ok_dgvm)
232	WRITE(numout,*) 'LPJ is activated: ',control%ok_dgvm
233
234	!!$ >> DS added for the externalised version
235
236	! 1. Number of PFTs defined by the user
237
238	!Config Key = NVM
239	!Config Desc = number of PFTs
240	!Config if = ANYTIME
241	!Config Def = 13
242	!Config Help = The number of vegetation types define by the user
243	!Config Units = NONE
244	!
245	CALL getin_p('NVM',nvm)
246	WRITE(numout,*)'the number of pfts is : ', nvm
247
248	! 2. Should we read the parameters in the run.def file ?
249
250	!Config Key = IMPOSE_PARAM
251	!Config Desc = Do you impose the values of the parameters?
252	!Config if = ANYTIME
253	!Config Def = y
254	!Config Help = This flag can deactivate the reading of some parameters.
255	! Useful if you want to use the standard values without commenting the run.def
256	!Config Units = NONE
257	!
258	CALL getin_p('IMPOSE_PARAM',impose_param)
259
260	! 3. Allocate and intialize the pft parameters
261
262	CALL pft_parameters_main(control)
263
264	! 4. Activation sub-models of ORCHIDEE
265
266	CALL activate_sub_models(control)
267
268	! 5. Vegetation configuration (impose_veg, land_use, lcchange...previously in slowproc)
269
270	CALL veget_config
271
272	! 6. Read the parameters in the run.def file
273
274	IF ( impose_param) THEN
275	CALL getin_pft_parameters
276	CALL getin_stomate_pft_parameters
277	CALL getin_co2_parameters
278	CALL getin_stomate_parameters
279	ENDIF
280	!-
281	IF (control%ok_dgvm) THEN
282	IF ( impose_param ) THEN
283	CALL getin_dgvm_parameters
284	ENDIF
285	ENDIF
286
287	!!$ >> DS
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(:,:) = 0.0
373	lat_g(:,:) = 0.0
374	lev(1) = 0.0
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	!!$ >> DS added for the externalised version
483	ALLOCATE (hist_PFTaxis(nvm),stat=ier)
484	a_er = a_er .OR. (ier.NE.0)
485	!!$ end add DS
486	ALLOCATE (indexveg(kjpindex*nvm), stat=ier)
487	a_er = a_er .OR. (ier.NE.0)
488	ALLOCATE (soiltype(kjpindex,nstm),stat=ier)
489	a_er = a_er .OR. (ier.NE.0)
490	ALLOCATE (veget_x(kjpindex,nvm),stat=ier)
491	a_er = a_er .OR. (ier.NE.0)
492	ALLOCATE (totfrac_nobio(kjpindex),stat=ier)
493	a_er = a_er .OR. (ier.NE.0)
494	ALLOCATE (frac_nobio(kjpindex,nnobio),stat=ier)
495	a_er = a_er .OR. (ier.NE.0)
496	ALLOCATE (veget_max_x(kjpindex,nvm),stat=ier)
497	a_er = a_er .OR. (ier.NE.0)
498	ALLOCATE (lai_x(kjpindex,nvm),stat=ier)
499	a_er = a_er .OR. (ier.NE.0)
500	ALLOCATE (veget_force_x(kjpindex,nvm),stat=ier)
501	a_er = a_er .OR. (ier.NE.0)
502	ALLOCATE (veget_max_force_x(kjpindex,nvm),stat=ier)
503	a_er = a_er .OR. (ier.NE.0)
504	ALLOCATE (lai_force_x(kjpindex,nvm),stat=ier)
505	a_er = a_er .OR. (ier.NE.0)
506	ALLOCATE (t2m(kjpindex),stat=ier)
507	a_er = a_er .OR. (ier.NE.0)
508	ALLOCATE (t2m_min(kjpindex),stat=ier)
509	a_er = a_er .OR. (ier.NE.0)
510	ALLOCATE (temp_sol(kjpindex),stat=ier)
511	a_er = a_er .OR. (ier.NE.0)
512	ALLOCATE (soiltemp(kjpindex,nbdl),stat=ier)
513	a_er = a_er .OR. (ier.NE.0)
514	ALLOCATE (soilhum(kjpindex,nbdl),stat=ier)
515	a_er = a_er .OR. (ier.NE.0)
516	ALLOCATE (humrel_x(kjpindex,nvm),stat=ier)
517	a_er = a_er .OR. (ier.NE.0)
518	ALLOCATE (litterhum(kjpindex),stat=ier)
519	a_er = a_er .OR. (ier.NE.0)
520	ALLOCATE (precip_rain(kjpindex),stat=ier)
521	a_er = a_er .OR. (ier.NE.0)
522	ALLOCATE (precip_snow(kjpindex),stat=ier)
523	a_er = a_er .OR. (ier.NE.0)
524	ALLOCATE (gpp_x(kjpindex,nvm),stat=ier)
525	a_er = a_er .OR. (ier.NE.0)
526	ALLOCATE (deadleaf_cover(kjpindex),stat=ier)
527	a_er = a_er .OR. (ier.NE.0)
528	ALLOCATE (assim_param_x(kjpindex,nvm,npco2),stat=ier)
529	a_er = a_er .OR. (ier.NE.0)
530	ALLOCATE (height_x(kjpindex,nvm),stat=ier)
531	a_er = a_er .OR. (ier.NE.0)
532	ALLOCATE (qsintmax_x(kjpindex,nvm),stat=ier)
533	a_er = a_er .OR. (ier.NE.0)
534	ALLOCATE (co2_flux(kjpindex,nvm),stat=ier)
535	a_er = a_er .OR. (ier.NE.0)
536	ALLOCATE (fco2_lu(kjpindex),stat=ier)
537	a_er = a_er .OR. (ier.NE.0)
538	IF (a_er) THEN
539	CALL ipslerr (3,'teststomate', &
540	& 'PROBLEM WITH ALLOCATION', &
541	& 'for local variables 1','')
542	ENDIF
543	!
544	! prepare forcing
545	!
546	max_totsize = 50
547	CALL getin_p ('STOMATE_FORCING_MEMSIZE',max_totsize)
548	max_totsize = max_totsize * 1000000
549
550	totsize_1step = SIZE(soiltype(:,3))*KIND(soiltype(:,3)) + &
551	SIZE(humrel_x)*KIND(humrel_x) + &
552	SIZE(litterhum)*KIND(litterhum) + &
553	SIZE(t2m)*KIND(t2m) + &
554	SIZE(t2m_min)*KIND(t2m_min) + &
555	SIZE(temp_sol)*KIND(temp_sol) + &
556	SIZE(soiltemp)*KIND(soiltemp) + &
557	SIZE(soilhum)*KIND(soilhum) + &
558	SIZE(precip_rain)*KIND(precip_rain) + &
559	SIZE(precip_snow)*KIND(precip_snow) + &
560	SIZE(gpp_x)*KIND(gpp_x) + &
561	SIZE(veget_force_x)*KIND(veget_force_x) + &
562	SIZE(veget_max_force_x)*KIND(veget_max_force_x) + &
563	SIZE(lai_force_x)*KIND(lai_force_x)
564	CALL reduce_sum(totsize_1step,totsize_tmp)
565	CALL bcast(totsize_tmp)
566	totsize_1step=totsize_tmp
567
568	! total number of forcing steps
569	IF ( nsft .NE. INT(one_year/(dt_force/one_day)) ) THEN
570	CALL ipslerr (3,'teststomate', &
571	& 'stomate: error with total number of forcing steps', &
572	& 'nsft','teststomate computation different with forcing file value.')
573	ENDIF
574	! number of forcing steps in memory
575	nsfm = MIN(nsft, &
576	& MAX(1,NINT( REAL(max_totsize,r_std) &
577	& /REAL(totsize_1step,r_std))))
578	!-
579	WRITE(numout,*) 'Offline forcing of Stomate:'
580	WRITE(numout,*) ' Total number of forcing steps:',nsft
581	WRITE(numout,*) ' Number of forcing steps in memory:',nsfm
582	!-
583	CALL init_forcing(kjpindex,nsfm,nsft)
584	!-
585	! ensure that we read all new forcing states
586	iisf = nsfm
587	! initialize that table that contains the indices
588	! of the forcing states that will be in memory
589	isf(:) = (/ (i,i=1,nsfm) /)
590
591	nf_written(:) = .FALSE.
592	nf_written(isf(:)) = .TRUE.
593
594	!-
595	! a time step for STOMATE corresponds to itau_step timesteps in SECHIBA
596	!-
597	itau_step = NINT(dt_force/dtradia)
598	IF (debug) WRITE(numout,*) "dtradia, dt_rest, dt_force, itau_step",dtradia, dt, dt_force, itau_step
599	!
600	CALL ioconf_startdate(date0)
601	CALL intsurf_time( itau_dep, date0, dtradia )
602	!-
603	! Length of integration
604	!-
605	WRITE(time_str,'(a)') '1Y'
606	CALL getin_p ('TIME_LENGTH', time_str)
607	! transform into itau
608	CALL tlen2itau(time_str, dt, date0, itau_len)
609	! itau_len*dtradia must be a multiple of dt_force
610	itau_len = NINT( MAX(1.,FLOAT(NINT(itau_len*dtradia/dt_force))) &
611	& *dt_force/dtradia)
612	!-
613	itau_fin = itau_dep+itau_len
614	!-
615	! set up STOMATE history file
616	!-
617	!Config Key = STOMATE_OUTPUT_FILE
618	!Config Desc = Name of file in which STOMATE's output is going
619	!Config to be written
620	!Config Def = stomate_history.nc
621	!Config Help = This file is going to be created by the model
622	!Config and will contain the output from the model.
623	!Config This file is a truly COADS compliant netCDF file.
624	!Config It will be generated by the hist software from
625	!Config the IOIPSL package.
626	!-
627	stom_histname='stomate_history.nc'
628	CALL getin_p ('STOMATE_OUTPUT_FILE', stom_histname)
629	WRITE(numout,*) 'STOMATE_OUTPUT_FILE', TRIM(stom_histname)
630	!-
631	!Config Key = STOMATE_HIST_DT
632	!Config Desc = STOMATE history time step (d)
633	!Config Def = 10.
634	!Config Help = Time step of the STOMATE history file
635	!-
636	hist_days_stom = 10.
637	CALL getin_p ('STOMATE_HIST_DT', hist_days_stom)
638	IF ( hist_days_stom == -1. ) THEN
639	hist_dt_stom = -1.
640	WRITE(numout,*) 'output frequency for STOMATE history file (d): one month.'
641	ELSE
642	hist_dt_stom = NINT( hist_days_stom ) * one_day
643	WRITE(numout,*) 'output frequency for STOMATE history file (d): ', &
644	hist_dt_stom/one_day
645	ENDIF
646	!-
647	!-
648	!- initialize
649	WRITE(numout,*) "before histbeg : ",date0,dt
650	IF ( rectilinear ) THEN
651	#ifdef CPP_PARA
652	CALL histbeg(stom_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
653	& itau_dep, date0, dt, hori_id, hist_id_stom, domain_id=orch_domain_id)
654	#else
655	CALL histbeg(stom_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
656	& itau_dep, date0, dt, hori_id, hist_id_stom)
657	#endif
658	ELSE
659	#ifdef CPP_PARA
660	CALL histbeg(stom_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
661	& itau_dep, date0, dt, hori_id, hist_id_stom, domain_id=orch_domain_id)
662	#else
663	CALL histbeg(stom_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
664	& itau_dep, date0, dt, hori_id, hist_id_stom)
665	#endif
666	ENDIF
667	!- define PFT axis
668	hist_PFTaxis = (/ ( REAL(i,r_std), i=1,nvm ) /)
669	!- declare this axis
670	CALL histvert (hist_id_stom, 'PFT', 'Plant functional type', &
671	& '1', nvm, hist_PFTaxis, hist_PFTaxis_id)
672	!- define Pool_10 axis
673	hist_pool_10axis = (/ ( REAL(i,r_std), i=1,10 ) /)
674	!- declare this axis
675	CALL histvert (hist_id_stom, 'P10', 'Pool 10 years', &
676	& '1', 10, hist_pool_10axis, hist_pool_10axis_id)
677
678	!- define Pool_100 axis
679	hist_pool_100axis = (/ ( REAL(i,r_std), i=1,100 ) /)
680	!- declare this axis
681	CALL histvert (hist_id_stom, 'P100', 'Pool 100 years', &
682	& '1', 100, hist_pool_100axis, hist_pool_100axis_id)
683
684	!- define Pool_11 axis
685	hist_pool_11axis = (/ ( REAL(i,r_std), i=1,11 ) /)
686	!- declare this axis
687	CALL histvert (hist_id_stom, 'P11', 'Pool 10 years + 1', &
688	& '1', 11, hist_pool_11axis, hist_pool_11axis_id)
689
690	!- define Pool_101 axis
691	hist_pool_101axis = (/ ( REAL(i,r_std), i=1,101 ) /)
692	!- declare this axis
693	CALL histvert (hist_id_stom, 'P101', 'Pool 100 years + 1', &
694	& '1', 101, hist_pool_101axis, hist_pool_101axis_id)
695
696	!- define STOMATE history file
697	CALL stom_define_history (hist_id_stom, nvm, iim, jjm, &
698	& dt, hist_dt_stom, hori_id, hist_PFTaxis_id, &
699	& hist_pool_10axis_id, hist_pool_100axis_id, &
700	& hist_pool_11axis_id, hist_pool_101axis_id)
701
702	!- end definition
703	CALL histend(hist_id_stom)
704	!-
705	!-
706	! STOMATE IPCC OUTPUTS IS ACTIVATED
707	!-
708	!Config Key = STOMATE_IPCC_OUTPUT_FILE
709	!Config Desc = Name of file in which STOMATE's output is going
710	!Config to be written
711	!Config Def = stomate_ipcc_history.nc
712	!Config Help = This file is going to be created by the model
713	!Config and will contain the output from the model.
714	!Config This file is a truly COADS compliant netCDF file.
715	!Config It will be generated by the hist software from
716	!Config the IOIPSL package.
717	!-
718	stom_ipcc_histname='stomate_ipcc_history.nc'
719	CALL getin_p('STOMATE_IPCC_OUTPUT_FILE', stom_ipcc_histname)
720	WRITE(numout,*) 'STOMATE_IPCC_OUTPUT_FILE', TRIM(stom_ipcc_histname)
721	!-
722	!Config Key = STOMATE_IPCC_HIST_DT
723	!Config Desc = STOMATE IPCC history time step (d)
724	!Config Def = 0.
725	!Config Help = Time step of the STOMATE IPCC history file
726	!-
727	hist_days_stom_ipcc = zero
728	CALL getin_p('STOMATE_IPCC_HIST_DT', hist_days_stom_ipcc)
729	IF ( hist_days_stom_ipcc == moins_un ) THEN
730	hist_dt_stom_ipcc = moins_un
731	WRITE(numout,*) 'output frequency for STOMATE IPCC history file (d): one month.'
732	ELSE
733	hist_dt_stom_ipcc = NINT( hist_days_stom_ipcc ) * one_day
734	WRITE(numout,*) 'output frequency for STOMATE IPCC history file (d): ', &
735	hist_dt_stom_ipcc/one_day
736	ENDIF
737
738	! test consistency between STOMATE_IPCC_HIST_DT and DT_SLOW parameters
739	dt_slow_ = one_day
740	CALL getin_p('DT_SLOW', dt_slow_)
741	IF ( hist_days_stom_ipcc > zero ) THEN
742	IF (dt_slow_ > hist_dt_stom_ipcc) THEN
743	WRITE(numout,*) "DT_SLOW = ",dt_slow_," , STOMATE_IPCC_HIST_DT = ",hist_dt_stom_ipcc
744	CALL ipslerr (3,'intsurf_history', &
745	& 'Problem with DT_SLOW > STOMATE_IPCC_HIST_DT','', &
746	& '(must be less or equal)')
747	ENDIF
748	ENDIF
749
750	IF ( hist_dt_stom_ipcc == 0 ) THEN
751	hist_id_stom_ipcc = -1
752	ELSE
753	!-
754	!- initialize
755	IF ( rectilinear ) THEN
756	#ifdef CPP_PARA
757	CALL histbeg(stom_ipcc_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
758	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc,domain_id=orch_domain_id)
759	#else
760	CALL histbeg(stom_ipcc_histname, iim, lon_rect, jjm, lat_rect, 1, iim, 1, jjm, &
761	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc)
762	#endif
763	ELSE
764	#ifdef CPP_PARA
765	CALL histbeg(stom_ipcc_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
766	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc,domain_id=orch_domain_id)
767	#else
768	CALL histbeg(stom_ipcc_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
769	& itau_dep, date0, dt, hori_id, hist_id_stom_ipcc)
770	#endif
771	ENDIF
772	!- declare this axis
773	CALL histvert (hist_id_stom_IPCC, 'PFT', 'Plant functional type', &
774	& '1', nvm, hist_PFTaxis, hist_IPCC_PFTaxis_id)
775
776	!- define STOMATE history file
777	CALL stom_IPCC_define_history (hist_id_stom_IPCC, nvm, iim, jjm, &
778	& dt, hist_dt_stom_ipcc, hori_id, hist_IPCC_PFTaxis_id)
779
780	!- end definition
781	CALL histend(hist_id_stom_IPCC)
782
783	ENDIF
784	!
785	CALL histwrite(hist_id_stom, 'Areas', itau_dep+itau_step, area, kjpindex, indices)
786	IF ( hist_id_stom_IPCC > 0 ) THEN
787	CALL histwrite(hist_id_stom_IPCC, 'Areas', itau_dep+itau_step, area, kjpindex, indices)
788	ENDIF
789	!
790	hist_id_sec = -1
791	hist_id_sec2 = -1
792	!-
793	! first call of slowproc to initialize variables
794	!-
795	itau = itau_dep
796	!
797	DO ji=1,kjpindex
798	DO jv=1,nvm
799	indexveg((jv-1)kjpindex + ji) = indices(ji) + (jv-1)kjpij
800	ENDDO
801	ENDDO
802	!-
803	!MM Problem here with dpu which depends on soil type
804	DO l = 1, nbdl-1
805	! first 2.0 is dpu
806	! second 2.0 is average
807	diaglev(l) = dpu_cste/(2*(nbdl-1) -1) ( ( 2(l-1) -1) + ( 2(l) -1) ) / 2.0
808	ENDDO
809	diaglev(nbdl) = dpu_cste
810	!
811
812	! >> DS : getin the value of slowproc parameters
813	!
814	!Config Key = CLAYFRACTION_DEFAULT
815	!Config Desc =
816	!Config If = OK_SECHIBA
817	!Config Def = 0.2
818	!Config Help =
819	!Config Units = NONE
820	CALL getin_p('CLAYFRACTION_DEFAULT',clayfraction_default)
821	!
822	!Config Key = MIN_VEGFRAC
823	!Config Desc = Minimal fraction of mesh a vegetation type can occupy
824	!Config If = OK_SECHIBA
825	!Config Def = 0.001
826	!Config Help =
827	!Config Units = NONE
828	CALL getin_p('MIN_VEGFRAC',min_vegfrac)
829	!
830	!Config Key = STEMPDIAG_BID
831	!Config Desc = only needed for an initial LAI if there is no restart file
832	!Config If = OK_SECHIBA
833	!Config Def = 280.
834	!Config Help =
835	!Config Units =
836	CALL getin_p('STEMPDIAG_BID',stempdiag_bid)
837	!
838	!Config Key = SOILTYPE_DEFAULT
839	!Config Desc = Default soil texture distribution in the following order : sand, loam and clay
840	!Config If = OK_SECHIBA
841	!Config Def = 0.0, 1.0, 0.0
842	!Config Help =
843	!Config Units = NONE
844	CALL getin_p('SOILTYPE_DEFAULT',soiltype_default)
845	! >> DS
846
847	CALL ioget_expval(val_exp)
848	ldrestart_read = .TRUE.
849	ldrestart_write = .FALSE.
850	!-
851	! read some variables we need from SECHIBA's restart file
852	!-
853	CALL slowproc_main &
854	& (itau, kjpij, kjpindex, dt_force, date0, &
855	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
856	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
857	& t2m, t2m_min, temp_sol, soiltemp, &
858	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
859	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
860	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
861	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
862	! correct date
863	day_counter = one_day - dt_force
864	date=1
865	!-
866	! time loop
867	!-
868	IF (debug) check_time=.TRUE.
869	CALL intsurf_time( itau_dep+itau_step, date0, dtradia )
870	l_first_intersurf=.FALSE.
871	!
872	DO itau = itau_dep+itau_step,itau_fin,itau_step
873	!
874	CALL intsurf_time( itau, date0, dtradia )
875	!
876	!-- next forcing state
877	iisf = iisf+1
878	IF (debug) WRITE(numout,*) "itau,iisf : ",itau,iisf
879	!---
880	IF (iisf .GT. nsfm) THEN
881	!-----
882	!---- we have to read new forcing states from the file
883	!-----
884	!---- determine blocks of forcing states that are contiguous in memory
885	!-----
886	CALL forcing_read(forcing_id,nsfm)
887
888	!--------------------------
889
890	!-----
891	!---- determine which forcing states must be read next time
892	!-----
893	isf(1) = isf(nsfm)+1
894	IF ( isf(1) .GT. nsft ) isf(1) = 1
895	DO iiisf = 2, nsfm
896	isf(iiisf) = isf(iiisf-1)+1
897	IF ( isf(iiisf) .GT. nsft ) isf(iiisf) = 1
898	ENDDO
899	nf_written(isf(:)) = .TRUE.
900	!---- start again at first forcing state
901	iisf = 1
902	ENDIF
903	! Bug here ! soiltype(:,3) != clay
904	! soiltype(:,3) = clay_fm(:,iisf)
905	humrel_x(:,:) = humrel_daily_fm(:,:,iisf)
906	litterhum(:) = litterhum_daily_fm(:,iisf)
907	t2m(:) = t2m_daily_fm(:,iisf)
908	t2m_min(:) = t2m_min_daily_fm(:,iisf)
909	temp_sol(:) = tsurf_daily_fm(:,iisf)
910	soiltemp(:,:) = tsoil_daily_fm(:,:,iisf)
911	soilhum(:,:) = soilhum_daily_fm(:,:,iisf)
912	precip_rain(:) = precip_fm(:,iisf)
913	gpp_x(:,:) = gpp_daily_fm(:,:,iisf)
914	veget_force_x(:,:) = veget_fm(:,:,iisf)
915	veget_max_force_x(:,:) = veget_max_fm(:,:,iisf)
916	lai_force_x(:,:) = lai_fm(:,:,iisf)
917	WHERE ( t2m(:) .LT. ZeroCelsius )
918	precip_snow(:) = precip_rain(:)
919	precip_rain(:) = 0.
920	ELSEWHERE
921	precip_snow(:) = 0.
922	ENDWHERE
923	!---
924	!-- scale GPP to new lai and veget_max
925	!---
926	WHERE ( lai_x(:,:) .EQ. 0.0 ) gpp_x(:,:) = 0.0
927	!-- scale GPP to new LAI
928	WHERE (lai_force_x(:,:) .GT. 0.0 )
929	gpp_x(:,:) = gpp_x(:,:)ATAN(2.lai_x(:,:)) &
930	& /ATAN( 2.*MAX(lai_force_x(:,:),0.01))
931	ENDWHERE
932	!-- scale GPP to new veget_max
933	WHERE (veget_max_force_x(:,:) .GT. 0.0 )
934	gpp_x(:,:) = gpp_x(:,:)*veget_max_x(:,:)/veget_max_force_x(:,:)
935	ENDWHERE
936	!---
937	!-- number crunching
938	!---
939	ldrestart_read = .FALSE.
940	ldrestart_write = .FALSE.
941	CALL slowproc_main &
942	& (itau, kjpij, kjpindex, dt_force, date0, &
943	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
944	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
945	& t2m, t2m_min, temp_sol, soiltemp, &
946	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
947	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
948	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
949	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
950	day_counter = one_day - dt_force
951	ENDDO ! end of the time loop
952	!-
953	! write restart files
954	!-
955	IF (is_root_prc) THEN
956	! first, read and write variables that are not managed otherwise
957	taboo_vars = &
958	& '$lat$ $lon$ $lev$ $veget_year$ '// &
959	& '$height$ $day_counter$ $veget$ $veget_max$ $frac_nobio$ '// &
960	& '$lai$ $soiltype_frac$ $clay_frac$ '// &
961	& '$nav_lon$ $nav_lat$ $nav_lev$ $time$ $time_steps$'
962	!-
963	CALL ioget_vname(rest_id_sec, nbvar, varnames)
964	!-
965	DO iv = 1, nbvar
966	!-- check if the variable is to be written here
967	IF (INDEX( taboo_vars,'$'//TRIM(varnames(iv))//'$') .EQ. 0 ) THEN
968	IF (debug) WRITE(numout,*) "restart var : ",TRIM(varnames(iv)),itau_dep,itau_fin
969
970	!---- get variable dimensions, especially 3rd dimension
971	CALL ioget_vdim &
972	& (rest_id_sec, varnames(iv), varnbdim_max, varnbdim, vardims)
973	l1d = ALL(vardims(1:varnbdim) .EQ. 1)
974	!---- read it
975	IF (l1d) THEN
976	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
977	1,1,1,itau_dep,.TRUE.,var_1d)
978	ELSE
979	ALLOCATE(var_3d(nbp_glo,vardims(3)),stat=ier)
980	IF (ier .NE. 0) STOP 'ALLOCATION PROBLEM'
981	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
982	nbp_glo,vardims(3),1,itau_dep,.TRUE.,var_3d, &
983	"gather",nbp_glo,indices_g)
984	ENDIF
985	!---- write it
986	IF (l1d) THEN
987	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
988	1,1,1,itau_fin,var_1d)
989	ELSE
990	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
991	nbp_glo,vardims(3),1,itau_fin,var_3d, &
992	'scatter',nbp_glo,indices_g)
993	DEALLOCATE (var_3d)
994	ENDIF
995	ENDIF
996	ENDDO
997	ENDIF
998	CALL barrier_para()
999
1000	! call slowproc to write restart files
1001	ldrestart_read = .FALSE.
1002	ldrestart_write = .TRUE.
1003	!-
1004	IF (debug) WRITE(numout,*) "Call slowproc for restart."
1005	CALL slowproc_main &
1006	& (itau_fin, kjpij, kjpindex, dt_force, date0, &
1007	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
1008	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
1009	& t2m, t2m_min, temp_sol, soiltemp, &
1010	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
1011	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
1012	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
1013	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_stom, hist_id_stom_IPCC, co2_flux, fco2_lu)
1014	!-
1015	! close files
1016	!-
1017	IF (is_root_prc) THEN
1018	CALL restclo
1019	IF ( debug ) WRITE(numout,*) 'REST CLOSED'
1020	ENDIF
1021	CALL histclo
1022
1023	IF (is_root_prc) &
1024	ier = NF90_CLOSE (forcing_id)
1025
1026	IF (is_root_prc) THEN
1027	CALL getin_dump()
1028	ENDIF
1029	#ifdef CPP_PARA
1030	CALL MPI_FINALIZE(ier)
1031	#endif
1032	WRITE(numout,*) "End of teststomate."
1033
1034	!---------------
1035	END PROGRAM teststomate
1036	!
1037	!===
1038	!

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

Download in other formats: