Integration of branch ORCHIDEE-DRIVER in the trunk

[jgipsl]

The branch ORCHIDEE-DRIVER with developments of Jan Polcher will be integrated in the trunk.

Main developments are :

• enhancement to run on different types of grids not only regular lat-lon grid
• added a new driver orchideedriver
• added interface to OASIS and new driver to run with OASIS in offline. The OASIS interface can be used to couple to WRF.

Validation is done for the standard driver dim2_driver. Only small numerical changes in results have been detected.

The validation of the new driver orchideedriver and the interface to OASIS will continue after integration in the trunk.

[jgipsl]

Revision 3446 of the branch ORCHIDEE-DRIVER has been integrated in the trunk rev [3447]

Note : the files in arch/ have not yet been added.

[nvuilsce]

With this development, ORCHIDEE does not work anymore on a single pixel.
Tests have been done with the following forcing files:

• /ccc/work/cont003/dsm/p86ipsl/IGCM/SRF/METEO/Princeton_GPCC/v1/1deg_3h/PGF_GPCC_1deg_3h_2000.nc on a restricted area (LIMIT_NORTH=4; LIMIT_SOUTH=41; LIMIT_WEST=-6; LIMIT_EAST=-5)
• /ccc/work/cont003/dsm/p86ipsl/IGCM/SRF/METEO/CRU-NCEP/v6.1/halfdeg/cruncep_halfdeg_2000.nc

on (LIMIT_NORTH=42; LIMIT_SOUTH=41.5; LIMIT_WEST=-6; LIMIT_EAST=-5.5)

• any fluxnet forcing, for instance:/ccc/work/cont003/dsm/p86ipsl/IGCM/SRF/METEO/FLUXNET/BR-Sa3_2000-2003.nc

It stops with the following message:

FATAL ERROR FROM ROUTINE haversine_laloseglen

--> The polygon here does not originate from a regular
--> latitude longitude grid.
-->

This is problematic ...

[nvuilsce]

I have solved the problem by modifying the routine haversine_reglatlontoploy. ORCHIDEE works now well on a single point. I don't know if this fix can be considered as a permanent solution and so commited in the trunk.

The initial block of code was:

DO i = 1, nbpt
!!$       ! index of latitude
!!$       jp = INT( (index_loc(i)-1) /iim ) + 1
!!$       ! index of longitude
!!$       ip = index_loc(i) - ( jp-1 ) * iim
          ip = iorig(i)
          jp = jorig(i)
          !
          ipm1 = ip-1
          ipp1 = ip+1
          jpm1 = jp-1
          jpp1 = jp+1
          !
          ! Compute the longitude increments
          !
          IF ( ipp1 <= iim ) THEN
             dlonp1 = (lon(ipp1,jp)-lon(ip,jp))/2.0
          ELSE 
             dlonp1 = (lon(ip,jp)-lon(ipm1,jp))/2.0
             IF ( global ) ipp1=1
          ENDIF
          IF ( ipm1 > 0 ) THEN
             dlonm1 = (lon(ip,jp)-lon(ipm1,jp))/2.0
          ELSE
             dlonm1 = (lon(ipp1,jp)-lon(ip,jp))/2.0
             IF ( global ) ipm1=iim
          ENDIF
          !
          ! Compute the latitude increments
          !
          IF ( jpp1 <= jjm ) THEN
             dlatp1 = (lat(ip,jpp1)-lat(ip,jp))/2.0
          ELSE
             dlatp1 = (lat(ip,jp)-lat(ip,jpm1))/2.0
          ENDIF
          IF ( jpm1 > 0 ) THEN
             dlatm1 = (lat(ip,jp)-lat(ip,jpm1))/2.0
          ELSE
             dlatm1 = (lat(ip,jpp1)-lat(ip,jp))/2.0
          ENDIF
          !
          !
          ! The longitude of all elements of the polygone
          !
          lonpoly(i,1) = lon(ip,jp)-dlonm1
          lonpoly(i,2) = lon(ip,jp)
          lonpoly(i,3) = lon(ip,jp)+dlonp1
          lonpoly(i,4) = lon(ip,jp)+dlonp1
          lonpoly(i,5) = lon(ip,jp)+dlonp1
          lonpoly(i,6) = lon(ip,jp)
          lonpoly(i,7) = lon(ip,jp)-dlonm1
          lonpoly(i,8) = lon(ip,jp)-dlonm1
          !
          ! The longitude of all elements of the polygone
          !
          latpoly(i,1) = lat(ip,jp)-dlatp1
          latpoly(i,2) = lat(ip,jp)-dlatp1
          latpoly(i,3) = lat(ip,jp)-dlatp1
          latpoly(i,4) = lat(ip,jp)
          latpoly(i,5) = lat(ip,jp)+dlatm1
          latpoly(i,6) = lat(ip,jp)+dlatm1
          latpoly(i,7) = lat(ip,jp)+dlatm1
          latpoly(i,8) = lat(ip,jp)
          !
          ! Keep the center of the gridbox
          !
          center(i,1) = lon(ip,jp)
          center(i,2) = lat(ip,jp)
          !
          ! Get the neighbours when they exist in the list of land points
          ! There are no neighbours over the North or South poles.
          !
          IF ( ipm1 > 0 .AND. jpm1 > 0 ) neighb_loc(i,1) = correspondance(ipm1,jpm1)
          IF ( jpm1 > 0 ) neighb_loc(i,2) = correspondance(ip,jpm1)
          IF ( ipp1 <= iim .AND. jpm1 > 0 ) neighb_loc(i,3) = correspondance(ipp1,jpm1)
          IF ( ipp1 <= iim ) neighb_loc(i,4) = correspondance(ipp1,jp)
          IF ( ipp1 <= iim .AND. jpp1 <= jjm ) neighb_loc(i,5) = correspondance(ipp1,jpp1)
          IF ( jpp1 <= jjm ) neighb_loc(i,6) = correspondance(ip,jpp1)
          IF ( ipm1 > 0 .AND. jpp1 <= jjm ) neighb_loc(i,7) = correspondance(ipm1,jpp1)
          IF ( ipm1 > 0 ) neighb_loc(i,8) = correspondance(ipm1,jp)
          !
       ENDDO

A IF statement has been added for accounting for a specific treatment for nbpt=1, as below:

    IF(nbpt == 1) THEN
       lonpoly(1,:) = lon(ip,jp)
       latpoly(1,:) = lat(ip,jp)
       center(i,1) = lon(ip,jp)
       center(i,2) = lat(ip,jp)
       neighb_loc(:,:) = 1.
    ELSE
       DO i = 1, nbpt
!!$       ! index of latitude
!!$       jp = INT( (index_loc(i)-1) /iim ) + 1
!!$       ! index of longitude
!!$       ip = index_loc(i) - ( jp-1 ) * iim
          ip = iorig(i)
          jp = jorig(i)
          !
          ...
          ...
          IF ( ipm1 > 0 ) neighb_loc(i,8) = correspondance(ipm1,jp)
          !
       ENDDO
  ENDIF 

[jpolcher]

I would like to propose the following solution.

In grid.f90 call a subroutine names haversine_singleploy which will be called in case nblan=iim=jjm=1.

Just as haversine_reglatlontoploy and haversine_regxytoploy it will generate a polygon around the point on which the model in run.

This will define an area for which the forcing data is representative. So this routine will need as exogenous data this area. This can be solved by adding to the run.def a parameter called for instance POINT_AREA.

This solution ensures that all other routines will then have a properly polygon to work on.

[nvuilsce]

The solution proposed by Jan has been tested. It works over a single point. Simulations over a single latudinal or longitudinal band are not handled but a message is printed out before stopping the simulation.
The solution has been commited in r3568.

[jgipsl]

Integration has been finalized with above commit in the trunk. Validation of the new drivers will continue in the trunk.