source: LMDZ6/trunk/libf/phylmd/surf_land_orchidee_mod.F90 @ 3551

Last change on this file since 3551 was 3551, checked in by oboucher, 5 years ago

remove old comment

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 29.1 KB
RevLine 
[781]1!
2MODULE surf_land_orchidee_mod
[1146]3#ifndef ORCHIDEE_NOOPENMP
[2571]4#ifndef ORCHIDEE_NOZ0H
[2952]5#ifndef ORCHIDEE_NOFREIN
[3435]6#ifndef ORCHIDEE_NOUNSTRUCT
[781]7!
[2410]8! This module controles the interface towards the model ORCHIDEE.
[781]9!
[2410]10! Compatibility with ORCHIDIEE :
[2952]11! The current version can be used with ORCHIDEE/trunk from revision 4465.
12! This interface is used if none of the cpp keys ORCHIDEE_NOOPENMP,
13! ORCHIDEE_NOZ0H or ORCHIDEE_NOFREIN is set.
[2410]14!
[781]15! Subroutines in this module : surf_land_orchidee
16!                              Init_orchidee_index
17!                              Get_orchidee_communicator
18!                              Init_neighbours
19
20  USE dimphy
[1067]21#ifdef CPP_VEGET
[781]22  USE intersurf     ! module d'ORCHIDEE
[1067]23#endif
[781]24  USE cpl_mod,      ONLY : cpl_send_land_fields
[996]25  USE surface_data, ONLY : type_ocean
[3435]26  USE geometry_mod, ONLY : dx, dy, boundslon, boundslat,longitude, latitude, cell_area,  ind_cell_glo
[781]27  USE mod_grid_phy_lmdz
[2351]28  USE mod_phys_lmdz_para, mpi_root_rank=>mpi_master
[3391]29  USE carbon_cycle_mod, ONLY : nbcf_in_orc, nbcf_out, fields_in, yfields_in, yfields_out, cfname_in, cfname_out
[3435]30  USE nrtype, ONLY : PI
31 
[781]32  IMPLICIT NONE
33
34  PRIVATE
35  PUBLIC  :: surf_land_orchidee
36
37CONTAINS
38!
39!****************************************************************************************
40
41  SUBROUTINE surf_land_orchidee(itime, dtime, date0, knon, &
[2410]42       knindex, rlon, rlat, yrmu0, pctsrf, &
[781]43       debut, lafin, &
[2240]44       plev,  u1_lay, v1_lay, gustiness, temp_air, spechum, epot_air, ccanopy, &
[781]45       tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
46       precip_rain, precip_snow, lwdown, swnet, swdown, &
[1146]47       ps, q2m, t2m, &
[781]48       evap, fluxsens, fluxlat, &             
[888]49       tsol_rad, tsurf_new, alb1_new, alb2_new, &
[2952]50       emis_new, z0m_new, z0h_new, qsurf, &
51       veget, lai, height )
[1279]52
53    USE mod_surf_para
54    USE mod_synchro_omp
[3391]55    USE carbon_cycle_mod
[1785]56    USE indice_sol_mod
[2311]57    USE print_control_mod, ONLY: lunout
[2344]58    USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
[2348]59#ifdef CPP_VEGET
60    USE time_phylmdz_mod, ONLY: itau_phy
61#endif
[781]62!   
63! Cette routine sert d'interface entre le modele atmospherique et le
64! modele de sol continental. Appel a sechiba
65!
66! L. Fairhead 02/2000
67!
68! input:
69!   itime        numero du pas de temps
70!   dtime        pas de temps de la physique (en s)
71!   nisurf       index de la surface a traiter (1 = sol continental)
72!   knon         nombre de points de la surface a traiter
73!   knindex      index des points de la surface a traiter
74!   rlon         longitudes de la grille entiere
75!   rlat         latitudes de la grille entiere
76!   pctsrf       tableau des fractions de surface de chaque maille
77!   debut        logical: 1er appel a la physique (lire les restart)
78!   lafin        logical: dernier appel a la physique (ecrire les restart)
79!                     (si false calcul simplifie des fluxs sur les continents)
80!   plev         hauteur de la premiere couche (Pa)     
81!   u1_lay       vitesse u 1ere couche
82!   v1_lay       vitesse v 1ere couche
83!   temp_air     temperature de l'air 1ere couche
84!   spechum      humidite specifique 1ere couche
85!   epot_air     temp pot de l'air
[1279]86!   ccanopy      concentration CO2 canopee, correspond au co2_send de
87!                carbon_cycle_mod ou valeur constant co2_ppm
[781]88!   tq_cdrag     cdrag
89!   petAcoef     coeff. A de la resolution de la CL pour t
90!   peqAcoef     coeff. A de la resolution de la CL pour q
91!   petBcoef     coeff. B de la resolution de la CL pour t
92!   peqBcoef     coeff. B de la resolution de la CL pour q
93!   precip_rain  precipitation liquide
94!   precip_snow  precipitation solide
95!   lwdown       flux IR descendant a la surface
96!   swnet        flux solaire net
97!   swdown       flux solaire entrant a la surface
98!   ps           pression au sol
99!   radsol       rayonnement net aus sol (LW + SW)
100!
101! output:
102!   evap         evaporation totale
103!   fluxsens     flux de chaleur sensible
104!   fluxlat      flux de chaleur latente
105!   tsol_rad     
106!   tsurf_new    temperature au sol
[888]107!   alb1_new     albedo in visible SW interval
108!   alb2_new     albedo in near IR interval
[781]109!   emis_new     emissivite
[2571]110!   z0m_new      surface roughness for momentum
111!   z0h_new      surface roughness for heat
[781]112!   qsurf        air moisture at surface
113!
[793]114    INCLUDE "YOMCST.h"
[2952]115    INCLUDE "dimpft.h"
[781]116!
117! Parametres d'entree
118!****************************************************************************************
119    INTEGER, INTENT(IN)                       :: itime
120    REAL, INTENT(IN)                          :: dtime
121    REAL, INTENT(IN)                          :: date0
122    INTEGER, INTENT(IN)                       :: knon
123    INTEGER, DIMENSION(klon), INTENT(IN)      :: knindex
124    LOGICAL, INTENT(IN)                       :: debut, lafin
125    REAL, DIMENSION(klon,nbsrf), INTENT(IN)   :: pctsrf
126    REAL, DIMENSION(klon), INTENT(IN)         :: rlon, rlat
[2410]127    REAL, DIMENSION(klon), INTENT(IN)         :: yrmu0 ! cosine of solar zenith angle
[781]128    REAL, DIMENSION(klon), INTENT(IN)         :: plev
[2240]129    REAL, DIMENSION(klon), INTENT(IN)         :: u1_lay, v1_lay, gustiness
[781]130    REAL, DIMENSION(klon), INTENT(IN)         :: temp_air, spechum
131    REAL, DIMENSION(klon), INTENT(IN)         :: epot_air, ccanopy
132    REAL, DIMENSION(klon), INTENT(IN)         :: tq_cdrag
133    REAL, DIMENSION(klon), INTENT(IN)         :: petAcoef, peqAcoef
134    REAL, DIMENSION(klon), INTENT(IN)         :: petBcoef, peqBcoef
135    REAL, DIMENSION(klon), INTENT(IN)         :: precip_rain, precip_snow
136    REAL, DIMENSION(klon), INTENT(IN)         :: lwdown, swnet, swdown, ps
[1146]137    REAL, DIMENSION(klon), INTENT(IN)         :: q2m, t2m
[781]138
139! Parametres de sortie
140!****************************************************************************************
141    REAL, DIMENSION(klon), INTENT(OUT)        :: evap, fluxsens, fluxlat, qsurf
[888]142    REAL, DIMENSION(klon), INTENT(OUT)        :: tsol_rad, tsurf_new
143    REAL, DIMENSION(klon), INTENT(OUT)        :: alb1_new, alb2_new
[2571]144    REAL, DIMENSION(klon), INTENT(OUT)        :: emis_new, z0m_new, z0h_new
[2952]145    REAL, DIMENSION(klon,nvm_lmdz), INTENT(OUT) :: veget
146    REAL, DIMENSION(klon,nvm_lmdz), INTENT(OUT) :: lai
147    REAL, DIMENSION(klon,nvm_lmdz), INTENT(OUT) :: height
[781]148
149! Local
150!****************************************************************************************
[3391]151    INTEGER                                   :: ij, jj, igrid, ireal, index, nb
[781]152    INTEGER                                   :: error
[1146]153    REAL, DIMENSION(klon)                     :: swdown_vrai
[781]154    CHARACTER (len = 20)                      :: modname = 'surf_land_orchidee'
155    CHARACTER (len = 80)                      :: abort_message
156    LOGICAL,SAVE                              :: check = .FALSE.
157    !$OMP THREADPRIVATE(check)
158
159! type de couplage dans sechiba
160!  character (len=10)   :: coupling = 'implicit'
161! drapeaux controlant les appels dans SECHIBA
162!  type(control_type), save   :: control_in
163! Preserved albedo
164    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: albedo_keep, zlev
165    !$OMP THREADPRIVATE(albedo_keep,zlev)
166! coordonnees geographiques
167    REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: lalo
168    !$OMP THREADPRIVATE(lalo)
[3435]169! boundaries of cells
170    REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE   :: bounds_lalo
171    !$OMP THREADPRIVATE(bounds_lalo)
[781]172! pts voisins
173    INTEGER,ALLOCATABLE, DIMENSION(:,:), SAVE :: neighbours
174    !$OMP THREADPRIVATE(neighbours)
175! fractions continents
176    REAL,ALLOCATABLE, DIMENSION(:), SAVE      :: contfrac
177    !$OMP THREADPRIVATE(contfrac)
178! resolution de la grille
179    REAL, ALLOCATABLE, DIMENSION (:,:), SAVE  :: resolution
180    !$OMP THREADPRIVATE(resolution)
181
182    REAL, ALLOCATABLE, DIMENSION (:,:), SAVE  :: lon_scat, lat_scat 
183    !$OMP THREADPRIVATE(lon_scat,lat_scat)
184
[3435]185! area of cells
186    REAL, ALLOCATABLE, DIMENSION (:), SAVE  :: area 
187    !$OMP THREADPRIVATE(area)
188
[781]189    LOGICAL, SAVE                             :: lrestart_read = .TRUE.
190    !$OMP THREADPRIVATE(lrestart_read)
191    LOGICAL, SAVE                             :: lrestart_write = .FALSE.
192    !$OMP THREADPRIVATE(lrestart_write)
193
194    REAL, DIMENSION(knon,2)                   :: albedo_out
195
196! Pb de nomenclature
197    REAL, DIMENSION(klon)                     :: petA_orc, peqA_orc
198    REAL, DIMENSION(klon)                     :: petB_orc, peqB_orc
199! Pb de correspondances de grilles
200    INTEGER, DIMENSION(:), SAVE, ALLOCATABLE  :: ig, jg
201    !$OMP THREADPRIVATE(ig,jg)
202    INTEGER :: indi, indj
203    INTEGER, SAVE, ALLOCATABLE,DIMENSION(:)   :: ktindex
204    !$OMP THREADPRIVATE(ktindex)
205
206! Essai cdrag
207    REAL, DIMENSION(klon)                     :: cdrag
208    INTEGER,SAVE                              :: offset
209    !$OMP THREADPRIVATE(offset)
210
211    REAL, DIMENSION(klon_glo)                 :: rlon_g,rlat_g
212    INTEGER, SAVE                             :: orch_comm
213    !$OMP THREADPRIVATE(orch_comm)
214
215    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: coastalflow
216    !$OMP THREADPRIVATE(coastalflow)
217    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: riverflow
218    !$OMP THREADPRIVATE(riverflow)
[987]219   
[3435]220    INTEGER :: orch_mpi_rank
221    INTEGER :: orch_mpi_size
[987]222    INTEGER :: orch_omp_rank
223    INTEGER :: orch_omp_size
[3435]224
225    REAL, ALLOCATABLE, DIMENSION(:)         :: longitude_glo
226    REAL, ALLOCATABLE, DIMENSION(:)         :: latitude_glo
227    REAL, ALLOCATABLE, DIMENSION(:,:)       :: boundslon_glo
228    REAL, ALLOCATABLE, DIMENSION(:,:)       :: boundslat_glo
229    INTEGER, ALLOCATABLE, DIMENSION(:)      :: ind_cell_glo_glo
230    INTEGER, ALLOCATABLE, SAVE,DIMENSION(:) :: ind_cell
231    !$OMP THREADPRIVATE(ind_cell)
232    INTEGER :: begin, end
[781]233!
234! Fin definition
235!****************************************************************************************
236
237    IF (check) WRITE(lunout,*)'Entree ', modname
238 
239! Initialisation
240 
241    IF (debut) THEN
[1067]242! Test of coherence between variable ok_veget and cpp key CPP_VEGET
243#ifndef CPP_VEGET
244       abort_message='Pb de coherence: ok_veget = .true. mais CPP_VEGET = .false.'
[2311]245       CALL abort_physic(modname,abort_message,1)
[1067]246#endif
247
[987]248       CALL Init_surf_para(knon)
[781]249       ALLOCATE(ktindex(knon))
250       IF ( .NOT. ALLOCATED(albedo_keep)) THEN
[987]251!ym          ALLOCATE(albedo_keep(klon))
252!ym bizarre que non alloué en knon precedement
253          ALLOCATE(albedo_keep(knon))
[781]254          ALLOCATE(zlev(knon))
255       ENDIF
256! Pb de correspondances de grilles
257       ALLOCATE(ig(klon))
258       ALLOCATE(jg(klon))
259       ig(1) = 1
260       jg(1) = 1
261       indi = 0
262       indj = 2
263       DO igrid = 2, klon - 1
264          indi = indi + 1
[2344]265          IF ( indi > nbp_lon) THEN
[781]266             indi = 1
267             indj = indj + 1
268          ENDIF
269          ig(igrid) = indi
270          jg(igrid) = indj
271       ENDDO
272       ig(klon) = 1
[2344]273       jg(klon) = nbp_lat
[781]274
[3435]275       IF ((.NOT. ALLOCATED(area))) THEN
276          ALLOCATE(area(knon), stat = error)
[781]277          IF (error /= 0) THEN
[3435]278             abort_message='Pb allocation area'
[2311]279             CALL abort_physic(modname,abort_message,1)
[781]280          ENDIF
281       ENDIF
282       DO igrid = 1, knon
[3435]283          area(igrid) = cell_area(knindex(igrid))
[781]284       ENDDO
[3435]285       
286       IF (grid_type==unstructured) THEN
[781]287
[3435]288
289         IF ((.NOT. ALLOCATED(lon_scat))) THEN
290            ALLOCATE(lon_scat(nbp_lon,nbp_lat), stat = error)
291            IF (error /= 0) THEN
292               abort_message='Pb allocation lon_scat'
293               CALL abort_physic(modname,abort_message,1)
294            ENDIF
295         ENDIF
296 
297         IF ((.NOT. ALLOCATED(lat_scat))) THEN
298            ALLOCATE(lat_scat(nbp_lon,nbp_lat), stat = error)
299            IF (error /= 0) THEN
300               abort_message='Pb allocation lat_scat'
301               CALL abort_physic(modname,abort_message,1)
302            ENDIF
303         ENDIF
304         CALL Gather(rlon,rlon_g)
305         CALL Gather(rlat,rlat_g)
306
307         IF (is_mpi_root) THEN
308            index = 1
309            DO jj = 2, nbp_lat-1
310               DO ij = 1, nbp_lon
311                  index = index + 1
312                  lon_scat(ij,jj) = rlon_g(index)
313                  lat_scat(ij,jj) = rlat_g(index)
314               ENDDO
315            ENDDO
316            lon_scat(:,1) = lon_scat(:,2)
317            lat_scat(:,1) = rlat_g(1)
318            lon_scat(:,nbp_lat) = lon_scat(:,2)
319            lat_scat(:,nbp_lat) = rlat_g(klon_glo)
320         ENDIF
321     
322         CALL bcast(lon_scat)
323         CALL bcast(lat_scat)
324               
325       ELSE IF (grid_type==regular_lonlat) THEN
326
327         IF ((.NOT. ALLOCATED(lalo))) THEN
328            ALLOCATE(lalo(knon,2), stat = error)
329            IF (error /= 0) THEN
330               abort_message='Pb allocation lalo'
331               CALL abort_physic(modname,abort_message,1)
332            ENDIF
333         ENDIF
[781]334       
[3435]335         IF ((.NOT. ALLOCATED(bounds_lalo))) THEN
336           ALLOCATE(bounds_lalo(knon,nvertex,2), stat = error)
337           IF (error /= 0) THEN
338             abort_message='Pb allocation lalo'
339             CALL abort_physic(modname,abort_message,1)
340           ENDIF
341         ENDIF
[781]342       
[3435]343         IF ((.NOT. ALLOCATED(lon_scat))) THEN
344            ALLOCATE(lon_scat(nbp_lon,nbp_lat), stat = error)
345            IF (error /= 0) THEN
346               abort_message='Pb allocation lon_scat'
347               CALL abort_physic(modname,abort_message,1)
348            ENDIF
349         ENDIF
350         IF ((.NOT. ALLOCATED(lat_scat))) THEN
351            ALLOCATE(lat_scat(nbp_lon,nbp_lat), stat = error)
352            IF (error /= 0) THEN
353               abort_message='Pb allocation lat_scat'
354               CALL abort_physic(modname,abort_message,1)
355            ENDIF
356         ENDIF
357         lon_scat = 0.
358         lat_scat = 0.
359         DO igrid = 1, knon
360            index = knindex(igrid)
361            lalo(igrid,2) = rlon(index)
362            lalo(igrid,1) = rlat(index)
363            bounds_lalo(igrid,:,2)=boundslon(index,:)*180./PI
364            bounds_lalo(igrid,:,1)=boundslat(index,:)*180./PI
365         ENDDO
[781]366
[3435]367       
368       
369         CALL Gather(rlon,rlon_g)
370         CALL Gather(rlat,rlat_g)
371
372         IF (is_mpi_root) THEN
373            index = 1
374            DO jj = 2, nbp_lat-1
375               DO ij = 1, nbp_lon
376                  index = index + 1
377                  lon_scat(ij,jj) = rlon_g(index)
378                  lat_scat(ij,jj) = rlat_g(index)
379               ENDDO
380            ENDDO
381            lon_scat(:,1) = lon_scat(:,2)
382            lat_scat(:,1) = rlat_g(1)
383            lon_scat(:,nbp_lat) = lon_scat(:,2)
384            lat_scat(:,nbp_lat) = rlat_g(klon_glo)
385         ENDIF
386   
387         CALL bcast(lon_scat)
388         CALL bcast(lat_scat)
389       
[781]390       ENDIF
391!
392! Allouer et initialiser le tableau des voisins et des fraction de continents
393!
394       IF (( .NOT. ALLOCATED(contfrac))) THEN
395          ALLOCATE(contfrac(knon), stat = error)
396          IF (error /= 0) THEN
397             abort_message='Pb allocation contfrac'
[2311]398             CALL abort_physic(modname,abort_message,1)
[781]399          ENDIF
400       ENDIF
401
402       DO igrid = 1, knon
403          ireal = knindex(igrid)
404          contfrac(igrid) = pctsrf(ireal,is_ter)
405       ENDDO
406
407
[3435]408       IF (grid_type==regular_lonlat) THEN
409 
410         IF ( (.NOT.ALLOCATED(neighbours))) THEN
411          ALLOCATE(neighbours(knon,8), stat = error)
412          IF (error /= 0) THEN
413             abort_message='Pb allocation neighbours'
414             CALL abort_physic(modname,abort_message,1)
415          ENDIF
416         ENDIF
417         neighbours = -1.
418         CALL Init_neighbours(knon,neighbours,knindex,pctsrf(:,is_ter))
[781]419
[3435]420       ELSE IF (grid_type==unstructured) THEN
421 
422         IF ( (.NOT.ALLOCATED(neighbours))) THEN
423          ALLOCATE(neighbours(knon,12), stat = error)
424          IF (error /= 0) THEN
425             abort_message='Pb allocation neighbours'
426             CALL abort_physic(modname,abort_message,1)
427          ENDIF
428         ENDIF
429         neighbours = -1.
430 
431       ENDIF
432         
433
[781]434!
435!  Allocation et calcul resolutions
436       IF ( (.NOT.ALLOCATED(resolution))) THEN
437          ALLOCATE(resolution(knon,2), stat = error)
438          IF (error /= 0) THEN
439             abort_message='Pb allocation resolution'
[2311]440             CALL abort_physic(modname,abort_message,1)
[781]441          ENDIF
442       ENDIF
[3435]443       
444       IF (grid_type==regular_lonlat) THEN
445         DO igrid = 1, knon
446            ij = knindex(igrid)
447            resolution(igrid,1) = dx(ij)
448           resolution(igrid,2) = dy(ij)
449         ENDDO
450       ENDIF
451       
[781]452       ALLOCATE(coastalflow(klon), stat = error)
453       IF (error /= 0) THEN
454          abort_message='Pb allocation coastalflow'
[2311]455          CALL abort_physic(modname,abort_message,1)
[781]456       ENDIF
457       
458       ALLOCATE(riverflow(klon), stat = error)
459       IF (error /= 0) THEN
460          abort_message='Pb allocation riverflow'
[2311]461          CALL abort_physic(modname,abort_message,1)
[781]462       ENDIF
[1279]463!
[1454]464! carbon_cycle_cpl not possible with this interface and version of ORHCHIDEE
[1279]465!
[3391]466! >> PC
467!       IF (carbon_cycle_cpl) THEN
468!          abort_message='carbon_cycle_cpl not yet possible with this interface of ORCHIDEE'
469!          CALL abort_physic(modname,abort_message,1)
470!       END IF
471! << PC
[1279]472       
[781]473    ENDIF                          ! (fin debut)
[1279]474 
[781]475!
476! Appel a la routine sols continentaux
477!
478    IF (lafin) lrestart_write = .TRUE.
479    IF (check) WRITE(lunout,*)'lafin ',lafin,lrestart_write
[987]480     
[781]481    petA_orc(1:knon) = petBcoef(1:knon) * dtime
482    petB_orc(1:knon) = petAcoef(1:knon)
483    peqA_orc(1:knon) = peqBcoef(1:knon) * dtime
484    peqB_orc(1:knon) = peqAcoef(1:knon)
485
486    cdrag = 0.
487    cdrag(1:knon) = tq_cdrag(1:knon)
488
489! zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/287.05*temp_air(1:knon))*9.80665)
[2011]490!    zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/RD*temp_air(1:knon))*RG)
491     zlev(1:knon) = plev(1:knon)*RD*temp_air(1:knon)/((ps(1:knon)*100.0)*RG)
[781]492
493
494! PF et PASB
495!   where(cdrag > 0.01)
496!     cdrag = 0.01
497!   endwhere
498!  write(*,*)'Cdrag = ',minval(cdrag),maxval(cdrag)
499
500 
501    IF (debut) THEN
[987]502       CALL Init_orchidee_index(knon,knindex,offset,ktindex)
[3435]503       CALL Get_orchidee_communicator(orch_comm,orch_mpi_size,orch_mpi_rank, orch_omp_size,orch_omp_rank)
504
505       IF (grid_type==unstructured) THEN
506         IF (knon==0) THEN
507           begin=1
508           end=0
509         ELSE
510           begin=offset+1
511           end=offset+ktindex(knon)
512         ENDIF
513       
514         IF (orch_mpi_rank==orch_mpi_size-1 .AND. orch_omp_rank==orch_omp_size-1) end=nbp_lon*nbp_lat
515         
516         ALLOCATE(lalo(end-begin+1,2))
517         ALLOCATE(bounds_lalo(end-begin+1,nvertex,2))
518         ALLOCATE(ind_cell(end-begin+1))
519         
520         ALLOCATE(longitude_glo(klon_glo))
521         CALL gather(longitude,longitude_glo)
522         CALL bcast(longitude_glo)
523         lalo(:,2)=longitude_glo(begin:end)*180./PI
524 
525         ALLOCATE(latitude_glo(klon_glo))
526         CALL gather(latitude,latitude_glo)
527         CALL bcast(latitude_glo)
528         lalo(:,1)=latitude_glo(begin:end)*180./PI
529
530         ALLOCATE(boundslon_glo(klon_glo,nvertex))
531         CALL gather(boundslon,boundslon_glo)
532         CALL bcast(boundslon_glo)
533         bounds_lalo(:,:,2)=boundslon_glo(begin:end,:)*180./PI
534 
535         ALLOCATE(boundslat_glo(klon_glo,nvertex))
536         CALL gather(boundslat,boundslat_glo)
537         CALL bcast(boundslat_glo)
538         bounds_lalo(:,:,1)=boundslat_glo(begin:end,:)*180./PI
539         
540         ALLOCATE(ind_cell_glo_glo(klon_glo))
541         CALL gather(ind_cell_glo,ind_cell_glo_glo)
542         CALL bcast(ind_cell_glo_glo)
543         ind_cell(:)=ind_cell_glo_glo(begin:end)
544         
545       ENDIF
[987]546       CALL Init_synchro_omp
[3504]547
548!$OMP BARRIER
[987]549       
[1023]550       IF (knon > 0) THEN
[1067]551#ifdef CPP_VEGET
[3435]552         CALL Init_intersurf(nbp_lon,nbp_lat,knon,ktindex,offset,orch_omp_size,orch_omp_rank,orch_comm,grid=grid_type)
[1067]553#endif
[987]554       ENDIF
[802]555
[3504]556       CALL Synchro_omp
557
[987]558       
[3435]559       IF (knon > 0) THEN
[987]560
[1067]561#ifdef CPP_VEGET
[3435]562
[2571]563         CALL intersurf_initialize_gathered (itime+itau_phy-1, nbp_lon, nbp_lat, knon, ktindex, dtime, &
564               lrestart_read, lrestart_write, lalo, contfrac, neighbours, resolution, date0, &
565               zlev,  u1_lay, v1_lay, spechum, temp_air, epot_air, &
[802]566               cdrag, petA_orc, peqA_orc, petB_orc, peqB_orc, &
567               precip_rain, precip_snow, lwdown, swnet, swdown, ps, &
568               evap, fluxsens, fluxlat, coastalflow, riverflow, &
[2952]569               tsol_rad, tsurf_new, qsurf, albedo_out, emis_new, z0m_new, &   
[3435]570               lon_scat, lat_scat, q2m(1:knon), t2m(1:knon), z0h_new(1:knon), nvm_orch, &
571               grid=grid_type, bounds_latlon=bounds_lalo, cell_area=area, ind_cell_glo=ind_cell, &
[3391]572               field_out_names=cfname_out, field_in_names=cfname_in(1:nbcf_in_orc))
[1067]573#endif         
[781]574       ENDIF
575
[987]576       CALL Synchro_omp
577
[781]578       albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
579
580    ENDIF
[987]581   
[781]582!  swdown_vrai(1:knon) = swnet(1:knon)/(1. - albedo_keep(1:knon))
583    swdown_vrai(1:knon) = swdown(1:knon)
[3465]584!$OMP BARRIER
[781]585
[987]586    IF (knon > 0) THEN
[1067]587#ifdef CPP_VEGET   
[2952]588       IF (nvm_orch .NE. nvm_lmdz ) THEN
589          abort_message='Pb de dimensiosn PFT: nvm_orch et nvm_lmdz differents.'
590          CALL abort_physic(modname,abort_message,1)
591       ENDIF
592
[2571]593       CALL intersurf_main_gathered (itime+itau_phy, nbp_lon, nbp_lat, knon, ktindex, dtime,  &
[802]594            lrestart_read, lrestart_write, lalo, &
595            contfrac, neighbours, resolution, date0, &
[781]596            zlev,  u1_lay(1:knon), v1_lay(1:knon), spechum(1:knon), temp_air(1:knon), epot_air(1:knon), ccanopy(1:knon), &
597            cdrag(1:knon), petA_orc(1:knon), peqA_orc(1:knon), petB_orc(1:knon), peqB_orc(1:knon), &
598            precip_rain(1:knon), precip_snow(1:knon), lwdown(1:knon), swnet(1:knon), swdown_vrai(1:knon), ps(1:knon), &
599            evap(1:knon), fluxsens(1:knon), fluxlat(1:knon), coastalflow(1:knon), riverflow(1:knon), &
[2571]600            tsol_rad(1:knon), tsurf_new(1:knon), qsurf(1:knon), albedo_out(1:knon,:), emis_new(1:knon), z0m_new(1:knon), &
[3435]601            lon_scat, lat_scat, q2m(1:knon), t2m(1:knon), z0h_new(1:knon),&
[2952]602            veget(1:knon,:),lai(1:knon,:),height(1:knon,:),&
[3391]603            fields_out=yfields_out(1:knon,1:nbcf_out),  &
604            fields_in=yfields_in(1:knon,1:nbcf_in_orc), &
[2952]605            coszang=yrmu0(1:knon))
[1067]606#endif       
[781]607    ENDIF
608
[987]609    CALL Synchro_omp
610   
[781]611    albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
612
613!* Send to coupler
614!
[996]615    IF (type_ocean=='couple') THEN
[781]616       CALL cpl_send_land_fields(itime, knon, knindex, &
617            riverflow, coastalflow)
618    ENDIF
619
[888]620    alb1_new(1:knon) = albedo_out(1:knon,1)
621    alb2_new(1:knon) = albedo_out(1:knon,2)
[781]622
623! Convention orchidee: positif vers le haut
624    fluxsens(1:knon) = -1. * fluxsens(1:knon)
625    fluxlat(1:knon)  = -1. * fluxlat(1:knon)
626   
627!  evap     = -1. * evap
628
629    IF (debut) lrestart_read = .FALSE.
630   
[987]631    IF (debut) CALL Finalize_surf_para
[1279]632
[3391]633! >> PC
634! Decompressing variables into LMDz for the module carbon_cycle_mod
635! nbcf_in can be zero, in which case the loop does not operate
636! fields_in can then used elsewhere in the model
637     
638     fields_in(:,:)=0.0
639
640     DO nb=1, nbcf_in_orc
641       DO igrid = 1, knon
642        ireal = knindex(igrid)
643        fields_in(ireal,nb)=yfields_in(igrid,nb)
644       ENDDO
645       WRITE(*,*) 'surf_land_orchidee_mod --- yfields_in :',cfname_in(nb)
646     ENDDO
647! >> PC
[987]648   
[781]649  END SUBROUTINE surf_land_orchidee
650!
651!****************************************************************************************
652!
[987]653  SUBROUTINE Init_orchidee_index(knon,knindex,offset,ktindex)
654  USE mod_surf_para
655  USE mod_grid_phy_lmdz
656 
657    INTEGER,INTENT(IN)    :: knon
658    INTEGER,INTENT(IN)    :: knindex(klon)   
659    INTEGER,INTENT(OUT)   :: offset
660    INTEGER,INTENT(OUT)   :: ktindex(klon)
[781]661   
[987]662    INTEGER               :: ktindex_glo(knon_glo)
663    INTEGER               :: offset_para(0:omp_size*mpi_size-1)
664    INTEGER               :: LastPoint
665    INTEGER               :: task
[781]666   
[987]667    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
[781]668   
[987]669    CALL gather_surf(ktindex(1:knon),ktindex_glo)
670   
671    IF (is_mpi_root .AND. is_omp_root) THEN
672      LastPoint=0
673      DO Task=0,mpi_size*omp_size-1
674        IF (knon_glo_para(Task)>0) THEN
675           offset_para(task)= LastPoint-MOD(LastPoint,nbp_lon)
676           LastPoint=ktindex_glo(knon_glo_end_para(task))
677        ENDIF
678      ENDDO
[781]679    ENDIF
680   
[987]681    CALL bcast(offset_para)
[781]682   
[987]683    offset=offset_para(omp_size*mpi_rank+omp_rank)
684   
685    ktindex(1:knon)=ktindex(1:knon)-offset
[781]686
[987]687  END SUBROUTINE Init_orchidee_index
688
[781]689!
[987]690!************************* ***************************************************************
[781]691!
[987]692
[3435]693  SUBROUTINE Get_orchidee_communicator(orch_comm, orch_mpi_size, orch_mpi_rank, orch_omp_size,orch_omp_rank)
[987]694  USE  mod_surf_para
695     
[1001]696#ifdef CPP_MPI
[781]697    INCLUDE 'mpif.h'
698#endif   
699
700    INTEGER,INTENT(OUT) :: orch_comm
[3435]701    INTEGER,INTENT(OUT) :: orch_mpi_size
702    INTEGER,INTENT(OUT) :: orch_mpi_rank
[987]703    INTEGER,INTENT(OUT) :: orch_omp_size
704    INTEGER,INTENT(OUT) :: orch_omp_rank
[781]705    INTEGER             :: color
[987]706    INTEGER             :: i,ierr
[802]707!
708! End definition
709!****************************************************************************************
[781]710   
[987]711    IF (is_omp_root) THEN         
712     
713      IF (knon_mpi==0) THEN
714         color = 0
715      ELSE
716         color = 1
717      ENDIF
718   
[1001]719#ifdef CPP_MPI   
[987]720      CALL MPI_COMM_SPLIT(COMM_LMDZ_PHY,color,mpi_rank,orch_comm,ierr)
[3435]721      CALL MPI_COMM_SIZE(orch_comm,orch_mpi_size,ierr)
722      CALL MPI_COMM_RANK(orch_comm,orch_mpi_rank,ierr)
[781]723#endif
[1023]724   
725    ENDIF
726    CALL bcast_omp(orch_comm)
727   
728    IF (knon_mpi /= 0) THEN
729      orch_omp_size=0
730      DO i=0,omp_size-1
731        IF (knon_omp_para(i) /=0) THEN
732          orch_omp_size=orch_omp_size+1
733          IF (i==omp_rank) orch_omp_rank=orch_omp_size-1
734        ENDIF
735      ENDDO
736    ENDIF
[781]737   
738  END SUBROUTINE Get_orchidee_communicator
739!
740!****************************************************************************************
741
[987]742
743  SUBROUTINE Init_neighbours(knon,neighbours,knindex,pctsrf)
744    USE mod_grid_phy_lmdz
745    USE mod_surf_para   
[1785]746    USE indice_sol_mod
[987]747
[1001]748#ifdef CPP_MPI
[781]749    INCLUDE 'mpif.h'
750#endif   
751
[802]752! Input arguments
753!****************************************************************************************
[781]754    INTEGER, INTENT(IN)                     :: knon
[987]755    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
[781]756    REAL, DIMENSION(klon), INTENT(IN)       :: pctsrf
757   
[802]758! Output arguments
759!****************************************************************************************
760    INTEGER, DIMENSION(knon,8), INTENT(OUT) :: neighbours
761
762! Local variables
763!****************************************************************************************
764    INTEGER                              :: i, igrid, jj, ij, iglob
765    INTEGER                              :: ierr, ireal, index
766    INTEGER, DIMENSION(8,3)              :: off_ini
767    INTEGER, DIMENSION(8)                :: offset 
[987]768    INTEGER, DIMENSION(nbp_lon,nbp_lat)  :: correspond
769    INTEGER, DIMENSION(knon_glo)         :: ktindex_glo
770    INTEGER, DIMENSION(knon_glo,8)       :: neighbours_glo
771    REAL, DIMENSION(klon_glo)            :: pctsrf_glo
772    INTEGER                              :: ktindex(klon)
[802]773!
774! End definition
775!****************************************************************************************
776
[987]777    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
[781]778   
[987]779    CALL gather_surf(ktindex(1:knon),ktindex_glo)
780    CALL gather(pctsrf,pctsrf_glo)
[781]781   
[987]782    IF (is_mpi_root .AND. is_omp_root) THEN
783      neighbours_glo(:,:)=-1
[781]784!  Initialisation des offset   
785!
786! offset bord ouest
[987]787       off_ini(1,1) = - nbp_lon   ; off_ini(2,1) = - nbp_lon + 1     ; off_ini(3,1) = 1
788       off_ini(4,1) = nbp_lon + 1 ; off_ini(5,1) = nbp_lon           ; off_ini(6,1) = 2 * nbp_lon - 1
789       off_ini(7,1) = nbp_lon -1  ; off_ini(8,1) = - 1
[781]790! offset point normal
[987]791       off_ini(1,2) = - nbp_lon   ; off_ini(2,2) = - nbp_lon + 1     ; off_ini(3,2) = 1
792       off_ini(4,2) = nbp_lon + 1 ; off_ini(5,2) = nbp_lon           ; off_ini(6,2) = nbp_lon - 1
793       off_ini(7,2) = -1          ; off_ini(8,2) = - nbp_lon - 1
[781]794! offset bord   est
[987]795       off_ini(1,3) = - nbp_lon   ; off_ini(2,3) = - 2 * nbp_lon + 1 ; off_ini(3,3) = - nbp_lon + 1
796       off_ini(4,3) =  1          ; off_ini(5,3) = nbp_lon           ; off_ini(6,3) = nbp_lon - 1
797       off_ini(7,3) = -1          ; off_ini(8,3) = - nbp_lon - 1
[781]798!
799! Attention aux poles
800!
[987]801       DO igrid = 1, knon_glo
802          index = ktindex_glo(igrid)
803          jj = INT((index - 1)/nbp_lon) + 1
804          ij = index - (jj - 1) * nbp_lon
[781]805          correspond(ij,jj) = igrid
806       ENDDO
[2126]807!sonia : Les mailles des voisines doivent etre toutes egales (pour couplage orchidee)
808       IF (knon_glo == 1) THEN
809         igrid = 1
810         DO i = 1,8
811           neighbours_glo(igrid, i) = igrid
812         ENDDO
813       ELSE
[781]814       
[987]815       DO igrid = 1, knon_glo
816          iglob = ktindex_glo(igrid)
817         
818          IF (MOD(iglob, nbp_lon) == 1) THEN
[781]819             offset = off_ini(:,1)
[987]820          ELSE IF(MOD(iglob, nbp_lon) == 0) THEN
[781]821             offset = off_ini(:,3)
822          ELSE
823             offset = off_ini(:,2)
824          ENDIF
[987]825         
[781]826          DO i = 1, 8
827             index = iglob + offset(i)
[987]828             ireal = (MIN(MAX(1, index - nbp_lon + 1), klon_glo))
829             IF (pctsrf_glo(ireal) > EPSFRA) THEN
830                jj = INT((index - 1)/nbp_lon) + 1
831                ij = index - (jj - 1) * nbp_lon
832                neighbours_glo(igrid, i) = correspond(ij, jj)
[781]833             ENDIF
834          ENDDO
835       ENDDO
[2126]836       ENDIF !fin knon_glo == 1
[781]837
838    ENDIF
839   
[987]840    DO i = 1, 8
841      CALL scatter_surf(neighbours_glo(:,i),neighbours(1:knon,i))
[781]842    ENDDO
843  END SUBROUTINE Init_neighbours
[987]844
[781]845!
846!****************************************************************************************
847!
[1146]848#endif
[2571]849#endif
[2952]850#endif
[3435]851#endif
[781]852END MODULE surf_land_orchidee_mod
Note: See TracBrowser for help on using the repository browser.