Merge-Hydro Branch

This version will merge the version of ORCHIDEE developed at LMD with the 11-layer hydrological scheme into the latest version of the trunk (1.9.5.2). Based on the differences between the LMD version and those from which it originates (either 1_9 or 1_9_4_1) (see the attachments), several notes have been written in order to prepare the merge. Here below are the links to these notes.

The hydrology group is : isabelle.gouttevin, jan.polcher, aducharne, aurelien.campoy, matthieu.guimberteau, gerhard.krinner, nathalie.de-noblet, catherine.ottle, pierre.brender, tao.wang, frederique.cheruy, marie-alice.foujols, patricia.cadule.
Everybody with a svn login can(since 12 april 2012) download this branch but only people in the group can commit changes.

How to download ORCHIDEE Hydrology branch for use in offline configuration

1. Get modipsl

   > svn co http://forge.ipsl.jussieu.fr/igcmg/svn/modipsl/trunk modipsl
   

2. Get the ORCHIDEE TOOLS-box. You will be prompted for username and password. Use personal svn login (firstname.lastname) or old sechiba cvs login.

   > cd modipsl/util
   > svn co svn://forge.ipsl.jussieu.fr/orchidee/trunk/TOOLS
   

3. Get the hydrology version. svn_login can be your personal svn login or sechiba

   > TOOLS/recup_my_ORCHIDEE svn_login branches/Hydrology
   

How to download ORCHIDEE Hydrology branch for use in coupled LMDZOR configuration

1. Get modipsl

   > svn co http://forge.ipsl.jussieu.fr/igcmg/svn/modipsl/trunk modipsl
   

2. Modifiy in modipsl/util/mod.def to get ORCHIDEE Hydrology sources instead of default.

   In modipsl/util/mod.def change following line :
   #-C- LMDZOR_v4  ORCHIDEE                    orchidee_1_9_5   2  .   modeles
   into
   #-C- LMDZOR_v4  branches/Hydrology/ORCHIDEE HEAD             15 .  modeles
   

3. Extract configuration LMDZOR_v4 and compile, for more information see ​wiki LMDZOR_v4

   cd modipsl/util
   ./model LMDZOR_v4       => you will be asked for passwd for sechiba
   ./ins_make
   cd ../config/LMDZOR_v4
   gmake
   

4. Some changes are needed to run this version using Choisnel 2 layers. No default set up is done yet for using CWRR.
   • Add in PARAM/orchidee.def :

     SOILCLASS_FILE = soils_param.nc
     TOPOGRAPHY_FILE =cartepente2d_15min.nc
     

   • Add in orchidee.card, in section [InitialStateFiles] following line (and add ", \" at line above):

     (${R_INIT}/SRF/${config_UserChoices_TagName}/cartepente2d_15min.nc, .)
     

   • NB! For the moment this version bug with default set up using RIVER_ROUTING=y. To have it runnig, change in orchidee.def to RIVER_ROUTING=n.

Notes to prepare the merge

Note on CWRR

Note on evapnu

Note on dpu

Note on flag couple

Note on logz0

Note on LAI

Note on Soil Map

Note on Interpolations

Martial notes for the merge

Testing the Merge-Hydro version

Matthieu Guimberteau: Off-line tests in the Amazon with bug corrections

Aurélien Campoy (Feb 2012) - test peformed off-line using 50 years of SAFRAN met forcing for the SIRTA site, without routing, without STOMATE : compares reasonably well with the CWRR version of ORCHIDEE before the merge

Jan Polcher (March 2012) - Off-line simulation over Europe using the WATCH forcing data set : simulation performed, to be analyzed

Nicolas Vuichard (March 2012) - off-line tests at FLUXNET sites using STOMATE : ongoing. See details

Fabienne Maignan (Feb-Jul 2012) - global off-line simulation

• forcing: ERA-Interim regular grid 0.72° 3-hourly
• computer: LSCE/obelix*
• 16 procs
• CWRR + STOMATE + ROUTING
  • segmentation fault SIGSEGV here:

     Projection arrays for Slope map                      : 
     nbvmax =           16
     We will work with      1000105  points of the fine grid
     Aggregate_2d : Slope map                     
     
     aggregate_2D nbvmax =           16 max used =           16
    

    The problem is in slowproc_slope after the call to aggregate_vec:

        DO ib = 1, nbpt
          idi=1
          !-
          !- Reinfiltration coefficient due to the slope: Calculation with parameteres maxlope_ro 
          !-
          slopecoef = zero
    !@BUG SOURCE      DO WHILE ( sub_area(ib,idi) > zero ) !@END BUG SOURCE 
    !@PROPOSED CORRECTION
           DO WHILE ( sub_area(ib,idi) > zero .AND. idi .LE. nbvmax) 
    !@END PROPOSED CORRECTION
             ip = sub_index(ib,idi,1)
             jp = sub_index(ib,idi,2)
             !
             slopecoef = slopecoef + MIN(slopemap(ip,jp)/slope_noreinf, un) * sub_area(ib,idi)
             idi = idi +1
          ENDDO
    

    idi may become greater than nbvmax, which causes the segmentation fault. The problem is resolved by the correction proposed in the above code.
  • stop in hydrol_soil_infilt:

    Error in the calculation of infilt tot -4.176114610371373E-008
     k, ji, jst, mc   1.17784270216726        2.16374183793683             2394
               1  0.277241865204730
    
    FATAL ERROR FROM ROUTINE hydrol_soil_infilt
     --> We will STOP after hydrol_soil_infilt.
     -->
     -->
    
    Fatal error from IOIPSL. STOP in ipslerr with code
    

    This is caused by slightly negative values in the ERA-Interim precipitations, as already mentionned by Nicolas Vuichard. Proposed (and validated) modification in intersurf_main_2d and intersurf_main_1d:

    !@ PROPOSED MODIFICATION
        WHERE(zprecip_rain(:) .LT. 0.)
    	zprecip_rain(:)=0.		
        ENDWHERE				
    !@ END PROPOSED MODIFICATION
        !
        IF (check_INPUTS) THEN
    ...
    

• test option CHECK_CWRR: OK (no stop, much longer ~*3)
• Comparison Choisnel against 196: on-going

Frédérique Cheruy - global simulation coupled to LMDz : to be done

Merge into the trunk (rev 941) (Didier Solyga)

This section lists the modifications needed to merge the DOC/Hydro branch into the trunk :

• hydrol : Move flag ok_throughfall_by_pft to pft_parameters :

       IF ( active_flags%hydrol_cwrr ) THEN
           
           !! 2.1 Read the flag ok_throughfall_by_pft to know if 
           !!      we have to use the parameter throughfall_by_pft
  
           !Config Key   = OK_THROUGHFALL_PFT
           !Config Desc  = Activate use of PERCENT_THROUGHFALL_PFT
           !Config If    = HYDROL_CWRR
           !Config Def   = FALSE
           !Config Help  = If NOT OFF_LINE_MODE it is always TRUE (coupled with a GCM)
           !Config Units = [FLAG]
           IF ( .NOT. OFF_LINE_MODE ) ok_throughfall_by_pft = .TRUE.
           CALL getin_p('OK_THROUGHFALL_PFT',ok_throughfall_by_pft)   
  
        END IF
  
  

The pft parameter throughfall_by_pft is initiliazed and read in pft_parameters.f90. Correct memory allocation and initialization for parameter throughfall_by_pft :

      IF ( .NOT.(active_flags%hydrol_cwrr) .OR. (active_flags%hydrol_cwrr .AND. ok_throughfall_by_pft) ) THEN
         ALLOCATE(throughfall_by_pft(nvm),stat=ier)
         l_error = l_error .OR. (ier /= 0)
         IF (l_error) THEN
            WRITE(numout,*) ' Memory allocation error for throughfall_by_pft. We stop. We need nvm words = ',nvm
            STOP 'pft_parameters_alloc'
         END IF
      END IF

      IF ( .NOT.(active_flags%hydrol_cwrr) .OR.  (active_flags%hydrol_cwrr .AND. ok_throughfall_by_pft) ) THEN
         throughfall_by_pft(:) = throughfall_by_mtc(pft_to_mtc(:))
      ENDIF

because throughfall_by_pft is still used with Choisnel hydrology but not automatically with CWRR. ===> Q: Inconsistency ?

• pft_parameters : humcste has different default values according the value of dpu_max. We know that if we use the 11-layers hydrology, dpu_max should be set to 2 and humcste the corresponding values

     If (active_flags%hydrol_cwrr ) THEN
        humcste(:) = humcste_cwrr(pft_to_mtc(:)) ! values for 2m soil depth
     ELSE
        humcste(:) = humcste_mtc(pft_to_mtc(:))  ! values for 4m soil depth 
     END IF
  

• constantes_soil : reintegrate the module constantes_soil. Add the corresponding "USE" in the code. The old hydrological parameters externalized which are obsolete now have been commented in constantes.f90. dpu_max is set at 4 meters by default (value used for AR5 + Choisnel)
• intersurf : add a consistency test for dpu_max value. dpu_max can't be set to 4 if hydrol_cwrr is activated.

   IF (control_flags%hydrol_cwrr .AND. (dpu_max /= 2.)) THEN
       WRITE (numout,*) "You can not use the 11-layers hydrology with dpu_max /= 2. We set it to 2."
       dpu_max = 2.
    END IF

• routing : Add documentation. the following lines have been adapted for the externalization :

    DO ig = 1, nbpt
       IF (MAXVAL(veget_max(ig,(nvm-3):nvm)) .GT. min_sechiba) THEN
          DO jv = nvm-3, nvm
             transpot_mean(ig) = transpot_mean(ig) + transpot(ig,jv) * veget_max(ig,jv)/ SUM(veget_max(ig,(nvm-3):nvm))
          ENDDO
       ELSE

The corresponding code now is :

   tot_vegfrac_nowoody(:) = zero
   DO jv  = 1, nvm
      IF (is_c3(jv) .OR. is_c4(jv)) THEN
         tot_vegfrac_nowoody(:) = tot_vegfrac_nowoody(:) + veget_max(:,jv) 
      END IF
   END DO

   DO ig = 1, nbpt
      IF ( tot_vegfrac_nowoody(ig) .GT. min_sechiba ) THEN
         DO jv = 1,nvm
            IF ( is_c3(jv) .OR. is_c4(jv) ) THEN
               transpot_mean(ig) = transpot_mean(ig) + transpot(ig,jv) * veget_max(ig,jv)/tot_vegfrac_nowoody(ig)  
            END IF
         END DO }}}