Ignore:
Location:
tags/ORCHIDEE_1_9_5_2
Files:
11 added
54 edited

Legend:

Unmodified
Added
Removed

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_parallel/transfert_para.f90

    r119 r405  
    2828  END INTERFACE 
    2929 
    30   INTERFACE gather_s 
    31     MODULE PROCEDURE gather_is, & 
    32                      gather_rs, & 
    33                      gather_ls 
    34   END INTERFACE 
     30!!$  INTERFACE gather_s 
     31!!$    MODULE PROCEDURE gather_is, & 
     32!!$                     gather_rs, & 
     33!!$                  gather_ls 
     34!!$  END INTERFACE 
    3535   
    3636  INTERFACE gather 
    … …  
    196196  IMPLICIT NONE 
    197197    LOGICAL,INTENT(INOUT) :: Var 
    198     
    199 #ifndef CPP_PARA 
    200     RETURN 
    201 #else 
    202     CALL bcast_lgen(Var,1) 
     198    LOGICAL,DIMENSION(1) :: Var1 
     199#ifndef CPP_PARA 
     200    RETURN 
     201#else 
     202    IF (is_root_prc) & 
     203         Var1(1)=Var 
     204    CALL bcast_lgen(Var1,1) 
     205    Var=Var1(1) 
    203206#endif 
    204207  END SUBROUTINE bcast_l 
    … …  
    548551!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
    549552 
    550   SUBROUTINE gather_is(VarIn, VarOut) 
    551     USE data_para 
    552     USE timer 
    553  
    554     IMPLICIT NONE 
    555    
    556 #ifdef CPP_PARA 
    557     INCLUDE 'mpif.h' 
    558 #endif 
    559      
    560     INTEGER,INTENT(IN) :: VarIn 
    561     INTEGER,INTENT(OUT),DIMENSION(:) :: VarOut 
    562    
    563 #ifdef CPP_PARA 
    564     INTEGER :: nb,i,index_para,rank 
    565     INTEGER :: ierr 
    566     LOGICAL :: flag=.FALSE. 
    567     LOGICAL, PARAMETER :: check=.FALSE. 
    568 #endif 
    569  
    570 #ifndef CPP_PARA 
    571     VarOut(:)=VarIn 
    572     RETURN 
    573 #else 
    574  
    575     IF (timer_state(timer_mpi)==running) THEN 
    576       flag=.TRUE. 
    577     ELSE 
    578       flag=.FALSE. 
    579     ENDIF 
    580      
    581     IF (flag) CALL suspend_timer(timer_mpi) 
    582  
    583     IF (check) & 
    584          WRITE(numout,*) "gather_rgen VarIn=",VarIn     
    585  
    586 #ifdef CPP_PARA 
    587     CALL MPI_GATHER(VarIn,1,MPI_INT_ORCH,VarOut,1,MPI_INT_ORCH,root_prc,MPI_COMM_ORCH,ierr) 
    588 #endif 
    589  
    590     IF (check) & 
    591          WRITE(numout,*) "gather_rgen VarOut=",VarOut 
    592     IF (flag) CALL resume_timer(timer_mpi) 
    593 #endif 
    594   END SUBROUTINE gather_is 
    595  
    596   SUBROUTINE gather_rs(VarIn, VarOut) 
    597     USE data_para 
    598     USE timer 
    599  
    600     IMPLICIT NONE 
    601    
    602 #ifdef CPP_PARA 
    603     INCLUDE 'mpif.h' 
    604 #endif 
    605  
    606     REAL,INTENT(IN) :: VarIn 
    607     REAL,INTENT(OUT),DIMENSION(:) :: VarOut 
    608    
    609 #ifdef CPP_PARA 
    610     INTEGER :: nb,i,index_para,rank 
    611     INTEGER :: ierr 
    612     LOGICAL :: flag=.FALSE. 
    613     LOGICAL, PARAMETER :: check=.FALSE. 
    614 #endif 
    615  
    616 #ifndef CPP_PARA 
    617     VarOut(:)=VarIn 
    618     RETURN 
    619 #else 
    620  
    621     IF (timer_state(timer_mpi)==running) THEN 
    622       flag=.TRUE. 
    623     ELSE 
    624       flag=.FALSE. 
    625     ENDIF 
    626      
    627     IF (flag) CALL suspend_timer(timer_mpi) 
    628  
    629     IF (check) & 
    630          WRITE(numout,*) "gather_rgen VarIn=",VarIn     
    631  
    632 #ifdef CPP_PARA 
    633     CALL MPI_GATHER(VarIn,1,MPI_REAL_ORCH,VarOut,1,MPI_REAL_ORCH,root_prc,MPI_COMM_ORCH,ierr) 
    634 #endif 
    635  
    636     IF (check) & 
    637          WRITE(numout,*) "gather_rgen VarOut=",VarOut 
    638  
    639     IF (flag) CALL resume_timer(timer_mpi) 
    640 #endif 
    641   END SUBROUTINE gather_rs 
    642  
    643   SUBROUTINE gather_ls(VarIn, VarOut) 
    644     USE data_para 
    645     USE timer 
    646  
    647     IMPLICIT NONE 
    648    
    649 #ifdef CPP_PARA 
    650     INCLUDE 'mpif.h' 
    651 #endif 
    652      
    653     LOGICAL,INTENT(IN) :: VarIn 
    654     LOGICAL,INTENT(OUT),DIMENSION(:) :: VarOut 
    655    
    656 #ifdef CPP_PARA 
    657     INTEGER :: nb,i,index_para,rank 
    658     INTEGER :: ierr 
    659     LOGICAL :: flag=.FALSE. 
    660     LOGICAL, PARAMETER :: check=.FALSE. 
    661 #endif 
    662  
    663 #ifndef CPP_PARA 
    664     VarOut(:)=VarIn 
    665     RETURN 
    666 #else 
    667  
    668     IF (timer_state(timer_mpi)==running) THEN 
    669       flag=.TRUE. 
    670     ELSE 
    671       flag=.FALSE. 
    672     ENDIF 
    673      
    674     IF (flag) CALL suspend_timer(timer_mpi) 
    675  
    676     IF (check) & 
    677          WRITE(numout,*) "gather_rgen VarIn=",VarIn     
    678  
    679 #ifdef CPP_PARA 
    680     CALL MPI_GATHER(VarIn,1,MPI_LOGICAL,VarOut,1,MPI_LOGICAL,root_prc,MPI_COMM_ORCH,ierr) 
    681 #endif 
    682  
    683     IF (check) & 
    684          WRITE(numout,*) "gather_rgen VarOut=",VarOut 
    685     IF (flag) CALL resume_timer(timer_mpi) 
    686 #endif 
    687   END SUBROUTINE gather_ls 
     553!!$  SUBROUTINE gather_is(VarIn, VarOut) 
     554!!$    USE data_para 
     555!!$    USE timer 
     556!!$ 
     557!!$    IMPLICIT NONE 
     558!!$   
     559!!$#ifdef CPP_PARA 
     560!!$    INCLUDE 'mpif.h' 
     561!!$#endif 
     562!!$     
     563!!$    INTEGER,INTENT(IN) :: VarIn 
     564!!$    INTEGER,INTENT(OUT),DIMENSION(:) :: VarOut 
     565!!$   
     566!!$#ifdef CPP_PARA 
     567!!$    INTEGER :: nb,i,index_para,rank 
     568!!$    INTEGER :: ierr 
     569!!$    LOGICAL :: flag=.FALSE. 
     570!!$    LOGICAL, PARAMETER :: check=.FALSE. 
     571!!$#endif 
     572!!$ 
     573!!$#ifndef CPP_PARA 
     574!!$    VarOut(:)=VarIn 
     575!!$    RETURN 
     576!!$#else 
     577!!$ 
     578!!$    IF (timer_state(timer_mpi)==running) THEN 
     579!!$      flag=.TRUE. 
     580!!$    ELSE 
     581!!$      flag=.FALSE. 
     582!!$    ENDIF 
     583!!$     
     584!!$    IF (flag) CALL suspend_timer(timer_mpi) 
     585!!$ 
     586!!$    IF (check) & 
     587!!$         WRITE(numout,*) "gather_rgen VarIn=",VarIn     
     588!!$ 
     589!!$#ifdef CPP_PARA 
     590!!$    CALL MPI_GATHER(VarIn,1,MPI_INT_ORCH,VarOut,1,MPI_INT_ORCH,root_prc,MPI_COMM_ORCH,ierr) 
     591!!$#endif 
     592!!$ 
     593!!$    IF (check) & 
     594!!$         WRITE(numout,*) "gather_rgen VarOut=",VarOut 
     595!!$    IF (flag) CALL resume_timer(timer_mpi) 
     596!!$#endif 
     597!!$  END SUBROUTINE gather_is 
     598!!$ 
     599!!$  SUBROUTINE gather_rs(VarIn, VarOut) 
     600!!$    USE data_para 
     601!!$    USE timer 
     602!!$ 
     603!!$    IMPLICIT NONE 
     604!!$   
     605!!$#ifdef CPP_PARA 
     606!!$    INCLUDE 'mpif.h' 
     607!!$#endif 
     608!!$ 
     609!!$    REAL,INTENT(IN) :: VarIn 
     610!!$    REAL,INTENT(OUT),DIMENSION(:) :: VarOut 
     611!!$   
     612!!$#ifdef CPP_PARA 
     613!!$    INTEGER :: nb,i,index_para,rank 
     614!!$    INTEGER :: ierr 
     615!!$    LOGICAL :: flag=.FALSE. 
     616!!$    LOGICAL, PARAMETER :: check=.FALSE. 
     617!!$#endif 
     618!!$ 
     619!!$#ifndef CPP_PARA 
     620!!$    VarOut(:)=VarIn 
     621!!$    RETURN 
     622!!$#else 
     623!!$ 
     624!!$    IF (timer_state(timer_mpi)==running) THEN 
     625!!$      flag=.TRUE. 
     626!!$    ELSE 
     627!!$      flag=.FALSE. 
     628!!$    ENDIF 
     629!!$     
     630!!$    IF (flag) CALL suspend_timer(timer_mpi) 
     631!!$ 
     632!!$    IF (check) & 
     633!!$         WRITE(numout,*) "gather_rgen VarIn=",VarIn     
     634!!$ 
     635!!$#ifdef CPP_PARA 
     636!!$    CALL MPI_GATHER(VarIn,1,MPI_REAL_ORCH,VarOut,1,MPI_REAL_ORCH,root_prc,MPI_COMM_ORCH,ierr) 
     637!!$#endif 
     638!!$ 
     639!!$    IF (check) & 
     640!!$         WRITE(numout,*) "gather_rgen VarOut=",VarOut 
     641!!$ 
     642!!$    IF (flag) CALL resume_timer(timer_mpi) 
     643!!$#endif 
     644!!$  END SUBROUTINE gather_rs 
     645!!$ 
     646!!$  SUBROUTINE gather_ls(VarIn, VarOut) 
     647!!$    USE data_para 
     648!!$    USE timer 
     649!!$ 
     650!!$    IMPLICIT NONE 
     651!!$   
     652!!$#ifdef CPP_PARA 
     653!!$    INCLUDE 'mpif.h' 
     654!!$#endif 
     655!!$     
     656!!$    LOGICAL,INTENT(IN) :: VarIn 
     657!!$    LOGICAL,INTENT(OUT),DIMENSION(:) :: VarOut 
     658!!$   
     659!!$#ifdef CPP_PARA 
     660!!$    INTEGER :: nb,i,index_para,rank 
     661!!$    INTEGER :: ierr 
     662!!$    LOGICAL :: flag=.FALSE. 
     663!!$    LOGICAL, PARAMETER :: check=.FALSE. 
     664!!$#endif 
     665!!$ 
     666!!$#ifndef CPP_PARA 
     667!!$    VarOut(:)=VarIn 
     668!!$    RETURN 
     669!!$#else 
     670!!$ 
     671!!$    IF (timer_state(timer_mpi)==running) THEN 
     672!!$      flag=.TRUE. 
     673!!$    ELSE 
     674!!$      flag=.FALSE. 
     675!!$    ENDIF 
     676!!$     
     677!!$    IF (flag) CALL suspend_timer(timer_mpi) 
     678!!$ 
     679!!$    IF (check) & 
     680!!$         WRITE(numout,*) "gather_rgen VarIn=",VarIn     
     681!!$ 
     682!!$#ifdef CPP_PARA 
     683!!$    CALL MPI_GATHER(VarIn,1,MPI_LOGICAL,VarOut,1,MPI_LOGICAL,root_prc,MPI_COMM_ORCH,ierr) 
     684!!$#endif 
     685!!$ 
     686!!$    IF (check) & 
     687!!$         WRITE(numout,*) "gather_rgen VarOut=",VarOut 
     688!!$    IF (flag) CALL resume_timer(timer_mpi) 
     689!!$#endif 
     690!!$  END SUBROUTINE gather_ls 
    688691 
    689692!!!!! --> Les entiers 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_sechiba/intersurf.f90

    r119 r405  
    3737 
    3838  PRIVATE 
    39   PUBLIC :: intersurf_main, stom_define_history, intsurf_time 
     39  PUBLIC :: intersurf_main, stom_define_history, stom_IPCC_define_history, intsurf_time 
    4040 
    4141  INTERFACE intersurf_main 
    … …  
    6464  REAL(r_std) :: julian0 
    6565  ! 
    66   LOGICAL :: check_INPUTS = .FALSE.         !! (very) long print of INPUTs in intersurf  
     66  LOGICAL, PARAMETER :: check_INPUTS = .FALSE.         !! (very) long print of INPUTs in intersurf  
    6767  LOGICAL, SAVE :: OFF_LINE_MODE = .FALSE.  
     68  LOGICAL, SAVE :: check_time = .FALSE. 
     69  PUBLIC check_time, l_first_intersurf 
    6870  ! 
    6971CONTAINS 
    … …  
    159161    REAL(r_std),DIMENSION (kjpindex)                      :: dcoastal      !! Work array to keep coastalflow 
    160162    REAL(r_std),DIMENSION (kjpindex)                      :: driver        !! Work array to keep riverflow 
     163    REAL(r_std),DIMENSION (kjpindex)                      :: znetco2       !! Work array to keep netco2flux 
     164    REAL(r_std),DIMENSION (kjpindex)                      :: zcarblu       !! Work array to keep fco2_land_use 
    161165    REAL(r_std),DIMENSION (kjpindex)                      :: ztsol_rad     !! Work array to keep tsol_rad 
    162166    REAL(r_std),DIMENSION (kjpindex)                      :: zvevapp       !! Work array to keep vevapp 
    … …  
    354358       & zprecip_rain ,zprecip_snow,  zlwdown, zswnet, zswdown, zpb, & 
    355359! Output : Fluxes 
    356        & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, & 
     360       & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, & 
    357361! Surface temperatures and surface properties 
    358362       & ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0, & 
    … …  
    698702    REAL(r_std),DIMENSION (kjpindex)                      :: dcoastal      !! Work array to keep coastal flow 
    699703    REAL(r_std),DIMENSION (kjpindex)                      :: driver        !! Work array to keep river out flow 
     704    REAL(r_std),DIMENSION (kjpindex)                      :: znetco2       !! Work array to keep netco2flux 
     705    REAL(r_std),DIMENSION (kjpindex)                      :: zcarblu       !! Work array to keep fco2_land_use 
    700706    REAL(r_std),DIMENSION (kjpindex)                      :: ztsol_rad     !! Work array to keep tsol_rad 
    701707    REAL(r_std),DIMENSION (kjpindex)                      :: zvevapp       !! Work array to keep vevapp 
    … …  
    871877       & zprecip_rain ,zprecip_snow,  zlwdown, zswnet, zswdown, zpb, & 
    872878! Output : Fluxes 
    873        & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, & 
     879       & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, & 
    874880! Surface temperatures and surface properties 
    875881       & ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0, & 
    … …  
    12081214    REAL(r_std),DIMENSION (kjpindex)                      :: dcoastal      !! Work array to keep coastal flow 
    12091215    REAL(r_std),DIMENSION (kjpindex)                      :: driver        !! Work array to keep river out flow 
     1216    REAL(r_std),DIMENSION (kjpindex)                      :: znetco2       !! Work array to keep netco2flux 
     1217    REAL(r_std),DIMENSION (kjpindex)                      :: zcarblu       !! Work array to keep fco2_land_use 
    12101218    REAL(r_std),DIMENSION (kjpindex)                      :: ztsol_rad     !! Work array to keep tsol_rad 
    12111219    REAL(r_std),DIMENSION (kjpindex)                      :: zvevapp       !! Work array to keep vevapp 
    … …  
    15681576       & zprecip_rain ,zprecip_snow,  lwdown, swnet, swdown, pb, & 
    15691577! Output : Fluxes 
    1570        & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, & 
     1578       & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, & 
    15711579! Surface temperatures and surface properties 
    15721580       & ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0, & 
    … …  
    18471855     & tsol_rad, temp_sol_new, qsurf, albedo, emis, z0, lon_scat_g, lat_scat_g, & 
    18481856! Ajout Nathalie - passage q2m/t2m pour calcul Rveget 
    1849      & q2m, t2m)   
     1857     & q2m, t2m, & 
     1858! Add emission/deposit fields 
     1859     & field_out_names, fields_out, field_in_names, fields_in)   
    18501860#else 
    18511861  SUBROUTINE intersurf_gathered_2m (kjit, iim_glo, jjm_glo, kjpindex, kindex, xrdt, & 
    … …  
    18631873     & tsol_rad, temp_sol_new, qsurf, albedo, emis, z0, lon_scat_g, lat_scat_g, & 
    18641874! Ajout Nathalie - passage q2m/t2m pour calcul Rveget 
    1865      & q2m, t2m)   
     1875     & q2m, t2m, & 
     1876! Add emission/deposit fields 
     1877     & field_out_names, fields_out, field_in_names, fields_in) 
    18661878#endif 
    18671879    ! routines called : sechiba_main 
    … …  
    19221934    REAL(r_std),DIMENSION (kjpindex), INTENT(out)         :: fluxlat       !! Latent chaleur flux 
    19231935    REAL(r_std),DIMENSION (kjpindex), INTENT(out)         :: emis          !! Emissivity 
     1936    ! 
     1937    ! Optional arguments 
     1938    ! 
     1939    ! Names and fields for emission variables : to be transport by GCM to chemistry model. 
     1940    CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_out_names 
     1941    REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: fields_out 
     1942    ! 
     1943    ! Names and fields for deposit variables : to be transport from chemistry model by GCM to ORCHIDEE. 
     1944    CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_in_names 
     1945    REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(IN) :: fields_in 
     1946    ! 
    19241947    ! LOCAL declaration 
    19251948    ! work arrays to scatter and/or gather information just before/after sechiba_main call's 
    … …  
    19341957    REAL(r_std),DIMENSION (kjpindex)                      :: dcoastal      !! Work array to keep coastal flow 
    19351958    REAL(r_std),DIMENSION (kjpindex)                      :: driver        !! Work array to keep river out flow 
     1959    REAL(r_std),DIMENSION (kjpindex)                      :: znetco2       !! Work array to keep netco2flux 
     1960    REAL(r_std),DIMENSION (kjpindex)                      :: zcarblu       !! Work array to keep fco2_land_use 
    19361961    REAL(r_std),DIMENSION (kjpindex)                      :: ztsol_rad     !! Work array to keep tsol_rad 
    19371962    REAL(r_std),DIMENSION (kjpindex)                      :: zvevapp       !! Work array to keep vevapp 
    … …  
    19451970    ! Optional arguments 
    19461971    ! 
    1947     REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN), OPTIONAL :: lon_scat_g, lat_scat_g !! The scattered values for longitude  
     1972    REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN) :: lon_scat_g, lat_scat_g !! The scattered values for longitude  
    19481973    ! 
    19491974    INTEGER(i_std)                          :: iim,jjm                                  !! local sizes 
    … …  
    19772002    LOGICAL, SAVE                                         :: fatmco2       !! Flag to force the value of atmospheric CO2 for vegetation. 
    19782003    REAL(r_std), SAVE                                     :: atmco2        !! atmospheric CO2  
     2004    ! 
     2005    ! Number of fields to give (nb_fields_out) or get from (nb_fields_in) GCM : 
     2006    INTEGER(i_std), SAVE                                  :: nb_fields_out, nb_fields_in 
     2007    ! Id of fields to give (nb_fields_out) or get from (nb_fields_in) GCM : 
     2008    INTEGER(i_std)                                        :: i_fields_out, i_fields_in 
    19792009    ! 
    19802010    CALL ipslnlf(old_number=old_fileout) 
    … …  
    20602090       !  we have to do the work here. 
    20612091       ! 
    2062        IF ( PRESENT(lon_scat_g) .AND. PRESENT(lat_scat_g)) THEN 
     2092       IF ( .TRUE. ) THEN 
    20632093           
    20642094          lon_scat(:,:)=zero 
    … …  
    20782108             lat_g(:,:) = lat_scat_g(:,:) 
    20792109          ENDIF 
    2080  
    2081        ELSE IF ( PRESENT(lon_scat_g) .OR. PRESENT(lat_scat_g)) THEN 
    2082  
    2083           WRITE(numout,*) 'You need to provide the longitude AND latitude on the' 
    2084           WRITE(numout,*) 'gathered grid in order to start ORCHIDEE.' 
    2085           STOP 'intersurf_gathered' 
    20862110 
    20872111       ELSE 
    … …  
    21952219       ENDIF 
    21962220       ! 
     2221 
     2222       ! Prepare fieds out/in for interface with GCM. 
     2223       IF (PRESENT(field_out_names)) THEN 
     2224          nb_fields_out=SIZE(field_out_names) 
     2225       ELSE 
     2226          nb_fields_out=0 
     2227       ENDIF 
     2228       IF (PRESENT(field_in_names)) THEN 
     2229          nb_fields_in=SIZE(field_in_names) 
     2230       ELSE 
     2231          nb_fields_in=0 
     2232       ENDIF 
     2233 
    21972234       IF ( check ) WRITE(numout,*) 'End of Initialisation of intersurf' 
    21982235       ! 
    … …  
    22482285       WRITE(numout,*) "Fraction of continent in the grid = ",zcontfrac 
    22492286    ENDIF 
     2287 
     2288 
     2289    ! Fields for deposit variables : to be transport from chemistry model by GCM to ORCHIDEE. 
     2290    WRITE(numout,*) "Get fields from atmosphere." 
     2291 
     2292    DO i_fields_in=1,nb_fields_in 
     2293       WRITE(numout,*) i_fields_in," Champ = ",TRIM(field_in_names(i_fields_in))  
     2294       SELECT CASE(TRIM(field_in_names(i_fields_in))) 
     2295       CASE DEFAULT  
     2296          CALL ipslerr (3,'intsurf_gathered_2m', & 
     2297            &          'You ask in GCM an unknown field '//TRIM(field_in_names(i_fields_in))//& 
     2298            &          ' to give to ORCHIDEE for this specific version.',& 
     2299            &          'This model won''t be able to continue.', & 
     2300            &          '(check your tracer parameters in GCM)') 
     2301       END SELECT 
     2302    ENDDO 
     2303     
    22502304    ! 
    22512305    ! 2. modification of co2 
    … …  
    22982352       & zprecip_rain ,zprecip_snow,  lwdown, swnet, swdown, pb, & 
    22992353! Output : Fluxes 
    2300        & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, & 
     2354       & zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, & 
    23012355! Surface temperatures and surface properties 
    23022356       & ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0, & 
    … …  
    25522606    ENDDO 
    25532607    ! 
     2608    WRITE(numout,*) "Give fields to atmosphere." 
     2609     
     2610    ! Fields for emission variables : to be transport by GCM to chemistry model. 
     2611    DO i_fields_out=1,nb_fields_out 
     2612       SELECT CASE(TRIM(field_out_names(i_fields_out))) 
     2613       CASE("fCO2_land")  
     2614          fields_out(:,i_fields_out)=znetco2(:) 
     2615       CASE("fCO2_land_use") 
     2616          fields_out(:,i_fields_out)=zcarblu(:) 
     2617       CASE DEFAULT  
     2618          CALL ipslerr (3,'intsurf_gathered_2m', & 
     2619            &          'You ask from GCM an unknown field '//TRIM(field_out_names(i_fields_out))//& 
     2620            &          ' to ORCHIDEE for this specific version.',& 
     2621            &          'This model won''t be able to continue.', & 
     2622            &          '(check your tracer parameters in GCM)') 
     2623       END SELECT 
     2624    ENDDO 
     2625    ! 
    25542626    IF ( lrestart_write .AND. ok_watchout .AND. is_root_prc ) THEN 
    25552627       CALL watchout_close() 
    … …  
    25772649    REAL(r_std), INTENT(in)                     :: dt        !! Time step 
    25782650    ! 
    2579     ! LOCAL 
    2580     LOGICAL     :: check=.FALSE. 
    25812651 
    25822652    IF (l_first_intersurf) THEN 
    … …  
    25902660       ENDIF 
    25912661 
    2592        IF (check) THEN 
     2662       IF (check_time) THEN 
    25932663          write(numout,*) "calendar_str =",calendar_str 
    25942664          write(numout,*) "one_year=",one_year,", one_day=",one_day 
    … …  
    25982668 
    25992669    ! 
    2600     IF (check) & 
     2670    IF (check_time) & 
    26012671         WRITE(numout,*) "---"  
    26022672    ! Dans diffuco (ie date0 == date0_shift !!)  
    … …  
    26122682!!$       julian_diff = in_julian 
    26132683!!$       month_len = ioget_mon_len (year,month) 
    2614 !!$       IF (check) THEN 
     2684!!$       IF (check_time) THEN 
    26152685!!$          write(numout,*) "in_julian, jur, julian_diff=",in_julian, jur, julian_diff 
    26162686!!$          write(numout,*) 'DATE ymds', year, month,'(',month_len,'d)', day, sec, '-- stp --', istp 
    … …  
    26242694       sec = NINT((julian_diff-REAL(INT(julian_diff)))*one_day) 
    26252695       month_len = ioget_mon_len (year,month) 
    2626        IF (check) THEN 
     2696       IF (check_time) THEN 
    26272697          write(numout,*) "2 in_julian, julian0, julian_diff=",in_julian, julian0, julian_diff 
    26282698          write(numout,*) '2 DATE ymds', year, month,'(',month_len,'d)', day, sec, '-- stp --', istp 
    … …  
    26342704!!$       julian_diff = in_julian 
    26352705!!$       month_len = ioget_mon_len (year,month) 
    2636 !!$       IF (check) THEN 
     2706!!$       IF (check_time) THEN 
    26372707!!$          write(numout,*) "in_julian=",in_julian, jur, julian_diff 
    26382708!!$          write(numout,*) 'DATE ymds', year, month,'(',month_len,'d)', day, sec, '-- stp --', istp 
    … …  
    26462716!!$       sec = NINT((julian_diff-REAL(INT(julian_diff)))*one_day) 
    26472717!!$       month_len = ioget_mon_len (year,month) 
    2648 !!$       IF (check) THEN 
     2718!!$       IF (check_time) THEN 
    26492719!!$          write(numout,*) "2 in_julian, jur, julian_diff=",in_julian, jur, julian_diff 
    26502720!!$          write(numout,*) '2 DATE ymds', year, month,'(',month_len,'d)', day, sec, '-- stp --', istp 
    … …  
    26522722 
    26532723 
    2654 !!$       IF (check) & 
     2724!!$       IF (check_time) & 
    26552725!!$            WRITE(numout,*) "-" 
    26562726 
    … …  
    26632733       julian_diff = in_julian 
    26642734       month_len = ioget_mon_len (year,month) 
    2665        IF (check) THEN 
     2735       IF (check_time) THEN 
    26662736          write(numout,*) "in_julian=",in_julian, julian0, julian_diff 
    26672737          write(numout,*) 'DATE ymds', year, month,'(',month_len,'d)', day, sec, '-- stp --', istp 
    26682738       ENDIF 
    26692739    ENDIF 
    2670 !!$    IF (check) & 
     2740!!$    IF (check_time) & 
    26712741!!$         WRITE(numout,*) "---"  
    26722742 
    … …  
    26942764    long_print = .FALSE. 
    26952765    CALL getin_p('LONGPRINT',long_print) 
     2766    ! 
     2767    !Config Key  = CHECKTIME 
     2768    !Config Desc = ORCHIDEE will print messages on time 
     2769    !Config Def  = n 
     2770    !Config Help = This flag permits to print debug messages on the time. 
     2771    ! 
     2772    check_time = .FALSE. 
     2773    CALL getin_p('CHECKTIME',check_time) 
    26962774    ! 
    26972775    ! 
    … …  
    27862864    CALL getin_p('STOMATE_OK_DGVM',control_flags%ok_dgvm) 
    27872865 
    2788     IF ( control_flags%ok_dgvm ) THEN 
    2789        WRITE(numout,*) 'You try to use LPJ ',control_flags%ok_dgvm, & 
    2790             ' with this version. ' 
    2791        WRITE(numout,*) 'It is not possible because it has to be modified ', & 
    2792             ' to give correct values.' 
    2793        CALL ipslerr (3,'intsurf_config', & 
    2794          &          'Use of STOMATE_OK_DGVM not allowed with this version.',& 
    2795          &          'ORCHIDEE will stop.', & 
    2796          &          'Please disable DGVM to use this version of ORCHIDEE.') 
    2797     ENDIF 
    27982866    ! 
    27992867    ! control initialisation with sechiba 
    … …  
    44484516            & hist_pool_10axis_id, hist_pool_100axis_id, & 
    44494517            & hist_pool_11axis_id, hist_pool_101axis_id) 
    4450 ! deforestation axis added as arguments 
    44514518 
    44524519       !- end definition 
    … …  
    48004867         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(3), dt, hist_dt) 
    48014868 
     4869    ! Adaptation to climate 
     4870    CALL histdef (hist_id_stom, & 
     4871         &               TRIM("ADAPTATION          "), & 
     4872         &               TRIM("Adaptation to climate (DGVM)                      "), & 
     4873         &               TRIM("-                   "), iim,jjm, hist_hori_id, & 
     4874         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(3), dt, hist_dt) 
     4875     
     4876    ! Probability from regenerative 
     4877    CALL histdef (hist_id_stom, & 
     4878         &               TRIM("REGENERATION        "), & 
     4879         &               TRIM("Probability from regenerative (DGVM)               "), & 
     4880         &               TRIM("-                   "), iim,jjm, hist_hori_id, & 
     4881         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(3), dt, hist_dt) 
     4882 
     4883    ! crown area of individuals (m**2) 
     4884    CALL histdef (hist_id_stom, & 
     4885         &               TRIM("CN_IND              "), & 
     4886         &               TRIM("crown area of individuals                         "), & 
     4887         &               TRIM("m^2                 "), iim,jjm, hist_hori_id, & 
     4888         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(3), dt, hist_dt) 
     4889 
     4890    ! woodmass of individuals (gC) 
     4891    CALL histdef (hist_id_stom, & 
     4892         &               TRIM("WOODMASS_IND        "), & 
     4893         &               TRIM("Woodmass of individuals                           "), & 
     4894         &               TRIM("gC/pft              "), iim,jjm, hist_hori_id, & 
     4895         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(3), dt, hist_dt) 
     4896 
    48024897    ! total living biomass 
    48034898    CALL histdef (hist_id_stom, & 
    … …  
    50305125         &               TRIM("1/day               "), iim,jjm, hist_hori_id, & 
    50315126         &               nvm,1,nvm, hist_PFTaxis_id,32, ave(6), dt, hist_dt) 
     5127 
     5128    ! Establish tree 
     5129    CALL histdef (hist_id_stom, & 
     5130         &               TRIM("ESTABTREE           "), & 
     5131         &               TRIM("Rate of tree establishement                       "), & 
     5132         &               TRIM("1/day               "), iim,jjm, hist_hori_id, & 
     5133         &               1,1,1, -99,32, ave(6), dt, hist_dt) 
     5134 
     5135    ! Establish grass 
     5136    CALL histdef (hist_id_stom, & 
     5137         &               TRIM("ESTABGRASS          "), & 
     5138         &               TRIM("Rate of grass establishement                      "), & 
     5139         &               TRIM("1/day               "), iim,jjm, hist_hori_id, & 
     5140         &               1,1,1, -99,32, ave(6), dt, hist_dt) 
    50325141 
    50335142    ! Fraction of plants that dies (light competition)   
    … …  
    52685377         &               TRIM("Carbon in Products of Land Use Change"), & 
    52695378         &               TRIM("kg C m-2"), iim,jjm, hist_hori_id, & 
     5379         &               1,1,1, -99,32, ave(1), dt, hist_dt) 
     5380    ! Carbon Mass Variation 
     5381    CALL histdef (hist_id_stom_IPCC, & 
     5382         &               TRIM("cMassVariation"), & 
     5383         &               TRIM("Terrestrial Carbon Mass Variation"), & 
     5384         &               TRIM("kg C m-2 s-1"), iim,jjm, hist_hori_id, & 
    52705385         &               1,1,1, -99,32, ave(1), dt, hist_dt) 
    52715386    ! Leaf Area Fraction 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_sechiba/routing.f90

    r119 r405  
    575575       CALL ioconf_setatt('LONG_NAME','Time counter for the routing scheme') 
    576576       CALL restget (rest_id, var_name, 1, 1, 1, kjit, .TRUE., tmp_day) 
    577        time_counter = tmp_day(1)  
     577       IF (tmp_day(1) == val_exp) THEN 
     578          time_counter = zero 
     579       ELSE 
     580          time_counter = tmp_day(1)  
     581       ENDIF 
    578582       CALL setvar (time_counter, val_exp, 'NO_KEYWORD', zero) 
    579583    ENDIF 
    … …  
    678682    CALL ioconf_setatt('LONG_NAME','Water in the lake reservoir') 
    679683    CALL restget_p (rest_id, var_name, nbp_glo, 1, 1, kjit, .TRUE., lake_reservoir, "gather", nbp_glo, index_g) 
    680     CALL setvar (lake_reservoir, val_exp, 'NO_KEYWORD', zero) 
     684    CALL setvar_p (lake_reservoir, val_exp, 'NO_KEYWORD', zero) 
    681685    ! 
    682686    ! Map of irrigated areas 
    … …  
    10401044!ym mais n'est pas la plus efficace 
    10411045 
    1042     IF (is_root_prc) & 
    1043          ALLOCATE( fast_flow_g(nbp_glo, nbasmax), slow_flow_g(nbp_glo, nbasmax), & 
    1044           stream_flow_g(nbp_glo, nbasmax), floods_g(nbp_glo, nbasmax), wdelay_g(nbp_glo, nbasmax) ) 
     1046    IF (is_root_prc)  THEN 
     1047       ALLOCATE( fast_flow_g(nbp_glo, nbasmax), slow_flow_g(nbp_glo, nbasmax), & 
     1048            stream_flow_g(nbp_glo, nbasmax), floods_g(nbp_glo, nbasmax),  & 
     1049            wdelay_g(nbp_glo, nbasmax) ) 
     1050    ELSE 
     1051       ALLOCATE( fast_flow_g(1,1), slow_flow_g(1,1), & 
     1052            stream_flow_g(1, 1), floods_g(1,1),  & 
     1053            wdelay_g(1,1) ) 
     1054    ENDIF 
    10451055     
    10461056        
    … …  
    10641074    ENDIF 
    10651075 
    1066     IF (is_root_prc) & 
    1067          DEALLOCATE( fast_flow_g, slow_flow_g, stream_flow_g, floods_g, wdelay_g ) 
     1076    DEALLOCATE( fast_flow_g, slow_flow_g, stream_flow_g, floods_g, wdelay_g ) 
    10681077    
    10691078    CALL scatter(transport_glo,transport) 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_sechiba/sechiba.f90

    r119 r405  
    187187    & precip_rain, precip_snow, lwdown, swnet, swdown, pb, & 
    188188         ! Output : Fluxes 
    189     & vevapp, fluxsens, fluxlat, coastalflow, riverflow, & 
     189    & vevapp, fluxsens, fluxlat, coastalflow, riverflow, netco2flux, fco2_lu, & 
    190190         ! Surface temperatures and surface properties 
    191191    & tsol_rad, temp_sol_new, qsurf_out, albedo_out, emis_out, z0_out, & 
    … …  
    250250    REAL(r_std),DIMENSION (kjpindex), INTENT (out)           :: fluxlat          !! Latent chaleur flux 
    251251    REAL(r_std),DIMENSION (kjpindex), INTENT (out)           :: emis_out         !! Emissivity 
    252  
     252    REAL(r_std),DIMENSION (kjpindex), INTENT (out)           :: netco2flux       !! Sum CO2 flux over PFTs (gC/m**2 of average ground/s) 
     253    REAL(r_std),DIMENSION (kjpindex), INTENT (out)           :: fco2_lu          !! Land Cover Change CO2 flux (gC/m**2 of average ground/s) 
     254 
     255    ! local declaration 
     256    INTEGER(i_std)                                :: jv 
    253257    REAL(r_std), ALLOCATABLE, DIMENSION (:)                  :: runoff1,drainage1, soilcap1,soilflx1 
    254258    REAL(r_std), ALLOCATABLE, DIMENSION (:,:)                :: shumdiag1 
    … …  
    318322            lai, height, veget, frac_nobio, veget_max, totfrac_nobio, qsintmax, & 
    319323            rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    320             co2_flux) 
     324            co2_flux, fco2_lu) 
     325       netco2flux(:) = zero 
     326       DO jv = 2,nvm 
     327          netco2flux(:) = netco2flux(:) + co2_flux(:,jv)*veget_max(:,jv) 
     328       ENDDO 
    321329       !  
    322330       ! computes initialisation of diffusion coeff 
    … …  
    566574         lai, height, veget, frac_nobio, veget_max, totfrac_nobio, qsintmax, & 
    567575         rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    568          co2_flux) 
    569  
     576         co2_flux, fco2_lu) 
     577    ! 
     578    ! Compute global CO2 flux 
     579    ! 
     580    netco2flux(:) = zero 
     581    DO jv = 2,nvm 
     582       netco2flux(:) = netco2flux(:) + co2_flux(:,jv)*veget_max(:,jv) 
     583    ENDDO 
    570584    ! 
    571585    ! call swap from new computed variables   
    … …  
    809823            lai, height, veget, frac_nobio, veget_max, totfrac_nobio, qsintmax, & 
    810824            rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    811             co2_flux) 
    812  
     825            co2_flux, fco2_lu) 
     826       netco2flux(:) = zero 
     827       DO jv = 2,nvm 
     828          netco2flux(:) = netco2flux(:) + co2_flux(:,jv)*veget_max(:,jv) 
     829       ENDDO 
    813830 
    814831       var_name= 'shumdiag'   

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_sechiba/slowproc.f90

    r119 r405  
    8181       lai, height, veget, frac_nobio, veget_max, totfrac_nobio, qsintmax, & 
    8282       rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    83        co2_flux) 
     83       co2_flux, fco2_lu) 
    8484 
    8585 
    … …  
    120120    ! output fields 
    121121    REAL(r_std), DIMENSION (kjpindex,nvm), INTENT(out)      :: co2_flux         !! CO2 flux in gC/m**2 of average ground/second 
     122    REAL(r_std),DIMENSION (kjpindex), INTENT (out)          :: fco2_lu          !! Land Cover Change CO2 flux (gC/m**2 of average ground/s) 
    122123    ! modified scalar 
    123124    ! modified fields 
    … …  
    193194               veget_nextyear, totfrac_nobio_nextyear, & 
    194195               hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    195                co2_flux,resp_maint,resp_hetero,resp_growth) 
     196               co2_flux, fco2_lu, resp_maint,resp_hetero,resp_growth) 
    196197          ! 
    197198       ENDIF 
    … …  
    289290               veget_nextyear, totfrac_nobio_nextyear, & 
    290291               hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    291                co2_flux,resp_maint,resp_hetero,resp_growth) 
     292               co2_flux, fco2_lu, resp_maint,resp_hetero,resp_growth) 
    292293       ENDIF 
    293294 
    … …  
    387388            veget_nextyear, totfrac_nobio_nextyear, & 
    388389            hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    389             co2_flux,resp_maint,resp_hetero,resp_growth) 
     390            co2_flux, fco2_lu, resp_maint,resp_hetero,resp_growth) 
    390391       IF ( control%ok_stomate .AND. control%ok_sechiba ) THEN 
    391392          CALL histwrite(hist_id, 'maint_resp', kjit, resp_maint, kjpindex*nvm, indexveg) 
    … …  
    775776       ! to be in sechiba when teststomate will have disapeared. 
    776777!MM Problem here with dpu which depends on soil type 
    777     DO jv = 1, nbdl-1 
     778    DO l = 1, nbdl-1 
    778779       ! first 2.0 is dpu  
    779780       ! second 2.0 is average 
    780        diaglev(jv) = dpu_cste/(2**(nbdl-1) -1) * ( ( 2**(jv-1) -1) + ( 2**(jv) -1) ) / 2.0 
     781       diaglev(l) = dpu_cste/(2**(nbdl-1) -1) * ( ( 2**(l-1) -1) + ( 2**(l) -1) ) / 2.0 
    781782    ENDDO 
    782783    diaglev(nbdl) = dpu_cste 
    … …  
    26792680          !    et PFT naturel / (somme des vegets - somme des vegets anthropiques) 
    26802681          !       est conservee. 
    2681           ! Sum veget_next = old (sum veget_next Naturel) + (sum veget_next Anthropic)  
    2682           !           = new (sum veget_next Naturel) + (sum veget_next Anthropic) 
    2683           !    a / (S-A) = e / (S-B) ; b/(S-A) = f/(S-B) 
     2682          ! Modification de Nathalie :  
     2683          ! Si les PFTs anthropique diminue, on les remplace plutôt par du sol nu. 
     2684          ! Le DGVM est chargé de ré-introduire les PFTs naturels. 
    26842685          IF (sumf > min_sechiba) THEN 
    26852686             sumvAnthro_old = zero 
    … …  
    26882689                IF ( .NOT. natural(jv) ) THEN 
    26892690                   veget_next(ib,jv) = veget_next(ib,jv) / sumf 
    2690                    sumvAnthro = sumvAnthro + veget_last(ib,jv) 
     2691                   sumvAnthro = sumvAnthro + veget_next(ib,jv) 
    26912692                   sumvAnthro_old = sumvAnthro_old + veget_last(ib,jv) 
    26922693                ENDIF 
    26932694             ENDDO 
    2694              ! conservation : 
    2695              rapport = ( sum_veg - sumvAnthro ) / ( sum_veg - sumvAnthro_old ) 
    2696              veget_next(ib,1) = veget_last(ib,1) * rapport 
    2697              DO jv = 2, nvm 
    2698                 IF ( .NOT. natural(jv) ) THEN 
    2699                    veget_next(ib,jv) = veget_last(ib,jv) * rapport 
    2700                 ENDIF 
    2701              ENDDO 
     2695 
     2696             IF ( sumvAnthro_old < sumvAnthro ) THEN 
     2697                ! deforestation 
     2698                ! conservation : 
     2699                rapport = ( sum_veg - sumvAnthro ) / ( sum_veg - sumvAnthro_old ) 
     2700                DO jv = 1, nvm 
     2701                   IF ( natural(jv) ) THEN 
     2702                      veget_next(ib,jv) = veget_last(ib,jv) * rapport 
     2703                   ENDIF 
     2704                ENDDO 
     2705             ELSE 
     2706                ! reforestation 
     2707                DO jv = 1, nvm 
     2708                   IF ( natural(jv) ) THEN 
     2709                      veget_next(ib,jv) = veget_last(ib,jv) 
     2710                   ENDIF 
     2711                ENDDO 
     2712                veget_next(ib,1) = veget_next(ib,1) + sumvAnthro_old - sumvAnthro 
     2713             ENDIF 
     2714 
    27022715             ! test 
    2703              IF ( ABS( SUM(veget_next(ib,:)) - sum_veg ) > EPSILON(un) ) THEN 
     2716             IF ( ABS( SUM(veget_next(ib,:)) - sum_veg ) > 10*EPSILON(un) ) THEN 
    27042717                WRITE(numout,*) "No conservation of sum of veget for point ",ib,",(",lalo(ib,1),",",lalo(ib,2),")"  
    27052718                WRITE(numout,*) "last sum of veget ",sum_veg," new sum of veget ",SUM(veget_next(ib,:))," error : ",& 
    27062719                     &                         SUM(veget_next(ib,:)) - sum_veg 
    2707                 WRITE(numout,*) "Anthropic modificaztions : last ",sumvAnthro_old," new ",sumvAnthro                 
     2720                WRITE(numout,*) "Anthropic modifications : last ",sumvAnthro_old," new ",sumvAnthro                 
    27082721                CALL ipslerr (3,'slowproc_update', & 
    27092722                     &          'No conservation of sum of veget_next', & 
    … …  
    28892902    ! 
    28902903    IF (MAXVAL(vegmap) .LT. nolson) THEN 
    2891       WRITE(*,*) 'WARNING -- WARNING' 
    2892       WRITE(*,*) 'The vegetation map has to few vegetation types.' 
    2893       WRITE(*,*) 'If you are lucky it will work but please check' 
     2904       WRITE(numout,*) 'WARNING -- WARNING' 
     2905       WRITE(numout,*) 'The vegetation map has to few vegetation types.' 
     2906       WRITE(numout,*) 'If you are lucky it will work but please check' 
    28942907    ELSE IF ( MAXVAL(vegmap) .GT. nolson) THEN 
    2895       WRITE(*,*) 'More vegetation types in file than the code can' 
    2896       WRITE(*,*) 'deal with.: ',  MAXVAL(vegmap),  nolson 
    2897       STOP 'slowproc_interpol' 
     2908       WRITE(numout,*) 'More vegetation types in file than the code can' 
     2909       WRITE(numout,*) 'deal with.: ',  MAXVAL(vegmap),  nolson 
     2910       STOP 'slowproc_interpol' 
    28982911    ENDIF 
    28992912    ! 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_constraints.f90

    r119 r405  
    147147          IF ( tree(j) .AND. ( pheno_crit%pheno_model(j) .NE. 'none' ) ) THEN 
    148148 
    149              WHERE ( when_growthinit(:,j) .GT. too_long*one_year ) 
     149             WHERE ( when_growthinit(:,j) .GT. too_long*one_year .AND. when_growthinit(:,j).LT. large_value) 
    150150                adapted(:,j) = zero 
    151151             ENDWHERE 
    … …  
    199199    ENDDO 
    200200 
     201    CALL histwrite (hist_id_stomate, 'ADAPTATION', itime, & 
     202         adapted, npts*nvm, horipft_index) 
     203    CALL histwrite (hist_id_stomate, 'REGENERATION', itime, & 
     204         regenerate, npts*nvm, horipft_index) 
     205 
    201206    IF (bavard.GE.4) WRITE(numout,*) 'Leaving constraints' 
    202207 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_cover.f90

    r119 r405  
    2323 
    2424  SUBROUTINE cover (npts, cn_ind, ind, biomass, & 
    25        veget_max, veget_max_old, veget, lai, litter, carbon) 
     25       veget_max, veget_max_old, veget, lai, litter, carbon, turnover_daily, bm_to_litter) 
    2626 
    2727    ! 
    … …  
    3737    ! density of individuals (1/(m**2 of ground)) 
    3838    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)          :: ind 
     39    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground at beginning of time step 
    3940    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)          :: veget_max_old 
    4041 
    … …  
    4445    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground 
    4546    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)       :: veget_max 
     47    ! Turnover rates (gC/(m**2 of ground)/day) 
     48    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)          :: turnover_daily 
     49    ! conversion of biomass to litter (g/m**2 / day 
     50    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)          :: bm_to_litter 
    4651 
    4752    ! 0.3 output 
    … …  
    5055    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)       :: veget 
    5156    ! leaf area index OF AN INDIVIDUAL PLANT 
    52     REAL(r_std), DIMENSION(npts,nvm), INTENT(in)         :: lai 
     57    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)         :: lai 
    5358 
    5459    ! metabolic and structural litter, above and below ground (gC/(m**2 of ground)) 
    … …  
    6065 
    6166    ! index 
    62     INTEGER(i_std)                                         :: i,j 
     67    INTEGER(i_std)                                         :: i,j,k,m 
    6368 
    6469    ! Litter dilution (gC/m²) 
    … …  
    6873 
    6974    ! conversion vectors 
    70     REAL(r_std),DIMENSION(nvm)                                         :: delta_veg 
     75    REAL(r_std),DIMENSION(nvm)                                         :: delta_veg,reduct 
    7176    ! vecteur de conversion 
    72     REAL(r_std)                                                        :: delta_veg_sum 
     77    REAL(r_std)                                                        :: delta_veg_sum,diff,sr 
     78    REAL(r_std), DIMENSION(npts)                                       :: frac_nat,sum_vegettree,sum_vegetgrass 
     79    REAL(r_std), DIMENSION(npts)                                       :: sum_veget_natveg 
    7380 
    7481    ! ========================================================================= 
    … …  
    8188    IF ( control%ok_dgvm ) THEN 
    8289 
    83        veget_max(:,ibare_sechiba) = 1. 
     90       ! some initialisations 
     91       frac_nat(:) = un 
     92       sum_veget_natveg(:) = zero 
     93       sum_vegettree(:) = zero 
     94       sum_vegetgrass(:) = zero 
     95 
     96       veget_max(:,ibare_sechiba) = un 
    8497 
    8598       DO j = 2,nvm 
    … …  
    88101 
    89102             veget_max(:,j) = ind(:,j) * cn_ind(:,j) 
    90  
    91           ENDIF 
    92  
     103             sum_veget_natveg(:) = sum_veget_natveg(:) + veget_max(:,j) 
     104 
     105          ELSE 
     106             !fraction occupied by agriculture needs to be substracted for the DGVM 
     107             !this is used below to constrain veget for natural vegetation, see below 
     108             frac_nat(:) = frac_nat(:) - veget_max(:,j) 
     109 
     110          ENDIF 
     111 
     112       ENDDO 
     113 
     114       DO i = 1, npts  
     115 
     116          IF (sum_veget_natveg(i) .GT. frac_nat(i) .AND. frac_nat(i) .GT. min_stomate) THEN 
     117 
     118             DO j = 2,nvm 
     119                IF( natural(j) ) THEN 
     120                   veget_max(i,j) =  veget_max(i,j) * frac_nat(i) / sum_veget_natveg(i) 
     121                ENDIF 
     122             ENDDO 
     123 
     124          ENDIF 
     125       ENDDO 
     126 
     127       DO j = 2,nvm 
    93128          veget_max(:,ibare_sechiba) = veget_max(:,ibare_sechiba) - veget_max(:,j) 
    94  
    95        ENDDO 
    96  
     129       ENDDO 
    97130       veget_max(:,ibare_sechiba) = MAX( veget_max(:,ibare_sechiba), zero ) 
    98131 
     132       ! 1.3 calculate carbon fluxes between PFTs to maintain mass balance 
     133       ! 
     134 
     135       DO i = 1, npts          
     136          ! Generation of the conversion vector 
     137 
     138          delta_veg(:) = veget_max(i,:)-veget_max_old(i,:) 
     139          delta_veg_sum = SUM(delta_veg,MASK=delta_veg.LT.zero) 
     140 
     141          dilu_lit(i,:,:) = zero 
     142          dilu_soil_carbon(i,:) = zero 
     143          DO j=1, nvm 
     144             IF ( delta_veg(j) < -min_stomate ) THEN  
     145                dilu_lit(i,:,:)=  dilu_lit(i,:,:) + delta_veg(j)*litter(i,:,j,:) / delta_veg_sum 
     146                dilu_soil_carbon(i,:)=  dilu_soil_carbon(i,:) + delta_veg(j) * carbon(i,:,j) / delta_veg_sum 
     147             ENDIF 
     148          ENDDO 
     149 
     150          DO j=1, nvm 
     151             IF ( delta_veg(j) > min_stomate) THEN 
     152 
     153                ! Dilution of reservoirs 
     154 
     155                ! Litter 
     156                litter(i,:,j,:)=(litter(i,:,j,:) * veget_max_old(i,j) + dilu_lit(i,:,:) * delta_veg(j)) / veget_max(i,j) 
     157 
     158                ! Soil carbon 
     159                carbon(i,:,j)=(carbon(i,:,j) * veget_max_old(i,j) + dilu_soil_carbon(i,:) * delta_veg(j)) / veget_max(i,j) 
     160 
     161             ENDIF 
     162 
     163             IF(j.GE.2.AND.veget_max_old(i,j).GT.min_stomate.AND.veget_max(i,j).GT.min_stomate) THEN 
     164 
     165                ! Correct biomass densities (i.e. also litter fall) to conserve mass  
     166                ! since it's defined on veget_max 
     167 
     168                biomass(i,j,:)=biomass(i,j,:)*veget_max_old(i,j)/veget_max(i,j) 
     169                turnover_daily(i,j,:)=turnover_daily(i,j,:)*veget_max_old(i,j)/veget_max(i,j) 
     170                bm_to_litter(i,j,:)=bm_to_litter(i,j,:)*veget_max_old(i,j)/veget_max(i,j) 
     171 
     172             ENDIF 
     173 
     174          ENDDO 
     175       ENDDO 
    99176    ENDIF 
    100  
    101     DO i = 1, npts          
    102        ! Generation of the conversion vector 
    103  
    104        delta_veg(:) = veget_max(i,:)-veget_max_old(i,:) 
    105        delta_veg_sum = SUM(delta_veg,MASK=delta_veg.LT.zero) 
    106  
    107        dilu_lit(i,:,:) = zero 
    108        dilu_soil_carbon(i,:) = zero 
    109        DO j=1, nvm 
    110           IF ( delta_veg(j) < -min_stomate ) THEN  
    111              dilu_lit(i,:,:)=  dilu_lit(i,:,:) - delta_veg(j)*litter(i,:,j,:) / delta_veg_sum 
    112              dilu_soil_carbon(i,:)=  dilu_soil_carbon(i,:) - delta_veg(j) * carbon(i,:,j) / delta_veg_sum 
    113           ENDIF 
    114        ENDDO 
    115  
    116        DO j=1, nvm 
    117           IF ( delta_veg(j) > min_stomate) THEN 
    118  
    119              ! Dilution of reservoirs 
    120  
    121              ! Litter 
    122              litter(i,:,j,:)=(litter(i,:,j,:) * veget_max_old(i,j) + dilu_lit(i,:,:) * delta_veg(j)) / veget_max(i,j) 
    123  
    124              ! Soil carbon 
    125              carbon(i,:,j)=(carbon(i,:,j) * veget_max_old(i,j) + dilu_soil_carbon(i,:) * delta_veg(j)) / veget_max(i,j) 
    126  
    127           ENDIF 
    128           !SZ correct biomass to conserve mass since it's defined on veget_max 
    129           IF(j.GE.2.AND.veget_max_old(i,j).GT.min_stomate.AND.veget_max(i,j).GT.min_stomate) THEN 
    130              biomass(i,j,:)=biomass(i,j,:)*veget_max_old(i,j)/veget_max(i,j) 
    131           ENDIF 
    132  
    133        ENDDO 
    134     ENDDO 
    135177 
    136178    ! 
    … …  
    140182    ! 
    141183    !MM in Soenke code but not in merge version ; must keep that ?? 
     184!NV, MM : we keep those comments for compatibility with CMIP5 computations. 
     185!! They have to be uncommented avec CMIP5 versions in the trunk ! 
    142186!!$    DO j = 2,nvm 
    143187!!$       lai(:,j) = biomass(:,j,ileaf,icarbon)*sla(j) 
    … …  
    153197             veget(i,j) = veget_max(i,j) 
    154198          ELSE 
    155              veget(i,j) = veget_max(i,j) * ( un - exp( - lai(i,j) * ext_coeff(j) ) ) 
     199             IF ( control%ok_dgvm ) THEN 
     200!!$SZneed to check this - this formulation will cause 100% veget, otherwise there will always  
     201!!$ be some percent bare ground 
     202                veget(i,j) = ind(i,j) * cn_ind(i,j)  * ( un - EXP( - lai(i,j) * ext_coeff(j) ) ) 
     203             ELSE 
     204                veget(i,j) = veget_max(i,j) * ( un - EXP( - lai(i,j) * ext_coeff(j) ) ) 
     205             ENDIF 
     206          ENDIF 
     207 
     208          ! check sums of fpc for natural vegetation (see correction below!) in dynamic mode 
     209          IF ( control%ok_dgvm ) THEN 
     210 
     211             IF(natural(j))THEN 
     212                IF(tree(j)) THEN 
     213                   sum_vegettree(i)=sum_vegettree(i)+veget(i,j) 
     214                ELSE  
     215                   sum_vegetgrass(i)=sum_vegetgrass(i)+veget(i,j) 
     216                ENDIF 
     217             ENDIF 
     218 
    156219          ENDIF 
    157220       ENDDO 
    158221    ENDDO 
    159     ! 
     222 
     223 
     224    ! 3.1 correct gridscale fpc for dynamic vegetation 
     225!!$SZ, this part should be obsolete now that veget_max is forced to 1.0 
     226!!$ nevertheless maintained just for savety. Whoever wants to test 
     227!!$ whether this works without is invited to do so. 
     228 
     229    ! in the DGVM mode, we can arrive at a sum of veget slighly exceeding 1.0, 
     230    ! because mainly of grass dynamics... 
     231    ! In this case, we devide the fpar over natural vegetation first such that  
     232    ! grasses are shadowed by trees, and in the theoretically impossible case that 
     233    ! this is not sufficient, reduce proportionally all veget's. 
     234    ! 
     235    IF ( control%ok_dgvm ) THEN 
     236 
     237       DO i = 1,npts 
     238 
     239          diff=sum_vegettree(i)+sum_vegetgrass(i)-frac_nat(i) 
     240          reduct(:) = 0. 
     241          ! ordinary case, the reason too much grasses  
     242          ! reduce grass veget to match the maximum 
     243          IF (diff .GT. 0. ) THEN 
     244 
     245             IF (sum_vegetgrass(i).GT.min_stomate) THEN 
     246                sr=0. 
     247                DO j=2,nvm 
     248                   IF(natural(j).AND..NOT.tree(j)) THEN 
     249                      reduct(j)=-MIN(diff,sum_vegetgrass(i))*veget(i,j)/sum_vegetgrass(i) 
     250                      sr=sr+reduct(j) 
     251                   ENDIF 
     252                ENDDO 
     253                diff=diff+sr 
     254             ENDIF 
     255 
     256          ENDIF 
     257 
     258          ! this is theoretically impossible, since trees can only occupy 95%, 
     259          ! but better be save than sorry 
     260          IF (diff .GT. min_stomate ) THEN 
     261 
     262             IF (sum_vegettree(i).GT.min_stomate) THEN 
     263                sr=0. 
     264                DO j=2,nvm 
     265                   IF(natural(j).AND.tree(j)) THEN 
     266                      reduct(j)=-MIN(diff,sum_vegettree(i))*veget(i,j)/sum_vegettree(i) 
     267                      sr=sr+reduct(j) 
     268                   ENDIF 
     269                ENDDO 
     270                diff=diff+sr  
     271             ENDIF 
     272 
     273          ENDIF 
     274 
     275!!$          ! tell user if the problem could not be resolved 
     276!!$          ! in theory the model should stop here! 
     277!!$          IF (diff .GT. min_stomate ) THEN 
     278!!$ 
     279!!$             write(numout,*) 'ATT, DGVM!: veget exceeds bareground without vegetation left' 
     280!!$             write(numout,*) 'ATT, DGVM!: is this a bug? cell: ',i 
     281!!$             write(numout,*) 'ATT, DGVM!: veget ',veget(i,:) 
     282!!$ 
     283!!$          ENDIF 
     284 
     285          ! finally, implement the reduction. (reduc is negative!) 
     286          veget(i,:)=veget(i,:)+reduct(:) 
     287 
     288       ENDDO 
     289 
     290    ENDIF 
     291 
    160292    veget(:,ibare_sechiba) = un 
    161293    DO j = 2,nvm 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_crown.f90

    r119 r405  
    66  !--------------------------------------------------------------------- 
    77  !- calculate individual crown area from stem mass. 
     8  !- SZ, I've put the woodmass calculation out of this routine 
     9  !      because after the very first establishment, woodmass 
     10  !      could not be calculated here as veget_max = zero and  
     11  !      d_ind not known... 
    812  !--------------------------------------------------------------------- 
    913  USE ioipsl 
    … …  
    2327  !- 
    2428  SUBROUTINE crown & 
    25        &  (npts, PFTpresent, ind, biomass, veget_max, cn_ind, height) 
     29       &  (npts, PFTpresent, ind, biomass, woodmass_ind, veget_max, cn_ind, height) 
    2630    !--------------------------------------------------------------------- 
    2731    ! 0 declarations 
    … …  
    3741    ! biomass (gC/(m**2 of ground)) 
    3842    REAL(r_std),DIMENSION(npts,nvm,nparts),INTENT(in) :: biomass 
     43    ! woodmass of the individual, needed to calculate crownarea in lpj_crownarea 
     44    REAL(r_std),DIMENSION(npts,nvm),INTENT(in) :: woodmass_ind 
    3945    !- 
    4046    ! 0.2 modified fields 
    … …  
    5864    ! wood mass of an individual 
    5965    !- 
    60     REAL(r_std),DIMENSION(npts) :: woodmass 
     66!!$    REAL(r_std),DIMENSION(npts) :: woodmass 
    6167    !- 
    6268    ! index 
    … …  
    7480    ! 1.1 check if DGVM activated 
    7581    !- 
    76     IF (.NOT.control%ok_dgvm) THEN 
     82    IF (.NOT.control%ok_dgvm .AND. lpj_gap_const_mort) THEN 
    7783       STOP 'crown: not to be called with static vegetation.' 
    7884    ENDIF 
    … …  
    9399          IF (natural(j)) THEN 
    94100             !------ 2.1.1 natural 
    95              WHERE (PFTpresent(:,j) .AND.ind(:,j).GT.min_stomate) 
    96                 !-------- 2.1.1.1 calculate individual wood mass 
    97                 woodmass(:) = & 
    98                      &         (biomass(:,j,isapabove)+biomass(:,j,isapbelow) & 
    99                      &         +biomass(:,j,iheartabove)+biomass(:,j,iheartbelow))/ind(:,j) 
     101             !WHERE (PFTpresent(:,j) .AND.ind(:,j).GT.min_stomate) 
     102             WHERE (PFTpresent(:,j) .AND.woodmass_ind(:,j).GT.min_stomate) 
     103!!$SZ note that woodmass_ind needs to be defined on the individual, hence 
     104!!$ biomass*veget_max/ind, not as stated here, correction MERGE 
     105!!$!-------- 2.1.1.1 calculate individual wood mass 
     106!!$          woodmass(:) = & 
     107!!$ &         (biomass(:,j,isapabove)+biomass(:,j,isapbelow) & 
     108!!$ &         +biomass(:,j,iheartabove)+biomass(:,j,iheartbelow))/ind(:,j) 
    100109                !-------- 2.1.1.2 stem diameter (pipe model) 
    101                 dia(:) = (woodmass(:)/(pipe_density*pi/4.*pipe_tune2)) & 
     110!!$          dia(:) = (woodmass(:)/(pipe_density*pi/4.*pipe_tune2)) & 
     111                dia(:) = (woodmass_ind(:,j)/(pipe_density*pi/4.*pipe_tune2)) & 
    102112                     &                **(1./(2.+pipe_tune3)) 
    103113                !-------- 2.1.1.3 height 
    104114                height(:,j) = pipe_tune2*(dia(:)**pipe_tune3) 
    105                 WHERE (height(:,j) > height_presc_12(j)) 
    106                    dia(:) = (height_presc_12(j)/pipe_tune2)**(1./pipe_tune3) 
    107                    height(:,j) = height_presc_12(j) 
    108                 ENDWHERE 
     115!!$SZ: The constraint on height has nothing to do with LPJ (for that purpose there's dia_max 
     116!!$ cannot see why this is necessary - it also blurrs the output, hence I leave it commented 
     117!!$                WHERE (height(:,j) > height_presc_12(j)) 
     118!!$                   dia(:) = (height_presc_12(j)/pipe_tune2)**(1./pipe_tune3) 
     119!!$                   height(:,j) = height_presc_12(j) 
     120!!$                ENDWHERE 
    109121                !-------- 2.1.1.4 crown area: for large truncs, crown area cannot 
    110122                !--------         exceed a certain value, prescribed through maxdia. 
    … …  
    128140       !       ind and cn_ind are 0 if not present 
    129141       !--- 
    130        !SZ isn't this physically inconsistent with the assumptions of sechiba?? 
    131        ! the actual LPJ formulation calculates fpc from maximum LAI, and fpar from actual LAI=veget 
    132        IF (natural(j).AND.control%ok_dgvm) THEN 
    133           veget_max(:,j) = ind(:,j) * cn_ind(:,j) 
    134        ENDIF 
     142!!$SZ: since now all state variables are defined on veget_max it is very 
     143!!$ dangerous to change this several times in stomate_lpj, as then GPP, turnover and allocated  
     144!!$ biomass are not defined on the same space! Hence, veget_max is now kept constant 
     145!!$ and updated at the end of stomate_lpj in lpj_cover.f90 
     146!!$ Eventually, this routine should only be called once at the beginning and the end of stomate_lpj 
     147!!$ or prefereably cn_ind made a saved state variable! 
     148!!$    IF (natural(j).AND.control%ok_dgvm) THEN 
     149!!$      veget_max(:,j) = ind(:,j) * cn_ind(:,j) 
     150!!$    ENDIF 
    135151    ENDDO 
    136152    !------------------- 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_establish.f90

    r119 r405  
    3333       neighbours, resolution, need_adjacent, herbivores, & 
    3434       precip_annual, gdd0, lm_lastyearmax, & 
    35        cn_ind, lai, avail_tree, avail_grass, & 
     35       cn_ind, lai, avail_tree, avail_grass,  npp_longterm, & 
    3636       leaf_age, leaf_frac, & 
    37        ind, biomass, age, everywhere, co2_to_bm,veget_max) 
    38  
     37       ind, biomass, age, everywhere, co2_to_bm,veget_max, woodmass_ind) 
    3938    ! 
    4039    ! 0 declarations 
    … …  
    7473    ! space availability for grasses 
    7574    REAL(r_std), DIMENSION(npts), INTENT(in)                    :: avail_grass 
     75    ! longterm NPP, for each PFT (gC/(m**2 of ground)) 
     76    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                :: npp_longterm 
    7677    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground  
    7778    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: veget_max 
    … …  
    9495    !NV passage 2D 
    9596    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                 :: co2_to_bm 
     97    ! woodmass of the individual, needed to calculate crownarea in lpj_crownarea 
     98    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                 :: woodmass_ind 
    9699 
    97100    ! 0.3 local 
    … …  
    111114    ! total natural fpc 
    112115    REAL(r_std), DIMENSION(npts)                                :: sumfpc 
     116    ! total fraction occupied by natural vegetation 
     117    REAL(r_std), DIMENSION(npts)                                :: fracnat 
    113118    ! total woody fpc 
    114119    REAL(r_std), DIMENSION(npts)                                :: sumfpc_wood 
    … …  
    129134    ! woodmass of an individual 
    130135    REAL(r_std), DIMENSION(npts)                                :: woodmass 
     136    ! carbon mass in youngest leaf age class (gC/m**2 PFT) 
     137    REAL(r_std), DIMENSION(npts)                                :: leaf_mass_young 
    131138    ! ratio of hw(above) to total hw, sm(above) to total sm 
    132139    REAL(r_std), DIMENSION(npts)                                :: sm_at 
    133140    ! reduction factor for establishment if many trees or grasses are present 
    134141    REAL(r_std), DIMENSION(npts)                                :: factor 
     142    ! Total carbon mass for all pools 
     143    REAL(r_std), DIMENSION(npts)                                :: total_bm_c 
     144    ! Total sappling biomass for all pools 
     145    REAL(r_std), DIMENSION(npts)                                :: total_bm_sapl 
    135146    ! from how many sides is the grid box invaded 
    136147    INTEGER(i_std)                                              :: nfrontx 
    137148    INTEGER(i_std)                                              :: nfronty 
    138149    ! daily establishment rate is large compared to present number of individuals 
    139     LOGICAL, DIMENSION(npts)                                   :: many_new 
     150    !LOGICAL, DIMENSION(npts)                                   :: many_new 
     151    ! flow due to new individuals 
     152    !   veget_max after establishment, to get a proper estimate of carbon and nitrogen  
     153    REAL(r_std), DIMENSION(npts)                                 :: vn 
     154    !   lai on each PFT surface  
     155    REAL(r_std), DIMENSION(npts)                                 :: lai_ind 
     156 
    140157    ! indices 
    141158    INTEGER(i_std)                                              :: i,j,k,m 
    … …  
    161178    ENDIF 
    162179 
    163     ! 
    164     ! 2 recalculate fpc 
    165     ! 
    166  
    167     ! 
    168     ! 2.1 Only natural part of the grid cell 
    169     ! 
    170  
    171     DO j = 2,nvm 
    172  
    173        IF ( natural(j) ) THEN 
    174           DO i = 1, npts 
    175              IF (lai(i,j) == val_exp) THEN                
    176                 fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) 
    177              ELSE 
    178                 fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * ( un - exp( -lai(i,j) * ext_coeff(j) ) ) 
    179              ENDIF 
    180           ENDDO 
    181        ELSE 
    182  
    183           fpc_nat(:,j) = zero 
    184  
    185        ENDIF 
    186  
    187     ENDDO 
    188  
    189     ! 
    190     ! 2.2 total natural fpc on grid 
    191     ! 
    192  
    193     sumfpc(:) = SUM( fpc_nat(:,:), DIM=2 ) 
    194  
    195     ! 
    196     ! 2.3 total woody fpc on grid and number of regenerative tree pfts 
    197     ! 
    198  
    199     sumfpc_wood(:) = zero 
    200     spacefight_tree(:) = zero 
    201  
    202     DO j = 2,nvm 
    203  
    204        IF ( tree(j) .AND. natural(j) ) THEN 
    205  
    206           ! total woody fpc 
    207  
    208           WHERE ( PFTpresent(:,j) ) 
    209              sumfpc_wood(:) = sumfpc_wood(:) + fpc_nat(:,j) 
    210           ENDWHERE 
    211  
    212           ! how many trees are competing? Count a PFT fully only if it is present 
    213           !   on the whole grid box. 
    214  
    215           WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) ) 
    216              spacefight_tree(:) = spacefight_tree(:) + everywhere(:,j) 
    217           ENDWHERE 
    218  
    219        ENDIF 
    220  
    221     ENDDO 
    222  
    223     ! 
    224     ! 2.4 number of natural grasses 
    225     ! 
    226  
    227     spacefight_grass(:) = zero 
    228  
    229     DO j = 2,nvm 
    230  
    231        IF ( .NOT. tree(j) .AND. natural(j) ) THEN 
    232  
    233           ! how many grasses are competing? Count a PFT fully only if it is present 
    234           !   on the whole grid box. 
    235  
    236           WHERE ( PFTpresent(:,j) ) 
    237              spacefight_grass(:) = spacefight_grass(:) + everywhere(:,j) 
    238           ENDWHERE 
    239  
    240        ENDIF 
    241  
    242     ENDDO 
    243  
    244     ! 
    245     ! 3 establishment rate 
    246     ! 
    247  
    248     ! 
    249     ! 3.1 maximum establishment rate, based on climate only 
    250     ! 
    251  
    252     WHERE ( ( precip_annual(:) .GE. precip_crit ) .AND. ( gdd0(:) .GE. gdd_crit ) ) 
    253  
    254        estab_rate_max_climate_tree(:) = estab_max_tree 
    255        estab_rate_max_climate_grass(:) = estab_max_grass 
    256  
    257     ELSEWHERE 
    258  
    259        estab_rate_max_climate_tree(:) = zero 
    260        estab_rate_max_climate_grass(:) = zero 
    261  
    262     ENDWHERE 
    263  
    264     ! 
    265     ! 3.2 reduce maximum tree establishment rate if many trees present. 
    266     !     In the original DGVM, this is done using a step function which yields a 
    267     !     reduction by factor 4 if sumfpc_wood(i) .GT.  fpc_crit - 0.05. 
    268     !     This can lead to small oscillations (without consequences however). 
    269     !     Here, a steady linear transition is used between fpc_crit-0.075 and 
    270     !     fpc_crit-0.025. 
    271     ! 
    272  
    273     factor(:) = un - 15. * ( sumfpc_wood(:) - (fpc_crit-.075) ) 
    274     factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:) ) ) 
    275  
    276     estab_rate_max_tree(:) = estab_rate_max_climate_tree(:) * factor(:) 
    277  
    278     ! 
    279     ! 3.3 Modulate grass establishment rate. 
    280     !     If canopy is not closed (fpc < fpc_crit-0.05), normal establishment. 
    281     !     If canopy is closed, establishment is reduced by a factor 4. 
    282     !     Factor is linear between these two bounds. 
    283     !     This is different from the original DGVM where a step function is 
    284     !     used at fpc_crit-0.05 (This can lead to small oscillations, 
    285     !     without consequences however). 
    286     ! 
    287  
    288     factor(:) = un - 15. * ( sumfpc(:) - (fpc_crit-.05) ) 
    289     factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:) ) ) 
    290  
    291     estab_rate_max_grass(:) = estab_rate_max_climate_grass(:) * factor(:) 
    292  
    293     ! 
    294     ! 4 do establishment for natural PFTs 
    295     ! 
    296  
    297     d_ind(:,:) = zero 
    298  
    299     DO j = 2,nvm 
    300  
    301        ! only for natural PFTs 
    302  
    303        IF ( natural(j) ) THEN 
    304  
    305           ! 
    306           ! 4.1 PFT expansion across the grid box. Not to be confused with areal 
    307           !     coverage. 
    308           ! 
    309  
    310           IF ( treat_expansion ) THEN 
    311  
    312              ! only treat plants that are regenerative and present and still can expand 
    313  
    314              DO i = 1, npts 
    315  
    316                 IF ( PFTpresent(i,j) .AND. & 
    317                      ( everywhere(i,j) .LT. un ) .AND. & 
    318                      ( regenerate(i,j) .GT. regenerate_crit ) ) THEN 
    319  
    320                    ! from how many sides is the grid box invaded (separate x and y directions 
    321                    ! because resolution may be strongly anisotropic) 
    322                    ! 
    323                    ! For the moment we only look into 4 direction but that can be extanded (JP)  
    324                    ! 
    325                    nfrontx = 0 
    326                    IF ( neighbours(i,3) .GT. 0 ) THEN 
    327                       IF ( everywhere(neighbours(i,3),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1 
    328                    ENDIF 
    329                    IF ( neighbours(i,7) .GT. 0 ) THEN 
    330                       IF ( everywhere(neighbours(i,7),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1 
    331                    ENDIF 
    332  
    333                    nfronty = 0 
    334                    IF ( neighbours(i,1) .GT. 0 ) THEN 
    335                       IF ( everywhere(neighbours(i,1),j) .GT. 1.-min_stomate ) nfronty = nfronty+1 
    336                    ENDIF 
    337                    IF ( neighbours(i,5) .GT. 0 ) THEN 
    338                       IF ( everywhere(neighbours(i,5),j) .GT. 1.-min_stomate ) nfronty = nfronty+1 
    339                    ENDIF 
    340  
    341                    everywhere(i,j) = & 
    342                         everywhere(i,j) + migrate(j) * dt/one_year * & 
    343                         ( nfrontx / resolution(i,1) + nfronty / resolution(i,2) ) 
    344  
    345                    IF ( .NOT. need_adjacent(i,j) ) THEN 
    346  
    347                       ! in that case, we also assume that the PFT expands from places within 
    348                       ! the grid box (e.g., oasis). 
     180 
     181    IF (control%ok_dgvm) THEN 
     182       ! 
     183       ! 2 recalculate fpc 
     184       ! 
     185 
     186       ! 
     187       ! 2.1 Only natural part of the grid cell 
     188       ! 
     189 
     190       fracnat(:) = 1. 
     191       do j = 2,nvm 
     192          IF ( .NOT. natural(j) ) THEN 
     193             fracnat(:) = fracnat(:) - veget_max(:,j) 
     194          ENDIF 
     195       ENDDO 
     196 
     197       ! 
     198       ! 2.2 total natural fpc on grid 
     199       ! 
     200       sumfpc(:) = zero 
     201       DO j = 2,nvm 
     202 
     203          IF ( natural(j) ) THEN 
     204             WHERE(fracnat(:).GT.min_stomate) 
     205                WHERE (lai(:,j) == val_exp)  
     206                   fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) 
     207                ELSEWHERE 
     208                   fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) &  
     209                        * ( 1. - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ) 
     210                ENDWHERE 
     211             ENDWHERE 
     212 
     213             WHERE ( PFTpresent(:,j) ) 
     214                sumfpc(:) = sumfpc(:) + fpc_nat(:,j) 
     215             ENDWHERE 
     216          ELSE 
     217 
     218             fpc_nat(:,j) = 0.0 
     219 
     220          ENDIF 
     221 
     222       ENDDO 
     223 
     224       ! 
     225       ! 2.3 total woody fpc on grid and number of regenerative tree pfts 
     226       ! 
     227 
     228       sumfpc_wood(:) = zero 
     229       spacefight_tree(:) = zero 
     230 
     231       DO j = 2,nvm 
     232 
     233          IF ( tree(j) .AND. natural(j) ) THEN 
     234 
     235             ! total woody fpc 
     236 
     237             WHERE ( PFTpresent(:,j) ) 
     238                sumfpc_wood(:) = sumfpc_wood(:) + fpc_nat(:,j) 
     239             ENDWHERE 
     240 
     241             ! how many trees are competing? Count a PFT fully only if it is present 
     242             !   on the whole grid box. 
     243 
     244             WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) ) 
     245                spacefight_tree(:) = spacefight_tree(:) + everywhere(:,j) 
     246             ENDWHERE 
     247 
     248          ENDIF 
     249 
     250       ENDDO 
     251 
     252       ! 
     253       ! 2.4 number of natural grasses 
     254       ! 
     255 
     256       spacefight_grass(:) = zero 
     257 
     258       DO j = 2,nvm 
     259 
     260          IF ( .NOT. tree(j) .AND. natural(j) ) THEN 
     261 
     262             ! how many grasses are competing? Count a PFT fully only if it is present 
     263             !   on the whole grid box. 
     264 
     265             WHERE ( PFTpresent(:,j) ) 
     266                spacefight_grass(:) = spacefight_grass(:) + everywhere(:,j) 
     267             ENDWHERE 
     268 
     269          ENDIF 
     270 
     271       ENDDO 
     272 
     273       ! 
     274       ! 3 establishment rate 
     275       ! 
     276 
     277       ! 
     278       ! 3.1 maximum establishment rate, based on climate only 
     279       ! 
     280 
     281       WHERE ( ( precip_annual(:) .GE. precip_crit ) .AND. ( gdd0(:) .GE. gdd_crit ) ) 
     282 
     283          estab_rate_max_climate_tree(:) = estab_max_tree 
     284          estab_rate_max_climate_grass(:) = estab_max_grass 
     285 
     286       ELSEWHERE 
     287 
     288          estab_rate_max_climate_tree(:) = zero 
     289          estab_rate_max_climate_grass(:) = zero 
     290 
     291       ENDWHERE 
     292 
     293       ! 
     294       ! 3.2 reduce maximum tree establishment rate if many trees present. 
     295       !     In the original DGVM, this is done using a step function which yields a 
     296       !     reduction by factor 4 if sumfpc_wood(i) .GT.  fpc_crit - 0.05. 
     297       !     This can lead to small oscillations (without consequences however). 
     298       !     Here, a steady linear transition is used between fpc_crit-0.075 and 
     299       !     fpc_crit-0.025. 
     300       ! 
     301 
     302       !       factor(:) = 1. - 15. * ( sumfpc_wood(:) - (fpc_crit-.075) ) 
     303       !       factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:) ) ) 
     304 
     305       !SZ modified according to Smith et al. 2001, 080806 
     306       factor(:)=(1.0-exp(-5.0*(1.0-sumfpc_wood(:))))*(1.0-sumfpc_wood(:)) 
     307 
     308       estab_rate_max_tree(:) = estab_rate_max_climate_tree(:) * factor(:) 
     309 
     310       ! 
     311       ! 3.3 Modulate grass establishment rate. 
     312       !     If canopy is not closed (fpc < fpc_crit-0.05), normal establishment. 
     313       !     If canopy is closed, establishment is reduced by a factor 4. 
     314       !     Factor is linear between these two bounds. 
     315       !     This is different from the original DGVM where a step function is 
     316       !     used at fpc_crit-0.05 (This can lead to small oscillations, 
     317       !     without consequences however). 
     318       ! 
     319 
     320       !       factor(:) = 1. - 15. * ( sumfpc(:) - (fpc_crit-.05) ) 
     321       !       factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:) ) ) 
     322       !       estab_rate_max_grass(:) = estab_rate_max_climate_grass(:) * factor(:) 
     323 
     324       !SZ modified to true LPJ formulation, grasses are only allowed in the 
     325       !fpc fraction not occupied by trees..., 080806 
     326!NVmodif       estab_rate_max_grass(:)=MAX(0.98-sumfpc(:),zero) 
     327       estab_rate_max_grass(:)=MAX(MIN(estab_rate_max_climate_grass(:),0.98-sumfpc(:)),zero) 
     328 
     329       ! SZ: longterm grass NPP for competition between C4 and C3 grasses 
     330       !     to avoid equal veget_max, the idea is that more reestablishment 
     331       !     is possible for the more productive PFT 
     332       factor(:)=min_stomate 
     333       DO j = 2,nvm 
     334          IF ( natural(j) .AND. .NOT.tree(j)) &  
     335               factor(:)=factor(:)+npp_longterm(:,j) * & 
     336               lm_lastyearmax(:,j) * sla(j) 
     337       ENDDO 
     338       ! 
     339       ! 
     340       ! 
     341       ! 4 do establishment for natural PFTs 
     342       ! 
     343 
     344       d_ind(:,:) = zero 
     345 
     346       DO j = 2,nvm 
     347 
     348          ! only for natural PFTs 
     349 
     350          IF ( natural(j) ) THEN 
     351 
     352             ! 
     353             ! 4.1 PFT expansion across the grid box. Not to be confused with areal 
     354             !     coverage. 
     355             ! 
     356 
     357             IF ( treat_expansion ) THEN 
     358 
     359                ! only treat plants that are regenerative and present and still can expand 
     360 
     361                DO i = 1, npts 
     362 
     363                   IF ( PFTpresent(i,j) .AND. & 
     364                        ( everywhere(i,j) .LT. un ) .AND. & 
     365                        ( regenerate(i,j) .GT. regenerate_crit ) ) THEN 
     366 
     367                      ! from how many sides is the grid box invaded (separate x and y directions 
     368                      ! because resolution may be strongly anisotropic) 
     369                      ! 
     370                      ! For the moment we only look into 4 direction but that can be extanded (JP)  
     371                      ! 
     372                      nfrontx = 0 
     373                      IF ( neighbours(i,3) .GT. 0 ) THEN 
     374                         IF ( everywhere(neighbours(i,3),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1 
     375                      ENDIF 
     376                      IF ( neighbours(i,7) .GT. 0 ) THEN 
     377                         IF ( everywhere(neighbours(i,7),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1 
     378                      ENDIF 
     379 
     380                      nfronty = 0 
     381                      IF ( neighbours(i,1) .GT. 0 ) THEN 
     382                         IF ( everywhere(neighbours(i,1),j) .GT. 1.-min_stomate ) nfronty = nfronty+1 
     383                      ENDIF 
     384                      IF ( neighbours(i,5) .GT. 0 ) THEN 
     385                         IF ( everywhere(neighbours(i,5),j) .GT. 1.-min_stomate ) nfronty = nfronty+1 
     386                      ENDIF 
    349387 
    350388                      everywhere(i,j) = & 
    351389                           everywhere(i,j) + migrate(j) * dt/one_year * & 
    352                            2. * SQRT( pi*everywhere(i,j)/(resolution(i,1)*resolution(i,2)) ) 
     390                           ( nfrontx / resolution(i,1) + nfronty / resolution(i,2) ) 
     391 
     392                      IF ( .NOT. need_adjacent(i,j) ) THEN 
     393 
     394                         ! in that case, we also assume that the PFT expands from places within 
     395                         ! the grid box (e.g., oasis). 
     396 
     397                         everywhere(i,j) = & 
     398                              everywhere(i,j) + migrate(j) * dt/one_year * & 
     399                              2. * SQRT( pi*everywhere(i,j)/(resolution(i,1)*resolution(i,2)) ) 
     400 
     401                      ENDIF 
     402 
     403                      everywhere(i,j) = MIN( everywhere(i,j), 1._r_std ) 
    353404 
    354405                   ENDIF 
    355406 
    356                    everywhere(i,j) = MIN( everywhere(i,j), 1._r_std ) 
    357  
    358                 ENDIF 
    359  
    360              ENDDO 
    361  
    362           ENDIF  ! treat expansion? 
    363  
    364           ! 
    365           ! 4.2 establishment rate 
    366           !     - Is lower if the PFT is only present in a small part of the grid box 
    367           !       (after its introduction), therefore multiplied by "everywhere". 
    368           !     - Is divided by the number of PFTs that compete ("spacefight"). 
    369           !     - Is modulated by space availability (avail_tree, avail_grass). 
    370           ! 
    371  
    372           IF ( tree(j) ) THEN 
    373  
    374              ! 4.2.1 present and regenerative trees 
    375  
    376              WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) ) 
    377  
    378  
    379                 d_ind(:,j) = estab_rate_max_tree(:)*everywhere(:,j)/spacefight_tree(:) * & 
    380                      avail_tree(:) * dt/one_year 
    381  
    382              ENDWHERE 
    383  
    384           ELSE 
    385  
    386              ! 4.2.2 present and regenerative grasses 
    387  
    388              WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) ) 
    389  
    390                 d_ind(:,j) = estab_rate_max_grass(:)*everywhere(:,j)/spacefight_grass(:) * & 
    391                      avail_grass(:) * dt/one_year 
    392  
    393              ENDWHERE 
    394  
    395           ENDIF  ! tree/grass 
     407                ENDDO 
     408 
     409             ENDIF  ! treat expansion? 
     410 
     411             ! 
     412             ! 4.2 establishment rate 
     413             !     - Is lower if the PFT is only present in a small part of the grid box 
     414             !       (after its introduction), therefore multiplied by "everywhere". 
     415             !     - Is divided by the number of PFTs that compete ("spacefight"). 
     416             !     - Is modulated by space availability (avail_tree, avail_grass). 
     417             ! 
     418 
     419             IF ( tree(j) ) THEN 
     420 
     421                ! 4.2.1 present and regenerative trees 
     422 
     423                WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) ) 
     424 
     425 
     426                   d_ind(:,j) = estab_rate_max_tree(:)*everywhere(:,j)/spacefight_tree(:) * & 
     427                        avail_tree(:) * dt/one_year 
     428 
     429                ENDWHERE 
     430 
     431             ELSE 
     432 
     433                ! 4.2.2 present and regenerative grasses 
     434 
     435                WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit )  &  
     436                     .AND.factor(:).GT.min_stomate .AND. spacefight_grass(:).GT. min_stomate)  
     437 
     438                   d_ind(:,j) = estab_rate_max_grass(:)*everywhere(:,j)/spacefight_grass(:) * & 
     439                        MAX(min_stomate,npp_longterm(:,j)*lm_lastyearmax(:,j)*sla(j)/factor(:)) * fracnat(:) * dt/one_year 
     440 
     441                ENDWHERE 
     442 
     443             ENDIF  ! tree/grass 
     444 
     445          ENDIF ! if natural 
     446       ENDDO ! PFTs 
     447 
     448    ELSE ! lpj establishment in static case, SZ 080806, account for real LPJ dynamics in  
     449       ! prescribed vegetation, i.e. population dynamics within a given area of the  
     450       ! grid cell 
     451 
     452       d_ind(:,:) = 0.0 
     453 
     454       DO j = 2,nvm 
     455 
     456          ! only for natural PFTs 
     457 
     458          WHERE(ind(:,j)*cn_ind(:,j).GT.min_stomate) 
     459             lai_ind(:)=sla(j) * lm_lastyearmax(:,j)/(ind(:,j)*cn_ind(:,j)) 
     460          ELSEWHERE 
     461             lai_ind(:)=0.0 
     462          ENDWHERE 
     463 
     464          IF ( natural(j) .AND. tree(j)) THEN  
     465 
     466             fpc_nat(:,j) =  MIN(1.0,cn_ind(:,j) * ind(:,j) * &  
     467                  MAX( ( 1._r_std - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover ) ) 
     468             !fpc_nat(:,j) = max(fpc_nat(:,j),1.-exp(-0.5*sla(j) * lm_lastyearmax(:,j))) 
     469 
     470 
     471             WHERE (veget_max(:,j).GT.min_stomate.AND.ind(:,j).LE.2.) 
     472 
     473                ! only establish into growing stands, ind can become very 
     474                ! large in the static mode because LAI is very low in poor  
     475                ! growing conditions, favouring continuous establishment. To avoid this 
     476                ! a maximum IND is set. BLARPP: This should be replaced by a  
     477                ! better stand density criteria 
     478                ! 
     479                factor(:)=(1.0-exp(-5.0*(1.0-fpc_nat(:,j))))*(1.0-fpc_nat(:,j)) 
     480 
     481                estab_rate_max_tree(:) = estab_max_tree * factor(:)  
     482                ! 
     483                ! 4 do establishment for natural PFTs 
     484                ! 
     485                d_ind(:,j) = MAX( 0.0, estab_rate_max_tree(:) * dt/one_year) 
     486 
     487             ENDWHERE 
     488 
     489             !SZ: quickfix: to simulate even aged stand, uncomment the following lines... 
     490             !where (ind(:,j) .LE. min_stomate) 
     491             !d_ind(:,j) = 0.1 !MAX( 0.0, estab_rate_max_tree(:) * dt/one_year) 
     492 
     493             WHERE (veget_max(:,j).GT.min_stomate.AND.ind(:,j).EQ.0.0) 
     494                d_ind(:,j) = ind_0*10. 
     495                !          elsewhere 
     496                !d_ind(:,j) =0.0 
     497             endwhere 
     498 
     499          ELSEIF ( natural(j) .AND. .NOT.tree(j)) THEN  
     500 
     501             WHERE (veget_max(:,j).GT.min_stomate) 
     502 
     503                fpc_nat(:,j) =  cn_ind(:,j) * ind(:,j) * &  
     504                     MAX( ( 1._r_std - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover ) 
     505 
     506                d_ind(:,j) = MAX(0.0 , (1.0-fpc_nat(:,j)) * dt/one_year ) 
     507 
     508             ENDWHERE 
     509 
     510             WHERE (veget_max(:,j).GT.min_stomate.AND.ind(:,j).EQ.0.0) 
     511                d_ind(:,j) = ind_0*10. 
     512             ENDWHERE 
     513 
     514          ENDIF 
     515 
     516       ENDDO 
     517 
     518    ENDIF ! DGVM OR NOT 
     519 
     520    DO j = 2,nvm 
     521 
     522       ! only for natural PFTs 
     523 
     524       IF ( natural(j) ) THEN 
    396525 
    397526          ! 
    … …  
    409538          ! 
    410539          ! 4.4 be sure that ind*cn_ind does not exceed 1 
    411           ! 
    412  
    413           WHERE ( ( d_ind(:,j) .GT. zero ) .AND. & 
    414                ( (ind(:,j)+d_ind(:,j))*cn_ind(:,j) .GT. un ) ) 
    415  
    416              d_ind(:,j) = MAX( 1._r_std / cn_ind(:,j) - ind(:,j), 0._r_std ) 
    417  
    418           ENDWHERE 
     540          !SZ This control is now moved to lpj_cover.f90 
     541          !SZ 
     542 
     543          !The aim is to control for sum(veget)=1., irrespective of ind*cnd (crowns can overlap as long as 
     544          ! there is enough light 
     545          ! 
     546          !SZ: This could be part of the dynamic vegetation problem of Orchidee 
     547          !in conjunction with the wrong formulation of establishment response  
     548          !to tree fpc above... 
     549          !          WHERE ( ( d_ind(:,j) .GT. zero ) .AND. & 
     550          !                  ( (ind(:,j)+d_ind(:,j))*cn_ind(:,j) .GT. un ) ) 
     551          ! 
     552          !            d_ind(:,j) = MAX( 1._stnd / cn_ind(:,j) - ind(:,j), zero ) 
     553          ! 
     554          !          ENDWHERE 
    419555 
    420556          ! 
    … …  
    428564 
    429565          ! compare establishment rate and present number of inidivuals 
    430           many_new(:) = ( d_ind(:,j) .GT. 0.1 * ind(:,j) ) 
     566          !many_new(:) = ( d_ind(:,j) .GT. 0.1 * ind(:,j) ) 
    431567 
    432568          ! gives a better vectorization of the VPP 
    433569 
    434           IF ( ANY( many_new(:) ) ) THEN 
    435  
    436              DO k = 1, nparts 
    437  
    438                 WHERE ( many_new(:) ) 
    439  
    440                    bm_new(:) = d_ind(:,j) * bm_sapl(j,k) / veget_max (:,j) 
    441  
    442                    biomass(:,j,k) = biomass(:,j,k) + bm_new(:) 
    443  
    444                    !NV passage 2D 
    445                    co2_to_bm(:,j) = co2_to_bm(:,j) + bm_new(:) / dt 
    446  
    447                 ENDWHERE 
    448  
     570          !IF ( ANY( many_new(:) ) ) THEN 
     571 
     572          ! save old leaf mass to calculate leaf age 
     573          leaf_mass_young(:) = leaf_frac(:,j,1) * biomass(:,j,ileaf) 
     574          ! total biomass of existing PFT to limit biomass added from establishment 
     575          total_bm_c(:) = zero 
     576 
     577          DO k = 1, nparts 
     578             total_bm_c(:)=total_bm_c(:)+biomass(:,j,k) 
     579          ENDDO 
     580          IF(control%ok_dgvm) THEN 
     581             vn(:)=veget_max(:,j) 
     582          ELSE 
     583             vn(:)=1.0 
     584          ENDIF 
     585          total_bm_sapl(:)=zero 
     586          DO k = 1, nparts 
     587             WHERE(d_ind(:,j).GT.min_stomate.AND.vn(:).GT.min_stomate) 
     588 
     589                total_bm_sapl(:) = total_bm_sapl(:) + &  
     590                     bm_sapl(j,k) * d_ind(:,j) / vn(:) 
     591             ENDWHERE 
     592          ENDDO 
     593 
     594          IF(control%ok_dgvm) THEN 
     595             ! SZ calculate new woodmass_ind and veget_max after establishment (needed for correct scaling!) 
     596             ! essential correction for MERGE! 
     597             IF(tree(j))THEN 
     598                DO i=1,npts 
     599                   IF((d_ind(i,j)+ind(i,j)).GT.min_stomate) THEN 
     600 
     601                      IF((total_bm_c(i).LE.min_stomate) .OR. (veget_max(i,j) .LE. min_stomate)) THEN 
     602 
     603                         ! new wood mass of PFT 
     604                         woodmass_ind(i,j) = & 
     605                              & (((biomass(i,j,isapabove)+biomass(i,j,isapbelow) & 
     606                              & +biomass(i,j,iheartabove)+biomass(i,j,iheartbelow))*veget_max(i,j)) & 
     607                              & +(bm_sapl(j,isapabove)+bm_sapl(j,isapbelow) & 
     608                              & +bm_sapl(j,iheartabove)+bm_sapl(j,iheartbelow))*d_ind(i,j))/(ind(i,j)+d_ind(i,j)) 
     609 
     610                      ELSE  
     611                         ! new biomass is added to the labile pool, hence there is no change in CA associated with establishment 
     612                         woodmass_ind(i,j) = & 
     613                              & (biomass(i,j,isapabove)+biomass(i,j,isapbelow) & 
     614                              & +biomass(i,j,iheartabove)+biomass(i,j,iheartbelow))*veget_max(i,j) & 
     615                              & /(ind(i,j)+d_ind(i,j)) 
     616 
     617                      ENDIF 
     618 
     619                      ! new diameter of PFT 
     620                      dia(i) = (woodmass_ind(i,j)/(pipe_density*pi/4.*pipe_tune2)) & 
     621                           &                **(1./(2.+pipe_tune3)) 
     622                      vn(i)=(ind(i,j)+d_ind(i,j))*pipe_tune1*MIN(dia(i),maxdia(j))**1.6 
     623 
     624                   ENDIF 
     625                ENDDO 
     626             ELSE ! for grasses, cnd=1, so the above calculation cancels 
     627                vn(:)=ind(:,j)+d_ind(:,j) 
     628             ENDIF 
     629          ELSE ! static 
     630             DO i=1,npts 
     631                IF(tree(j).AND.(d_ind(i,j)+ind(i,j)).GT.min_stomate) THEN 
     632                   IF(total_bm_c(i).LE.min_stomate) THEN 
     633 
     634                      ! new wood mass of PFT 
     635                      woodmass_ind(i,j) = & 
     636                           & (((biomass(i,j,isapabove)+biomass(i,j,isapbelow) & 
     637                           & +biomass(i,j,iheartabove)+biomass(i,j,iheartbelow))) & 
     638                           & +(bm_sapl(j,isapabove)+bm_sapl(j,isapbelow) & 
     639                           & +bm_sapl(j,iheartabove)+bm_sapl(j,iheartbelow))*d_ind(i,j))/(ind(i,j)+d_ind(i,j)) 
     640 
     641                   ELSE ! new biomass is added to the labile pool, hence there is no change in CA associated with establishment 
     642 
     643                      woodmass_ind(i,j) = & 
     644                           & (biomass(i,j,isapabove)+biomass(i,j,isapbelow) & 
     645                           & +biomass(i,j,iheartabove)+biomass(i,j,iheartbelow)) & 
     646                           & /(ind(i,j)+d_ind(i,j)) 
     647 
     648                   ENDIF 
     649                ENDIF 
    449650             ENDDO 
    450651 
    451              ! reset leaf ages. Should do a real calculation like in the npp routine,  
    452              ! but this case is rare and not worth messing around. 
    453  
    454              WHERE ( many_new(:) ) 
    455                 leaf_age(:,j,1) = zero 
    456                 leaf_frac(:,j,1) = un 
    457              ENDWHERE 
    458  
    459              DO m = 2, nleafages 
    460  
    461                 WHERE ( many_new(:) ) 
    462                    leaf_age(:,j,m) = zero 
    463                    leaf_frac(:,j,m) = zero 
    464                 ENDWHERE 
    465  
    466              ENDDO 
    467  
    468           ENDIF   ! establishment rate is large 
    469  
    470           WHERE ( d_ind(:,j) .GT. zero ) 
    471  
    472              ! 4.5.2 age decreases 
     652             vn(:)=1.0 ! cannot change in static!, and veget_max implicit in d_ind 
     653 
     654          ENDIF 
     655          ! total biomass of PFT added by establishment defined over veget_max ... 
     656          total_bm_sapl(:)=zero 
     657          DO k = 1, nparts 
     658             WHERE(d_ind(:,j).GT.min_stomate.AND.total_bm_c(:).GT.min_stomate.AND.vn(:).GT.min_stomate) 
     659 
     660                total_bm_sapl(:) = total_bm_sapl(:) + &  
     661                     bm_sapl(j,k) * d_ind(:,j) / vn(:) 
     662             ENDWHERE 
     663          ENDDO 
     664 
     665          DO k = 1, nparts 
     666 
     667             bm_new(:)=zero 
     668 
     669             ! first ever establishment, C flows 
     670             WHERE( d_ind(:,j).GT.min_stomate .AND. & 
     671                  total_bm_c(:).LE.min_stomate .AND. & 
     672                  vn(:).GT.min_stomate) 
     673                ! WHERE ( many_new(:) ) 
     674 
     675                !bm_new(:) = d_ind(:,j) * bm_sapl(j,k) / veget_max (:,j) 
     676                bm_new(:) = d_ind(:,j) * bm_sapl(j,k) / vn(:) 
     677 
     678                biomass(:,j,k) = biomass(:,j,k) + bm_new(:) 
     679 
     680                co2_to_bm(:,j) = co2_to_bm(:,j) + bm_new(:) / dt 
     681 
     682             ENDWHERE 
     683 
     684             ! establishment into existing population, C flows 
     685             WHERE(d_ind(:,j).GT.min_stomate.AND.total_bm_c(:).GT.min_stomate) 
     686 
     687                bm_new(:) = total_bm_sapl(:) * biomass(:,j,k) / total_bm_c(:) 
     688 
     689                biomass(:,j,k) = biomass(:,j,k) + bm_new(:) 
     690                co2_to_bm(:,j) = co2_to_bm(:,j) + bm_new(:) / dt 
     691 
     692             ENDWHERE 
     693          ENDDO 
     694 
     695          ! reset leaf ages. Should do a real calculation like in the npp routine,  
     696          ! but this case is rare and not worth messing around. 
     697          ! SZ 080806, added real calculation now, because otherwise leaf_age/leaf_frac 
     698          ! are not initialised for the calculation of vmax, and hence no growth at all. 
     699          ! logic follows that of stomate_npp.f90, just that it's been adjusted for the code here 
     700          ! 
     701          ! 4.5.2 Decrease leaf age in youngest class if new leaf biomass is higher than old one. 
     702          ! 
     703 
     704!!$          WHERE ( many_new(:) ) 
     705!!$             leaf_age(:,j,1) = zero 
     706!!$             leaf_frac(:,j,1) = un 
     707!!$          ENDWHERE 
     708!!$ 
     709!!$          DO m = 2, nleafages 
     710!!$ 
     711!!$             WHERE ( many_new(:) ) 
     712!!$                leaf_age(:,j,m) = zero 
     713!!$                leaf_frac(:,j,m) = zero 
     714!!$             ENDWHERE 
     715!!$ 
     716!!$          ENDDO 
     717 
     718          WHERE ( d_ind(:,j) * bm_sapl(j,ileaf) .GT. min_stomate )  
     719 
     720             leaf_age(:,j,1) = leaf_age(:,j,1) * leaf_mass_young(:) / & 
     721                  ( leaf_mass_young(:) + d_ind(:,j) * bm_sapl(j,ileaf) ) 
     722 
     723          ENDWHERE 
     724 
     725          leaf_mass_young(:) = leaf_mass_young(:) + d_ind(:,j) * bm_sapl(j,ileaf) 
     726 
     727          ! 
     728          ! new age class fractions (fraction in youngest class increases) 
     729          ! 
     730 
     731          ! youngest class: new mass in youngest class divided by total new mass 
     732 
     733          WHERE ( biomass(:,j,ileaf) .GT. min_stomate ) 
     734 
     735             leaf_frac(:,j,1) = leaf_mass_young(:) / biomass(:,j,ileaf) 
     736 
     737          ENDWHERE 
     738 
     739          ! other classes: old mass in leaf age class divided by new mass 
     740 
     741          DO m = 2, nleafages 
     742 
     743             WHERE ( biomass(:,j,ileaf) .GT. min_stomate ) 
     744 
     745                leaf_frac(:,j,m) = leaf_frac(:,j,m) * &  
     746                     ( biomass(:,j,ileaf) + d_ind(:,j) * bm_sapl(j,ileaf) ) /  biomass(:,j,ileaf) 
     747 
     748             ENDWHERE 
     749 
     750          ENDDO 
     751 
     752          !ENDIF   ! establishment rate is large 
     753 
     754          WHERE ( d_ind(:,j) .GT. min_stomate ) 
     755 
     756             ! 4.5.3 age decreases 
    473757 
    474758             age(:,j) = age(:,j) * ind(:,j) / ( ind(:,j) + d_ind(:,j) ) 
    475759 
    476              ! 4.5.3 new number of individuals 
     760             ! 4.5.4 new number of individuals 
    477761 
    478762             ind(:,j) = ind(:,j) + d_ind(:,j) 
    … …  
    484768          ! 
    485769 
     770          !SZ to clarify with Gerhard Krinner: This is theoretically inconsistent because  
     771          ! the allocation to sapwood and leaves do not follow the LPJ logic in stomate_alloc.f90 
     772          ! hence imposing this here not only solves for the uneveness of age (mixing new and average individual) 
     773          ! but also corrects for the discrepancy between SLAVE and LPJ logic of allocation, thus leads to excess heartwood 
     774          ! and thus carbon accumulation! 
     775          ! should be removed. 
     776 
    486777          IF ( tree(j) ) THEN 
    487778 
    488              sm2(:) = zero 
    489  
    490              WHERE ( d_ind(:,j) .GT. zero )  
    491  
    492                 ! ratio of above / total sap parts 
    493                 sm_at(:) = biomass(:,j,isapabove) / & 
    494                      ( biomass(:,j,isapabove) + biomass(:,j,isapbelow) ) 
    495  
    496                 ! woodmass of an individual 
    497  
    498                 woodmass(:) = & 
    499                      ( biomass(:,j,isapabove) + biomass(:,j,isapbelow) + & 
    500                      biomass(:,j,iheartabove) + biomass(:,j,iheartbelow) ) / ind(:,j) 
    501  
    502                 ! crown area (m**2) depends on stem diameter (pipe model) 
    503                 dia(:) = ( woodmass(:) / ( pipe_density * pi/4. * pipe_tune2 ) ) & 
    504                      ** ( un / ( 2. + pipe_tune3 ) ) 
    505  
    506                 b1(:) = pipe_k1 / ( sla(j) * pipe_density*pipe_tune2 * dia(:)**pipe_tune3 ) * & 
    507                      ind(:,j) 
    508                 sm2(:) = lm_lastyearmax(:,j) / b1(:) 
    509  
    510              ENDWHERE 
    511  
    512              WHERE ( ( d_ind(:,j) .GT. zero ) .AND. & 
     779!!$             sm2(:) = 0.0 
     780!!$             WHERE ( d_ind(:,j) .GT. 0.0 )  
     781!!$ 
     782!!$                ! ratio of above / total sap parts 
     783!!$                sm_at(:) = biomass(:,j,isapabove) / & 
     784!!$                     ( biomass(:,j,isapabove) + biomass(:,j,isapbelow) ) 
     785!!$ 
     786!!$                ! woodmass of an individual 
     787!!$ 
     788!!$                woodmass(:) = & 
     789!!$                     ( biomass(:,j,isapabove) + biomass(:,j,isapbelow) + & 
     790!!$                     biomass(:,j,iheartabove) + biomass(:,j,iheartbelow) ) / ind(:,j) 
     791!!$ 
     792!!$                ! crown area (m**2) depends on stem diameter (pipe model) 
     793!!$                dia(:) = ( woodmass(:) / ( pipe_density * pi/4. * pipe_tune2 ) ) & 
     794!!$                     ** ( 1. / ( 2. + pipe_tune3 ) ) 
     795!!$ 
     796!!$                b1(:) = pipe_k1 / ( sla(j) * pipe_density*pipe_tune2 * dia(:)**pipe_tune3 ) * & 
     797!!$                     ind(:,j) 
     798!!$                sm2(:) = lm_lastyearmax(:,j) / b1(:) 
     799!!$ 
     800!!$             ENDWHERE 
     801 
     802             sm2(:)=biomass(:,j,isapabove) + biomass(:,j,isapbelow) 
     803 
     804             WHERE ( ( d_ind(:,j) .GT. min_stomate ) .AND. & 
    513805                  ( biomass(:,j,isapabove) + biomass(:,j,isapbelow) ) .GT. sm2(:) ) 
    514806 
    … …  
    536828 
    537829    CALL histwrite (hist_id_stomate, 'IND_ESTAB', itime, d_ind, npts*nvm, horipft_index) 
     830    CALL histwrite (hist_id_stomate, 'ESTABTREE', itime, estab_rate_max_tree, npts, hori_index) 
     831    CALL histwrite (hist_id_stomate, 'ESTABGRASS', itime, estab_rate_max_grass, npts, hori_index) 
    538832 
    539833    IF (bavard.GE.4) WRITE(numout,*) 'Leaving establish' 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_fire.f90

    r119 r405  
    421421       !       individuals. 
    422422 
    423        IF ( control%ok_dgvm .AND. tree(j) ) THEN 
     423       IF ( (control%ok_dgvm .OR. .NOT.lpj_gap_const_mort) .AND. tree(j) ) THEN 
    424424 
    425425          ! fraction of plants that dies each day. 
    … …  
    472472       !       into CO2) 
    473473 
    474        residue(:) = litter(:,istructural,j,iabove) * firefrac(:,j) * & 
    475             struc_residual(:) 
    476        !MM in SZ        residue(:) = firefrac(:,j) * struc_residual(:) 
     474!NV,MM : We add this test to keep coherence with CMIP5 computations without DGVM. 
     475!        It has to be removed in trunk version after CMIP5. 
     476       IF (control%ok_dgvm .OR. .NOT.lpj_gap_const_mort) THEN 
     477          residue(:) = firefrac(:,j) * struc_residual(:) 
     478       ELSE 
     479          residue(:) = litter(:,istructural,j,iabove) * firefrac(:,j) * & 
     480               struc_residual(:) 
     481       ENDIF 
    477482 
    478483       ! 5.2.4 determine fraction of black carbon in the residue. 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_gap.f90

    r119 r405  
    3838  SUBROUTINE gap (npts, dt, & 
    3939       npp_longterm, turnover_longterm, lm_lastyearmax, & 
    40        PFTpresent, biomass, ind, bm_to_litter) 
     40       PFTpresent, biomass, ind, bm_to_litter, mortality) 
    4141 
    4242    ! 
    … …  
    6767    ! biomass taken away (gC/(m**2 of ground)) 
    6868    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)  :: bm_to_litter 
     69    ! mortality (fraction of trees that is dying per time step), per day in history file 
     70    REAL(r_std), DIMENSION(npts,nvm),INTENT(out)             :: mortality 
    6971 
    7072    ! 0.3 local 
    7173 
    72     ! which kind of mortality 
    73     LOGICAL, SAVE                                           :: constant_mortality 
    7474    ! biomass increase 
    7575    REAL(r_std), DIMENSION(npts)                             :: delta_biomass 
     76    ! biomass increase 
     77    REAL(r_std), DIMENSION(npts)                             :: dmortality 
    7678    ! vigour 
    7779    REAL(r_std), DIMENSION(npts)                             :: vigour 
    7880    ! natural availability, based on vigour 
    7981    REAL(r_std), DIMENSION(npts)                             :: availability 
    80     ! mortality (fraction of trees that is dying per time step), per day in history file 
    81     REAL(r_std), DIMENSION(npts,nvm)                        :: mortality 
    8282    ! indices 
    83     INTEGER(i_std)                                           :: j,k 
     83    INTEGER(i_std)                                           :: j,k,m 
     84    REAL(r_std) :: ref_greff 
    8485 
    8586    ! ========================================================================= 
    … …  
    8990       firstcall = .FALSE. 
    9091 
    91        !Config  Key  = LPJ_GAP_CONST_MORT 
    92        !Config  Desc = constant tree mortality 
    93        !Config  Def  = y 
    94        !Config  Help = If yes, then a constant mortality is applied to trees.  
    95        !Config         Otherwise, mortality is a function of the trees'  
    96        !Config         vigour (as in LPJ). 
    97  
    98        constant_mortality = .TRUE. 
    99        CALL getin_p('LPJ_GAP_CONST_MORT', constant_mortality)      
    100        WRITE(numout,*) 'gap: constant mortality:', constant_mortality 
    101  
    10292    ENDIF 
    10393 
    104     IF (bavard.GE.3) WRITE(numout,*) 'Entering gap' 
     94    IF (bavard.GE.3) WRITE(numout,*) 'Entering gap',lpj_gap_const_mort 
    10595 
    10696    mortality(:,:) = zero 
    10797 
     98    ref_greff =  0.035 
     99 
    108100    DO j = 2,nvm 
    109101 
    … …  
    116108          ! 
    117109 
    118           IF ( .NOT. constant_mortality ) THEN 
     110          IF ( .NOT.  lpj_gap_const_mort ) THEN 
    119111 
    120112             ! 
    … …  
    124116             WHERE ( PFTpresent(:,j) .AND. ( lm_lastyearmax(:,j) .GT. min_stomate ) ) 
    125117 
     118!SZ 080806, changed to LPJ formulation according to Smith et al., 2001  
     119 
    126120                ! how much did the tree grow per year? 
    127121 
    128                 delta_biomass(:) = & 
    129                      MAX( npp_longterm(:,j) - ( turnover_longterm(:,j,ileaf) + & 
    130                      turnover_longterm(:,j,iroot) + turnover_longterm(:,j,ifruit) ), & 
    131                      0._r_std ) 
     122!!$                delta_biomass(:) = & 
     123!!$                     MAX( npp_longterm(:,j) - ( turnover_longterm(:,j,ileaf) + & 
     124!!$                     turnover_longterm(:,j,iroot) + turnover_longterm(:,j,ifruit) ), & 
     125!!$                     0._r_std ) 
     126 
     127            ! note that npp_longterm is now actually longterm growth efficiency (NPP/LAI) 
     128            ! to be fair to deciduous trees 
     129             delta_biomass(:) = MAX( npp_longterm(:,j) - ( turnover_longterm(:,j,ileaf) + & 
     130                  turnover_longterm(:,j,iroot) + turnover_longterm(:,j,ifruit) + &  
     131                  turnover_longterm(:,j,isapabove) + turnover_longterm(:,j,isapbelow) ) ,zero) 
    132132 
    133133                ! scale this to the leaf surface of the tree 
    134  
    135                 vigour(:) = delta_biomass(:) / (lm_lastyearmax(:,j)*sla(j)) / 70. 
     134!!$                vigour(:) = delta_biomass(:) / (lm_lastyearmax(:,j)*sla(j)) / 70. 
     135             vigour(:) = delta_biomass(:) / (lm_lastyearmax(:,j)*sla(j)) 
    136136 
    137137             ELSEWHERE 
    … …  
    146146                ! low vigour. 
    147147 
    148                 availability(:) = 0.02 / ( 1.+vigour(:)/0.17 ) 
     148!SZ 080806, changed to LPJ formulation according to Smith et al., 2001  
     149! tuned maximal mortality to 0.05 to get realistic range of avergage age to get ~100 years at GREFF=100 
     150! for the range of modelled annual NPP 
     151!!$                availability(:) = min_avail / ( 1.+vigour(:)/0.17 ) 
     152                availability(:) = 0.1 / ( 1.+ref_greff*vigour(:) ) 
    149153 
    150154                ! Mortality (fraction per time step). 
    … …  
    157161                ! approximation ok as availability < 0.02 << 1 
    158162 
    159                 mortality(:,j) = availability(:) * dt/one_year 
    160  
     163                mortality(:,j) = MAX(min_avail,availability(:))  * dt/one_year   
     164!!$                mortality(:,j) = availability(:) * dt/one_year 
     165                 
    161166             ENDWHERE 
    162167 
    … …  
    198203             WHERE ( PFTpresent(:,j) ) 
    199204 
    200                 bm_to_litter(:,j,k) = bm_to_litter(:,j,k) + mortality(:,j) * biomass(:,j,k) 
    201  
    202                 biomass(:,j,k) = biomass(:,j,k) * ( un - mortality(:,j) ) 
     205                dmortality(:) =  mortality(:,j) * biomass(:,j,k) 
     206                bm_to_litter(:,j,k) = bm_to_litter(:,j,k) + dmortality(:) 
     207                 
     208                biomass(:,j,k) = biomass(:,j,k) - dmortality(:) 
    203209 
    204210             ENDWHERE 
    … …  
    210216          ! 
    211217 
    212           IF ( control%ok_dgvm ) THEN 
     218!SZ 080806, allow changing density in static case when mortality is dynamic 
     219          IF ( control%ok_dgvm .OR. .NOT.lpj_gap_const_mort) THEN 
    213220 
    214221             WHERE ( PFTpresent(:,j) ) 
    … …  
    219226 
    220227          ENDIF 
    221  
     228       ELSE  
     229 
     230          IF ( .NOT.control%ok_dgvm .AND. .NOT.lpj_gap_const_mort) THEN 
     231 
     232             WHERE ( PFTpresent(:,j) .AND. ( npp_longterm(:,j) .LE. 10. ) ) 
     233 
     234                mortality(:,j) = 1. 
     235 
     236             ENDWHERE 
     237             DO k = 1, nparts 
     238 
     239                WHERE ( PFTpresent(:,j) ) 
     240 
     241                   dmortality(:) =  mortality(:,j) * biomass(:,j,k) 
     242                    
     243                   bm_to_litter(:,j,k) = bm_to_litter(:,j,k) + dmortality(:) 
     244                    
     245                   biomass(:,j,k) = biomass(:,j,k) - dmortality(:) 
     246 
     247                ENDWHERE 
     248             ENDDO 
     249              
     250          ENDIF 
     251           
    222252       ENDIF       ! only trees 
    223253 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_kill.f90

    r119 r405  
    2424  SUBROUTINE kill (npts, whichroutine, lm_lastyearmax, & 
    2525       ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
    26        lai, age, leaf_age, leaf_frac, & 
     26       lai, age, leaf_age, leaf_frac, npp_longterm, & 
    2727       when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    2828 
    … …  
    7171    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground 
    7272    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: veget_max 
     73    ! "long term" net primary productivity (gC/(m**2 of ground)/year) 
     74    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: npp_longterm  
    7375    ! conversion of biomass to litter (gC/(m**2 of ground)) / day 
    7476    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)    :: bm_to_litter 
    … …  
    9799          ! the "was_killed" business is necessary for a more efficient code on the VPP 
    98100 
    99           WHERE ( PFTpresent(:,j) .AND. & 
    100                ( ( ind(:,j) .LT. min_stomate ) .OR. & 
    101                ( lm_lastyearmax(:,j) .LT. min_stomate ) ) ) 
    102  
    103              was_killed(:) = .TRUE. 
    104  
    105           ENDWHERE 
     101          IF ( control%ok_dgvm ) THEN 
     102             WHERE ( PFTpresent(:,j) .AND. & 
     103                  ( ( ind(:,j) .LT. min_stomate ) .OR. & 
     104                  ( lm_lastyearmax(:,j) .LT. min_stomate ) ) ) 
     105 
     106                was_killed(:) = .TRUE. 
     107 
     108             ENDWHERE 
     109 
     110          ELSE 
     111             WHERE ( PFTpresent(:,j) .AND. &  
     112                  (biomass(:,j,icarbres) .LE.zero .OR. &  
     113                  biomass(:,j,iroot).LT.-min_stomate .OR. biomass(:,j,ileaf).LT.-min_stomate ).AND. &  
     114                  ind(:,j).GT. zero) 
     115 
     116                was_killed(:) = .TRUE. 
     117 
     118             ENDWHERE 
     119 
     120             IF(.NOT.tree(j).AND..NOT.lpj_gap_const_mort)THEN 
     121                WHERE ( was_killed(:) ) 
     122 
     123                   npp_longterm(:,j)=500. 
     124 
     125                ENDWHERE 
     126             ENDIF 
     127 
     128          ENDIF 
    106129 
    107130          IF ( ANY( was_killed(:) ) ) THEN 
    108131 
    109132             WHERE ( was_killed(:) ) 
    110  
    111                 ind(:,j) = zero 
    112133 
    113134                bm_to_litter(:,j,ileaf) = bm_to_litter(:,j,ileaf) + biomass(:,j,ileaf) 
    … …  
    131152                biomass(:,j,icarbres) = zero 
    132153 
    133                 PFTpresent(:,j) = .FALSE. 
     154             ENDWHERE   ! number of individuals very low 
     155 
     156             IF (control%ok_dgvm) THEN 
     157 
     158                WHERE ( was_killed(:) ) 
     159                   PFTpresent(:,j) = .FALSE. 
     160 
     161                   veget_max(:,j) = zero 
     162                    
     163                   RIP_time(:,j) = zero 
     164 
     165                ENDWHERE  ! number of individuals very low 
     166 
     167             ENDIF 
     168 
     169             WHERE ( was_killed(:) ) 
     170 
     171                ind(:,j) = zero 
    134172 
    135173                cn_ind(:,j) = zero 
    … …  
    140178                age(:,j) = zero 
    141179 
    142                 when_growthinit(:,j) = undef 
     180                ! SZ: why undef ??? this causes a delay in reestablishment 
     181                !when_growthinit(:,j) = undef 
     182                when_growthinit(:,j) = large_value  
    143183 
    144184                everywhere(:,j) = zero 
    145185 
    146186                veget(:,j) = zero 
    147  
    148                 veget_max(:,j) = zero 
    149  
    150                 RIP_time(:,j) = zero 
    151187 
    152188             ENDWHERE   ! number of individuals very low 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_light.f90

    r119 r405  
    1414! Exclude agricultural pfts from competition 
    1515! 
     16! SZ: added light competition for the static case if the mortality is not  
     17!     assumed to be constant. 
     18! other modifs: 
     19! -1      FPC is now always calculated from lm_lastyearmax*sla, since the aim of this DGVM is  
     20!         to represent community ecology effects; seasonal variations in establishment related to phenology 
     21!         may be relevant, but beyond the scope of a 1st generation DGVM  
     22! -2      problem, if agriculture is present, fpc can never reach 1.0 since natural veget_max < 1.0. To 
     23!         correct for this, ind must be recalculated to correspond to the natural density... 
     24!         since ind is 1/m2 grid cell, this can be achived by dividing ind by the agricultural fraction 
     25 
     26! 
    1627! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/lpj_light.f90,v 1.8 2009/01/06 15:01:25 ssipsl Exp $ 
    1728! IPSL (2006) 
    … …  
    4354 
    4455  SUBROUTINE light (npts, dt, & 
    45        PFTpresent, cn_ind, lai, maxfpc_lastyear, & 
    46        ind, biomass, veget_lastlight, bm_to_litter) 
     56       veget_max, fpc_max, PFTpresent, cn_ind, lai, maxfpc_lastyear, & 
     57       lm_lastyearmax, ind, biomass, veget_lastlight, bm_to_litter, mortality) 
    4758 
    4859    ! 
    … …  
    6475    ! last year's maximum fpc for each natural PFT, on ground 
    6576    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)             :: maxfpc_lastyear 
     77    ! last year's maximum leafmass for each natural PFT, on ground 
     78    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)             :: lm_lastyearmax 
     79    ! last year's maximum fpc for each natural PFT, on ground 
     80    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)             :: veget_max 
     81    ! last year's maximum fpc for each natural PFT, on ground 
     82    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)             :: fpc_max 
    6683 
    6784    ! 0.2 modified fields 
    … …  
    7592    ! biomass taken away (gC/m**2) 
    7693    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)   :: bm_to_litter 
     94    ! fraction of individuals that died this time step 
     95    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)          :: mortality 
    7796 
    7897    ! 0.3 local 
    … …  
    86105    LOGICAL, PARAMETER                                       :: annual_increase = .TRUE. 
    87106    ! index 
    88     INTEGER(i_std)                                            :: i,j 
     107    INTEGER(i_std)                                            :: i,j,k,m 
    89108    ! total natural fpc 
    90109    REAL(r_std), DIMENSION(npts)                              :: sumfpc 
     110    ! fraction of natural vegetation at grid cell level 
     111    REAL(r_std), DIMENSION(npts)                              :: fracnat 
    91112    ! total natural woody fpc 
    92113    REAL(r_std)                                               :: sumfpc_wood 
    … …  
    107128    ! Fraction of plants that survive 
    108129    REAL(r_std), DIMENSION(nvm)                              :: survive 
     130    ! FPC for static mode 
     131    REAL(r_std), DIMENSION(npts)                              :: fpc_real 
     132    ! FPC mortality for static mode 
     133    REAL(r_std), DIMENSION(npts)                              :: lai_ind 
    109134    ! number of grass PFTs present in the grid box 
    110     INTEGER(i_std)                                            :: num_grass 
     135    !    INTEGER(i_std)                                            :: num_grass 
    111136    ! New total grass fpc 
    112137    REAL(r_std)                                               :: sumfpc_grass2 
    113138    ! fraction of plants that dies each day (1/day) 
    114139    REAL(r_std), DIMENSION(npts,nvm)                         :: light_death 
     140    ! Relative change of number of individuals for trees 
     141    REAL(r_std)                                               :: fpc_dec 
    115142 
    116143    ! ========================================================================= 
    … …  
    146173    ENDIF 
    147174 
    148     ! 
    149     ! 2 fpc characteristics 
    150     ! 
    151  
    152     ! 
    153     ! 2.1 calculate fpc on natural part of grid cell. 
    154     ! 
    155  
    156     DO j = 2, nvm 
    157  
    158        IF ( natural(j) ) THEN 
    159  
    160           ! 2.1.1 natural PFTs 
    161  
    162           IF ( tree(j) ) THEN 
    163  
    164              ! 2.1.1.1 trees: minimum cover due to stems, branches etc. 
    165  
    166              DO i = 1, npts 
    167                 IF (lai(i,j) == val_exp) THEN 
    168                    fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) 
    169                 ELSE 
    170                    fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & 
    171                         MAX( ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) 
    172                 ENDIF 
    173              ENDDO 
     175    IF (control%ok_dgvm) THEN 
     176       ! 
     177       ! 2 fpc characteristics 
     178       ! 
     179 
     180       ! 2.0 Only natural part of the grid cell: 
     181       ! calculate fraction of natural and agricultural (1-fracnat) surface 
     182 
     183       fracnat(:) = 1. 
     184       DO j = 2,nvm 
     185          IF ( .NOT. natural(j) ) THEN 
     186             fracnat(:) = fracnat(:) - veget_max(:,j) 
     187          ENDIF 
     188       ENDDO 
     189       ! 
     190       ! 2.1 calculate fpc on natural part of grid cell. 
     191       ! 
     192       fpc_nat(:,:)=zero 
     193       fpc_nat(:,ibare_sechiba)=un 
     194 
     195       DO j = 2, nvm 
     196 
     197          IF ( natural(j) ) THEN 
     198 
     199             ! 2.1.1 natural PFTs 
     200 
     201             IF ( tree(j) ) THEN 
     202 
     203                ! 2.1.1.1 trees: minimum cover due to stems, branches etc. 
     204 
     205                !          DO i = 1, npts 
     206                !             IF (lai(i,j) == val_exp) THEN 
     207                !                fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) 
     208                !             ELSE 
     209                !                fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & 
     210                !                     MAX( ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) 
     211                !             ENDIF 
     212                !          ENDDO 
     213 
     214                !NV : modif from SZ version : fpc is based on veget_max, not veget. 
     215                WHERE(fracnat(:).GE.min_stomate) 
     216                   !            WHERE(LAI(:,j) == val_exp) 
     217                   !               fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) 
     218                   !            ELSEWHERE 
     219                   !               fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * & 
     220                   !                    MAX( ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ), min_cover ) 
     221                   !            ENDWHERE 
     222                   fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) 
     223                ENDWHERE 
     224 
     225             ELSE 
     226 
     227                !NV : modif from SZ version : fpc is based on veget_max, not veget. 
     228                WHERE(fracnat(:).GE.min_stomate) 
     229                   !            WHERE(LAI(:,j) == val_exp) 
     230                   !               fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) 
     231                   !            ELSEWHERE 
     232                   !               fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * & 
     233                   !                    ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ) 
     234                   !            ENDWHERE 
     235                   fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) 
     236                ENDWHERE 
     237 
     238!!$                ! 2.1.1.2 bare ground  
     239!!$                IF (j == ibare_sechiba) THEN 
     240!!$                   fpc_nat(:,j) = cn_ind(:,j) * ind(:,j)  
     241!!$ 
     242!!$                   ! 2.1.1.3 grasses 
     243!!$                ELSE 
     244!!$                   DO i = 1, npts 
     245!!$                      IF (lai(i,j) == val_exp) THEN 
     246!!$                         fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) 
     247!!$                      ELSE 
     248!!$                         fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & 
     249!!$                              ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ) 
     250!!$                      ENDIF 
     251!!$                   ENDDO 
     252!!$                ENDIF 
     253 
     254             ENDIF  ! tree/grass 
    174255 
    175256          ELSE 
    176257 
    177              ! 2.1.1.2 bare ground  
    178              IF (j == ibare_sechiba) THEN 
    179                 fpc_nat(:,j) = cn_ind(:,j) * ind(:,j)  
    180  
    181                 ! 2.1.1.3 grasses 
     258             ! 2.1.2 agricultural PFTs: not present on natural part 
     259 
     260             fpc_nat(:,j) = zero 
     261 
     262          ENDIF    ! natural/agricultural 
     263 
     264       ENDDO 
     265 
     266       ! 
     267       ! 2.2 sum natural fpc for every grid point 
     268       ! 
     269 
     270       sumfpc(:) = zero 
     271       DO j = 2,nvm 
     272          !SZ bug correction MERGE: need to subtract agricultural area! 
     273          sumfpc(:) = sumfpc(:) + fpc_nat(:,j) 
     274       ENDDO 
     275 
     276       ! 
     277       ! 3 Light competition 
     278       ! 
     279 
     280       light_death(:,:) = zero 
     281 
     282       DO i = 1, npts ! SZ why this loop and not a vector statement ? 
     283 
     284          ! Only if vegetation cover is dense 
     285 
     286          IF ( sumfpc(i) .GT. fpc_crit ) THEN 
     287 
     288             ! fpc change for each pft 
     289             ! There are two possibilities: either we compare today's fpc with the fpc after the last 
     290             ! time step, or we compare it to last year's maximum fpc of that PFT. In the first case, 
     291             ! the fpc increase will be strong for seasonal PFTs at the beginning of the growing season. 
     292             ! As for trees, the cutback is proportional to this increase, this means that seasonal trees 
     293             ! will be disadvantaged compared to evergreen trees. In the original LPJ model, with its  
     294             ! annual time step, the second method was used (this corresponds to annual_increase=.TRUE.) 
     295 
     296             IF ( annual_increase ) THEN 
     297                deltafpc(:) = MAX( (fpc_nat(i,:)-maxfpc_lastyear(i,:)), 0._r_std ) 
    182298             ELSE 
    183                 DO i = 1, npts 
    184                    IF (lai(i,j) == val_exp) THEN 
    185                       fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) 
     299                deltafpc(:) = MAX( (fpc_nat(i,:)-veget_lastlight(i,:)), 0._r_std ) 
     300             ENDIF 
     301 
     302             ! default: survive 
     303 
     304             survive(:) = 1.0 
     305 
     306             ! 
     307             ! 3.1 determine some characteristics of the fpc distribution 
     308             ! 
     309 
     310             sumfpc_wood = zero 
     311             sumdelta_fpc_wood = zero 
     312             maxfpc_wood = zero 
     313             optpft_wood = 0 
     314             sumfpc_grass = zero 
     315             !        num_grass = 0 
     316 
     317             DO j = 2,nvm 
     318 
     319                ! only natural pfts 
     320 
     321                IF ( natural(j) ) THEN 
     322 
     323                   IF ( tree(j) ) THEN 
     324 
     325                      ! trees 
     326 
     327                      ! total woody fpc 
     328 
     329                      sumfpc_wood = sumfpc_wood + fpc_nat(i,j) 
     330 
     331                      ! how much did the woody fpc increase 
     332 
     333                      sumdelta_fpc_wood = sumdelta_fpc_wood + deltafpc(j) 
     334 
     335                      ! which woody pft is preponderant 
     336 
     337                      IF ( fpc_nat(i,j) .GT. maxfpc_wood ) THEN 
     338 
     339                         optpft_wood = j 
     340 
     341                         maxfpc_wood = fpc_nat(i,j) 
     342 
     343                      ENDIF 
     344 
    186345                   ELSE 
    187                       fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & 
    188                            ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ) 
    189                    ENDIF 
    190                 ENDDO 
    191              ENDIF 
    192           ENDIF  ! tree/grass 
    193  
    194        ELSE 
    195  
    196           ! 2.1.2 agricultural PFTs: not present on natural part 
    197  
    198           fpc_nat(:,j) = zero 
    199  
    200        ENDIF    ! natural/agricultural 
    201  
    202     ENDDO 
    203  
    204     ! 
    205     ! 2.2 sum natural fpc for every grid point 
    206     ! 
    207  
    208     sumfpc(:) = zero 
    209     DO j = 2,nvm 
    210        !SZ bug correction MERGE: need to subtract agricultural area! 
    211        sumfpc(:) = sumfpc(:) + fpc_nat(:,j) 
    212     ENDDO 
    213  
    214     ! 
    215     ! 3 Light competition 
    216     ! 
    217  
    218     light_death(:,:) = zero 
    219  
    220     DO i = 1, npts ! SZ why this loop and not a vector statement ? 
    221  
    222        ! Only if vegetation cover is dense 
    223  
    224        IF ( sumfpc(i) .GT. fpc_crit ) THEN 
    225  
    226           ! fpc change for each pft 
    227           ! There are two possibilities: either we compare today's fpc with the fpc after the last 
    228           ! time step, or we compare it to last year's maximum fpc of that PFT. In the first case, 
    229           ! the fpc increase will be strong for seasonal PFTs at the beginning of the growing season. 
    230           ! As for trees, the cutback is proportional to this increase, this means that seasonal trees 
    231           ! will be disadvantaged compared to evergreen trees. In the original LPJ model, with its  
    232           ! annual time step, the second method was used (this corresponds to annual_increase=.TRUE.) 
    233  
    234           IF ( annual_increase ) THEN 
    235              deltafpc(:) = MAX( (fpc_nat(i,:)-maxfpc_lastyear(i,:)), 0._r_std ) 
    236           ELSE 
    237              deltafpc(:) = MAX( (fpc_nat(i,:)-veget_lastlight(i,:)), 0._r_std ) 
    238           ENDIF 
    239  
    240           ! default: survive 
    241  
    242           survive(:) = 1.0 
    243  
    244           ! 
    245           ! 3.1 determine some characteristics of the fpc distribution 
    246           ! 
    247  
    248           sumfpc_wood = zero 
    249           sumdelta_fpc_wood = zero 
    250           maxfpc_wood = zero 
    251           optpft_wood = 0 
    252           sumfpc_grass = zero 
    253           num_grass = 0 
    254  
    255           DO j = 2,nvm 
    256  
    257              ! only natural pfts 
    258  
    259              IF ( natural(j) ) THEN 
    260  
    261                 IF ( tree(j) ) THEN 
    262  
    263                    ! trees 
    264  
    265                    ! total woody fpc 
    266  
    267                    sumfpc_wood = sumfpc_wood + fpc_nat(i,j) 
    268  
    269                    ! how much did the woody fpc increase 
    270  
    271                    sumdelta_fpc_wood = sumdelta_fpc_wood + deltafpc(j) 
    272  
    273                    ! which woody pft is preponderant 
    274  
    275                    IF ( fpc_nat(i,j) .GT. maxfpc_wood ) THEN 
    276  
    277                       optpft_wood = j 
    278  
    279                       maxfpc_wood = fpc_nat(i,j) 
    280  
    281                    ENDIF 
    282  
    283                 ELSE 
    284  
    285                    ! grasses 
    286  
    287                    ! total (natural) grass fpc 
    288  
    289                    sumfpc_grass = sumfpc_grass + fpc_nat(i,j) 
    290  
    291                    ! number of grass PFTs present in the grid box 
    292  
    293                    IF ( PFTpresent(i,j) ) THEN 
    294                       num_grass = num_grass + 1 
    295                    ENDIF 
    296  
    297                 ENDIF   ! tree or grass 
    298  
    299              ENDIF   ! natural 
    300  
    301           ENDDO     ! loop over pfts 
    302  
    303           ! 
    304           ! 3.2 light competition: assume wood outcompetes grass 
    305           ! 
    306  
    307           IF (sumfpc_wood .GE. fpc_crit ) THEN 
     346 
     347                      ! grasses 
     348 
     349                      ! total (natural) grass fpc 
     350 
     351                      sumfpc_grass = sumfpc_grass + fpc_nat(i,j) 
     352 
     353                      ! number of grass PFTs present in the grid box 
     354 
     355                      ! IF ( PFTpresent(i,j) ) THEN 
     356                      !    num_grass = num_grass + 1 
     357                      ! ENDIF 
     358 
     359                   ENDIF   ! tree or grass 
     360 
     361                ENDIF   ! natural 
     362 
     363             ENDDO     ! loop over pfts 
     364 
     365             ! 
     366             ! 3.2 light competition: assume wood outcompetes grass 
     367             ! 
     368             !SZ 
     369!!$             IF (sumfpc_wood .GE. fpc_crit ) THEN 
    308370 
    309371             ! 
    … …  
    326388                      ! 
    327389 
    328                       IF ( maxfpc_wood .GE. fpc_crit ) THEN 
    329  
    330                          ! 3.2.1.1.1 one single woody pft is overwhelming 
    331  
    332                          IF ( j .eq. optpft_wood ) THEN 
    333  
    334                             ! reduction for this dominant pft 
    335  
    336                             reduct = un - fpc_crit / fpc_nat(i,j) 
    337  
    338                          ELSE 
    339  
    340                             ! strongly reduce all other woody pfts 
    341                             !   (original DGVM: tree_mercy = 0.0 ) 
    342  
    343                             reduct = un - tree_mercy 
    344  
    345                          ENDIF   ! pft = dominant woody pft 
     390                      ! no single woody pft is overwhelming 
     391                      ! (original DGVM: tree_mercy = 0.0 ) 
     392                      ! The reduction rate is proportional to the ratio deltafpc/fpc. 
     393 
     394                      IF (sumfpc_wood .GE. fpc_crit .AND. fpc_nat(i,j) .GT. min_stomate .AND. &  
     395                           sumdelta_fpc_wood .GT. min_stomate) THEN 
     396 
     397                         ! reduct = MIN( ( ( deltafpc(j)/sumdelta_fpc_wood * & 
     398                         !     (sumfpc_wood-fpc_crit) ) / fpc_nat(i,j) ), & 
     399                         !     ( 1._r_std - tree_mercy ) ) 
     400                         reduct = un - MIN((fpc_nat(i,j)-(sumfpc_wood-fpc_crit) &  
     401                              * deltafpc(j)/sumdelta_fpc_wood)/fpc_nat(i,j), un ) 
    346402 
    347403                      ELSE 
    348404 
    349                          ! 3.2.1.1.2 no single woody pft is overwhelming 
    350                          !           (original DGVM: tree_mercy = 0.0 ) 
    351                          !           The reduction rate is proportional to the ratio deltafpc/fpc. 
    352  
    353                          IF ( fpc_nat(i,j) .GE. min_stomate ) THEN 
    354  
    355                             reduct = MIN( ( ( deltafpc(j)/sumdelta_fpc_wood * & 
    356                                  (sumfpc_wood-fpc_crit) ) / fpc_nat(i,j) ), & 
    357                                  ( 1._r_std - tree_mercy ) ) 
    358  
    359                          ELSE 
    360  
    361                             ! tree fpc didn't icrease or it started from nothing 
    362  
    363                             reduct = zero 
    364  
    365                          ENDIF 
    366  
    367                       ENDIF   ! maxfpc_wood > fpc_crit 
     405                         ! tree fpc didn't icrease or it started from nothing 
     406 
     407                         reduct = zero 
     408 
     409                      ENDIF 
    368410 
    369411                      survive(j) = un - reduct 
    … …  
    379421                      ! 
    380422 
    381                       survive(j) = ( grass_mercy / REAL( num_grass,r_std ) ) / ind(i,j) 
    382  
    383                       survive(j) = MIN( 1._r_std, survive(j) ) 
    384  
     423                      ! survive(j) = ( grass_mercy / REAL( num_grass,r_std ) ) / ind(i,j) 
     424 
     425                      ! survive(j) = MIN( 1._r_std, survive(j) ) 
     426 
     427                      IF(sumfpc_grass .GE. 1.0-MIN(fpc_crit,sumfpc_wood).AND. &  
     428                           sumfpc_grass.GE.min_stomate) THEN 
     429 
     430                         fpc_dec=(sumfpc_grass-1.+MIN(fpc_crit,sumfpc_wood))*fpc_nat(i,j)/sumfpc_grass 
     431 
     432                         reduct=fpc_dec 
     433                      ELSE  
     434                         reduct = zero 
     435                      ENDIF 
     436                      survive(j) = ( un -  reduct )  
    385437                   ENDIF   ! tree or grass 
    386438 
    … …  
    389441             ENDDO       ! loop over pfts 
    390442 
     443             !SZ 
     444!!$          ELSE 
     445!!$ 
     446!!$             ! 
     447!!$             ! 3.2.2 not too much wood so that grasses can subsist 
     448!!$             ! 
     449!!$ 
     450!!$             ! new total grass fpc 
     451!!$             sumfpc_grass2 = fpc_crit - sumfpc_wood 
     452!!$ 
     453!!$             DO j = 2,nvm 
     454!!$ 
     455!!$                ! only present and natural PFTs compete 
     456!!$ 
     457!!$                IF ( PFTpresent(i,j) .AND. natural(j) ) THEN 
     458!!$ 
     459!!$                   IF ( tree(j) ) THEN 
     460!!$ 
     461!!$                      ! no change for trees 
     462!!$ 
     463!!$                      survive(j) = 1.0 
     464!!$ 
     465!!$                   ELSE 
     466!!$ 
     467!!$                      ! grass: fractional loss is the same for all grasses 
     468!!$ 
     469!!$                      IF ( sumfpc_grass .GT. min_stomate ) THEN 
     470!!$                         survive(j) = sumfpc_grass2 / sumfpc_grass 
     471!!$                      ELSE 
     472!!$                         survive(j)=  zero 
     473!!$                      ENDIF 
     474!!$ 
     475!!$                   ENDIF 
     476!!$ 
     477!!$                ENDIF    ! pft there and natural 
     478!!$ 
     479!!$             ENDDO       ! loop over pfts 
     480!!$ 
     481!!$          ENDIF    ! sumfpc_wood > fpc_crit 
     482 
     483             ! 
     484             ! 3.3 update output variables 
     485             ! 
     486 
     487             DO j = 2,nvm 
     488 
     489                IF ( PFTpresent(i,j) .AND. natural(j) ) THEN 
     490 
     491                   bm_to_litter(i,j,:) = bm_to_litter(i,j,:) + & 
     492                        biomass(i,j,:) * ( un - survive(j) ) 
     493 
     494                   biomass(i,j,:) = biomass(i,j,:) * survive(j) 
     495 
     496                   IF ( control%ok_dgvm ) THEN 
     497                      ind(i,j) = ind(i,j) * survive(j) 
     498                   ENDIF 
     499 
     500                   ! fraction of plants that dies each day.  
     501                   ! exact formulation: light_death(i,j) = un - survive(j) / dt 
     502                   light_death(i,j) = ( un - survive(j) ) / dt 
     503 
     504                ENDIF      ! pft there and natural 
     505 
     506             ENDDO        ! loop over pfts 
     507 
     508          ENDIF      ! sumfpc > fpc_crit 
     509 
     510       ENDDO        ! loop over grid points 
     511 
     512       ! 
     513       ! 4 recalculate fpc on natural part of grid cell (for next light competition) 
     514       ! 
     515 
     516       DO j = 2,nvm 
     517 
     518          IF ( natural(j) ) THEN 
     519 
     520             ! 
     521             ! 4.1 natural PFTs 
     522             ! 
     523 
     524             IF ( tree(j) ) THEN 
     525 
     526                ! 4.1.1 trees: minimum cover due to stems, branches etc. 
     527 
     528                DO i = 1, npts 
     529                   !NVMODIF          
     530                   !    IF (lai(i,j) == val_exp) THEN 
     531                   !                veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     532                   !             ELSE 
     533                   !                veget_lastlight(i,j) = & 
     534                   !                     cn_ind(i,j) * ind(i,j) * & 
     535                   !                     MAX( ( un - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) 
     536                   !             ENDIF 
     537                   !!                veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     538                   IF (lai(i,j) == val_exp) THEN 
     539                      veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     540                   ELSE 
     541                      veget_lastlight(i,j) = & 
     542                           cn_ind(i,j) * ind(i,j) * & 
     543                           MAX( ( un - EXP( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) ), min_cover ) 
     544                   ENDIF 
     545                ENDDO 
     546 
     547             ELSE 
     548 
     549                ! 4.1.2 grasses 
     550                DO i = 1, npts 
     551                   !NVMODIF          
     552                   !            IF (lai(i,j) == val_exp) THEN 
     553                   !                veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     554                   !             ELSE 
     555                   !                veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * & 
     556                   !                     ( un - exp( -lai(i,j) * ext_coeff(j) ) ) 
     557                   !             ENDIF 
     558                   !!veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     559                   IF (lai(i,j) == val_exp) THEN 
     560                      veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
     561                   ELSE 
     562                      veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * & 
     563                           ( un - exp( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) ) 
     564                   ENDIF 
     565                ENDDO 
     566             ENDIF    ! tree/grass 
     567 
    391568          ELSE 
    392569 
    393570             ! 
    394              ! 3.2.2 not too much wood so that grasses can subsist 
    395              ! 
    396  
    397              ! new total grass fpc 
    398              sumfpc_grass2 = fpc_crit - sumfpc_wood 
    399  
    400              DO j = 2,nvm 
    401  
    402                 ! only present and natural PFTs compete 
    403  
    404                 IF ( PFTpresent(i,j) .AND. natural(j) ) THEN 
    405  
    406                    IF ( tree(j) ) THEN 
    407  
    408                       ! no change for trees 
    409  
    410                       survive(j) = 1.0 
    411  
    412                    ELSE 
    413  
    414                       ! grass: fractional loss is the same for all grasses 
    415  
    416                       IF ( sumfpc_grass .GT. min_stomate ) THEN 
    417                          survive(j) = sumfpc_grass2 / sumfpc_grass 
    418                       ELSE 
    419                          survive(j)=  zero 
    420                       ENDIF 
    421  
    422                    ENDIF 
    423  
    424                 ENDIF    ! pft there and natural 
    425  
    426              ENDDO       ! loop over pfts 
    427  
    428           ENDIF    ! sumfpc_wood > fpc_crit 
    429  
    430           ! 
    431           ! 3.3 update output variables 
    432           ! 
    433  
    434           DO j = 2,nvm 
    435  
    436              IF ( PFTpresent(i,j) .AND. natural(j) ) THEN 
    437  
    438                 bm_to_litter(i,j,:) = bm_to_litter(i,j,:) + & 
    439                      biomass(i,j,:) * ( un - survive(j) ) 
    440  
    441                 biomass(i,j,:) = biomass(i,j,:) * survive(j) 
    442  
    443                 IF ( control%ok_dgvm ) THEN 
    444                    ind(i,j) = ind(i,j) * survive(j) 
    445                 ENDIF 
    446  
    447                 ! fraction of plants that dies each day.  
    448                 ! exact formulation: light_death(i,j) = un - survive(j) ** (1/dt) 
    449                 light_death(i,j) = ( un - survive(j) ) / dt 
    450  
    451              ENDIF      ! pft there and natural 
    452  
    453           ENDDO        ! loop over pfts 
    454  
    455        ENDIF      ! sumfpc > fpc_crit 
    456  
    457     ENDDO        ! loop over grid points 
    458  
    459     ! 
    460     ! 4 recalculate fpc on natural part of grid cell (for next light competition) 
    461     ! 
    462  
    463     DO j = 2,nvm 
    464  
    465        IF ( natural(j) ) THEN 
    466  
    467           ! 
    468           ! 4.1 natural PFTs 
    469           ! 
    470  
    471           IF ( tree(j) ) THEN 
    472  
    473              ! 4.1.1 trees: minimum cover due to stems, branches etc. 
    474  
    475              DO i = 1, npts 
    476                 IF (lai(i,j) == val_exp) THEN 
    477                    veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
    478                 ELSE 
    479                    veget_lastlight(i,j) = & 
    480                         cn_ind(i,j) * ind(i,j) * & 
    481                         MAX( ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) 
    482                 ENDIF 
     571             ! 4.2 agricultural PFTs: not present on natural part 
     572             ! 
     573 
     574             veget_lastlight(:,j) = zero 
     575 
     576          ENDIF      ! natural/agricultural 
     577 
     578       ENDDO 
     579 
     580    ELSE ! static 
     581 
     582       light_death(:,:)=0.0 
     583 
     584       DO j = 2, nvm 
     585 
     586          IF ( natural(j) ) THEN 
     587 
     588             ! 2.1.1 natural PFTs, in the one PFT only case there needs to be no special case for grasses, 
     589             ! neither a redistribution of mortality (delta fpc) 
     590 
     591             WHERE( ind(:,j)*cn_ind(:,j) .GT. min_stomate )  
     592                lai_ind(:)=sla(j) * lm_lastyearmax(:,j) / ( ind(:,j) * cn_ind(:,j) ) 
     593             ELSEWHERE 
     594                lai_ind(:)=zero 
     595             ENDWHERE 
     596 
     597             fpc_nat(:,j) =  cn_ind(:,j) * ind(:,j) * &  
     598                  MAX( ( 1._r_std - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover ) 
     599 
     600             WHERE(fpc_nat(:,j).GT.fpc_max(:,j)) 
     601 
     602                light_death(:,j)=MIN(1.0,1.0-fpc_max(:,j)/fpc_nat(:,j))  
     603 
     604             ENDWHERE 
     605 
     606             DO k=1,nparts 
     607 
     608                bm_to_litter(:,j,k)=bm_to_litter(:,j,k)+light_death(:,j)*biomass(:,j,k) 
     609                biomass(:,j,k)=biomass(:,j,k)-light_death(:,j)*biomass(:,j,k) 
     610 
    483611             ENDDO 
    484  
    485           ELSE 
    486  
    487              ! 4.1.2 grasses 
    488              DO i = 1, npts 
    489                 IF (lai(i,j) == val_exp) THEN 
    490                    veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)  
    491                 ELSE 
    492                    veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * & 
    493                         ( un - exp( -lai(i,j) * ext_coeff(j) ) ) 
    494                 ENDIF 
    495              ENDDO 
    496           ENDIF    ! tree/grass 
    497  
    498        ELSE 
    499  
    500           ! 
    501           ! 4.2 agricultural PFTs: not present on natural part 
    502           ! 
    503  
    504           veget_lastlight(:,j) = zero 
    505  
    506        ENDIF      ! natural/agricultural 
    507  
    508     ENDDO 
     612             ind(:,j)=ind(:,j)-light_death(:,j)*ind(:,j) 
     613             ! if (j==10) print *,'ind10bis=',ind(:,j),light_death(:,j)*ind(:,j) 
     614          ENDIF 
     615       ENDDO 
     616 
     617       light_death(:,:)=light_death(:,:)/dt 
     618 
     619    ENDIF 
    509620 
    510621    ! 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_pftinout.f90

    r119 r405  
    3232  SUBROUTINE pftinout (npts, dt, adapted, regenerate, & 
    3333       neighbours, veget, veget_max, & 
    34        biomass, ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, & 
     34       biomass, ind, cn_ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, & 
    3535       PFTpresent, everywhere, when_growthinit, need_adjacent, RIP_time, & 
    3636       co2_to_bm, & 
    … …  
    6565    ! density of individuals 1/m**2 
    6666    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: ind 
     67    ! crownarea of individuals m**2 
     68    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: cn_ind 
    6769    ! mean age (years) 
    6870    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: age 
    … …  
    105107    REAL(r_std), DIMENSION(npts)                               :: avail 
    106108    ! indices 
    107     INTEGER(i_std)                                             :: i,j 
     109    INTEGER(i_std)                                             :: i,j,m 
    108110    ! total woody vegetation cover 
    109111    REAL(r_std), DIMENSION(npts)                               :: sumfrac_wood 
    … …  
    112114    ! we can introduce this PFT 
    113115    LOGICAL, DIMENSION(npts)                                  :: can_introduce 
     116    ! no real need for dimension(ntps) except for vectorisation 
     117    REAL(r_std), DIMENSION(npts)                               :: fracnat 
    114118 
    115119    ! ========================================================================= 
    … …  
    133137    ! 
    134138 
    135     ! need to know total woody vegetation fraction 
    136  
     139    ! 2.1 Only natural part of the grid cell 
     140    ! 
     141    !SZ bug correction MERGE: need to subtract agricultural area! 
     142    ! fraction of agricultural surface 
     143    fracnat(:) = 1. 
     144    do j = 2,nvm 
     145       IF ( .NOT. natural(j) ) THEN 
     146          fracnat(:) = fracnat(:) - veget_max(:,j) 
     147       ENDIF 
     148    ENDDO 
     149 
     150    ! 
     151    ! 2.2 total woody fpc on grid 
     152    ! 
    137153    sumfrac_wood(:) = zero 
    138154 
    139155    DO j = 2,nvm 
    140  
    141        IF ( tree(j) ) THEN 
    142  
    143           sumfrac_wood(:) = sumfrac_wood(:) + veget(:,j) 
    144  
     156       !SZ problem here: agriculture, not convinced that this representation of LPJ is correct 
     157       !if agriculture is present, ind must be recalculated to correspond to the natural density... 
     158       ! since ind is per grid cell, can be achived by discounting for agricultura fraction 
     159       IF ( natural(j).AND.tree(j) ) THEN 
     160          WHERE(fracnat(:).GT.min_stomate) 
     161                sumfrac_wood(:) = sumfrac_wood(:) + cn_ind(:,j) * ind(:,j) / fracnat(:) &  
     162                     * ( 1. - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ) 
     163                !lai changed to lm_last 
     164          ENDWHERE 
    145165       ENDIF 
    146  
    147166    ENDDO 
    148167 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate.f90

    r119 r405  
    3030  IMPLICIT NONE 
    3131  PRIVATE 
    32   PUBLIC stomate_main,stomate_clear 
     32  PUBLIC stomate_main,stomate_clear,init_forcing,forcing_read 
    3333  ! 
    3434  INTEGER,PARAMETER :: r_typ =nf90_real4 
    … …  
    231231  REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:)  :: harvest_above_monthly, cflux_prod_monthly 
    232232 
     233  ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground  
     234  REAL(r_std), ALLOCATABLE,SAVE,DIMENSION(:,:)              :: fpc_max 
     235 
    233236  ! Date and EndOfYear, intialize and update in slowproc 
    234237  ! (Now managed in slowproc for land_use) 
    … …  
    263266  REAL(r_std),ALLOCATABLE,DIMENSION(:,:)     :: precip_fm 
    264267  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: gpp_daily_fm 
    265   REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:,:) :: resp_maint_part_fm 
    266268  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: veget_fm 
    267269  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: veget_max_fm 
    268270  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: lai_fm 
     271  PUBLIC clay_fm, humrel_daily_fm, litterhum_daily_fm, t2m_daily_fm, t2m_min_daily_fm, tsurf_daily_fm, tsoil_daily_fm, & 
     272       soilhum_daily_fm, precip_fm, gpp_daily_fm, veget_fm, veget_max_fm, lai_fm 
    269273 
    270274  REAL(r_std),ALLOCATABLE,DIMENSION(:,:)     :: clay_fm_g 
    … …  
    278282  REAL(r_std),ALLOCATABLE,DIMENSION(:,:)     :: precip_fm_g 
    279283  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: gpp_daily_fm_g 
    280   REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:,:) :: resp_maint_part_fm_g 
    281284  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: veget_fm_g 
    282285  REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:)   :: veget_max_fm_g 
    … …  
    286289  LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:)      :: nf_written 
    287290  INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: nf_cumul 
     291  PUBLIC isf, nf_written 
     292   
    288293  ! first call 
    289294  LOGICAL,SAVE :: l_first_stomate = .TRUE. 
    … …  
    312317  ! harvest above ground biomass for agriculture 
    313318  REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:)                            :: harvest_above 
     319 
     320  ! Carbon Mass total 
     321  REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:)                            :: carb_mass_total 
    314322 
    315323CONTAINS 
    … …  
    327335       &  veget_max_new, totfrac_nobio_new, & 
    328336       &  hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, & 
    329        &  co2_flux,resp_maint,resp_hetero,resp_growth) 
     337       &  co2_flux, fco2_lu, resp_maint,resp_hetero,resp_growth) 
    330338    !--------------------------------------------------------------------- 
    331339    ! 
    … …  
    417425    !NV champs 2D  
    418426    REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out)      :: co2_flux 
     427    REAL(r_std),DIMENSION(kjpindex),INTENT(out)      :: fco2_lu 
    419428    ! autotrophic respiration in gC/m**2 of surface/dt 
    420429    REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out)  :: resp_maint 
    … …  
    490499    ! for forcing file: "daily" gpp 
    491500    REAL(r_std),DIMENSION(kjpindex,nvm)           :: gpp_daily_x 
    492     ! for forcing file: "daily" auto resp 
    493     REAL(r_std),DIMENSION(kjpindex,nvm,nparts)    :: resp_maint_part_x 
    494501    ! total "vegetation" cover 
    495502    REAL(r_std),DIMENSION(kjpindex)               :: cvegtot 
    … …  
    511518    INTEGER(i_std),SAVE        :: nparan            ! Number of time steps per year for carbon spinup 
    512519    INTEGER(i_std),SAVE        :: nbyear 
    513     INTEGER(i_std),PARAMETER   :: nparanmax=36      ! Number max of time steps per year for carbon spinup 
     520    INTEGER(i_std),PARAMETER   :: nparanmax=366     ! Number max of time steps per year for carbon spinup 
    514521    REAL(r_std)                 :: sf_time 
    515     INTEGER(i_std),SAVE        :: iatt=1 
     522    INTEGER(i_std),SAVE        :: iatt 
    516523    INTEGER(i_std),SAVE        :: iatt_old=1 
    517524    INTEGER(i_std)             :: max_totsize, totsize_1step,totsize_tmp 
    … …  
    591598            rest_id_stom, hist_id_stom, hist_id_stom_IPCC) 
    592599 
    593        co2_flux_monthly(:,:) = zero 
    594600       ! 
    595601       ! 1.2 read PFT data 
    … …  
    600606       ! 
    601607       ! 1.3.1 read STOMATE's start file 
     608       ! 
     609       co2_flux(:,:) = zero 
     610       fco2_lu(:) = zero 
    602611       ! 
    603612       CALL readstart & 
    … …  
    629638            &         carbon, black_carbon, lignin_struc,turnover_time,& 
    630639            &         prod10,prod100,flux10, flux100, & 
    631             &         convflux, cflux_prod10, cflux_prod100, bm_to_litter) 
     640            &         convflux, cflux_prod10, cflux_prod100, bm_to_litter, carb_mass_total) 
    632641 
    633642       ! 1.4 read the boundary conditions 
    … …  
    726735                  &     +SIZE(precip_daily)*KIND(precip_daily) & 
    727736                  &     +SIZE(gpp_daily_x)*KIND(gpp_daily_x) & 
    728                   &     +SIZE(resp_maint_part_x)*KIND(resp_maint_part_x) & 
    729737                  &     +SIZE(veget)*KIND(veget) & 
    730738                  &     +SIZE(veget_max)*KIND(veget_max) & 
    … …  
    813821                ier = NF90_DEF_VAR (forcing_id,'lai', & 
    814822                     &                            r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid) 
    815                 ier = NF90_DEF_VAR (forcing_id,'resp_maint_part', & 
    816                      &                       r_typ,(/ d_id(1),d_id(3),d_id(7),d_id(6) /),vid) 
    817823                ier = NF90_ENDDEF (forcing_id) 
    818824                !- 
    … …  
    867873             !Config  Key  = FORCESOIL_STEP_PER_YEAR 
    868874             !Config  Desc = Number of time steps per year for carbon spinup. 
    869              !Config  Def  = 12 
     875             !Config  Def  = 365 
    870876             !Config  Help = Number of time steps per year for carbon spinup. 
    871              nparan = 12 
     877             nparan = 365 
    872878             CALL getin_p('FORCESOIL_STEP_PER_YEAR', nparan) 
    873879 
    … …  
    10081014            &          carbon, black_carbon, lignin_struc,turnover_time,& 
    10091015            &          prod10,prod100,flux10, flux100, & 
    1010             &          convflux, cflux_prod10, cflux_prod100, bm_to_litter) 
     1016            &          convflux, cflux_prod10, cflux_prod100, bm_to_litter,carb_mass_total) 
    10111017 
    10121018       IF (ldforcing_write .AND. TRIM(forcing_name) /= 'NONE' ) THEN   
    … …  
    13361342               &             t_photo_min, t_photo_opt, t_photo_max,bm_to_litter,& 
    13371343               &             prod10, prod100, flux10, flux100, veget_cov_max_new,& 
    1338                &             convflux, cflux_prod10, cflux_prod100, harvest_above, lcchange) 
    1339  
     1344               &             convflux, cflux_prod10, cflux_prod100, harvest_above, carb_mass_total, lcchange,& 
     1345               &             fpc_max) 
     1346 
     1347          ! 
     1348          ! fco2_lu --> luccarb 
     1349          fco2_lu(:) = convflux(:) & 
     1350               &             + cflux_prod10(:)  & 
     1351               &             + cflux_prod100(:) & 
     1352               &             + harvest_above(:) 
    13401353          ! 
    13411354          ! 6.4 output: transform from dimension nvm to dimension nvm 
    … …  
    13971410             ENDDO 
    13981411             sf_time = MODULO(REAL(date,r_std)-1,one_year*REAL(nbyear,r_std)) 
    1399              iatt=FLOOR(sf_time/dt_forcesoil)+1 
    1400              IF ((iatt < 1) .OR. (iatt > nparan*nbyear)) THEN 
    1401                 WRITE(numout,*) 'Error with iatt=',iatt 
    1402                 CALL ipslerr (3,'stomate', & 
    1403                      &          'Error with iatt.', '', & 
    1404                      &          '(Problem with dt_forcesoil ?)') 
    1405              ENDIF 
     1412             iatt=FLOOR(sf_time/dt_forcesoil) 
     1413             IF (iatt == 0) iatt = iatt_old + 1 
    14061414 
    14071415             IF ((iatt<iatt_old) .and. (.not. cumul_Cforcing)) THEN 
    14081416                nforce(:)=0 
    1409                 soilcarbon_input(:,:,:,:) = 0 
    1410                 control_moist(:,:,:) = 0 
    1411                 control_temp(:,:,:) = 0 
    1412                 npp_equil(:,:) = 0 
     1417                soilcarbon_input(:,:,:,:) = zero 
     1418                control_moist(:,:,:) = zero 
     1419                control_temp(:,:,:) = zero 
     1420                npp_equil(:,:) = zero 
    14131421             ENDIF 
    14141422             iatt_old=iatt 
    … …  
    14371445 
    14381446          gpp_daily_x(:,:) = zero 
    1439           resp_maint_part_x(:,:,:) = zero 
    14401447          !gpp needs to be multiplied by coverage for forcing (see above) 
    14411448          DO j = 2, nvm              
    14421449             gpp_daily_x(:,j) = gpp_daily_x(:,j) + & 
    14431450                  &                              gpp_daily(:,j) * dt_slow / one_day * veget_cov_max(:,j) 
    1444              resp_maint_part_x(:,j,:) = resp_maint_part_x(:,j,:) + & 
    1445                   &                              resp_maint_part(:,j,:) * dt_slow / one_day 
    14461451          ENDDO 
    14471452          ! 
    … …  
    14791484             gpp_daily_fm(:,:,iisf) = & 
    14801485                  &                (xn*gpp_daily_fm(:,:,iisf) + gpp_daily_x(:,:))/(xn+1.) 
    1481              resp_maint_part_fm(:,:,:,iisf) = & 
    1482                   &                ( xn*resp_maint_part_fm(:,:,:,iisf) & 
    1483                   &         +resp_maint_part_x(:,:,:) )/(xn+1.) 
    14841486             veget_fm(:,:,iisf) = & 
    14851487                  &                (xn*veget_fm(:,:,iisf) + veget(:,:) )/(xn+1.) 
    … …  
    14911493             clay_fm(:,iisf) = clay(:) 
    14921494             humrel_daily_fm(:,:,iisf) = humrel_daily(:,:) 
    1493              litterhum_daily_fm(:,iisf) = +litterhum_daily(:) 
     1495             litterhum_daily_fm(:,iisf) = litterhum_daily(:) 
    14941496             t2m_daily_fm(:,iisf) = t2m_daily(:) 
    14951497             t2m_min_daily_fm(:,iisf) =t2m_min_daily(:) 
    … …  
    14991501             precip_fm(:,iisf) = precip_daily(:) 
    15001502             gpp_daily_fm(:,:,iisf) =gpp_daily_x(:,:) 
    1501              resp_maint_part_fm(:,:,:,iisf) = resp_maint_part_x(:,:,:) 
    15021503             veget_fm(:,:,iisf) = veget(:,:) 
    15031504             veget_max_fm(:,:,iisf) =veget_max(:,:) 
    … …  
    17161717    ! allocation error 
    17171718    LOGICAL                                     :: l_error 
    1718     ! Global world fraction of vegetation type map 
    1719     REAL(r_std),DIMENSION(360,180,nvm)           :: veget_ori_on_disk 
    17201719    INTEGER(i_std)                              :: ier 
    17211720    ! indices 
    … …  
    19861985    ALLOCATE (harvest_above(kjpindex), stat=ier) 
    19871986    l_error = l_error .OR. (ier.NE.0) 
     1987    ALLOCATE (carb_mass_total(kjpindex), stat=ier) 
     1988    l_error = l_error .OR. (ier.NE.0) 
    19881989    ALLOCATE (soilcarbon_input_daily(kjpindex,ncarb,nvm), stat=ier) 
    19891990    l_error = l_error .OR. (ier.NE.0) 
    … …  
    19931994    l_error = l_error .OR. (ier.NE.0) 
    19941995    ! 
     1996    ALLOCATE (fpc_max(kjpindex,nvm), stat=ier) 
     1997    l_error = l_error .OR. (ier.NE.0) 
     1998    ! 
    19951999    IF (l_error) THEN 
    19962000       STOP 'stomate_init: error in memory allocation' 
    … …  
    20662070    WRITE(numout,*) & 
    20672071         &  'expansion across a grid cell is treated: ',treat_expansion 
     2072 
     2073    !Config Key  = LPJ_GAP_CONST_MORT 
     2074    !Config Desc = prescribe mortality if not using DGVM? 
     2075    !Config Def  = y 
     2076    !Config Help = set to TRUE if constant mortality is to be activated 
     2077    !              ignored if DGVM=true! 
     2078    ! 
     2079    lpj_gap_const_mort=.TRUE. 
     2080    CALL getin('LPJ_GAP_CONST_MORT', lpj_gap_const_mort) 
     2081    WRITE(numout,*) 'LPJ GAP: constant mortality:', lpj_gap_const_mort 
    20682082 
    20692083    !Config  Key  = HARVEST_AGRI 
    … …  
    20982112    cflux_prod10(:) = zero 
    20992113    cflux_prod100(:)= zero 
     2114 
     2115    fpc_max(:,:)=zero 
    21002116    !-------------------------- 
    21012117  END SUBROUTINE stomate_init 
    … …  
    22032219    IF (ALLOCATED(precip_fm)) DEALLOCATE(precip_fm) 
    22042220    IF (ALLOCATED(gpp_daily_fm))  DEALLOCATE(gpp_daily_fm) 
    2205     IF (ALLOCATED(resp_maint_part_fm))  DEALLOCATE(resp_maint_part_fm) 
    22062221    IF (ALLOCATED(veget_fm)) DEALLOCATE(veget_fm) 
    22072222    IF (ALLOCATED(veget_max_fm)) DEALLOCATE(veget_max_fm) 
    … …  
    22192234       IF (ALLOCATED(precip_fm_g)) DEALLOCATE(precip_fm_g) 
    22202235       IF (ALLOCATED(gpp_daily_fm_g))  DEALLOCATE(gpp_daily_fm_g) 
    2221        IF (ALLOCATED(resp_maint_part_fm_g))  DEALLOCATE(resp_maint_part_fm_g) 
    22222236       IF (ALLOCATED(veget_fm_g)) DEALLOCATE(veget_fm_g) 
    22232237       IF (ALLOCATED(veget_max_fm_g)) DEALLOCATE(veget_max_fm_g) 
    … …  
    22472261    IF ( ALLOCATED (control_temp_daily)) DEALLOCATE (control_temp_daily) 
    22482262    IF ( ALLOCATED (control_moist_daily)) DEALLOCATE (control_moist_daily) 
     2263 
     2264    IF ( ALLOCATED (fpc_max)) DEALLOCATE (fpc_max) 
    22492265 
    22502266    ! 2. reset l_first 
    … …  
    24592475    ALLOCATE(gpp_daily_fm(kjpindex,nvm,nsfm),stat=ier) 
    24602476    l_error = l_error .OR. (ier /= 0) 
    2461     ALLOCATE(resp_maint_part_fm(kjpindex,nvm,nparts,nsfm),stat=ier) 
    2462     l_error = l_error .OR. (ier /= 0) 
    24632477    ALLOCATE(veget_fm(kjpindex,nvm,nsfm),stat=ier) 
    24642478    l_error = l_error .OR. (ier /= 0) 
    … …  
    24732487    ALLOCATE(nf_cumul(nsft),stat=ier) 
    24742488    l_error = l_error .OR. (ier /= 0) 
     2489    IF (l_error) THEN 
     2490       WRITE(numout,*) 'Problem with memory allocation: forcing variables' 
     2491       STOP 'init_forcing' 
     2492    ENDIF 
    24752493 
    24762494    IF (is_root_prc) THEN 
    … …  
    24952513       ALLOCATE(gpp_daily_fm_g(nbp_glo,nvm,nsfm),stat=ier) 
    24962514       l_error = l_error .OR. (ier /= 0) 
    2497        ALLOCATE(resp_maint_part_fm_g(nbp_glo,nvm,nparts,nsfm),stat=ier) 
    2498        l_error = l_error .OR. (ier /= 0) 
    24992515       ALLOCATE(veget_fm_g(nbp_glo,nvm,nsfm),stat=ier) 
    25002516       l_error = l_error .OR. (ier /= 0) 
    … …  
    25032519       ALLOCATE(lai_fm_g(nbp_glo,nvm,nsfm),stat=ier) 
    25042520       l_error = l_error .OR. (ier /= 0) 
     2521       IF (l_error) THEN 
     2522          WRITE(numout,*) 'Problem with memory allocation: forcing variables' 
     2523          STOP 'init_forcing' 
     2524       ENDIF 
     2525    ELSE 
     2526       ALLOCATE(clay_fm_g(0,nsfm),stat=ier) 
     2527       ALLOCATE(humrel_daily_fm_g(0,nvm,nsfm),stat=ier) 
     2528       ALLOCATE(litterhum_daily_fm_g(0,nsfm),stat=ier) 
     2529       ALLOCATE(t2m_daily_fm_g(0,nsfm),stat=ier) 
     2530       ALLOCATE(t2m_min_daily_fm_g(0,nsfm),stat=ier) 
     2531       ALLOCATE(tsurf_daily_fm_g(0,nsfm),stat=ier) 
     2532       ALLOCATE(tsoil_daily_fm_g(0,nbdl,nsfm),stat=ier) 
     2533       ALLOCATE(soilhum_daily_fm_g(0,nbdl,nsfm),stat=ier) 
     2534       ALLOCATE(precip_fm_g(0,nsfm),stat=ier) 
     2535       ALLOCATE(gpp_daily_fm_g(0,nvm,nsfm),stat=ier) 
     2536       ALLOCATE(veget_fm_g(0,nvm,nsfm),stat=ier) 
     2537       ALLOCATE(veget_max_fm_g(0,nvm,nsfm),stat=ier) 
     2538       ALLOCATE(lai_fm_g(0,nvm,nsfm),stat=ier) 
    25052539    ENDIF 
    25062540    ! 
    … …  
    25282562    precip_fm(:,:) = zero 
    25292563    gpp_daily_fm(:,:,:) = zero 
    2530     resp_maint_part_fm(:,:,:,:)=zero 
    25312564    veget_fm(:,:,:) = zero 
    25322565    veget_max_fm(:,:,:) = zero 
    … …  
    25802613    CALL gather(precip_fm,precip_fm_g) 
    25812614    CALL gather(gpp_daily_fm,gpp_daily_fm_g) 
    2582     CALL gather(resp_maint_part_fm,resp_maint_part_fm_g) 
    25832615    CALL gather(veget_fm,veget_fm_g) 
    25842616    CALL gather(veget_max_fm,veget_max_fm_g) 
    … …  
    26672699                  &            gpp_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    26682700                  &            start=start(1:ndim), count=count_force(1:ndim)) 
    2669              ndim = 4; 
    2670              start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    2671              count_force(1:ndim)=SHAPE(resp_maint_part_fm_g) 
    2672              count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    2673              ier = NF90_INQ_VARID (forcing_id,'resp_maint_part',vid) 
    2674              ier = NF90_PUT_VAR (forcing_id,vid, & 
    2675                   &            resp_maint_part_fm_g(:,:,:,ifirst(iblocks):ilast(iblocks)), & 
    2676                   &            start=start(1:ndim), count=count_force(1:ndim)) 
    26772701             ndim = 3; 
    26782702             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    27152739    INTEGER(i_std)                :: iisf, iblocks, nblocks 
    27162740    INTEGER(i_std)                :: ier 
     2741    LOGICAL    :: a_er 
    27172742    INTEGER(i_std),DIMENSION(0:2) :: ifirst, ilast 
    27182743    INTEGER(i_std),PARAMETER      :: ndm = 10 
    27192744    INTEGER(i_std),DIMENSION(ndm) :: start, count_force 
    27202745    INTEGER(i_std)                :: ndim, vid 
     2746 
     2747    LOGICAL, PARAMETER :: check=.FALSE. 
     2748 
     2749    IF (check) WRITE(numout,*) "forcing_read " 
    27212750    !--------------------------------------------------------------------- 
    27222751    ! 
    … …  
    27362765          precip_fm(:,iisf) = zero 
    27372766          gpp_daily_fm(:,:,iisf) = zero 
    2738           resp_maint_part_fm(:,:,:,iisf) = zero 
    27392767          veget_fm(:,:,iisf) = zero 
    27402768          veget_max_fm(:,:,iisf) = zero 
    … …  
    27652793       ENDIF 
    27662794    ENDDO 
     2795    IF (check) WRITE(numout,*) "forcing_read nblocks, ifirst, ilast",nblocks, ifirst, ilast 
    27672796    ! 
    27682797    IF (is_root_prc) THEN 
    27692798       DO iblocks = 1, nblocks 
     2799          IF (check) WRITE(numout,*) "forcing_read iblocks, ifirst(iblocks), ilast(iblocks)",iblocks, & 
     2800               ifirst(iblocks), ilast(iblocks) 
    27702801          IF (ifirst(iblocks) /= ilast(iblocks)) THEN 
     2802             a_er=.FALSE. 
    27712803             ndim = 2; 
    27722804             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    27742806             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    27752807             ier = NF90_INQ_VARID (forcing_id,'clay',vid) 
     2808             a_er = a_er.OR.(ier.NE.0) 
    27762809             ier = NF90_GET_VAR (forcing_id, vid, & 
    27772810                  &            clay_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    27782811                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2812             a_er = a_er.OR.(ier.NE.0) 
     2813!--------- 
    27792814             ndim = 3; 
    27802815             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    27822817             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    27832818             ier = NF90_INQ_VARID (forcing_id,'humrel',vid) 
     2819             a_er = a_er.OR.(ier.NE.0) 
    27842820             ier = NF90_GET_VAR (forcing_id, vid, & 
    27852821                  &            humrel_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    27862822                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2823             a_er = a_er.OR.(ier.NE.0) 
     2824!--------- 
    27872825             ndim = 2; 
    27882826             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    27902828             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    27912829             ier = NF90_INQ_VARID (forcing_id,'litterhum',vid) 
     2830             a_er = a_er.OR.(ier.NE.0) 
    27922831             ier = NF90_GET_VAR (forcing_id, vid, & 
    27932832                  &              litterhum_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    27942833                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2834             a_er = a_er.OR.(ier.NE.0) 
     2835!--------- 
    27952836             ndim = 2; 
    27962837             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    27982839             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    27992840             ier = NF90_INQ_VARID (forcing_id,'t2m',vid) 
     2841             a_er = a_er.OR.(ier.NE.0) 
    28002842             ier = NF90_GET_VAR (forcing_id, vid, & 
    28012843                  &              t2m_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    28022844                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2845             a_er = a_er.OR.(ier.NE.0) 
     2846!--------- 
    28032847             ndim = 2; 
    28042848             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28062850             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28072851             ier = NF90_INQ_VARID (forcing_id,'t2m_min',vid) 
     2852             a_er = a_er.OR.(ier.NE.0) 
    28082853             ier = NF90_GET_VAR (forcing_id, vid, & 
    28092854                  &              t2m_min_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    28102855                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2856             a_er = a_er.OR.(ier.NE.0) 
     2857!--------- 
    28112858             ndim = 2; 
    28122859             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28142861             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28152862             ier = NF90_INQ_VARID (forcing_id,'tsurf',vid) 
     2863             a_er = a_er.OR.(ier.NE.0) 
    28162864             ier = NF90_GET_VAR (forcing_id, vid, & 
    28172865                  &              tsurf_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    28182866                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2867             a_er = a_er.OR.(ier.NE.0) 
     2868!--------- 
    28192869             ndim = 3; 
    28202870             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28222872             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28232873             ier = NF90_INQ_VARID (forcing_id,'tsoil',vid) 
     2874             a_er = a_er.OR.(ier.NE.0) 
    28242875             ier = NF90_GET_VAR (forcing_id, vid, & 
    28252876                  &              tsoil_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28262877                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2878             a_er = a_er.OR.(ier.NE.0) 
     2879!--------- 
    28272880             ndim = 3; 
    28282881             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28302883             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28312884             ier = NF90_INQ_VARID (forcing_id,'soilhum',vid) 
     2885             a_er = a_er.OR.(ier.NE.0) 
    28322886             ier = NF90_GET_VAR (forcing_id, vid, & 
    28332887                  &              soilhum_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28342888                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2889             a_er = a_er.OR.(ier.NE.0) 
     2890!--------- 
    28352891             ndim = 2; 
    28362892             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28382894             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28392895             ier = NF90_INQ_VARID (forcing_id,'precip',vid) 
     2896             a_er = a_er.OR.(ier.NE.0) 
    28402897             ier = NF90_GET_VAR (forcing_id, vid, & 
    28412898                  &              precip_fm_g(:,ifirst(iblocks):ilast(iblocks)), & 
    28422899                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2900             a_er = a_er.OR.(ier.NE.0) 
     2901!--------- 
    28432902             ndim = 3; 
    28442903             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28462905             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28472906             ier = NF90_INQ_VARID (forcing_id,'gpp',vid) 
     2907             a_er = a_er.OR.(ier.NE.0) 
    28482908             ier = NF90_GET_VAR (forcing_id, vid, & 
    28492909                  &            gpp_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28502910                  &            start=start(1:ndim), count=count_force(1:ndim)) 
    2851              ndim = 4; 
    2852              start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    2853              count_force(1:ndim)=SHAPE(resp_maint_part_fm_g) 
    2854              count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    2855              ier = NF90_INQ_VARID (forcing_id,'resp_maint_part',vid) 
    2856              ier = NF90_GET_VAR (forcing_id,vid, & 
    2857                   &       resp_maint_part_fm_g(:,:,:,ifirst(iblocks):ilast(iblocks)), & 
    2858                   &            start=start(1:ndim), count=count_force(1:ndim)) 
     2911             a_er = a_er.OR.(ier.NE.0) 
     2912!--------- 
    28592913             ndim = 3; 
    28602914             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28622916             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28632917             ier = NF90_INQ_VARID (forcing_id,'veget',vid) 
     2918             a_er = a_er.OR.(ier.NE.0) 
    28642919             ier = NF90_GET_VAR (forcing_id, vid, & 
    28652920                  &            veget_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28662921                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2922             a_er = a_er.OR.(ier.NE.0) 
     2923!--------- 
    28672924             ndim = 3; 
    28682925             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28702927             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28712928             ier = NF90_INQ_VARID (forcing_id,'veget_max',vid) 
     2929             a_er = a_er.OR.(ier.NE.0) 
    28722930             ier = NF90_GET_VAR (forcing_id, vid, & 
    28732931                  &            veget_max_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28742932                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2933             a_er = a_er.OR.(ier.NE.0) 
     2934!--------- 
    28752935             ndim = 3; 
    28762936             start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks)); 
    … …  
    28782938             count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1 
    28792939             ier = NF90_INQ_VARID (forcing_id,'lai',vid) 
     2940             a_er = a_er.OR.(ier.NE.0) 
    28802941             ier = NF90_GET_VAR (forcing_id, vid, & 
    28812942                  &            lai_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), & 
    28822943                  &            start=start(1:ndim), count=count_force(1:ndim)) 
     2944             a_er = a_er.OR.(ier.NE.0) 
     2945             IF (a_er) THEN 
     2946                CALL ipslerr (3,'forcing_read', & 
     2947                     &        'PROBLEM when read forcing file', & 
     2948                     &        '','') 
     2949             ENDIF 
    28832950          ENDIF 
    28842951       ENDDO 
    … …  
    28942961    CALL scatter(precip_fm_g,precip_fm) 
    28952962    CALL scatter(gpp_daily_fm_g,gpp_daily_fm) 
    2896     CALL scatter(resp_maint_part_fm_g,resp_maint_part_fm) 
    28972963    CALL scatter(veget_fm_g,veget_fm) 
    28982964    CALL scatter(veget_max_fm_g,veget_max_fm) 
    2899     CALL scatter(lai_fm_g,lai_fm_g) 
     2965    CALL scatter(lai_fm_g,lai_fm) 
    29002966    !-------------------------- 
    29012967  END SUBROUTINE forcing_read 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate_constants.f90

    r119 r405  
    160160! initial density of individuals 
    161161  REAL(r_std),PARAMETER :: ind_0 = 0.02 
     162  ! min npp to test competition between grass 
     163  REAL(r_std), PARAMETER :: npp_min = 100. 
    162164!- 
    163165! Do we treat PFT expansion across a grid point after introduction? 
    164166! default = .FALSE. 
    165167  LOGICAL,SAVE :: treat_expansion = .FALSE. 
     168! Do we treat calculate constant mortality if vegetation is static? 
     169! default = .TRUE. 
     170  LOGICAL, SAVE :: lpj_gap_const_mort = .TRUE. 
    166171!- 
    167172! herbivores? 
    … …  
    193198! fraction of GPP which is lost as growth respiration 
    194199  REAL(r_std),PARAMETER :: frac_growthresp = 0.28 
     200!- 
     201! minimum availability to calculate mortality 
     202  REAL(r_std),PARAMETER :: min_avail = 0.02 
    195203!- 
    196204! description of the PFT 
    … …  
    498506! critical tmin, tabulated (C) 
    499507  tmin_crit_tab(2:nvm) =    & 
    500  & (/     0.0,     0.0,   -45.0,   -10.0,   -45.0,   -60.0, & 
    501  &      -60.0,   undef,   undef,   undef,   undef,   undef /) 
     508 & (/     0.0,     0.0,   -30.0,   -14.0,   -30.0,   -45.0, & 
     509 &      -45.0,   undef,   undef,   undef,   undef,   undef /) 
    502510! critical tcm, tabulated (C) 
    503511  tcm_crit_tab(2:nvm) =     & 
    504  & (/   undef,   undef,     5.0,    15.5,    15.5,    -2.0, & 
    505  &        5.0,    -2.0,   undef,   undef,   undef,   undef /) 
     512 & (/   undef,   undef,     5.0,    15.5,    15.5,    -8.0, & 
     513 &        -8.0,    -8.0,   undef,   undef,   undef,   undef /) 
    506514! critical gdd, tabulated (C), constant c of aT^2+bT+c 
    507515  gdd_crit1_tab(2:nvm) =    & 
    … …  
    552560 &            1.,      1.,      1.,      1.,      1.,      1.      /) 
    553561! Maximum rate of carboxylation 
     562  !Config Key  = vcmax_opt 
     563  !Config Desc = Maximum rate of carboxylation 
     564  !Config Def  = undef, 65., 65., 35., 45., 55., 35., 45., 35., 70., 70., 70., 70. 
     565  !Config Help =  
     566  ! 
    554567!Shilong 
    555568  vcmax_opt(:) =     & 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate_io.f90

    r119 r405  
    5555       &  carbon, black_carbon, lignin_struc,turnover_time, & 
    5656       &  prod10,prod100,flux10, flux100, & 
    57        &  convflux, cflux_prod10, cflux_prod100, bm_to_litter) 
     57       &  convflux, cflux_prod10, cflux_prod100, bm_to_litter, carb_mass_total) 
    5858    !--------------------------------------------------------------------- 
    5959    !- read start file 
    … …  
    275275    REAL(r_std), DIMENSION(npts), INTENT(out)                            :: cflux_prod100 
    276276    REAL(r_std),DIMENSION(npts,nvm,nparts),INTENT(out)                   :: bm_to_litter 
     277    REAL(r_std),DIMENSION(npts),INTENT(out)                              :: carb_mass_total 
    277278    !--------------------------------------------------------------------- 
    278279    IF (bavard >= 3) WRITE(numout,*) 'Entering readstart' 
    … …  
    342343       date = NINT(date_real) 
    343344    ENDIF 
    344     CALL bcast(date_real) 
     345    CALL bcast(date) 
    345346    !- 
    346347    ! 3 daily meteorological variables 
    … …  
    940941    ENDDO 
    941942 
     943    carb_mass_total(:) = val_exp 
     944    var_name = 'carb_mass_total' 
     945    CALL restget_p (rest_id_stomate, var_name, nbp_glo, 1     , 1, itime, & 
     946         &              .TRUE., carb_mass_total, 'gather', nbp_glo, index_g) 
     947    IF (ALL(carb_mass_total(:) == val_exp)) carb_mass_total(:) = zero 
    942948    !- 
    943949 
    … …  
    971977       &  carbon, black_carbon, lignin_struc, turnover_time, & 
    972978       &  prod10,prod100 ,flux10, flux100, & 
    973        &  convflux, cflux_prod10, cflux_prod100, bm_to_litter) 
     979       &  convflux, cflux_prod10, cflux_prod100, bm_to_litter, carb_mass_total) 
    974980 
    975981    !--------------------------------------------------------------------- 
    … …  
    11791185    REAL(r_std), DIMENSION(npts), INTENT(in)                            :: cflux_prod100 
    11801186    REAL(r_std),DIMENSION(npts,nvm,nparts),INTENT(in)                   :: bm_to_litter 
     1187    REAL(r_std),DIMENSION(npts),INTENT(in)                              :: carb_mass_total 
    11811188    !--------------------------------------------------------------------- 
    11821189    IF (bavard >= 3) WRITE(numout,*) 'Entering writerestart' 
    … …  
    16431650            &                bm_to_litter(:,:,k), 'scatter', nbp_glo, index_g) 
    16441651    ENDDO 
     1652    var_name = 'carb_mass_total' 
     1653    CALL restput_p (rest_id_stomate, var_name, nbp_glo, 1, 1, itime, & 
     1654         &              carb_mass_total, 'scatter', nbp_glo, index_g) 
    16451655    !- 
    16461656    IF (bavard >= 4) WRITE(numout,*) 'Leaving writerestart' 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate_lpj.f90

    r119 r405  
    9292       t_photo_min, t_photo_opt, t_photo_max,bm_to_litter, & 
    9393       prod10,prod100,flux10, flux100, veget_max_new, & 
    94        convflux,cflux_prod10,cflux_prod100, harvest_above, lcchange) 
     94       convflux,cflux_prod10,cflux_prod100, harvest_above, carb_mass_total, lcchange, & 
     95       fpc_max) 
    9596 
    9697    ! 
    … …  
    166167    ! maintenance respiration of different plant parts (gC/day/m**2 of ground) 
    167168    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)             :: resp_maint_part 
     169    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground  
     170    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: fpc_max 
    168171 
    169172    ! 0.2 modified fields 
    … …  
    292295    ! harvest above ground biomass for agriculture 
    293296    REAL(r_std), DIMENSION(npts), INTENT(inout)                       :: harvest_above 
     297    ! Carbon Mass total 
     298    REAL(r_std), DIMENSION(npts), INTENT(inout)                       :: carb_mass_total 
    294299 
    295300    ! land cover change flag 
    … …  
    319324    ! total soil carbon (gC/(m**2)) 
    320325    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_soil_carb 
     326    ! Carbon Mass variation 
     327    REAL(r_std), DIMENSION(npts)                                      :: carb_mass_variation 
    321328    ! crown area of individuals (m**2) 
    322329    REAL(r_std), DIMENSION(npts,nvm)                               :: cn_ind 
     330    ! woodmass of individuals (gC) 
     331    REAL(r_std), DIMENSION(npts,nvm)                               :: woodmass_ind 
    323332    ! fraction that goes into plant part 
    324333    REAL(r_std), DIMENSION(npts,nvm,nparts)                        :: f_alloc 
    … …  
    337346    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground 
    338347    REAL(r_std),DIMENSION(npts,nvm)                                :: veget_max_old 
     348 
     349    ! fraction of individual dying this time step 
     350    REAL(r_std), DIMENSION(npts,nvm)                               :: mortality 
    339351 
    340352    REAL(r_std), DIMENSION(npts)                                   :: vartmp 
    … …  
    367379    bm_to_litter(:,:,:) = zero 
    368380    cn_ind(:,:) = zero 
     381    woodmass_ind(:,:) = zero 
    369382    veget_max_old(:,:) = veget_max(:,:) 
    370383 
    371     ! 
    372     ! 1.3 Prescribe some vegetation characteristics if the vegetation is not dynamic 
     384    ! 1.3 Calculate some vegetation characteristics 
     385 
     386    ! 
     387    ! 1.3.1 Calculate some vegetation characteristics (cn_ind and height) from 
     388    !     state variables if running DGVM or dynamic mortality in static cover mode 
     389    ! 
     390    IF ( control%ok_dgvm .OR. .NOT.lpj_gap_const_mort) THEN 
     391       IF(control%ok_dgvm) THEN 
     392          WHERE (ind(:,:).GT.min_stomate) 
     393             woodmass_ind(:,:) = & 
     394                  ((biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     395                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow)) &  
     396                  *veget_max(:,:))/ind(:,:) 
     397          ENDWHERE 
     398       ELSE 
     399          WHERE (ind(:,:).GT.min_stomate) 
     400             woodmass_ind(:,:) = & 
     401                  (biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     402                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow))/ind(:,:) 
     403          ENDWHERE 
     404       ENDIF 
     405 
     406       CALL crown (npts,  PFTpresent, & 
     407            ind, biomass, woodmass_ind, & 
     408            veget_max, cn_ind, height) 
     409    ENDIF 
     410 
     411    ! 
     412    ! 1.3.2 Prescribe some vegetation characteristics if the vegetation is not dynamic 
    373413    !     IF the DGVM is not activated, the density of individuals and their crown 
    374414    !     areas don't matter, but they should be defined for the case we switch on 
    … …  
    389429 
    390430    CALL constraints (npts, dt_days, & 
    391          t2m_month, t2m_min_daily, when_growthinit, & 
     431         t2m_month, t2m_min_daily,when_growthinit, & 
    392432         adapted, regenerate) 
    393433 
    … …  
    404444       CALL pftinout (npts, dt_days, adapted, regenerate, & 
    405445            neighbours, veget, veget_max, & 
    406             biomass, ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, & 
     446            biomass, ind, cn_ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, & 
    407447            PFTpresent, everywhere, when_growthinit, need_adjacent, RIP_time, & 
    408448            co2_to_bm, & 
    … …  
    417457       CALL kill (npts, 'pftinout  ', lm_lastyearmax, & 
    418458            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
    419             lai, age, leaf_age, leaf_frac, & 
     459            lai, age, leaf_age, leaf_frac, npp_longterm, & 
    420460            when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    421461 
    … …  
    423463       ! 3.3 calculate new crown area and maximum vegetation cover 
    424464       ! 
     465       ! 
     466       ! unsure whether this is really required 
     467       ! - in theory this could ONLY be done at the END of stomate_lpj 
     468       ! 
     469 
     470       ! calculate woodmass of individual tree 
     471       WHERE ((ind(:,:).GT.min_stomate)) 
     472          WHERE  ( veget_max(:,:) .GT. min_stomate) 
     473             woodmass_ind(:,:) = & 
     474                  ((biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     475                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow))*veget_max(:,:))/ind(:,:) 
     476          ELSEWHERE 
     477             woodmass_ind(:,:) =(biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     478                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow))/ind(:,:) 
     479          ENDWHERE 
     480 
     481       ENDWHERE 
    425482 
    426483       CALL crown (npts, PFTpresent, & 
    427             ind, biomass, & 
     484            ind, biomass, woodmass_ind, & 
    428485            veget_max, cn_ind, height) 
    429486 
    … …  
    487544         resp_maint, resp_growth, npp_daily) 
    488545 
    489     IF ( control%ok_dgvm ) THEN 
     546    IF ( control%ok_dgvm .OR. .NOT.lpj_gap_const_mort) THEN 
     547       CALL kill (npts, 'npp       ', lm_lastyearmax,  & 
     548            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
     549            lai, age, leaf_age, leaf_frac, npp_longterm, & 
     550            when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    490551 
    491552       ! new provisional crown area and maximum vegetation cover after growth 
     553       IF(control%ok_dgvm) THEN 
     554          WHERE (ind(:,:).GT.min_stomate) 
     555             woodmass_ind(:,:) = & 
     556                  ((biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     557                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow)) &  
     558                  *veget_max(:,:))/ind(:,:) 
     559          ENDWHERE 
     560       ELSE 
     561          WHERE (ind(:,:).GT.min_stomate) 
     562             woodmass_ind(:,:) = & 
     563                  (biomass(:,:,isapabove)+biomass(:,:,isapbelow) & 
     564                  +biomass(:,:,iheartabove)+biomass(:,:,iheartbelow))/ind(:,:) 
     565          ENDWHERE 
     566       ENDIF 
    492567 
    493568       CALL crown (npts, PFTpresent, & 
    494             ind, biomass, & 
     569            ind, biomass, woodmass_ind,& 
    495570            veget_max, cn_ind, height) 
    496571 
    … …  
    513588       CALL kill (npts, 'fire      ', lm_lastyearmax, & 
    514589            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
    515             lai, age, leaf_age, leaf_frac, & 
     590            lai, age, leaf_age, leaf_frac, npp_longterm, & 
    516591            when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    517592 
    … …  
    524599    CALL gap (npts, dt_days, & 
    525600         npp_longterm, turnover_longterm, lm_lastyearmax, & 
    526          PFTpresent, biomass, ind, bm_to_litter) 
     601         PFTpresent, biomass, ind, bm_to_litter, mortality) 
    527602 
    528603    IF ( control%ok_dgvm ) THEN 
    … …  
    532607       CALL kill (npts, 'gap       ', lm_lastyearmax, & 
    533608            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
    534             lai, age, leaf_age, leaf_frac, & 
     609            lai, age, leaf_age, leaf_frac, npp_longterm, & 
    535610            when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    536611 
    … …  
    570645 
    571646       CALL light (npts, dt_days, & 
    572             PFTpresent, cn_ind, lai, maxfpc_lastyear, & 
    573             ind, biomass, veget_lastlight, bm_to_litter) 
     647            veget_max, fpc_max, PFTpresent, cn_ind, lai, maxfpc_lastyear, & 
     648            lm_lastyearmax, ind, biomass, veget_lastlight, bm_to_litter, mortality) 
    574649 
    575650       ! 
    … …  
    579654       CALL kill (npts, 'light     ', lm_lastyearmax, & 
    580655            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, & 
    581             lai, age, leaf_age, leaf_frac, & 
     656            lai, age, leaf_age, leaf_frac, npp_longterm, & 
    582657            when_growthinit, everywhere, veget, veget_max, bm_to_litter) 
    583658 
    … …  
    588663    ! 
    589664 
    590     IF ( control%ok_dgvm ) THEN 
     665    IF ( control%ok_dgvm .OR. .NOT.lpj_gap_const_mort ) THEN 
    591666 
    592667       ! 
    … …  
    597672            neighbours, resolution, need_adjacent, herbivores, & 
    598673            precip_lastyear, gdd0_lastyear, lm_lastyearmax, & 
    599             cn_ind, lai, avail_tree, avail_grass, & 
     674            cn_ind, lai, avail_tree, avail_grass, npp_longterm, & 
    600675            leaf_age, leaf_frac, & 
    601             ind, biomass, age, everywhere, co2_to_bm, veget_max) 
     676            ind, biomass, age, everywhere, co2_to_bm, veget_max, woodmass_ind) 
    602677 
    603678       ! 
    … …  
    606681 
    607682       CALL crown (npts, PFTpresent, & 
    608             ind, biomass, & 
     683            ind, biomass, woodmass_ind, & 
    609684            veget_max, cn_ind, height) 
    610685 
    … …  
    617692    CALL cover (npts, cn_ind, ind, biomass, & 
    618693         veget_max, veget_max_old, veget, & 
    619          lai, litter, carbon) 
     694         lai, litter, carbon, turnover_daily, bm_to_litter) 
    620695 
    621696    ! 
    … …  
    647722       ENDIF 
    648723    ENDIF 
    649 !MM déplacement pour initialisation correcte des grandeurs cumulées : 
     724    !MM déplacement pour initialisation correcte des grandeurs cumulées : 
    650725    cflux_prod_total(:) = convflux(:) + cflux_prod10(:) + cflux_prod100(:) 
    651726    prod10_total(:)=SUM(prod10,dim=2) 
    … …  
    681756         &             bm_to_litter(:,:,iheartabove) + bm_to_litter(:,:,iroot) + & 
    682757         &             bm_to_litter(:,:,ifruit) + bm_to_litter(:,:,icarbres) 
     758 
     759    carb_mass_variation(:)=-carb_mass_total(:) 
     760    carb_mass_total(:)=SUM((tot_live_biomass+tot_litter_carb+tot_soil_carb)*veget_max,dim=2) + & 
     761         &                 (prod10_total + prod100_total) 
     762    carb_mass_variation(:)=carb_mass_total(:)+carb_mass_variation(:) 
    683763 
    684764    ! 
    … …  
    759839    CALL histwrite (hist_id_stomate, 'IND', itime, & 
    760840         ind, npts*nvm, horipft_index) 
     841    CALL histwrite (hist_id_stomate, 'CN_IND', itime, & 
     842         cn_ind, npts*nvm, horipft_index) 
     843    CALL histwrite (hist_id_stomate, 'WOODMASS_IND', itime, & 
     844         woodmass_ind, npts*nvm, horipft_index) 
    761845    CALL histwrite (hist_id_stomate, 'TOTAL_M', itime, & 
    762846         tot_live_biomass, npts*nvm, horipft_index) 
    … …  
    832916       vartmp(:)=SUM(tot_live_biomass*veget_max,dim=2)/1e3*contfrac 
    833917       CALL histwrite (hist_id_stomate_IPCC, "cVeg", itime, & 
    834          vartmp, npts, hori_index) 
     918            vartmp, npts, hori_index) 
    835919       vartmp(:)=SUM(tot_litter_carb*veget_max,dim=2)/1e3*contfrac 
    836920       CALL histwrite (hist_id_stomate_IPCC, "cLitter", itime, & 
    837          vartmp, npts, hori_index) 
     921            vartmp, npts, hori_index) 
    838922       vartmp(:)=SUM(tot_soil_carb*veget_max,dim=2)/1e3*contfrac 
    839923       CALL histwrite (hist_id_stomate_IPCC, "cSoil", itime, & 
    840          vartmp, npts, hori_index) 
     924            vartmp, npts, hori_index) 
    841925       vartmp(:)=(prod10_total + prod100_total)/1e3 
    842926       CALL histwrite (hist_id_stomate_IPCC, "cProduct", itime, & 
    843          vartmp, npts, hori_index) 
     927            vartmp, npts, hori_index) 
     928       vartmp(:)=carb_mass_variation/1e3/one_day*contfrac 
     929       CALL histwrite (hist_id_stomate_IPCC, "cMassVariation", itime, & 
     930            vartmp, npts, hori_index) 
     931 
    844932       vartmp(:)=SUM(lai*veget_max,dim=2)*contfrac 
    845933       CALL histwrite (hist_id_stomate_IPCC, "lai", itime, & 
    846          vartmp, npts, hori_index) 
     934            vartmp, npts, hori_index) 
    847935       vartmp(:)=SUM(gpp_daily*veget_max,dim=2)/1e3/one_day*contfrac 
    848936       CALL histwrite (hist_id_stomate_IPCC, "gpp", itime, & 
    849          vartmp, npts, hori_index) 
     937            vartmp, npts, hori_index) 
    850938       vartmp(:)=SUM((resp_maint+resp_growth)*veget_max,dim=2)/1e3/one_day*contfrac 
    851939       CALL histwrite (hist_id_stomate_IPCC, "ra", itime, & 
    852          vartmp, npts, hori_index) 
     940            vartmp, npts, hori_index) 
    853941       vartmp(:)=SUM(npp_daily*veget_max,dim=2)/1e3/one_day*contfrac 
    854942       CALL histwrite (hist_id_stomate_IPCC, "npp", itime, & 
    855          vartmp, npts, hori_index) 
     943            vartmp, npts, hori_index) 
    856944       vartmp(:)=SUM(resp_hetero*veget_max,dim=2)/1e3/one_day*contfrac 
    857945       CALL histwrite (hist_id_stomate_IPCC, "rh", itime, & 
    858          vartmp, npts, hori_index) 
     946            vartmp, npts, hori_index) 
    859947       vartmp(:)=SUM(co2_fire*veget_max,dim=2)/1e3/one_day*contfrac 
    860948       CALL histwrite (hist_id_stomate_IPCC, "fFire", itime, & 
    861          vartmp, npts, hori_index) 
     949            vartmp, npts, hori_index) 
    862950       vartmp(:)=harvest_above/1e3/one_day*contfrac 
    863951       CALL histwrite (hist_id_stomate_IPCC, "fHarvest", itime, & 
    864          vartmp, npts, hori_index) 
     952            vartmp, npts, hori_index) 
    865953       vartmp(:)=cflux_prod_total/1e3/one_day*contfrac 
    866954       CALL histwrite (hist_id_stomate_IPCC, "fLuc", itime, & 
    867          vartmp, npts, hori_index) 
     955            vartmp, npts, hori_index) 
    868956       vartmp(:)=(SUM((gpp_daily-(resp_maint+resp_growth+resp_hetero)-co2_fire) & 
    869957            &        *veget_max,dim=2)-cflux_prod_total-harvest_above)/1e3/one_day*contfrac 
    870958       CALL histwrite (hist_id_stomate_IPCC, "nbp", itime, & 
    871          vartmp, npts, hori_index) 
     959            vartmp, npts, hori_index) 
    872960       vartmp(:)=SUM(tot_bm_to_litter*veget_max,dim=2)/1e3/one_day*contfrac 
    873961       CALL histwrite (hist_id_stomate_IPCC, "fVegLitter", itime, & 
    874          vartmp, npts, hori_index) 
     962            vartmp, npts, hori_index) 
    875963       vartmp(:)=SUM(SUM(soilcarbon_input,dim=2)*veget_max,dim=2)/1e3/one_day*contfrac 
    876964       CALL histwrite (hist_id_stomate_IPCC, "fLitterSoil", itime, & 
    877          vartmp, npts, hori_index) 
     965            vartmp, npts, hori_index) 
    878966       vartmp(:)=SUM(biomass(:,:,ileaf)*veget_max,dim=2)/1e3*contfrac 
    879967       CALL histwrite (hist_id_stomate_IPCC, "cLeaf", itime, & 
    880          vartmp, npts, hori_index) 
     968            vartmp, npts, hori_index) 
    881969       vartmp(:)=SUM((biomass(:,:,isapabove)+biomass(:,:,iheartabove))*veget_max,dim=2)/1e3*contfrac 
    882970       CALL histwrite (hist_id_stomate_IPCC, "cWood", itime, & 
    883          vartmp, npts, hori_index) 
     971            vartmp, npts, hori_index) 
    884972       vartmp(:)=SUM(( biomass(:,:,iroot) + biomass(:,:,isapbelow) + biomass(:,:,iheartbelow) ) & 
    885973            &        *veget_max,dim=2)/1e3*contfrac 
    886974       CALL histwrite (hist_id_stomate_IPCC, "cRoot", itime, & 
    887          vartmp, npts, hori_index) 
     975            vartmp, npts, hori_index) 
    888976       vartmp(:)=SUM(( biomass(:,:,icarbres) + biomass(:,:,ifruit))*veget_max,dim=2)/1e3*contfrac 
    889977       CALL histwrite (hist_id_stomate_IPCC, "cMisc", itime, & 
    890          vartmp, npts, hori_index) 
     978            vartmp, npts, hori_index) 
    891979       vartmp(:)=SUM((litter(:,istructural,:,iabove)+litter(:,imetabolic,:,iabove))*veget_max,dim=2)/1e3*contfrac 
    892980       CALL histwrite (hist_id_stomate_IPCC, "cLitterAbove", itime, & 
    893          vartmp, npts, hori_index) 
     981            vartmp, npts, hori_index) 
    894982       vartmp(:)=SUM((litter(:,istructural,:,ibelow)+litter(:,imetabolic,:,ibelow))*veget_max,dim=2)/1e3*contfrac 
    895983       CALL histwrite (hist_id_stomate_IPCC, "cLitterBelow", itime, & 
    896          vartmp, npts, hori_index) 
     984            vartmp, npts, hori_index) 
    897985       vartmp(:)=SUM(carbon(:,iactive,:)*veget_max,dim=2)/1e3*contfrac 
    898986       CALL histwrite (hist_id_stomate_IPCC, "cSoilFast", itime, & 
    899          vartmp, npts, hori_index) 
     987            vartmp, npts, hori_index) 
    900988       vartmp(:)=SUM(carbon(:,islow,:)*veget_max,dim=2)/1e3*contfrac 
    901989       CALL histwrite (hist_id_stomate_IPCC, "cSoilMedium", itime, & 
    902          vartmp, npts, hori_index) 
     990            vartmp, npts, hori_index) 
    903991       vartmp(:)=SUM(carbon(:,ipassive,:)*veget_max,dim=2)/1e3*contfrac 
    904992       CALL histwrite (hist_id_stomate_IPCC, "cSoilSlow", itime, & 
    905          vartmp, npts, hori_index) 
     993            vartmp, npts, hori_index) 
    906994       DO j=1,nvm 
    907995          histvar(:,j)=veget_max(:,j)*contfrac(:)*100 
    908996       ENDDO 
    909997       CALL histwrite (hist_id_stomate_IPCC, "landCoverFrac", itime, & 
    910          histvar, npts*nvm, horipft_index) 
     998            histvar, npts*nvm, horipft_index) 
    911999       vartmp(:)=(veget_max(:,3)+veget_max(:,6)+veget_max(:,8)+veget_max(:,9))*contfrac*100 
    9121000       CALL histwrite (hist_id_stomate_IPCC, "treeFracPrimDec", itime, & 
    913           vartmp, npts, hori_index) 
     1001            vartmp, npts, hori_index) 
    9141002       vartmp(:)=(veget_max(:,2)+veget_max(:,4)+veget_max(:,5)+veget_max(:,7))*contfrac*100 
    9151003       CALL histwrite (hist_id_stomate_IPCC, "treeFracPrimEver", itime, & 
    916          vartmp, npts, hori_index) 
     1004            vartmp, npts, hori_index) 
    9171005       vartmp(:)=(veget_max(:,10)+veget_max(:,12))*contfrac*100 
    9181006       CALL histwrite (hist_id_stomate_IPCC, "c3PftFrac", itime, & 
    919          vartmp, npts, hori_index) 
     1007            vartmp, npts, hori_index) 
    9201008       vartmp(:)=(veget_max(:,11)+veget_max(:,13))*contfrac*100 
    9211009       CALL histwrite (hist_id_stomate_IPCC, "c4PftFrac", itime, & 
    922          vartmp, npts, hori_index) 
     1010            vartmp, npts, hori_index) 
    9231011       vartmp(:)=SUM(resp_growth*veget_max,dim=2)/1e3/one_day*contfrac 
    9241012       CALL histwrite (hist_id_stomate_IPCC, "rGrowth", itime, & 
    925          vartmp, npts, hori_index) 
     1013            vartmp, npts, hori_index) 
    9261014       vartmp(:)=SUM(resp_maint*veget_max,dim=2)/1e3/one_day*contfrac 
    9271015       CALL histwrite (hist_id_stomate_IPCC, "rMaint", itime, & 
    928          vartmp, npts, hori_index) 
     1016            vartmp, npts, hori_index) 
    9291017       vartmp(:)=SUM(bm_alloc(:,:,ileaf)*veget_max,dim=2)/1e3/one_day*contfrac 
    9301018       CALL histwrite (hist_id_stomate_IPCC, "nppLeaf", itime, & 
    931          vartmp, npts, hori_index) 
     1019            vartmp, npts, hori_index) 
    9321020       vartmp(:)=SUM(bm_alloc(:,:,isapabove)*veget_max,dim=2)/1e3/one_day*contfrac 
    9331021       CALL histwrite (hist_id_stomate_IPCC, "nppWood", itime, & 
    934          vartmp, npts, hori_index) 
     1022            vartmp, npts, hori_index) 
    9351023       vartmp(:)=SUM(( bm_alloc(:,:,isapbelow) + bm_alloc(:,:,iroot) )*veget_max,dim=2)/1e3/one_day*contfrac 
    9361024       CALL histwrite (hist_id_stomate_IPCC, "nppRoot", itime, & 
    937          vartmp, npts, hori_index) 
     1025            vartmp, npts, hori_index) 
    9381026 
    9391027       CALL histwrite (hist_id_stomate_IPCC, 'RESOLUTION_X', itime, & 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate_prescribe.f90

    r119 r405  
    8989      ! only when the DGVM is not activated or agricultural PFT. 
    9090 
    91       IF ( ( .NOT. control%ok_dgvm ) .OR. ( .NOT. natural(j) ) ) THEN 
     91      IF ( ( .NOT. control%ok_dgvm .AND. lpj_gap_const_mort ) .OR. ( .NOT. natural(j) ) ) THEN 
    9292 
    9393        ! 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/stomate_season.f90

    r119 r405  
    163163    ! rapport maximal GPP/GGP_max pour dormance 
    164164    REAL(r_std), PARAMETER                                  :: gppfrac_dormance = 0.2 
    165 ! 
    166 !NVADD 
    167      ! minimum gpp considered as not "lowgpp" 
     165    ! 
     166    !NVADD 
     167    ! minimum gpp considered as not "lowgpp" 
    168168    REAL(r_std), PARAMETER                                  :: min_gpp_allowed = 0.3 
    169      ! tau (year) for "climatologic variables 
     169    ! tau (year) for "climatologic variables 
    170170    REAL(r_std), PARAMETER                                  :: tau_climatology = 20 
    171 !ENDNVADD 
     171    !ENDNVADD 
    172172    ! maximum ncd (d) (to avoid floating point underflows) 
    173173    REAL(r_std)                                             :: ncd_max  
    … …  
    186186    ! herbivore consumption (gC/m**2/day) 
    187187    REAL(r_std), DIMENSION(npts)                            :: consumption 
     188    ! fraction of each gridcell occupied by natural vegetation 
     189    REAL(r_std), DIMENSION(npts)                            :: fracnat 
    188190 
    189191    ! ========================================================================= 
    … …  
    226228 
    227229       ! 1.2.1.1 "monthly" 
    228 !MM PAS PARALLELISE!! 
     230       !MM PAS PARALLELISE!! 
    229231       IF ( ABS( SUM( moiavail_month(:,2:nvm) ) ) .LT. min_stomate ) THEN 
    230232 
    … …  
    278280 
    279281       ! 1.2.3 "monthly" soil temperatures 
    280 !MM PAS PARALLELISE!! 
     282       !MM PAS PARALLELISE!! 
    281283       IF ( ABS( SUM( tsoil_month(:,:) ) ) .LT. min_stomate ) THEN 
    282284 
    … …  
    465467    !         detect a beginning of the growing season by declaring it dormant 
    466468    ! 
    467 !NVMODIF 
     469    !NVMODIF 
    468470    DO j = 2,nvm 
    469471       WHERE ( ( gpp_week(:,j) .LT. min_gpp_allowed ) .OR. &  
    … …  
    471473            ( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. & 
    472474            ( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) ) 
    473 !       WHERE ( ( gpp_week(:,j) .EQ. zero ) .OR. &  
    474 !            ( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. & 
    475 !            ( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. & 
    476 !            ( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) ) 
    477            
     475          !       WHERE ( ( gpp_week(:,j) .EQ. zero ) .OR. &  
     476          !            ( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. & 
     477          !            ( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. & 
     478          !            ( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) ) 
     479 
    478480          time_lowgpp(:,j) = time_lowgpp(:,j) + dt 
    479            
     481 
    480482       ELSEWHERE 
    481            
     483 
    482484          time_lowgpp(:,j) = zero 
    483485 
    … …  
    817819    ! 
    818820 
     821    IF(control%ok_dgvm ) THEN 
     822 
     823       fracnat(:) = un 
     824       DO j = 2,nvm 
     825          IF ( .NOT. natural(j) ) THEN 
     826             fracnat(:) = fracnat(:) - veget_max(:,j) 
     827          ENDIF 
     828       ENDDO 
     829 
     830    ENDIF 
     831 
    819832    IF ( control%ok_stomate ) THEN 
    820  
    821        DO j = 2,nvm 
    822           WHERE ( biomass(:,j,ileaf) .GT. lm_thisyearmax(:,j) ) 
    823              lm_thisyearmax(:,j) = biomass(:,j,ileaf) 
    824           ENDWHERE 
    825        ENDDO 
    826  
     833       IF(control%ok_dgvm ) THEN 
     834          DO j=2,nvm 
     835 
     836             IF ( natural(j) .AND. control%ok_dgvm ) THEN 
     837 
     838                WHERE ( fracnat(:) .GT. min_stomate .AND. biomass(:,j,ileaf).GT. lm_lastyearmax(:,j)*0.75 ) 
     839                   maxfpc_lastyear(:,j) = ( maxfpc_lastyear(:,j) * ( one_year/leaflife_tab(j)- dt ) + & 
     840                        veget(:,j) / fracnat(:) * dt ) / (one_year/leaflife_tab(j)) 
     841                ENDWHERE 
     842                maxfpc_thisyear(:,j) = maxfpc_lastyear(:,j) ! just to initialise value 
     843 
     844             ENDIF 
     845 
     846!NV : correct initialization 
     847!!$             WHERE(biomass(:,j,ileaf).GT. lm_lastyearmax(:,j)*0.75) 
     848!!$                lm_lastyearmax(:,j) = ( lm_lastyearmax(:,j) * ( one_year/leaflife_tab(j)- dt ) + & 
     849!!$                     biomass(:,j,ileaf) * dt ) / (one_year/leaflife_tab(j)) 
     850!!$             ENDWHERE 
     851!!$             lm_thisyearmax(:,j)=lm_lastyearmax(:,j) ! just to initialise value 
     852             WHERE (lm_thisyearmax(:,j) .GT. min_stomate) 
     853                WHERE(biomass(:,j,ileaf).GT. lm_thisyearmax(:,j)*0.75) 
     854                   lm_thisyearmax(:,j) = ( lm_thisyearmax(:,j) * ( one_year/leaflife_tab(j)- dt ) + & 
     855                        biomass(:,j,ileaf) * dt ) / (one_year/leaflife_tab(j)) 
     856                ENDWHERE 
     857             ELSEWHERE 
     858                lm_thisyearmax(:,j) =biomass(:,j,ileaf) 
     859             ENDWHERE 
     860 
     861          ENDDO 
     862 
     863       ELSE 
     864 
     865          DO j = 2,nvm 
     866             WHERE ( biomass(:,j,ileaf) .GT. lm_thisyearmax(:,j) ) 
     867                lm_thisyearmax(:,j) = biomass(:,j,ileaf) 
     868             ENDWHERE 
     869          ENDDO 
     870 
     871       ENDIF 
    827872    ELSE 
    828873 
    … …  
    852897       ! 21.1 replace old values 
    853898       ! 
    854 !NVMODIF 
    855       maxmoiavail_lastyear(:,:) = (maxmoiavail_lastyear(:,:)*(tau_climatology-1)+ maxmoiavail_thisyear(:,:))/tau_climatology 
    856       minmoiavail_lastyear(:,:) = (minmoiavail_lastyear(:,:)*(tau_climatology-1)+ minmoiavail_thisyear(:,:))/tau_climatology 
    857       maxgppweek_lastyear(:,:) =( maxgppweek_lastyear(:,:)*(tau_climatology-1)+ maxgppweek_thisyear(:,:))/tau_climatology 
    858 !       maxmoiavail_lastyear(:,:) = maxmoiavail_thisyear(:,:) 
    859 !       minmoiavail_lastyear(:,:) = minmoiavail_thisyear(:,:) 
    860 !       maxgppweek_lastyear(:,:) = maxgppweek_thisyear(:,:) 
     899       !NVMODIF 
     900       maxmoiavail_lastyear(:,:) = (maxmoiavail_lastyear(:,:)*(tau_climatology-1)+ maxmoiavail_thisyear(:,:))/tau_climatology 
     901       minmoiavail_lastyear(:,:) = (minmoiavail_lastyear(:,:)*(tau_climatology-1)+ minmoiavail_thisyear(:,:))/tau_climatology 
     902       maxgppweek_lastyear(:,:) =( maxgppweek_lastyear(:,:)*(tau_climatology-1)+ maxgppweek_thisyear(:,:))/tau_climatology 
     903       !       maxmoiavail_lastyear(:,:) = maxmoiavail_thisyear(:,:) 
     904       !       minmoiavail_lastyear(:,:) = minmoiavail_thisyear(:,:) 
     905       !       maxgppweek_lastyear(:,:) = maxgppweek_thisyear(:,:) 
    861906 
    862907       gdd0_lastyear(:) = gdd0_thisyear(:) 
    … …  
    910955       !        fpc_crit. 
    911956 
    912        ! calculate the sum of maxfpc_lastyear 
    913        sumfpc_nat(:) = zero 
    914        DO j = 2,nvm 
    915           sumfpc_nat(:) = sumfpc_nat(:) + maxfpc_lastyear(:,j) 
    916        ENDDO 
    917  
    918        ! scale so that the new sum is fpc_crit 
    919        DO j = 2,nvm  
    920           WHERE ( sumfpc_nat(:) .GT. fpc_crit ) 
    921              maxfpc_lastyear(:,j) = maxfpc_lastyear(:,j) * (fpc_crit/sumfpc_nat(:)) 
    922           ENDWHERE 
    923        ENDDO 
     957!!$       ! calculate the sum of maxfpc_lastyear 
     958!!$       sumfpc_nat(:) = zero 
     959!!$       DO j = 2,nvm 
     960!!$          sumfpc_nat(:) = sumfpc_nat(:) + maxfpc_lastyear(:,j) 
     961!!$       ENDDO 
     962!!$ 
     963!!$       ! scale so that the new sum is fpc_crit 
     964!!$       DO j = 2,nvm  
     965!!$          WHERE ( sumfpc_nat(:) .GT. fpc_crit ) 
     966!!$             maxfpc_lastyear(:,j) = maxfpc_lastyear(:,j) * (fpc_crit/sumfpc_nat(:)) 
     967!!$          ENDWHERE 
     968!!$       ENDDO 
    924969 
    925970    ENDIF  ! EndOfYear 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE_OL/FLUXNET/Job_FLUXNET_validation

    r119 r405  
    212212 
    213213    IGCM_card_WriteOption ${New_SUBMIT_DIR}/COMP/spinup.card UserChoices DRIVER_NORESTART y 
     214    IGCM_card_WriteOption ${New_SUBMIT_DIR}/COMP/spinup.card UserChoices DRIVER_TIMELENGTH n 
    214215 
    215216    eval IGCM_card_WriteOption ${New_SUBMIT_DIR}/COMP/spinup.card UserChoices duree_nostomate $( correct_duree ${fluxnet_SPINUP_duree_nostomate} ${TIME_YEAR} ) 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE_OL/FLUXNET/PARAM/sechiba.def

    r119 r405  
    503503 
    504504# Total depth of soil reservoir 
    505 HYDROL_SOIL_DEPTH = 2. 
     505HYDROL_SOIL_DEPTH = 4. 
    506506# default = 2. 
    507507 
    … …  
    510510# For 4 meters soil depth, you may use those ones : 
    511511# 5., .4, .4, 1., .8, .8, 1., 1., .8, 4., 1., 4., 1. 
    512 HYDROL_HUMCSTE = 5., .8, .8, 1., .8, .8, 1., 1., .8, 4., 4., 4., 4. 
     512HYDROL_HUMCSTE = 5., .4, .4, 1., .8, .8, 1., 1., .8, 4., 1., 4., 1. 
    513513# default =  5., .8, .8, 1., .8, .8, 1., 1., .8, 4., 4., 4., 4. 
    514514 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE_OL/FLUXNET/fluxnet.card

    r119 r405  
    11[FLUXNET] 
    22# - Fluxnet files path 
    3 FluxnetPath=${R_BC}/OOL/${config_UserChoices_TagName}/FLUXNET/BC 
     3FluxnetPath=/home/orchidee01/vuichard/Input_Fluxnet 
    44 
    55# - Number of PFTs : 
    66NbPFTs= 13 
    77 
    8 # - Information on the sites to be treated : 
    9 #   * Number of physical parameters on each sites per PFTs 
    10 NbSitesParam= 2 
    11 # 4 first parameters are Name, Forcing file, Begin date, Number of years in forcing file 
    12  
    13 #   * ORCHIDEE name for physical parameters on each sites 
    14 #   PFT (IMPOSE_VEG), \ 
    15 #   initial LAI (IMPOSE_VEG) 
    16 NameSitesParam= ( SECHIBA_VEGMAX, SECHIBA_LAI ) 
    17 # by Default :  
    18 # 1) first line is for PFT 
    19 # 2) second line is for LAI default for SLOWPROC lai model with :  
    20 #    llaimax = 0.,  8.,  8.,  4., 4.5, 4.5,  4., 4.5,  4.,  2.,  2.,  2.,  2.) 
    21  
    22 #   * Name of component for each physical parameter described in NameSitesParam 
    23 #     (in SECHIBA, STOMATE, DRIVER) 
    24 CompSitesParam= ( SECHIBA, SECHIBA ) 
    25  
     8#**** Information on the sites to be treated ************************* 
     9#  Number of parameters to modify for each site 
     10NbSitesParam= 1 
     11#  Name of the parameters to modify for on each site 
     12NameSitesParam= ( SECHIBA_VEGMAX ) 
     13#  Name of the component for each parameter described in NameSitesParam (either, SECHIBA, STOMATE, or DRIVER) 
     14CompSitesParam= ( SECHIBA ) 
    2615 
    2716# Sites descriptions 
    28 #       Abbrv,  Filename ,      Inital year (for gregorian calendar) , Length (Y),  \ 
    29 #param 1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,  13 
    30 Sites= ( GU, GU.nc   ,          1996,                         3     , \ 
    31      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.7, 0.0, 0.3, 0.0, 0.0, 0.0, \ 
    32      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0) \ 
    33 \ 
    34        ( FL, FL.nc   ,          1996,                         3     , \ 
    35      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    36      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    37 \ 
    38        ( HY, HY.nc   ,          1996,                         5     , \ 
    39      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.9, 0.0, 0.0, 0.1, 0.0, 0.0, 0.0, \ 
    40      0.0, 8.0, 8.0, 4.0, 4.5, 4.5, 3.0, 2.5, 4.0, 3.0, 2.0, 2.0, 2.0) \ 
    41 \ 
    42        ( NB, NB.nc   ,          1994,                         5     , \ 
    43      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.7, 0.3, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    44      0.0, 0.0, 0.0, 0.0, 0.0, 4.5, 4.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    45 \ 
    46        ( NO, NO.nc   ,          1996,                         3     , \ 
    47      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    48      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 4.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    49 \ 
    50        ( HV, HV.nc   ,          1992,                         8     , \ 
    51      0.0, 0.0, 0.0, 0.3, 0.0, 0.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    52      0.0, 0.0, 0.0, 2.8, 0.0, 2.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    53 \ 
    54        ( SO, SO.nc   ,          1997,                         4     , \ 
    55      0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    56      0.0, 0.0, 0.0, 0.0, 0.0, 2.9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    57 \ 
    58        ( VI, VI.nc   ,          1996,                         3     , \ 
    59      0.0, 0.0, 0.0, 0.3, 0.0, 0.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    60      0.0, 0.0, 0.0, 2.5, 0.0, 2.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    61 \ 
    62        ( WB, WB.nc   ,          1995,                         3     , \ 
    63      0.0, 0.0, 0.0, 0.2, 0.0, 0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    64      0.0, 0.0, 0.0, 3.0, 0.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    65 \ 
    66        ( AB, AB.nc   ,          1997,                         3     , \ 
    67      0.0, 0.0, 0.0, 0.9, 0.0, 0.1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    68      0.0, 0.0, 0.0, 7.5, 0.0, 0.0, 0.0, 0.0, 0.0, 7.5, 0.0, 0.0, 0.0) \ 
    69 \ 
    70        ( BR, BR.nc   ,          1996,                         4     , \ 
    71      0.0, 0.0, 0.0, 0.6, 0.0, 0.4, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    72      0.0, 0.0, 0.0, 2.5, 0.0, 2.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    73 \ 
    74        ( LO, LO.nc   ,          1996,                         5     , \ 
    75      0.0, 0.0, 0.0, 0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.2, 0.0, 0.0, 0.0, \ 
    76      1.5, 8.0, 8.0, 3.0, 1.6, 5.5, 3.0, 2.5, 4.0, 3.2, 2.9, 5.0, 2.0) \ 
    77 \ 
    78        ( ME, ME.nc   ,          1996,                         2     , \ 
    79      0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    80      0.0, 0.0, 0.0, 3.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    81 \ 
    82        ( TH, TH.nc   ,          1996,                         5     , \ 
    83      0.0, 0.0, 0.0, 0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.2, 0.0, 0.0, 0.0, \ 
    84      1.5, 8.0, 8.0, 6.0, 1.6, 5.5, 3.0, 2.5, 4.0, 6.0, 2.9, 5.0, 2.0) \ 
    85 \ 
    86        ( WE, WE.nc   ,          1996,                         4     , \ 
    87      0.0, 0.0, 0.0, 0.8, 0.0, 0.0, 0.0, 0.0, 0.0, 0.2, 0.0, 0.0, 0.0, \ 
    88      0.0, 0.0, 0.0, 6.0, 0.0, 0.0, 0.0, 0.0, 0.0, 6.0, 0.0, 0.0, 0.0) \ 
    89 \ 
    90        ( MA, MA.nc   ,          1996,                         1     , \ 
    91      0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    92      0.0, 5.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    93 \ 
    94        ( LW, LW.nc   ,          1997,                         2     , \ 
    95      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, \ 
    96      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.5, 0.0, 0.0, 0.0) \ 
    97 \ 
    98        ( SH, SH.nc   ,          1997,                         1     , \ 
    99      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, \ 
    100      0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0) \ 
     17# 4 first parameters are Name, Forcing file, Initial Year, Number of years in forcing file 
     18# following parameters are NameSitesParam 
     19Sites= ( NL-Loo, NL-Loo.nc,  1996,                                 11, \ 
     20     0,   0,   0, 0.8,   0,   0,   0,   0,   0, 0.2,   0,   0,   0.) \ 
     21 \ 
     22       ( DE-Hai, DE-Hai.nc,  2000,                                  7, \ 
     23     0,   0,   0,   0,    0, 0.8,   0,   0,   0, 0.2,   0,   0,   0) \ 
     24 \ 
     25       ( BW-Ma1, BW-Ma1.nc,  1999,                                  3, \ 
     26    0.1, 0.2,   0,   0,   0,   0,   0,   0,   0, 0.7,   0,   0,   0) \ 
     27 \ 
     28       ( FI-Sod, FI-Sod.nc,  2000,                                  7, \ 
     29      0,   0,   0,   0,   0,   0, 0.8,   0,   0, 0.2,   0,   0,   0) \ 
     30 \ 
     31       ( BR-Sa1, BR-Sa1.nc,  2002,                                  3, \ 
     32      0, 0.8,   0,   0,   0,   0,   0,   0,   0, 0.2,   0,   0,   0) \ 
     33 \ 
     34       ( RU-Zot, RU-Zot.nc,  2002,                                  3, \ 
     35     0,   0,   0,   0,   0,   0, 0.8,   0,   0, 0.2,   0,   0,   0) \ 
     36 \ 
     37        ( BR-Ma2, BR-Ma2.nc,  2002,                                  5, \ 
     38      0, 0.8,   0,   0,   0,   0,   0,   0,   0, 0.2,   0,   0,   0) 
    10139 
    102 #??? 
    103 #        ( BX, BX.nc   ,     2     , \ 
    104 #      0.0, 0.0, 0.0, 0.7, 0.0, 0.0, 0.0, 0.0, 0.0, 0.3, 0.0, 0.0, 0.0, \ 
    105 #      0.0, 0.0, 0.0, 2.9, 0.0, 0.0, 0.0, 0.0, 0.0, 2.9, 0.0, 0.0, 0.0) \ 
    106 # \ 
    107  
    108 #\ 
    109 #       ( ??, ??.nc   ,     2000,     0     , \ 
    110 #     0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, \ 
    111 #                                   0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) \ 
    112 #\ 
    113  
    114  
    115 # The following tables of parameters for SECHIBA 
    116 # are in the following order : 
     40# To fill the VEGMAX for each site 
     41# here below is the standard PFT list 
    11742# 
    11843#    1 - Bare soil 
    … …  
    13156 
    13257[SPINUP] 
    133     # SPINUP configuration :  
    134     # ---------------------- 
    135 # !! Step of time in N Years !! 
    136 # !! The spinup will change if the fluxnet file contains more than one year !! 
    137 #    ( N = Number of years contain in fluxnet forcing file )  
    138 #     each Year * N 
    139  
    14058# Initialisation for spin-up : 
    14159# orchidee with sechiba alone (!!! if ok_stomate == n !!!) 
    … …  
    14361# orchidee with stomate 
    14462duree_inistomate=1 
    145 # teststomate (only if duree_nostomate or duree_inistomate > 0) 
     63# teststomate (only if duree_inistomate > 0) 
    14664duree_offlineini=0 
    14765 
    148 # Loop configuration for spin-up : 
     66# Loop over ORCHIDEE runs (used for spin-up) 
    14967# The whole job is restarted n_iter times 
    15068n_iter=1 
    … …  
    15674duree_carbonsol=10000 
    15775 
    158 # Finalization for spin-up : 
    159 # all orchidee 
    160 duree_final=20 
     76# Final run (full ORCHIDEE) 
    16177# This last parameter must be non-zero. 
     78duree_final=200 
    16279 
    16380 
    16481    # POST configuration :  
    16582    # -------------------- 
    166 # ATLAS fix parameters : 
    16783# Atlas Name :  
    168 AtlasCfg=atlas_FLUXNET.cfg 
    169 #atlas_FLUXNET.cfg 
    170 #atlas_FLUXNET_soenke.cfg 
     84AtlasCfg=atlas_FLUXNET_LATHUILE.cfg 
    17185 
    172 # observation_file 
    173 observation_file_path='${R_BC}/OOL/${config_UserChoices_TagName}/FLUXNET/BC/${Site}.nc' 
    174 #'${R_BC}/OOL/${config_UserChoices_TagName}/FLUXNET/OLD/${Site}.nc' 
    175 #'${R_BC}/OOL/${config_UserChoices_TagName}/FLUXNET/OBS/${Site}_obs_gapfilled.nc' 
     86# Observation_file 
     87observation_file_path='/home/orchidee01/vuichard/Input_Fluxnet/${Site}.nc' 
    17688 
    177 # old history file 
    178 reference_file_path='/dmnfs/cont003/p86manci/VALID_OL/OK_STOMATE/${Site}_sechiba_hist.nc' 
    179 # 3 choices : SECHIBA, OK_CO2, OK_STOMATE 
    180 #'/dmnfs/cont003/p86manci/VALID_OL/SECHIBA/${Site}_sechiba_hist.nc' 
     89# History file of former ORCHIDEE runs (Reference) to compare with the current simulations 
     90reference_file_path='/home/orchidee01/vuichard/ORCHIDEE_1951/IGCM_OUT/OL2/Fluxnet_Vuichard/${Site}_sechiba_hist.nc' 
    18191 
    18292# Modulo for SpinUp years 
    … …  
    18797 
    18898[UserChoices] 
    189  
    190 # 
    191 ###-- STOMATE flag 
    192 # 
     99# stomate activated or not ? 
    193100ok_stomate=y 
    194 # 
    195 ###-- OK_CO2 flag 
    196 # 
     101# Photosynthesis activated or not ? 
    197102ok_co2=y 
    198  
    199 # 
    200 ###-- NEW HYDROL CWRR flag 
    201 # 
     103# New hydrology (deRosnay) activated or not ? 
    202104ok_newhydrol=n 
    203105 
    204 # 
    205 ## DEBUG mode for SPINUP  
     106# DEBUG mode for SPINUP  
    206107# 
    207108# This mode keep all SPINUP directory in ARCHIVE 
    208109# If disable, all ARCHIVE is automaticly cleaned. 
    209 #  
    210110DEBUG_SPIN=n 
    211111# If you don't want to keep old spinup steps, but last one 
    212 CONSERVE=y 
     112CONSERVE=n  
    213113 
    214114[SubJobParams] 
    215 # You can specify here any parameters to be modified in sechiba.def, stomate.def or driver.def for SpinUp Subjobs. 
    216 # NEW : due to split of orchidee.def in component specific parameter files, 
    217 #       you must add here a prefix for the specific parameter file. 
    218 driver_DEBUG_INFO=n 
    219 sechiba_LONGPRINT=n 
     115# You can specify here any parameters to be modified in sechiba.def, stomate.def or driver.def 
     116# due to split of orchidee.def in component specific parameter files, 
     117# you must add here a prefix for the specific parameter file. 
    220118stomate_BAVARD=0 
    221119sechiba_ALMA_OUTPUT=y 
    222 driver_ALLOW_WEATHERGEN=n 
    223120sechiba_SECHIBA_reset_time=y 
    224 ## To begin with half water stress 
    225 #sechiba_HYDROL_HUMR=0.5 
    226 # FLUXNET files have hour frequency values. 
    227121driver_SPLIT_DT=1 

  • tags/ORCHIDEE_1_9_5_2/ORCHIDEE_OL/FLUXNET/fluxnet_LATHUILE.card

    r119 r405  
    5757\ 
    5858       ( RU-Zot, RU-Zot.nc,  2002,                        &nbs