wiki:DevelopmentActivities/ORCHIDEE-MICT-IMBALANCE-P/Evaluation

Version 24 (modified by ajornet, 5 years ago) (diff)

--

Solved Issues

Spitfire with negative values at fuel_XXhr variables by A.Jornet 15/12/2015

Supposedly solved at r3021 and r3022. It seems that it is not the case.

There are negative values at diff_frac in stomate_litter.

It only fails when DGVM is active.

I disscussed it with Chao, he thinks it comes from somewhere else. I also disscussed with Dan about DGVM. She did not find anything.

It remains to check if this is a problem related only to precission. Due to the use of min_stomate.

It is actually the source the makes fuel_XXhr variables become negative (day 41). It is difficult to debug.

Biomass is very small Xe-09 or Xe-10 at the beginning (day 3). This value is skipped due to min_stomate (1e-08). As the simulation evolves, this makes a difference at diff_frac (there is always a very small difference). So, diff_frac is more than 4 times bigger than the current fuel_10hr itself. So it becomes negative.

diff_frac(:) = litter(:,imetabolic,m,l,icarbon) - fuel_nlitt_total_pft(:,m,imetabolic,icarbon)
WHERE(diff_frac(:).GT.zero .OR. -1.0*diff_frac(:).GT.zero)
                   ! simply divided even fraction to each fuel class, that is 1/4 of buffer
                   fuel_1hr(:,m,imetabolic,icarbon) = fuel_1hr(:,m,imetabolic,icarbon) + 0.25 * diff_frac(:)
                   fuel_10hr(:,m,imetabolic,icarbon) = fuel_10hr(:,m,imetabolic,icarbon) + 0.25 * diff_frac(:)
                   fuel_100hr(:,m,imetabolic,icarbon) = fuel_100hr(:,m,imetabolic,icarbon) + 0.25 * diff_frac(:)
                   fuel_1000hr(:,m,imetabolic,icarbon) = fuel_1000hr(:,m,imetabolic,icarbon) + 0.25 * diff_frac(:)
                   fuel_nlitt_total_pft(:,m,imetabolic,icarbon) = fuel_1hr(:,m,imetabolic,icarbon) + &
                        fuel_10hr(:,m,imetabolic,icarbon) + fuel_100hr(:,m,imetabolic,icarbon) + &
                        fuel_1000hr(:,m,imetabolic,icarbon)
ENDWHERE

I replaced all min_stomate to zero in the stomate_litter. The simulation finishes. This makes lignin_struc get a value of ~0.7 with very small litter Xe-09 instead of 0. Due to min_stomate.

stomate_litter -> littercalc:

     !! 2.2.7 new lignin content: add old lignin and lignin increase, divide by
     !!       total structural litter (above/below)
     DO l = 1, nlevs !Loop over litter levels (above and below ground)
        DO m = 1,nvm !Loop over PFTs
           WHERE( litter(:,istructural,m,l,icarbon) .GT. min_stomate )
 
        !MM : Soenke modif
        ! Best vectorization ?
 !!$       lignin_struc(:,:,:) = &
 !!$            ( lignin_struc(:,:,:)*old_struc(:,:,:) + lignin_struc_inc(:,:,:) ) / &
 !!$            litter(:,istructural,:,:,icarbon)
 
              lignin_struc(:,m,l) = lignin_struc(:,m,l) * old_struc(:,m,l)
              lignin_struc(:,m,l) = lignin_struc(:,m,l) + lignin_struc_inc(:,m,l)
              lignin_struc(:,m,l) = lignin_struc(:,m,l) / litter(:,istructural,m,l,icarbon)
           ELSEWHERE
              lignin_struc(:,m,l) = zero
           ENDWHERE
 
          IF (m == 11) THEN
              WRITE(numout, *) 'lignin_struc = ', lignin_struc(:,m,l)
              IF (lignin_struc(1,m,l) > 0) THEN
                  WRITE (*,*) 'ooooooh'
              ENDIF
          ENDIF
 
        ENDDO
     ENDDO

Solved: [3075]

Bug report in DGVM by F. Maignan on 2015/12/04

Problem: There are no trees with DGVM (without SPITFIRE).

I'm looking at DGVM without SPITFIRE. In dynamic mode, we have first crops (PFTs 12 & 13) but natural PFTs appear later: grasses (PFTs 10 & 11) and after a few years trees (PFTs 2 to 10). At least that's how it behaves with the TRUNK version.

  • TRUNK 2800 + latest correction of the DGVM bug (#215)

Configuration: /home/satellites2/maignan/ORCHIDEE/TESTS/test_TRUNK_2800/modipsl/config/ORCHIDEE_OL/OOL_SEC_STO_DGVM
We see natural vegetation growing (check variable VEGET_MAX):
/home/satellites2/maignan/ORCHIDEE/TESTS/test_TRUNK_2800/ARCHIVE/IGCM_OUT/OL2/DEVT/secsto/SECHSTOM.2800/SBG/Output/MO/
Grasses are already well extended in 1903, same for trees in 1907.

  • ORCHIDEE-MICT branch r3039

Configuration: /home/surface3/maignan/ORCHIDEE/TESTS/test_MICT_branch/modipsl/config/ORCHIDEE_OL/OOL_SEC_STO_DGVM_NOFIRE
Trees appear but much later and much less than with the TRUNK version.
/home/surface3/maignan/ORCHIDEE/TESTS/test_MICT_branch/ARCHIVE/IGCM_OUT/OL2/DEVT/secsto/SECHSTOM.DGVM.NOFIRE/SBG/Output/MO
Grasses present in 1903, PFT 11 similar to TRUNK, PFT 10 quite different. No trees in 1907.

If you want to know more about PFTs, you have a list in pft_parameters:

      !! 5.1 Read the name of the PFTs given by the user
      !
      !Config Key   = PFT_NAME
      !Config Desc  = Name of a PFT
      !Config if    = OK_SECHIBA or OK_STOMATE
      !Config Def   = bare ground, tropical broad-leaved evergreen, tropical broad-leaved raingreen,
      !Config         temperate needleleaf evergreen, temperate broad-leaved evergreen temperate broad-leaved summergreen,
      !Config         boreal needleleaf evergreen, boreal broad-leaved summergreen, boreal needleleaf summergreen,
      !Config         C3 grass, C4 grass, C3 agriculture, C4 agriculture
      !Config Help  = the user can name the new PFTs he/she introducing for new species
      !Config Units = [-]
      CALL getin_p('PFT_NAME',pft_name)

Solved in r3057.

Bug report in stomate_wet_ch4_pt_ter_0.f90 by M. Guimberteau on 2015/08/27

The following parameters need to be defined and this part of the code should be added:

       ! define diffusion coefficients (after Penman)                                                                                                                          
       ! soil air                                                                                                                                                              
         diffsa=diffair*0.66*rpv
       ! soil water                                                                                                                                                            
         diffsw=0.0001*diffsa
       ! standing water                                                                                                                                                        
         diffwa=0.0001*diffair

Albert: changed r2906

Bug report in thermosoil by M. Guimberteau on 2015/08/24

The variable SoilTemp? should be in the almaoutput condition of the loop and related to ptn_pftmean and not ptn:

Before in the code:

       IF ( .NOT. almaoutput ) THEN
          CALL histwrite_p(hist_id, 'SoilTemp', kjit, ptn_pftmean, kjpindex*ngrnd, indexgrnd)                                                                                   
          CALL histwrite_p(hist2_id, 'ptn', kjit, ptn, kjpindex*ngrnd, indexgrnd)                                                                                               
       ELSE
          CALL histwrite_p(hist2_id, 'SoilTemp', kjit, ptn, kjpindex*ngrnd, indexgrnd)                                                                                          
          CALL histwrite_p(hist2_id, 'Qg', kjit, soilflx, kjpindex, index)
          CALL histwrite_p(hist2_id, 'DelSurfHeat', kjit, surfheat_incr, kjpindex, index)
          CALL histwrite_p(hist2_id, 'DelColdCont', kjit, coldcont_incr, kjpindex, index)
       ENDIF
    ENDIF

After corrections:

     IF ( hist2_id > 0 ) THEN
       IF ( .NOT. almaoutput ) THEN
          CALL histwrite_p(hist2_id, 'ptn_pftmean', kjit, ptn_pftmean, kjpindex*ngrnd, indexgrnd)
       ELSE
          CALL histwrite_p(hist2_id, 'SoilTemp', kjit, ptn_pftmean, kjpindex*ngrnd, indexgrnd)
          CALL histwrite_p(hist2_id, 'Qg', kjit, soilflx, kjpindex, index)
          CALL histwrite_p(hist2_id, 'DelSurfHeat', kjit, surfheat_incr, kjpindex, index)
          CALL histwrite_p(hist2_id, 'DelColdCont', kjit, coldcont_incr, kjpindex, index)
       ENDIF
    ENDIF

Albert: done in r2907.

Bugs report in intersurf and ioipslctrl by M. Guimberteau on 2015/08/24

The variable ptn_snow_pftmean does not exist in the code anymore. The call histdef should be removed:

           CALL histdef(hist_id, 'ptn_snow_pftmean', 'Snow temperature,PFT-mean', 'K', &                                                                                        
                & iim,jjm, hori_id, nsnow, 1, nsnow, snowax_id, 32,avescatter(6), dt,dw)

Albert: done in r2907.

The variable soiltemp should be put outside the if (ok_explicitsnow) condition

              CALL histdef(hist_id, 'soiltemp','Soil temperature profile','K', &
               & iim,jjm, hori_id, ngrnd, 1, ngrnd, solax_id, 32,avescatter(6),dt,dw)

Albert: done in r2907.

The histdef high frequency of the variable ptn_pftmean is missing in SECHIBA_HISTLEVEL2 = 7 level:

               CALL histdef(hist2_id, 'ptn_pftmean', 'Soil temperature, PFT-mean','K', &
                & iim,jjm, hori_id2, ngrnd, 1, ngrnd, solax_id, 32,avescatter2(7), dt,dw2)

Albert: done in r2907 and r2908.

Bug report in thermosoil.f90 / thermosoil_finalize by F. Maignan on 2015/09/10 & 2015/08/31

I've started a 2 degree V6 evaluation at LSCE here:
/home/surface3/maignan/ORCHIDEE/TESTS/test_MICT_V6/modipsl/modeles/ORCHIDEE/modipsl/config/ORCHIDEE_OL/OOL_SEC_STO_V6

The modified code is in /home/surface3/maignan/ORCHIDEE/TESTS/test_MICT_V6/modipsl

The gthick, gtemp and gpkappa variables were not saved in the restart file.
The following correction should be applied:

SUBROUTINE thermosoil_finalize(kjit, kjpindex, rest_id, temp_sol_new, gthick, gtemp, gpkappa, soilcap, soilflx)		!@FM_150910
    !! 0. Variable and parameter declaration 
    !! 0.1 Input variables 
    INTEGER(i_std), INTENT(in)                            :: kjit             !! Time step number (unitless)  
    INTEGER(i_std), INTENT(in)                            :: kjpindex         !! Domain size (unitless) 
    INTEGER(i_std),INTENT (in)                            :: rest_id          !! Restart file identifier(unitless) 
    REAL(r_std),DIMENSION (kjpindex), INTENT (in)         :: temp_sol_new     !! Surface temperature at the present time-step, 
    REAL(r_std),DIMENSION (kjpindex),INTENT(in)           :: gthick           !! First soil layer thickness		!@FM_150910
    REAL(r_std),DIMENSION (kjpindex),INTENT(in)           :: gtemp            !! First soil layer temperature		!@FM_150831
    REAL(r_std),DIMENSION (kjpindex),INTENT(in)           :: gpkappa          !! First soil layer thermal conductivity	!@FM_150910

...
        var_name= 'temp_sol_beg'
        CALL restput_p(rest_id, var_name, nbp_glo, 1, 1, kjit, temp_sol_beg, 'scatter', nbp_glo, index_g)

        CALL restput_p(rest_id, 'gthick', nbp_glo, 1, 1, kjit, gthick, 'scatter', nbp_glo, index_g)	!@FM_150910
        CALL restput_p(rest_id, 'gtemp', nbp_glo, 1, 1, kjit, gtemp, 'scatter', nbp_glo, index_g)	!@FM_150831
        CALL restput_p(rest_id, 'gpkappa', nbp_glo, 1, 1, kjit, gpkappa, 'scatter', nbp_glo, index_g)	!@FM_150910
...

with the corresponding call in sechiba_finalize:

    !! 5. Call soil thermodynamic to write restart files
!@FM_150831    CALL thermosoil_finalize(kjit, kjpindex, rest_id, temp_sol_new, soilcap, soilflx)
    CALL thermosoil_finalize(kjit, kjpindex, rest_id, temp_sol_new, gthick, gtemp, gpkappa, soilcap, soilflx)	!@FM_150910

The MICT thermosoil_finalize is different from the Trunk version so probably other variables have to be checked (done 2015/09/10).

Albert: done in r2909.

Bug report in stomate.f90 / call writerestart by F. Maignan on 2015/08/27

The following correction should be applied:

...
        &          litterpart, litter, dead_leaves, &
!@FM_150827        &          carbon, lignin_struc,turnover_time,&
        &          carbon, lignin_struc,&		!@FM_150827
        &          ni_acc,fuel_1hr,fuel_10hr,fuel_100hr,fuel_1000hr, &
        & 	   turnover_time, &			!@FM_150827
...

to have a call consistent with the interface of writerestart in stomate_io.f90, otherwise some variables are crushed in the restart file.

Albert: done in r2910.

Bug report in lpj_spitfire.f90 / spitfire by F. Maignan on 2015/08/27

Add:

ratio_C3_livegrass(:) = zero
ratio_C4_livegrass(:) = zero
char_dens_fuel_ave(:) = zero

as mentionned in Albert's mail (spitfire value too big, 2015/08/21).

Albert: done in r2910.

Bug report in stomate_permafrost_carbon.f90 / snow_interpol / ! 5. inter- or extrapolate by F. Maignan on 2015/08/24

IF ( dzio(ip,iv) .GT. 0. ) THEN

should be replaced with:

IF ( dzio(ip,iv) .GT. min_stomate ) THEN

dzio(ip,iv) can be very small and there is a subsequent division per zero and orchidee crashes.

Albert: done in r2910.

Bug report in stomate_Cforcing_name and stomate_forcing_name by F. Maignan on 2015/08/21

I've found a problem with the variable Cforcing_name, which is declared both as a global save variable of the stomate module and as a save variable of the stomate_main procedure. It is read in stomate_initialize using a getin_p but indeed the value is not correctly transmitted to the stomate_main procedure. There is also a stomate_Cforcing_name in constantes_var.

I corrected this using the constantes_var variable stomate_Cforcing_name and suppressing the two other Cforcing_name in the stomate module.

Then I did the same for forcing_name.

I don't know where this weird code originates from, I hope it's OK with forcesoil and teststomate, I'll check later.

Albert: done in r2910.

Bug report on LIGHTNING_FILE by F. Maignan on 2015/08/24

LIGHTNING_FILE should not be read if FIRE_DISABLE=y.
Albert: done in r2911.

To be corrected

  • 2015/08/21
    • GRAZING_MAP should not be read if grassland management is not activated.
    • If OK_EXPLICITSNOW is FALSE then snowdz is not initialized but used anyway (this is probably the case of other variables, to be checked).

Bug report by T. Wang

There is a small bug in our MICT version that might induce the failure of some point simulations. In subroutine enerbil_flux of enerbil.f90:

CHANGE

!! To get the extra energy used to melt the snowpack
IF (ok_explicitsnow .AND. temp_sol_new (ji) > tp_00 .AND. SUM(snowdz(ji,:)) .GT. zero .AND. snowcap(ji) .GT. zero) THEN

lwup_tmp(ji) = emis(ji) * c_stefan * temp_sol(ji)**4 + &
&     quatre * emis(ji) * c_stefan * temp_sol(ji)**3 * &
&     (tp_00 - temp_sol(ji))

To

REAL(r_std), PARAMETER :: EPS =1E-3_r_std        
!! To get the extra energy used to melt the snowpack
IF (ok_explicitsnow .AND. temp_sol_new (ji) > tp_00 .AND. SUM(snowdz(ji,:)) .GT. zero .AND. snowcap(ji) .GT. EPS) THEN

lwup_tmp(ji) = emis(ji) * c_stefan * temp_sol(ji)**4 + &
&     quatre * emis(ji) * c_stefan * temp_sol(ji)**3 * &
&     (tp_00 - temp_sol(ji))

Albert: done in r2912.

Bug report by C. Yue

Please put these changes in the src_stomate/Grassland_Management.f90. They're used to protect the division-by-zero from happening. But I don't think this correction is related negative fuel problem. But it could lead to Nan value in fuel if you activate the GM (which is not your case as I checked the run.def you sent before).

From

!spitfire
fuel_all_type(:,:,:) = fuel_1hr(:,:,:) + fuel_10hr(:,:,:) + &
                            fuel_100hr(:,:,:) + fuel_1000hr(:,:,:)
fuel_type_frac(:,:,:,1) = fuel_1hr(:,:,:)/fuel_all_type(:,:,:)
fuel_type_frac(:,:,:,2) = fuel_10hr(:,:,:)/fuel_all_type(:,:,:)
fuel_type_frac(:,:,:,3) = fuel_100hr(:,:,:)/fuel_all_type(:,:,:)
fuel_type_frac(:,:,:,4) = fuel_1000hr(:,:,:)/fuel_all_type(:,:,:)
!end spitfire   

to

!spitfire
fuel_type_frac(:,:,:,:) = zero
fuel_all_type(:,:,:) = fuel_1hr(:,:,:) + fuel_10hr(:,:,:) + &
                            fuel_100hr(:,:,:) + fuel_1000hr(:,:,:)
WHERE (fuel_all_type(:,:,:) .GT. min_stomate)
   fuel_type_frac(:,:,:,1) = fuel_1hr(:,:,:)/fuel_all_type(:,:,:)
   fuel_type_frac(:,:,:,2) = fuel_10hr(:,:,:)/fuel_all_type(:,:,:)
   fuel_type_frac(:,:,:,3) = fuel_100hr(:,:,:)/fuel_all_type(:,:,:)
   fuel_type_frac(:,:,:,4) = fuel_1000hr(:,:,:)/fuel_all_type(:,:,:)
ENDWHERE
!end spitfire   

Albert: done in r2913.