source: LMDZ6/trunk/libf/phylmd/surf_land_orchidee_nolic_mod.F90 @ 5831

Last change on this file since 5831 was 5823, checked in by fcheruy, 2 months ago

pour tourner en 1D avec Orchidee apres modif LMDZ pour separer dynamique et physique (rev. 2344)

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