= Fiche : flag couplé =

== Description ==
Plusieurs formulations (humrel, drysoil_frac, rsol, albedo) avaient été modifiées pour des applications en mode couplé afin de compenser des biais du modèle. Ces formulations proposées par Nathalie et Frédéric ne sont pas présentes dans la version LMD, version qui conserve donc les formulations originelles. 

== Choix ==
On décide de rajouter un flag spécifique pour le couplé (flag non actif par défaut) permettant d'activer les formulations proposées par Nathalie et Frédéric. Les modules impactées sont condveg.f90 pour alb_bare et hydrolc.f90 pour rsol, drysoilfrac, throughfall et humrel

== Mise en oeuvre ==
La version à modifier est la version d'origine (1.9 ou 1.9.4)

 * Traitement de alb_bare_model
  * dans condveg.f90, on remplace 
  {{{
      alb_bare_model=.FALSE.
      CALL getin_p('ALB_BARE_MODEL', alb_bare_model)
  }}}
  par
  {{{
      ! Declaration d'une variable du module en tete de module
      LOGICAL, SAVE                         :: couple = .FALSE.

+

      alb_bare_model = .TRUE.
      CALL getin_p('COUPLE', couple)	
      IF(couple) THEN
             alb_bare_model = .FALSE.
       ENDIF
       CALL getin_p('ALB_BARE_MODEL', alb_bare_model)
  }}}

 * Traitement de rsol
  * On remplace dans hydrolc_var_init (hydrolc.f90) 
  {{{
    !
    ! Compute the resistance to bare soil evaporation
    !
    rsol(:) = -un
    DO ji = 1, kjpindex
       IF (veget(ji,1) .GE. min_sechiba) THEN
          !
          ! Correction Nathalie - le 28 mars 2006 - sur conseils Fred Hourdin
          ! on modifie le rsol pour que la resistance croisse subitement si on s'approche
          ! du fond. En gros, rsol=hdry*rsol_cste pour hdry < 1m70
          !Ancienne formulation
          rsol(ji) = dss(ji,1) * rsol_cste
          !Nouvelle formulation Nath          
          !rsol(ji) =  ( drysoil_frac(ji) + 1./(10.*(dpu_cste - drysoil_frac(ji))+1.e-10)**2 ) * rsol_cste
       ENDIF
    ENDDO
  }}}
  par
  {{{
    ! Declaration d'une variable du module en tete de module
    LOGICAL, SAVE                                     ::  couple = .FALSE.

+

    ! Nouvelle Declaration locale
    REAL(r_std), DIMENSION (kjpindex) :: frsoil 

+  

    !
    ! Compute the resistance to bare soil evaporation
    !
    CALL getin_p('COUPLE', couple)	
    IF(couple) THEN
       frsoil(:)=( drysoil_frac(:) + 1./(10.*(dpu_cste – drysoil_frac(:))+1.e-10)**2 ) 
    ELSE
       frsoil(:)=dss(:,1) 
    ENDIF
    rsol(:) = -un
    DO ji = 1, kjpindex
       IF (veget(ji,1) .GE. min_sechiba) THEN
	 rsol(ji) = dss(ji,1) * rsol_cste
       ENDIF   
    ENDDO
  }}}
  * On remplace dans hydrolc_soil (hydrolc.f90)
  {{{
    ! The fraction of visibly dry soil (dry when dss(:,1) = 0.1 m)
    ! Correction Nathalie - le 28 Mars 2006 - re-ecriture drysoil_frac/hdry - Fred Hourdin
    !drysoil_frac(:) = MIN(MAX(dss(:,1),0.)*10._r_std, un)
    drysoil_frac(:) = a_subgrd(:,1)*dss(:,1) + (1.-a_subgrd(:,1))*dsp(:,1)
    !
    ! Compute the resistance to bare soil evaporation.
    !
    rsol(:) = -un
    DO ji = 1, kjpindex
       IF (veget(ji,1) .GE. min_sechiba) THEN
          !
          ! Correction Nathalie - le 28 mars 2006 - sur conseils Fred Hourdin
          ! on modifie le rsol pour que la resistance croisse subitement si on s'approche
          ! du fond. En gros, rsol=hdry*rsol_cste pour hdry < 1m70
          !rsol(ji) = dss(ji,1) * rsol_cste
          rsol(ji) =  ( drysoil_frac(ji) + 1./(10.*(dpu_cste - drysoil_frac(ji))+1.e-10)**2 ) * rsol_cste
       ENDIF
    ENDDO
  }}}
  par
  {{{
    ! Nouvelle Declaration locale
    REAL(r_std), DIMENSION (kjpindex) :: frsoil 

    ! The fraction of visibly dry soil (dry when dss(:,1) = 0.1 m)
    IF (couple) THEN
       ! Correction Nathalie - le 28 Mars 2006 - re-ecriture drysoil_frac/hdry - Fred Hourdin	 
       ! on modifie le rsol pour que la resistance croisse subitement si on s'approche
       ! du fond. En gros, rsol=hdry*rsol_cste pour hdry < 1m70
       drysoil_frac(:) = a_subgrd(:,1)*dss(:,1) + (1.-a_subgrd(:,1))*dsp(:,1)
       frsoil(:)=( drysoil_frac(:) + 1./(10.*(dpu_cste – drysoil_frac(:))+1.e-10)**2 ) 
    ELSE
       drysoil_frac(:) = MIN(MAX(dss(:,1),0.)*10._r_std, un)
       frsoil(:)=dss(:,1)
    ENDIF
    !
    ! Compute the resistance to bare soil evaporation.
    !
    rsol(:) = -un
    DO ji = 1, kjpindex
       IF (veget(ji,1) .GE. min_sechiba) THEN
          rsol(ji) = frsoil(ji) * rsol_cste
       ENDIF
    ENDDO
  }}}

 * Traitement de throughfall (pour calcul de precisol et qsintveg)
  * On remplace dans hydrolc_canop (hydrolc.f90)
  {{{
    IF ( firstcall ) THEN
       !Config  Key  = PERCENT_THROUGHFALL_PFT
       !Config  Desc = Percent by PFT of precip that is not intercepted by the canopy
       !Config  Def  = 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30.
       !Config  Help = During one rainfall event, PERCENT_THROUGHFALL_PFT% of the incident rainfall
       !Config         will get directly to the ground without being intercepted, for each PFT.
       
       throughfall_by_pft = (/ 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30. /)
       CALL getin_p('PERCENT_THROUGHFALL_PFT',throughfall_by_pft)
       throughfall_by_pft = throughfall_by_pft / 100.

       firstcall=.FALSE.
    ENDIF
  }}}
  par 
  {{{
    IF ( firstcall ) THEN
       !Config  Key  = PERCENT_THROUGHFALL_PFT
       !Config  Desc = Percent by PFT of precip that is not intercepted by the canopy
       !Config  Def  = 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30. 30.
       !Config  Help = During one rainfall event, PERCENT_THROUGHFALL_PFT% of the incident rainfall
       !Config         will get directly to the ground without being intercepted, for each PFT.
       
       throughfall_by_pft = (/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0. /)

       IF (couple) THEN
	throughfall_by_pft = (/ 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30., 30. /)
       ENDIF
       CALL getin_p('PERCENT_THROUGHFALL_PFT',throughfall_by_pft)
 
      throughfall_by_pft = throughfall_by_pft / 100.

       firstcall=.FALSE.
    ENDIF
  }}}

 * Traitement de humrel
  * On remplace dans hydrolc_soil (hydrolc.f90)
  {{{
            ELSE
               ! Corrections Nathalie - le 28 Mars 2006 - sur conseils Fred Hourdin
               !zhumrel_lo(ji) = EXP( - humcste(jv) * dpu_cste * (dsp(ji,jv)/dpu_cste) )
               !zhumrel_up(ji) = EXP( - humcste(jv) * dpu_cste * (dss(ji,jv)/dsg(ji,jv)) )
               !humrel(ji,jv) = MAX(zhumrel_lo(ji),zhumrel_up(ji))
               !
               ! As we need a slower variable for vegetation growth the stress is computed
               ! differently than in humrel.
               !
               zhumrel_lo(ji) = EXP( - humcste(jv) * dsp(ji,jv))
               zhumrel_up(ji) = EXP( - humcste(jv) * dss(ji,jv))
               ! Ajouts Nathalie - Fred - le 28 Mars 2006
               a_subgrd(ji,jv)=MIN(MAX(dsg(ji,jv)-dss(ji,jv),0.)/dsg_min,1.)
               humrel(ji,jv)=a_subgrd(ji,jv)*zhumrel_up(ji)+(1.-a_subgrd(ji,jv))*zhumrel_lo(ji)
               !
               vegstress(ji,jv) = zhumrel_lo(ji) + zhumrel_up(ji) - EXP( - humcste(jv) * dsg(ji,jv) ) 
               !
            ENDIF
  }}}
  par
  {{{ 
           ELSE
	   IF(couple)
	       ! Ajouts Nathalie - Fred - le 28 Mars 2006
                     zhumrel_lo(ji) = EXP( - humcste(jv) * dsp(ji,jv))
                     zhumrel_up(ji) = EXP( - humcste(jv) * dss(ji,jv))
                    a_subgrd(ji,jv)=MIN(MAX(dsg(ji,jv)-dss(ji,jv),0.)/dsg_min,1.)
                    humrel(ji,jv)=a_subgrd(ji,jv)*zhumrel_up(ji)+(1.-a_subgrd(ji,jv))*zhumrel_lo(ji)
	   ELSE
                    zhumrel_lo(ji) = EXP( - humcste(jv) * dpu_cste * (dsp(ji,jv)/dpu_cste) )
                    zhumrel_up(ji) = EXP( - humcste(jv) * dpu_cste * (dss(ji,jv)/dsg(ji,jv)) )
                    humrel(ji,jv) = MAX(zhumrel_lo(ji),zhumrel_up(ji))
	   ENDIF
             ENDIF
  }}}


== Tests ==