Context Navigation

source: tags/ORCHIDEE_2_1/ORCHIDEE/src_driver/forcesoil.f90 @ 6593

Last change on this file since 6593 was 4470, checked in by josefine.ghattas, 7 years ago
Change adresse for orchidee-help in the comments.
File size: 26.9 KB

Line
1	! =================================================================================================================================
2	! PROGRAM : forcesoil
3	!
4	! CONTACT : orchidee-help _at_ listes.ipsl.fr
5	!
6	! LICENCE : IPSL (2006). This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7	!
8	!>\BRIEF This subroutine runs the soilcarbon submodel using specific initial conditions
9	!! and driving variables in order to obtain soil carbon stocks closed to the steady-state values
10	!! quicker than when using the ''full'' ORCHIDEE.
11	!!
12	!!\n DESCRIPTION: None
13	!! This subroutine computes the soil carbon stocks by calling the soilcarbon routine at each time step. \n
14	!! The aim is to obtain soil carbon stocks closed to the steady-state values and ultimately to create
15	!! an updated stomate restart file for the stomate component. The state variables of the subsystem are the clay content
16	!! (fixed value) and the soil carbon stocks. Initial conditions for the state variables are read in an
17	!! input stomate restart file. Driving variables are Soil carbon input, Water and Temperature stresses on
18	!! Organic Matter decomposition. Driving variables are read from a specific forcing file produced by a former run of ORCHIDEE
19	!! (SECHIBA+STOMATE). \n
20	!! The FORCESOIL program first consists in reading a set of input files, allocating variables and
21	!! preparing output stomate restart file. \n
22	!! Then, a loop over time is performed in which the soilcarbon routine is called at each time step. \n
23	!! Last, final values of the soil carbon stocks are written into the output stomate restart file. \n
24	!! No flag is associated with the use of the FORCESOIL program. \n
25	!!
26	!! RECENT CHANGE(S): None
27	!!
28	!! REFERENCE(S) : None
29	!!
30	!! FLOWCHART : None
31	!!
32	!! SVN :
33	!! $HeadURL: $
34	!! $Date: $
35	!! $Revision: $
36	!! \n
37	!_ =================================================================================================================================
38
39	PROGRAM forcesoil
40
41	USE netcdf
42	!-
43	USE defprec
44	USE constantes
45	USE constantes_mtc
46	USE pft_parameters
47	USE stomate_data
48	USE ioipsl_para
49	USE mod_orchidee_para
50	USE stomate_soilcarbon
51
52	!-
53	IMPLICIT NONE
54	!-
55	CHARACTER(LEN=80) :: sto_restname_in,sto_restname_out
56	INTEGER(i_std) :: iim,jjm !! Indices (unitless)
57
58	INTEGER(i_std),PARAMETER :: llm = 1 !! Vertical Layers (requested by restini routine) (unitless)
59	INTEGER(i_std) :: kjpindex !! Domain size (unitless)
60
61	INTEGER(i_std) :: itau_dep,itau_len !! Time step read in the restart file (?)
62	!! and number of time steps of the simulation (unitless)
63	CHARACTER(LEN=30) :: time_str !! Length of the simulation (year)
64	REAL(r_std) :: dt_files !! time step between two successive itaus (?)
65	!! (requested by restini routine) (seconds)
66	REAL(r_std) :: date0 !! Time at which itau = 0 (requested by restini routine) (?)
67	INTEGER(i_std) :: rest_id_sto !! ID of the input restart file (unitless)
68	CHARACTER(LEN=20), SAVE :: thecalendar = 'noleap' !! Type of calendar defined in the input restart file
69	!! (unitless)
70	!-
71	CHARACTER(LEN=100) :: Cforcing_name !! Name of the forcing file (unitless)
72	INTEGER :: Cforcing_id !! ID of the forcing file (unitless)
73	INTEGER :: v_id !! ID of the variable 'Index' stored in the forcing file
74	!! (unitless)
75	REAL(r_std) :: dt_forcesoil !! Time step at which soilcarbon routine is called (days)
76	INTEGER :: nparan !! Number of values stored per year in the forcing file
77	!! (unitless)
78	INTEGER :: nbyear
79	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indices !! Grid Point Index used per processor (unitless)
80	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indices_g !! Grid Point Index for all processor (unitless)
81	REAL(r_std),DIMENSION(:),ALLOCATABLE :: x_indices_g !! Grid Point Index for all processor (unitless)
82	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lon, lat !! Longitude and Latitude of each grid point defined
83	!! in lat/lon (2D) (degrees)
84	REAL(r_std),DIMENSION(llm) :: lev !! Number of level (requested by restini routine) (unitless)
85
86
87	INTEGER :: i,m,iatt,iv,iyear !! counters (unitless)
88
89	CHARACTER(LEN=80) :: var_name
90	CHARACTER(LEN=800) :: taboo_vars !! string used for storing the name of the variables
91	!! of the stomate restart file that are not automatically
92	!! duplicated from input to output restart file (unitless)
93	REAL(r_std),DIMENSION(1) :: xtmp !! scalar read/written in restget/restput routines (unitless)
94	INTEGER(i_std),PARAMETER :: nbvarmax=300 !! maximum # of variables assumed in the stomate restart file
95	!! (unitless)
96	INTEGER(i_std) :: nbvar !! # of variables effectively present
97	!! in the stomate restart file (unitless)
98	CHARACTER(LEN=50),DIMENSION(nbvarmax) :: varnames !! list of the names of the variables stored
99	!! in the stomate restart file (unitless)
100	INTEGER(i_std) :: varnbdim !! # of dimensions of a given variable
101	!! of the stomate restart file
102	INTEGER(i_std),PARAMETER :: varnbdim_max=20 !! maximal # of dimensions assumed for any variable
103	!! of the stomate restart file
104	INTEGER,DIMENSION(varnbdim_max) :: vardims !! length of each dimension of a given variable
105	!! of the stomate restart file
106	LOGICAL :: l1d !! boolean : TRUE if all dimensions of a given variable
107	!! of the stomate restart file are of length 1 (ie scalar)
108	!! (unitless)
109	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: var_3d !! matrix read/written in restget/restput routines (unitless)
110	REAL(r_std) :: x_tmp !! temporary variable used to store return value
111	!! from nf90_get_att (unitless)
112	CHARACTER(LEN=10) :: part_str !! string suffix indicating the index of a PFT
113	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: clay_g !! clay fraction (nbpglo) (unitless)
114	REAL(r_std),DIMENSION(:,:,:,:),ALLOCATABLE :: soilcarbon_input_g !! soil carbon input (nbpglob,ncarb,nvm,time)
115	!! (\f$gC m^{-2} dt_forcesoil^{-1}\f$)
116	REAL(r_std),DIMENSION(:,:,:),ALLOCATABLE :: control_temp_g !! Temperature control (nbp_glo,above/below,time) on OM decomposition
117	!! (unitless)
118	REAL(r_std),DIMENSION(:,:,:),ALLOCATABLE :: control_moist_g !! Moisture control (nbp_glo,abo/below,time) on OM decomposition
119	!! ?? Should be defined per PFT as well (unitless)
120	REAL(r_std),DIMENSION(:,:,:),ALLOCATABLE :: carbon_g !! Soil carbon stocks (nbp_glo,ncarb,nvm) (\f$gC m^{-2}\f$)
121
122	REAL(r_std),ALLOCATABLE :: clay(:) !! clay fraction (nbp_loc) (unitless)
123	REAL(r_std),ALLOCATABLE :: soilcarbon_input(:,:,:,:) !! soil carbon input (nbp_loc,ncarb,nvm,time)
124	!! (\f$gC m^{-2} dt_forcesoil^{-1}\f$)
125	REAL(r_std),ALLOCATABLE :: control_temp(:,:,:) !! Temperature control (nbp_loc,above/below,time) on OM decomposition
126	!! (unitless)
127	REAL(r_std),ALLOCATABLE :: control_moist(:,:,:) !! Moisture control (nbp_loc,abo/below,time) on OM decomposition
128	!! ?? Should be defined per PFT as well (unitless)
129	REAL(r_std),ALLOCATABLE :: carbon(:,:,:) !! Soil carbon stocks (nbp_loc,ncarb,nvm) (\f$gC m^{-2}\f$)
130	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: resp_hetero_soil !! Heterotrophic respiration (\f$gC m^{-2} dt_forcesoil^{-1}\f$)
131	!! (requested by soilcarbon routine but not used here)
132
133	INTEGER(i_std) :: printlev_loc !! Local write level
134	INTEGER(i_std) :: ier,iret !! Used for hangling errors
135
136	CHARACTER(LEN=30) :: temp_name
137	LOGICAL :: l_error !! boolean for memory allocation
138	REAL(r_std), ALLOCATABLE, DIMENSION(:,:,:,:) :: matrixA !! Carbon fluxes matrix
139
140	!_ =================================================================================================================================
141
142	CALL Init_orchidee_para
143	CALL init_timer
144
145	! Set specific write level to forcesoil using PRINTLEV_forcesoil=[0-4] in run.def.
146	! The global printlev is used as default value.
147	printlev_loc=get_printlev('forcesoil')
148
149	!-
150	! Configure the number of PFTS
151	!-
152
153	! 1. Read the number of PFTs
154	!
155	!Config Key = NVM
156	!Config Desc = number of PFTs
157	!Config If = OK_SECHIBA or OK_STOMATE
158	!Config Def = 13
159	!Config Help = The number of vegetation types define by the user
160	!Config Units = [-]
161	CALL getin_p('NVM',nvm)
162
163	! 2. Allocation
164	l_error = .FALSE.
165	ALLOCATE(pft_to_mtc(nvm),stat=ier)
166	l_error = l_error .OR. (ier .NE. 0)
167	IF (l_error) THEN
168	STOP 'pft_to_mtc (forcesoil only) : error in memory allocation'
169	ENDIF
170
171	! 3. Initialisation of the correspondance table
172	pft_to_mtc(:) = undef_int
173
174	! 4.Reading of the conrrespondance table in the .def file
175	!
176	!Config Key = PFT_TO_MTC
177	!Config Desc = correspondance array linking a PFT to MTC
178	!Config if = OK_SECHIBA or OK_STOMATE
179	!Config Def = 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13
180	!Config Help =
181	!Config Units = [-]
182	CALL getin_p('PFT_TO_MTC',pft_to_mtc)
183
184	! 4.1 if nothing is found, we use the standard configuration
185	IF(nvm == nvmc ) THEN
186	IF(pft_to_mtc(1) == undef_int) THEN
187	WRITE(numout,*) 'Note to the user : we will use ORCHIDEE to its standard configuration'
188	pft_to_mtc(:) = (/ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13 /)
189	ENDIF
190	ELSE
191	IF(pft_to_mtc(1) == undef_int) THEN
192	WRITE(numout,*)' The array PFT_TO_MTC is empty : we stop'
193	ENDIF
194	ENDIF
195
196	! 4.2 What happened if pft_to_mtc(j) > nvmc (if the mtc doesn't exist)?
197	DO i = 1, nvm
198	IF(pft_to_mtc(i) > nvmc) THEN
199	WRITE(numout,*) "the MTC you chose doesn't exist"
200	STOP 'we stop reading pft_to_mtc'
201	ENDIF
202	ENDDO
203
204	! 4.3 Check if pft_to_mtc(1) = 1
205	IF(pft_to_mtc(1) /= 1) THEN
206	WRITE(numout,*) 'the first pft has to be the bare soil'
207	STOP 'we stop reading next values of pft_to_mtc'
208	ELSE
209	DO i = 2,nvm
210	IF(pft_to_mtc(i) == 1) THEN
211	WRITE(numout,*) 'only pft_to_mtc(1) has to be the bare soil'
212	STOP 'we stop reading pft_to_mtc'
213	ENDIF
214	ENDDO
215	ENDIF
216
217	! 5. Allocate and initialize natural ans is_c4
218
219	! 5.1 Memory allocation
220	l_error = .FALSE.
221	ALLOCATE(natural(nvm),stat=ier)
222	l_error = l_error .OR. (ier .NE. 0)
223	ALLOCATE(is_c4(nvm),stat=ier)
224
225	IF (l_error) THEN
226	STOP 'natural or is_c4 (forcesoil only) : error in memory allocation'
227	ENDIF
228
229	! 5.2 Initialisation
230	DO i = 1, nvm
231	natural(i) = natural_mtc(pft_to_mtc(i))
232	is_c4(i) = is_c4_mtc(pft_to_mtc(i))
233	ENDDO
234
235	!!-
236	!! 1. Initialisation stage
237	!! Reading a set of input files, allocating variables and preparing output restart file.
238	!!-
239	! Define restart file name
240	! for reading initial conditions (sto_restname_in var) and for writting final conditions (sto_restname_out var).
241	! User values are used if present in the .def file.
242	! If not present, default values (stomate_start.nc and stomate_rest_out.c) are used.
243	!-
244	IF (is_root_prc) THEN
245	sto_restname_in = 'stomate_start.nc'
246	CALL getin ('STOMATE_RESTART_FILEIN',sto_restname_in)
247	WRITE(numout,*) 'STOMATE INPUT RESTART_FILE: ',TRIM(sto_restname_in)
248	sto_restname_out = 'stomate_rest_out.nc'
249	CALL getin ('STOMATE_RESTART_FILEOUT',sto_restname_out)
250	WRITE(numout,*) 'STOMATE OUTPUT RESTART_FILE: ',TRIM(sto_restname_out)
251	!-
252	! Open the input file and Get some Dimension and Attributes ID's
253	!-
254	iret = NF90_OPEN (sto_restname_in, NF90_NOWRITE, rest_id_sto)
255	iret = NF90_INQUIRE_DIMENSION (rest_id_sto,1,len=iim_g)
256	iret = NF90_INQUIRE_DIMENSION (rest_id_sto,2,len=jjm_g)
257	iret = NF90_INQ_VARID (rest_id_sto, "time", iv)
258	iret = NF90_GET_ATT (rest_id_sto, iv, 'calendar',thecalendar)
259	iret = NF90_CLOSE (rest_id_sto)
260	i=INDEX(thecalendar,ACHAR(0))
261	IF ( i > 0 ) THEN
262	thecalendar(i:20)=' '
263	ENDIF
264	!-
265	! Allocate longitudes and latitudes
266	!-
267	ALLOCATE (lon(iim_g,jjm_g))
268	ALLOCATE (lat(iim_g,jjm_g))
269	lon(:,:) = zero
270	lat(:,:) = zero
271	lev(1) = zero
272	!-
273	CALL restini &
274	& (sto_restname_in, iim_g, jjm_g, lon, lat, llm, lev, &
275	& sto_restname_out, itau_dep, date0, dt_files, rest_id_sto)
276	ENDIF
277
278	CALL bcast(date0)
279	CALL bcast(thecalendar)
280	WRITE(numout,*) "calendar = ",thecalendar
281	!-
282	! calendar
283	!-
284	CALL ioconf_calendar (thecalendar)
285	CALL ioget_calendar (one_year,one_day)
286	CALL ioconf_startdate(date0)
287	!
288	!! For master process only
289	!
290	IF (is_root_prc) THEN
291	!-
292	! define forcing file's name (Cforcing_name var)
293	! User value is used if present in the .def file
294	! If not, default (NONE) is used
295	!-
296	Cforcing_name = 'NONE'
297	CALL getin ('STOMATE_CFORCING_NAME',Cforcing_name)
298	!-
299	! Open FORCESOIL's forcing file to read some basic info (dimensions, variable ID's)
300	! and allocate variables.
301	!-
302	iret = NF90_OPEN (TRIM(Cforcing_name),NF90_NOWRITE,Cforcing_id)
303	IF (iret /= NF90_NOERR) THEN
304	CALL ipslerr (3,'forcesoil', &
305	& 'Could not open file : ', &
306	& Cforcing_name,'(Do you have forget it ?)')
307	ENDIF
308	!-
309	! Total Domain size is stored in nbp_glo variable
310	!-
311	ier = NF90_GET_ATT (Cforcing_id,NF90_GLOBAL,'kjpindex',x_tmp)
312	nbp_glo = NINT(x_tmp)
313	!-
314	! Number of values stored per year in the forcing file is stored in nparan var.
315	!-
316	ier = NF90_GET_ATT (Cforcing_id,NF90_GLOBAL,'nparan',x_tmp)
317	nparan = NINT(x_tmp)
318	ier = NF90_GET_ATT (Cforcing_id,NF90_GLOBAL,'nbyear',x_tmp)
319	nbyear = NINT(x_tmp)
320	!-
321	ALLOCATE (indices_g(nbp_glo))
322	ALLOCATE (clay_g(nbp_glo))
323	!-
324	ALLOCATE (x_indices_g(nbp_glo),stat=ier)
325	ier = NF90_INQ_VARID (Cforcing_id,'index',v_id)
326	ier = NF90_GET_VAR (Cforcing_id,v_id,x_indices_g)
327	indices_g(:) = NINT(x_indices_g(:))
328	WRITE(numout,*) mpi_rank,"indices globaux : ",indices_g
329	DEALLOCATE (x_indices_g)
330	!-
331	ier = NF90_INQ_VARID (Cforcing_id,'clay',v_id)
332	ier = NF90_GET_VAR (Cforcing_id,v_id,clay_g)
333	!-
334	! time step of forcesoil program (in days)
335	!-
336	dt_forcesoil = one_year / FLOAT(nparan)
337	! Initialize dt_sechiba in module constantes_var. This is needed in carbonsoil.
338	dt_sechiba=dt_forcesoil*one_day
339	WRITE(numout,*) 'time step (d): ',dt_forcesoil
340	WRITE(numout,*) 'nparan: ',nparan
341	WRITE(numout,*) 'nbyear: ',nbyear
342	!-
343	! read and write the variables in the output restart file we do not modify within the Forcesoil program
344	! ie all variables stored in the input restart file except those stored in taboo_vars
345	!-
346	taboo_vars ='$lon$ $lat$ $lev$ $nav_lon$ $nav_lat$ $nav_lev$ $time$ $time_steps$ '// &
347	& '$day_counter$ $dt_days$ $date$ '
348	!-
349	DO m = 1,nvm
350	WRITE(part_str,'(I2)') m
351	IF (m < 10) part_str(1:1) = '0'
352	temp_name = '$carbon_'//part_str(1:LEN_TRIM(part_str))//'$'
353	taboo_vars = TRIM(taboo_vars)//' '//TRIM(temp_name)
354	ENDDO
355	!-
356	CALL ioget_vname(rest_id_sto, nbvar, varnames)
357	!-
358	! read and write some special variables (1D or variables that we need)
359	!-
360	var_name = 'day_counter'
361	CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
362	CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
363	!-
364	var_name = 'dt_days'
365	CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
366	CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
367	!-
368	var_name = 'date'
369	CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
370	CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
371	!-
372	DO iv=1,nbvar
373	!-- check if the variable is to be written here
374	IF (INDEX(taboo_vars,'$'//TRIM(varnames(iv))//'$') == 0 ) THEN
375	!---- get variable dimensions, especially 3rd dimension
376	CALL ioget_vdim &
377	& (rest_id_sto, varnames(iv), varnbdim_max, varnbdim, vardims)
378	l1d = ALL(vardims(1:varnbdim) == 1)
379	!---- read it
380	IF (l1d) THEN
381	CALL restget &
382	& (rest_id_sto, TRIM(varnames(iv)), 1, vardims(3), &
383	& 1, itau_dep, .TRUE., xtmp)
384	ELSE
385	ALLOCATE( var_3d(nbp_glo,vardims(3)), stat=ier)
386	IF (ier /= 0) STOP 'ALLOCATION PROBLEM'
387	!----
388	CALL restget &
389	& (rest_id_sto, TRIM(varnames(iv)), nbp_glo, vardims(3), &
390	& 1, itau_dep, .TRUE., var_3d, "gather", nbp_glo, indices_g)
391	ENDIF
392	!---- write it
393	IF (l1d) THEN
394	CALL restput &
395	& (rest_id_sto, TRIM(varnames(iv)), 1, vardims(3), &
396	& 1, itau_dep, xtmp)
397	ELSE
398	CALL restput &
399	& (rest_id_sto, TRIM(varnames(iv)), nbp_glo, vardims(3), &
400	& 1, itau_dep, var_3d, 'scatter', nbp_glo, indices_g)
401	!----
402	DEALLOCATE(var_3d)
403	ENDIF
404	ENDIF
405	ENDDO
406	!-
407	! read soil carbon stocks values stored in the input restart file
408	!-
409	ALLOCATE(carbon_g(nbp_glo,ncarb,nvm))
410	carbon_g(:,:,:) = val_exp
411	DO m = 1, nvm
412	WRITE (part_str, '(I2)') m
413	IF (m<10) part_str(1:1)='0'
414	var_name = 'carbon_'//part_str(1:LEN_TRIM(part_str))
415	CALL restget &
416	& (rest_id_sto, var_name, nbp_glo, ncarb , 1, itau_dep, &
417	& .TRUE., carbon_g(:,:,m), 'gather', nbp_glo, indices_g)
418	IF (ALL(carbon_g(:,:,m) == val_exp)) carbon_g(:,:,m) = zero
419	!-- do not write this variable: it will be modified.
420	ENDDO
421	WRITE(numout,*) "date0 : ",date0, itau_dep
422	!-
423	! Analytical spinup is set to false
424	!
425	spinup_analytic = .FALSE.
426
427	! Length of the run (in Years)
428	! User value is used if present in the .def file
429	! If not, default value (10000 Years) is used
430	!-
431	WRITE(time_str,'(a)') '10000Y'
432	CALL getin('TIME_LENGTH', time_str)
433	write(numout,*) 'Number of years for carbon spinup : ',time_str
434	! transform into itau
435	CALL tlen2itau(time_str, dt_forcesoil*one_day, date0, itau_len)
436	write(numout,*) 'Number of time steps to do: ',itau_len
437	!-
438	! read soil carbon inputs, water and temperature stresses on OM decomposition
439	! into the forcing file - We read an average year.
440	!-
441	ALLOCATE(soilcarbon_input_g(nbp_glo,ncarb,nvm,nparan*nbyear))
442	ALLOCATE(control_temp_g(nbp_glo,nlevs,nparan*nbyear))
443	ALLOCATE(control_moist_g(nbp_glo,nlevs,nparan*nbyear))
444	!-
445	ier = NF90_INQ_VARID (Cforcing_id,'soilcarbon_input',v_id)
446	ier = NF90_GET_VAR (Cforcing_id,v_id,soilcarbon_input_g)
447	ier = NF90_INQ_VARID (Cforcing_id, 'control_moist',v_id)
448	ier = NF90_GET_VAR (Cforcing_id,v_id,control_moist_g)
449	ier = NF90_INQ_VARID (Cforcing_id, 'control_temp',v_id)
450	ier = NF90_GET_VAR (Cforcing_id,v_id,control_temp_g)
451	!-
452	ier = NF90_CLOSE (Cforcing_id)
453	!-
454	ENDIF
455	CALL bcast(nparan)
456	CALL bcast(nbyear)
457	CALL bcast(dt_forcesoil)
458	CALL bcast(iim_g)
459	CALL bcast(jjm_g)
460	CALL bcast(nbp_glo)
461	CALL bcast(itau_dep)
462	CALL bcast(itau_len)
463	IF (.NOT. ALLOCATED(indices_g)) ALLOCATE (indices_g(nbp_glo))
464	CALL bcast(indices_g)
465
466	!
467	! We must initialize data_para :
468	CALL init_orchidee_data_para_driver(nbp_glo,indices_g)
469
470	kjpindex=nbp_loc
471	jjm=jj_nb
472	iim=iim_g
473	IF (printlev_loc>=3) WRITE(numout,*) "Local grid : ",kjpindex,iim,jjm
474
475	!---
476	!--- Create the index table
477	!---
478	!--- This job returns a LOCAL kindex.
479	!---
480	ALLOCATE (indices(kjpindex),stat=ier)
481	!
482	!! scattering to all processes in parallel mode
483	!
484	CALL scatter(indices_g,indices)
485	indices(1:kjpindex)=indices(1:kjpindex)-(jj_begin-1)*iim_g
486	IF (printlev_loc>=3) WRITE(numout,*) mpi_rank,"indices locaux = ",indices(1:kjpindex)
487	!-
488	! Allocation of the variables for a processor
489	!-
490	ALLOCATE(clay(kjpindex))
491	ALLOCATE(soilcarbon_input(kjpindex,ncarb,nvm,nparan*nbyear))
492	ALLOCATE(control_temp(kjpindex,nlevs,nparan*nbyear))
493	ALLOCATE(control_moist(kjpindex,nlevs,nparan*nbyear))
494	ALLOCATE(carbon(kjpindex,ncarb,nvm))
495	ALLOCATE(resp_hetero_soil(kjpindex,nvm))
496	ALLOCATE(matrixA(kjpindex,nvm,nbpools,nbpools))
497	DO i = 1,nbpools
498	matrixA(:,:,i,i) = un
499	ENDDO
500	iatt = 0
501
502	!-
503	! Initialization of the variables for a processor
504	!-
505	CALL Scatter(clay_g,clay)
506	CALL Scatter(soilcarbon_input_g,soilcarbon_input)
507	CALL Scatter(control_temp_g,control_temp)
508	CALL Scatter(control_moist_g,control_moist)
509	CALL Scatter(carbon_g,carbon)
510
511	!-
512	! Configuration of the parameters
513	!-
514
515	!Config Key = FRAC_CARB_AP
516	!Config Desc = frac carb coefficients from active pool: depends on clay content
517	!Config if = OK_STOMATE
518	!Config Def = 0.004
519	!Config Help = fraction of the active pool going to the passive pool
520	!Config Units = [-]
521	CALL getin_p('FRAC_CARB_AP',frac_carb_ap)
522	!
523	!Config Key = FRAC_CARB_SA
524	!Config Desc = frac_carb_coefficients from slow pool
525	!Config if = OK_STOMATE
526	!Config Def = 0.42
527	!Config Help = fraction of the slow pool going to the active pool
528	!Config Units = [-]
529	CALL getin_p('FRAC_CARB_SA',frac_carb_sa)
530	!
531	!Config Key = FRAC_CARB_SP
532	!Config Desc = frac_carb_coefficients from slow pool
533	!Config if = OK_STOMATE
534	!Config Def = 0.03
535	!Config Help = fraction of the slow pool going to the passive pool
536	!Config Units = [-]
537	CALL getin_p('FRAC_CARB_SP',frac_carb_sp)
538	!
539	!Config Key = FRAC_CARB_PA
540	!Config Desc = frac_carb_coefficients from passive pool
541	!Config if = OK_STOMATE
542	!Config Def = 0.45
543	!Config Help = fraction of the passive pool going to the passive pool
544	!Config Units = [-]
545	CALL getin_p('FRAC_CARB_PA',frac_carb_pa)
546	!
547	!Config Key = FRAC_CARB_PS
548	!Config Desc = frac_carb_coefficients from passive pool
549	!Config if = OK_STOMATE
550	!Config Def = 0.0
551	!Config Help = fraction of the passive pool going to the passive pool
552	!Config Units = [-]
553	CALL getin_p('FRAC_CARB_PS',frac_carb_ps)
554	!
555	!Config Key = ACTIVE_TO_PASS_CLAY_FRAC
556	!Config Desc =
557	!Config if = OK_STOMATE
558	!Config Def = .68
559	!Config Help =
560	!Config Units = [-]
561	CALL getin_p('ACTIVE_TO_PASS_CLAY_FRAC',active_to_pass_clay_frac)
562	!
563	!Config Key = CARBON_TAU_IACTIVE
564	!Config Desc = residence times in carbon pools
565	!Config if = OK_STOMATE
566	!Config Def = 0.149
567	!Config Help =
568	!Config Units = [days]
569	CALL getin_p('CARBON_TAU_IACTIVE',carbon_tau_iactive)
570	!
571	!Config Key = CARBON_TAU_ISLOW
572	!Config Desc = residence times in carbon pools
573	!Config if = OK_STOMATE
574	!Config Def = 5.48
575	!Config Help =
576	!Config Units = [days]
577	CALL getin_p('CARBON_TAU_ISLOW',carbon_tau_islow)
578	!
579	!Config Key = CARBON_TAU_IPASSIVE
580	!Config Desc = residence times in carbon pools
581	!Config if = OK_STOMATE
582	!Config Def = 241.
583	!Config Help =
584	!Config Units = [days]
585	CALL getin_p('CARBON_TAU_IPASSIVE',carbon_tau_ipassive)
586	!
587	!Config Key = FLUX_TOT_COEFF
588	!Config Desc =
589	!Config if = OK_STOMATE
590	!Config Def = 1.2, 1.4,.75
591	!Config Help =
592	!Config Units = [days]
593	CALL getin_p('FLUX_TOT_COEFF',flux_tot_coeff)
594
595	!!-
596	!! 2. Computational step
597	!! Loop over time - Call of soilcarbon routine at each time step
598	!! Updated soil carbon stocks are stored into carbon variable
599	!! We only keep the last value of carbon variable (no time dimension).
600	!!-
601	iyear=1
602	DO i=1,itau_len
603	iatt = iatt+1
604	IF (iatt > nparan*nbyear) THEN
605	IF (printlev_loc>=3) WRITE(numout,*) iyear
606	iatt = 1
607	iyear=iyear+1
608	ENDIF
609	CALL soilcarbon &
610	& (kjpindex, clay, &
611	& soilcarbon_input(:,:,:,iatt), &
612	& control_temp(:,:,iatt), control_moist(:,:,iatt), &
613	& carbon, resp_hetero_soil, &
614	& matrixA)
615	ENDDO
616	WRITE(numout,*) "End of soilcarbon LOOP."
617
618	!
619	!! Gathering of variables towards main processor in parallel mode
620	!
621	CALL Gather(carbon,carbon_g)
622	!!-
623	!! 3. write new carbon stocks into the ouput restart file
624	!!-
625	IF (is_root_prc) THEN
626	DO m=1,nvm
627	WRITE (part_str, '(I2)') m
628	IF (m<10) part_str(1:1)='0'
629	var_name = 'carbon_'//part_str(1:LEN_TRIM(part_str))
630	CALL restput &
631	& (rest_id_sto, var_name, nbp_glo, ncarb , 1, itau_dep, &
632	& carbon_g(:,:,m), 'scatter', nbp_glo, indices_g)
633	ENDDO
634	!-
635	CALL getin_dump
636	CALL restclo
637	ENDIF
638	#ifdef CPP_PARA
639	CALL MPI_FINALIZE(ier)
640	#endif
641	WRITE(numout,*) "End of forcesoil."
642	!--------------------
643	END PROGRAM forcesoil

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

Download in other formats: