Context Navigation

stomate_soilcarbon.f90 @ 7475

Last change on this file since 7475 was 2945, checked in by sebastiaan.luyssaert, 9 years ago
DEV: tested 1 year global. This code contains the latest version for anthropogenic tree species channges, several bug fixes to forest management as well as the code for the fully integrated multi-layer energy budget. This implies that the multi-layer energy budget makes use Pinty's albedo scheme, the rognostic canopy structure as well as a vertical profile for stomatal conductance. This is an intermediate version because species change code is not complete as some management changes have not been implemented yet. Further the multi-layer albedo code needs more work in terms of calculating average fluxes at the pixel rather than the PFT level
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 24.3 KB

Line
1	! =================================================================================================================================
2	! MODULE : stomate_soilcarbon
3	!
4	! CONTACT : orchidee-help _at_ ipsl.jussieu.fr
5	!
6	! LICENCE : IPSL (2006)
7	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8	!
9	!>\BRIEF Calculate soil dynamics largely following the Century model
10	!!
11	!!\n DESCRIPTION: None
12	!!
13	!! RECENT CHANGE(S): None
14	!!
15	!! SVN :
16	!! $HeadURL$
17	!! $Date$
18	!! $Revision$
19	!! \n
20	!_ ================================================================================================================================
21
22	MODULE stomate_soilcarbon
23
24	! modules used:
25
26	USE ioipsl_para
27	USE stomate_data
28	USE constantes
29
30	IMPLICIT NONE
31
32	! private & public routines
33
34	PRIVATE
35	PUBLIC soilcarbon,soilcarbon_clear
36
37	! Variables shared by all subroutines in this module
38
39	LOGICAL, SAVE :: firstcall = .TRUE. !! Is this the first call? (true/false)
40	!$OMP THREADPRIVATE(firstcall)
41
42	CONTAINS
43
44
45	!! ================================================================================================================================
46	!! SUBROUTINE : soilcarbon_clear
47	!!
48	!>\BRIEF Set the flag ::firstcall to .TRUE. and as such activate sections 1.1.2 and 1.2 of the subroutine soilcarbon
49	!! (see below).
50	!!
51	!_ ================================================================================================================================
52
53	SUBROUTINE soilcarbon_clear
54	firstcall=.TRUE.
55	END SUBROUTINE soilcarbon_clear
56
57
58	!! ================================================================================================================================
59	!! SUBROUTINE : soilcarbon
60	!!
61	!>\BRIEF Computes the soil respiration and carbon stocks, essentially
62	!! following Parton et al. (1987).
63	!!
64	!! DESCRIPTION : The soil is divided into 3 carbon pools, with different
65	!! characteristic turnover times : active (1-5 years), slow (20-40 years)
66	!! and passive (200-1500 years).\n
67	!! There are three types of carbon transferred in the soil:\n
68	!! - carbon input in active and slow pools from litter decomposition,\n
69	!! - carbon fluxes between the three pools,\n
70	!! - carbon losses from the pools to the atmosphere, i.e., soil respiration.\n
71	!!
72	!! The subroutine performs the following tasks:\n
73	!!
74	!! Section 1.\n
75	!! The flux fractions (f) between carbon pools are defined based on Parton et
76	!! al. (1987). The fractions are constants, except for the flux fraction from
77	!! the active pool to the slow pool, which depends on the clay content,\n
78	!! \latexonly
79	!! \input{soilcarbon_eq1.tex}
80	!! \endlatexonly\n
81	!! In addition, to each pool is assigned a constant turnover time.\n
82	!!
83	!! Section 2.\n
84	!! The carbon input, calculated in the stomate_litter module, is added to the
85	!! carbon stock of the different pools.\n
86	!!
87	!! Section 3.\n
88	!! First, the outgoing carbon flux of each pool is calculated. It is
89	!! proportional to the product of the carbon stock and the ratio between the
90	!! iteration time step and the residence time:\n
91	!! \latexonly
92	!! \input{soilcarbon_eq2.tex}
93	!! \endlatexonly
94	!! ,\n
95	!! Note that in the case of crops, the additional multiplicative factor
96	!! integrates the faster decomposition due to tillage (following Gervois et
97	!! al. (2008)).
98	!! In addition, the flux from the active pool depends on the clay content:\n
99	!! \latexonly
100	!! \input{soilcarbon_eq3.tex}
101	!! \endlatexonly
102	!! ,\n
103	!! Each pool is then cut from the carbon amount corresponding to each outgoing
104	!! flux:\n
105	!! \latexonly
106	!! \input{soilcarbon_eq4.tex}
107	!! \endlatexonly\n
108	!! Second, the flux fractions lost to the atmosphere is calculated in each pool
109	!! by subtracting from 1 the pool-to-pool flux fractions. The soil respiration
110	!! is then the summed contribution of all the pools,\n
111	!! \latexonly
112	!! \input{soilcarbon_eq5.tex}
113	!! \endlatexonly\n
114	!! Finally, each carbon pool accumulates the contribution of the other pools:
115	!! \latexonly
116	!! \input{soilcarbon_eq6.tex}
117	!! \endlatexonly
118	!!
119	!! Section 4.\n
120	!! If the flag SPINUP_ANALYTIC is set to true, the matrix A is updated following
121	!! Lardy (2011).
122	!!
123	!! RECENT CHANGE(S): None
124	!!
125	!! MAIN OUTPUTS VARIABLE(S): carbon, resp_hetero_soil
126	!!
127	!! REFERENCE(S) :
128	!! - Parton, W.J., D.S. Schimel, C.V. Cole, and D.S. Ojima. 1987. Analysis of
129	!! factors controlling soil organic matter levels in Great Plains grasslands.
130	!! Soil Sci. Soc. Am. J., 51, 1173-1179.
131	!! - Gervois, S., P. Ciais, N. de Noblet-Ducoudre, N. Brisson, N. Vuichard,
132	!! and N. Viovy (2008), Carbon and water balance of European croplands
133	!! throughout the 20th century, Global Biogeochem. Cycles, 22, GB2022,
134	!! doi:10.1029/2007GB003018.
135	!! - Lardy, R, et al., A new method to determine soil organic carbon equilibrium,
136	!! Environmental Modelling & Software (2011), doi:10.1016\|j.envsoft.2011.05.016
137	!!
138	!! FLOWCHART :
139	!! \latexonly
140	!! \includegraphics[scale=0.5]{soilcarbon_flowchart.jpg}
141	!! \endlatexonly
142	!! \n
143	!_ ================================================================================================================================
144
145	SUBROUTINE soilcarbon (npts, dt, clay, soilcarbon_input, &
146	control_temp, control_moist, carbon, resp_hetero_soil, &
147	MatrixA,veget_max)
148
149	!! 0. Variable and parameter declaration
150
151	!! 0.1 Input variables
152
153	INTEGER(i_std), INTENT(in) :: npts !! Domain size (unitless)
154	REAL(r_std), INTENT(in) :: dt !! Time step \f$(dtradia one_day^{-1})$\f
155	REAL(r_std), DIMENSION(npts), INTENT(in) :: clay !! Clay fraction (unitless, 0-1)
156	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max !! Coverage fraction of a PFT within the pixel
157	!! (unitless, 0-1)
158	REAL(r_std), DIMENSION(npts,nlevs), INTENT(in) :: control_temp !! Temperature control of heterotrophic respiration
159	!! (unitless: 0->1)
160	REAL(r_std), DIMENSION(npts,nlevs), INTENT(in) :: control_moist !! Moisture control of heterotrophic respiration
161	!! (unitless: 0.25->1)
162	REAL(r_std), DIMENSION(npts,ncarb,nvm), INTENT(in) :: soilcarbon_input !! Amount of carbon going into the carbon pools from
163	!! litter decomposition \f$(gC m^{-2} day^{-1})$\f
164
165	!! 0.2 Output variables
166
167	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: resp_hetero_soil !! Soil heterotrophic respiration
168	!! \f$(gC m^{-2} (dtradia one_day^{-1})^{-1})$\f
169
170	!! 0.3 Modified variables
171
172	REAL(r_std), DIMENSION(npts,ncarb,nvm), INTENT(inout) :: carbon !! Soil carbon pools: active, slow, or passive,
173	!! \f$(gC m^{2})$\f
174	REAL(r_std), DIMENSION(npts,nvm,nbpools,nbpools), &
175	INTENT(inout) :: MatrixA !! Matrix containing the fluxes between the carbon
176	!! pools per sechiba time step
177	!! @tex $(gC.m^2.day^{-1})$ @endtex
178
179	!! 0.4 Local variables
180	REAL(r_std), SAVE, DIMENSION(ncarb) :: carbon_tau !! Residence time in carbon pools (days)
181	!$OMP THREADPRIVATE(carbon_tau)
182	REAL(r_std), DIMENSION(npts,ncarb,ncarb) :: frac_carb !! Flux fractions between carbon pools
183	!! (second index=origin, third index=destination)
184	!! (unitless, 0-1)
185	REAL(r_std), DIMENSION(npts,ncarb) :: frac_resp !! Flux fractions from carbon pools to the atmosphere
186	!! (respiration) (unitless, 0-1)
187	REAL(r_std), DIMENSION(npts,ncarb,nelements) :: fluxtot !! Total flux out of carbon pools \f$(gC m^{2})$\f
188	REAL(r_std), DIMENSION(npts,ncarb,ncarb,nelements) :: flux !! Fluxes between carbon pools \f$(gC m^{2})$\f
189	CHARACTER(LEN=7), DIMENSION(ncarb) :: carbon_str !! Name of the carbon pools for informative outputs
190	!! (unitless)
191	INTEGER(i_std) :: iicarb,m,j !! Indices (unitless)
192	INTEGER(i_std) :: icarb, imbc !! Indices (unitless)
193	INTEGER(i_std) :: ipts, ivm !! Indices (unitless)
194	REAL(r_std), DIMENSION(npts,nvm,nmbcomp,nelements) :: check_intern !! Contains the components of the internal
195	!! mass balance chech for this routine
196	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
197	REAL(r_std), DIMENSION(npts,nvm,nelements) :: closure_intern !! Check closure of internal mass balance
198	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
199	REAL(r_std), DIMENSION(npts,nvm,nelements) :: pool_start !! Start and end pool of this routine
200	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
201	REAL(r_std), DIMENSION(npts,nvm,nelements) :: pool_end !! Start and end pool of this routine
202	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
203
204	!_ ================================================================================================================================
205
206	!! bavard is the level of diagnostic information, 0 (none) to 4 (full)
207	IF (bavard.GE.3) WRITE(numout,*) 'Entering soilcarbon'
208
209	!! 1. Initializations
210
211	!! 1.1 Initialize check for mass balance closure
212	! The mass balance is calculated at the end of this routine
213	! in section 14. Initial biomass and harvest pool and all other
214	! relevant pools were set to zero before calling this routine.
215	pool_start(:,:,:) = zero
216	DO icarb = 1,ncarb
217	pool_start(:,:,icarbon) = pool_start(:,:,icarbon) + &
218	( carbon(:,icarb,:) + (soilcarbon_input(:,icarb,:) * dt) ) * &
219	veget_max(:,:)
220	ENDDO
221
222	!! 1.2 Get soil "constants"
223	!! 1.2.1 Flux fractions between carbon pools: depend on clay content,
224	! recalculated each time.
225	! From active pool: depends on clay content
226	frac_carb(:,iactive,iactive) = zero
227	frac_carb(:,iactive,ipassive) = frac_carb_ap
228	frac_carb(:,iactive,islow) = un - frac_carb(:,iactive,ipassive) - &
229	(metabolic_ref_frac - active_to_pass_clay_frac*clay(:))
230
231	! 1.2.2 from slow pool
232	frac_carb(:,islow,islow) = zero
233	frac_carb(:,islow,iactive) = frac_carb_sa
234	frac_carb(:,islow,ipassive) = frac_carb_sp
235
236	! From passive pool
237	frac_carb(:,ipassive,ipassive) = zero
238	frac_carb(:,ipassive,iactive) = frac_carb_pa
239	frac_carb(:,ipassive,islow) = frac_carb_ps
240
241	IF ( firstcall ) THEN
242
243	!! 1.2.3 Residence times in carbon pools (days)
244	! 1.5 years. This is same as CENTURY. But, in Parton et
245	! al. (1987), it's weighted by moisture and temperature
246	! dependences.
247	carbon_tau(iactive) = carbon_tau_iactive * one_year
248
249	! 25 years. This is same as CENTURY. But, in Parton et
250	! al. (1987), it's weighted by moisture and temperature
251	! dependences.
252	carbon_tau(islow) = carbon_tau_islow * one_year
253
254	! 1000 years. This is same as CENTURY. But, in Parton et
255	! al. (1987), it's weighted by moisture and temperature
256	! dependences.
257	carbon_tau(ipassive) = carbon_tau_ipassive * one_year
258
259	!! 1.2.4 Messages : display the residence times
260	carbon_str(iactive) = 'active'
261	carbon_str(islow) = 'slow'
262	carbon_str(ipassive) = 'passive'
263
264	IF(ld_litter)THEN
265	WRITE(numout,*) 'soilcarbon:'
266
267	WRITE(numout,*) ' > minimal carbon residence time in carbon pools (d):'
268	DO icarb = 1, ncarb ! Loop over carbon pools
269	WRITE(numout,*) '(1, ::carbon_str(icarb)):',carbon_str(icarb), &
270	' : (1, ::carbon_tau(icarb)):',carbon_tau(icarb)
271	ENDDO
272
273	WRITE(numout,*) ' > flux fractions between carbon pools: depend on clay content'
274	ENDIF
275
276	firstcall = .FALSE.
277
278	ENDIF
279
280	!! 1.3 Set soil respiration to zero
281	resp_hetero_soil(:,:) = zero
282
283
284	!! 2. Update the carbon stocks with the different soil carbon input
285
286	carbon(:,:,:) = carbon(:,:,:) + soilcarbon_input(:,:,:) * dt
287
288
289	!! 3. Fluxes between carbon reservoirs, and to the atmosphere (respiration) \n
290
291	!! 3.1. Determine the respiration fraction : what's left after
292	! subtracting all the 'pool-to-pool' flux fractions
293	! Diagonal elements of frac_carb are zero
294	! VPP killer:
295	! frac_resp(:,:) = 1. - SUM( frac_carb(:,:,:), DIM=3 )
296	frac_resp(:,:) = un - frac_carb(:,:,iactive) - frac_carb(:,:,islow) - &
297	frac_carb(:,:,ipassive)
298
299	!! 3.2. Calculate fluxes
300	DO m = 1,nvm ! Loop over #Â PFTs
301
302	!! 3.2.1. Flux out of pools
303	DO icarb = 1, ncarb ! Loop over carbon pools from which the flux comes
304
305	! Determine total flux out of pool
306	! For crops multiply tillage factor of decomposition
307	! Not crop
308	IF ( natural(m) ) THEN
309
310	fluxtot(:,icarb,icarbon) = dt/carbon_tau(icarb) * carbon(:,icarb,m) * &
311	control_moist(:,ibelow) * control_temp(:,ibelow)
312
313	! C3 crop
314	ELSEIF ( (.NOT. natural(m)) .AND. (.NOT. is_c4(m)) ) THEN
315
316	fluxtot(:,icarb,icarbon) = dt/carbon_tau(icarb) * carbon(:,icarb,m) * &
317	control_moist(:,ibelow) * control_temp(:,ibelow) * &
318	flux_tot_coeff(1)
319
320	! C4 Crop
321	ELSEIF ( (.NOT. natural(m)) .AND. is_c4(m) ) THEN
322
323	fluxtot(:,icarb,icarbon) = dt/carbon_tau(icarb) * carbon(:,icarb,m) * &
324	control_moist(:,ibelow) * control_temp(:,ibelow) * &
325	flux_tot_coeff(2)
326
327	ENDIF ! natural(m)
328
329	! Carbon flux from active pools depends on clay content
330	IF ( icarb .EQ. iactive ) THEN
331
332	fluxtot(:,icarb,icarbon) = fluxtot(:,icarb,icarbon) * &
333	( un - flux_tot_coeff(3) * clay(:) )
334
335	ENDIF
336
337	! Update the loss in each carbon pool
338	carbon(:,icarb,m) = carbon(:,icarb,m) - fluxtot(:,icarb,icarbon)
339
340	! Fluxes towards the other pools (icarb -> iicarb)
341	! Loop over the carbon pools where the flux goes
342	DO iicarb = 1, ncarb
343
344	flux(:,icarb,iicarb,icarbon) = frac_carb(:,icarb,iicarb) * &
345	fluxtot(:,icarb,icarbon)
346
347	ENDDO
348
349	ENDDO ! End of loop over carbon pools
350
351	!! 3.2.2 respiration
352	! BE CAREFUL: Here resp_hetero_soil is divided by dt to have a value which corresponds to
353	! the sechiba time step but then in stomate.f90 resp_hetero_soil is multiplied by dt.
354	! Perhaps it could be simplified. Moreover, we should entirely adapt the routines to the
355	! dtradia/one_day time step and avoid some constructions that could create bugs during
356	! future developments.
357	resp_hetero_soil(:,m) = &
358	( frac_resp(:,iactive) * fluxtot(:,iactive,icarbon) + &
359	frac_resp(:,islow) * fluxtot(:,islow,icarbon) + &
360	frac_resp(:,ipassive) * fluxtot(:,ipassive,icarbon) ) / dt
361
362	!! 3.2.3 add fluxes to active, slow, and passive pools
363	DO icarb = 1, ncarb ! Loop over carbon pools
364
365	carbon(:,icarb,m) = carbon(:,icarb,m) + &
366	flux(:,iactive,icarb,icarbon) + flux(:,ipassive,icarb,icarbon) + &
367	flux(:,islow,icarb,icarbon)
368
369	ENDDO ! Loop over carbon pools
370
371	ENDDO ! End loop over PFTs
372
373
374	!! 4. Check mass balance closure
375
376	!! 4.1 Calculate components of the mass balance
377	pool_end(:,:,:) = zero
378	DO icarb = 1,ncarb
379	pool_end(:,:,icarbon) = pool_end(:,:,icarbon) + &
380	( carbon(:,icarb,:) * veget_max(:,:) )
381	ENDDO
382
383	!! 4.2 Calculate mass balance
384	! Note that soilcarbon is transfered to other pools but that the pool
385	! was not updated. We should not account for it in ::pool_end
386	check_intern(:,:,:,:) = zero
387	check_intern(:,:,iatm2land,icarbon) = zero
388	check_intern(:,:,iland2atm,icarbon) = -un * resp_hetero_soil(:,:) * &
389	veget_max(:,:) * dt
390	check_intern(:,:,ilat2out,icarbon) = zero
391	check_intern(:,:,ilat2in,icarbon) = -un * zero
392	check_intern(:,:,ipoolchange,icarbon) = -un * (pool_end(:,:,icarbon) - &
393	pool_start(:,:,icarbon))
394	closure_intern = zero
395	DO imbc = 1,nmbcomp
396	closure_intern(:,:,icarbon) = closure_intern(:,:,icarbon) + &
397	check_intern(:,:,imbc,icarbon)
398	ENDDO
399
400	!! 4.3 Write verdict
401	DO ipts=1,npts
402	DO ivm=1,nvm
403	IF(ABS(closure_intern(ipts,ivm,icarbon)) .LE. min_stomate)THEN
404	IF (ld_massbal) WRITE(numout,*) 'Mass balance closure in stomate_soilcarbon.f90'
405	ELSE
406	WRITE(numout,*) 'Error: mass balance is not closed in stomate_soilcarbon.f90'
407	WRITE(numout,*) ' ipts,ivm; ', ipts,ivm
408	WRITE(numout,*) ' Difference is, ', closure_intern(ipts,ivm,icarbon)
409	WRITE(numout,*) ' pool_end,pool_start: ', pool_end(ipts,ivm,icarbon), pool_start(ipts,ivm,icarbon)
410	WRITE(numout,*) ' resp_hetero_soil,veget_max: ', resp_hetero_soil(ipts,ivm),veget_max(ipts,ivm)
411	IF(ld_stop)THEN
412	CALL ipslerr_p (3,'soilcarbon', 'Mass balance error.','','')
413	ENDIF
414	ENDIF
415	ENDDO
416	ENDDO
417
418
419	!! 5. (Quasi-)Analytical Spin-up
420
421	!! 5.2 Finish to fill MatrixA with fluxes between soil pools
422
423	IF (spinup_analytic) THEN
424
425	DO m = 2,nvm
426
427	! flux leaving the active pool
428	MatrixA(:,m,iactive_pool,iactive_pool) = moins_un * &
429	dt/carbon_tau(iactive) * &
430	control_moist(:,ibelow) * control_temp(:,ibelow) * &
431	( 1. - flux_tot_coeff(3) * clay(:))
432
433	! flux received by the active pool from the slow pool
434	MatrixA(:,m,iactive_pool,islow_pool) = &
435	frac_carb(:,islow,iactive)dt/carbon_tau(islow) &
436	control_moist(:,ibelow) * control_temp(:,ibelow)
437
438	! flux received by the active pool from the passive pool
439	MatrixA(:,m,iactive_pool,ipassive_pool) = &
440	frac_carb(:,ipassive,iactive)dt/carbon_tau(ipassive) &
441	control_moist(:,ibelow) * control_temp(:,ibelow)
442
443	! flux received by the slow pool from the active pool
444	MatrixA(:,m,islow_pool,iactive_pool) = frac_carb(:,iactive,islow) *&
445	dt/carbon_tau(iactive) * &
446	control_moist(:,ibelow) * control_temp(:,ibelow) * &
447	( 1. - flux_tot_coeff(3) * clay(:) )
448
449	! flux leaving the slow pool
450	MatrixA(:,m,islow_pool,islow_pool) = moins_un * &
451	dt/carbon_tau(islow) * &
452	control_moist(:,ibelow) * control_temp(:,ibelow)
453
454	! flux received by the passive pool from the active pool
455	MatrixA(:,m,ipassive_pool,iactive_pool) = &
456	frac_carb(:,iactive,ipassive)* &
457	dt/carbon_tau(iactive) * &
458	control_moist(:,ibelow) * control_temp(:,ibelow) *&
459	( 1. - flux_tot_coeff(3) * clay(:) )
460
461	! flux received by the passive pool from the slow pool
462	MatrixA(:,m,ipassive_pool,islow_pool) = &
463	frac_carb(:,islow,ipassive) * &
464	dt/carbon_tau(islow) * &
465	control_moist(:,ibelow) * control_temp(:,ibelow)
466
467	! flux leaving the passive pool
468	MatrixA(:,m,ipassive_pool,ipassive_pool) = moins_un * &
469	dt/carbon_tau(ipassive) * &
470	control_moist(:,ibelow) * control_temp(:,ibelow)
471
472
473	IF ( (.NOT. natural(m)) .AND. (.NOT. is_c4(m)) ) THEN ! C3crop
474
475	! flux leaving the active pool
476	MatrixA(:,m,iactive_pool,iactive_pool) = &
477	MatrixA(:,m,iactive_pool,iactive_pool) * &
478	flux_tot_coeff(1)
479
480	! flux received by the active pool from the slow pool
481	MatrixA(:,m,iactive_pool,islow_pool)= &
482	MatrixA(:,m,iactive_pool,islow_pool) * &
483	flux_tot_coeff(1)
484
485	! flux received by the active pool from the passive pool
486	MatrixA(:,m,iactive_pool,ipassive_pool) = &
487	MatrixA(:,m,iactive_pool,ipassive_pool) * &
488	flux_tot_coeff(1)
489
490	! flux received by the slow pool from the active pool
491	MatrixA(:,m,islow_pool,iactive_pool) = &
492	MatrixA(:,m,islow_pool,iactive_pool) * &
493	flux_tot_coeff(1)
494
495	! flux leaving the slow pool
496	MatrixA(:,m,islow_pool,islow_pool) = &
497	MatrixA(:,m,islow_pool,islow_pool) * &
498	flux_tot_coeff(1)
499
500	! flux received by the passive pool from the active pool
501	MatrixA(:,m,ipassive_pool,iactive_pool) = &
502	MatrixA(:,m,ipassive_pool,iactive_pool) * &
503	flux_tot_coeff(1)
504
505	! flux received by the passive pool from the slow pool
506	MatrixA(:,m,ipassive_pool,islow_pool) = &
507	MatrixA(:,m,ipassive_pool,islow_pool) * &
508	flux_tot_coeff(1)
509
510	! flux leaving the passive pool
511	MatrixA(:,m,ipassive_pool,ipassive_pool) = &
512	MatrixA(:,m,ipassive_pool,ipassive_pool) *&
513	flux_tot_coeff(1)
514
515	ENDIF ! (.NOT. natural(m)) .AND. (.NOT. is_c4(m))
516
517
518	IF ( (.NOT. natural(m)) .AND. is_c4(m) ) THEN ! C4crop
519
520	! flux leaving the active pool
521	MatrixA(:,m,iactive_pool,iactive_pool) = &
522	MatrixA(:,m,iactive_pool,iactive_pool) * &
523	flux_tot_coeff(2)
524
525	! flux received by the active pool from the slow pool
526	MatrixA(:,m,iactive_pool,islow_pool)= &
527	MatrixA(:,m,iactive_pool,islow_pool) * &
528	flux_tot_coeff(2)
529
530	! flux received by the active pool from the passive pool
531	MatrixA(:,m,iactive_pool,ipassive_pool) = &
532	MatrixA(:,m,iactive_pool,ipassive_pool) * &
533	flux_tot_coeff(2)
534
535	! flux received by the slow pool from the active pool
536	MatrixA(:,m,islow_pool,iactive_pool) = &
537	MatrixA(:,m,islow_pool,iactive_pool) * &
538	flux_tot_coeff(2)
539
540	! flux leaving the slow pool
541	MatrixA(:,m,islow_pool,islow_pool) = &
542	MatrixA(:,m,islow_pool,islow_pool) * &
543	flux_tot_coeff(2)
544
545	! flux received by the passive pool from the active pool
546	MatrixA(:,m,ipassive_pool,iactive_pool) = &
547	MatrixA(:,m,ipassive_pool,iactive_pool) * &
548	flux_tot_coeff(2)
549
550	! flux received by the passive pool from the slow pool
551	MatrixA(:,m,ipassive_pool,islow_pool) = &
552	MatrixA(:,m,ipassive_pool,islow_pool) * &
553	flux_tot_coeff(2)
554
555	! flux leaving the passive pool
556	MatrixA(:,m,ipassive_pool,ipassive_pool) = &
557	MatrixA(:,m,ipassive_pool,ipassive_pool) * &
558	flux_tot_coeff(2)
559
560	ENDIF ! (.NOT. natural(m)) .AND. is_c4(m)
561
562	IF (bavard.GE.4) WRITE(numout,*)'Finish to fill MatrixA'
563
564	ENDDO ! Loop over # PFTS
565
566
567	! 5.2 Add Identity for each submatrix(7,7)
568
569	DO j = 1,nbpools
570	MatrixA(:,:,j,j) = MatrixA(:,:,j,j) + un
571	ENDDO
572
573	ENDIF ! (spinup_analytic)
574
575	IF (bavard.GE.4) WRITE(numout,*) 'Leaving soilcarbon'
576
577	END SUBROUTINE soilcarbon
578
579	END MODULE stomate_soilcarbon

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

Context Navigation

source: branches/publications/ORCHIDEE_CAN_r3069/src_stomate/stomate_soilcarbon.f90 @ 7475

Download in other formats: