source: LMDZ4/trunk/libf/phylmd/surf_land_orchidee_mod.F90 @ 1563

Last change on this file since 1563 was 1279, checked in by Laurent Fairhead, 15 years ago

Merged LMDZ4-dev branch changes r1241:1278 into the trunk
Running trunk and LMDZ4-dev in LMDZOR configuration on local
machine (sequential) and SX8 (4-proc) yields identical results
(restart and restartphy are identical binarily)
Log history from r1241 to r1278 is available by switching to
source:LMDZ4/branches/LMDZ4-dev-20091210

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 22.7 KB
Line 
1!
2MODULE surf_land_orchidee_mod
3#ifndef ORCHIDEE_NOOPENMP
4!
5! This module controles the interface towards the model ORCHIDEE
6!
7! Subroutines in this module : surf_land_orchidee
8!                              Init_orchidee_index
9!                              Get_orchidee_communicator
10!                              Init_neighbours
11
12  USE dimphy
13#ifdef CPP_VEGET
14  USE intersurf     ! module d'ORCHIDEE
15#endif
16  USE cpl_mod,      ONLY : cpl_send_land_fields
17  USE surface_data, ONLY : type_ocean
18  USE comgeomphy,   ONLY : cuphy, cvphy
19  USE mod_grid_phy_lmdz
20  USE mod_phys_lmdz_para, mpi_root_rank=>mpi_root
21
22  IMPLICIT NONE
23
24  PRIVATE
25  PUBLIC  :: surf_land_orchidee
26
27  LOGICAL, ALLOCATABLE, SAVE :: flag_omp(:)
28CONTAINS
29!
30!****************************************************************************************
31
32  SUBROUTINE surf_land_orchidee(itime, dtime, date0, knon, &
33       knindex, rlon, rlat, pctsrf, &
34       debut, lafin, &
35       plev,  u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &
36       tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
37       precip_rain, precip_snow, lwdown, swnet, swdown, &
38       ps, q2m, t2m, &
39       evap, fluxsens, fluxlat, &             
40       tsol_rad, tsurf_new, alb1_new, alb2_new, &
41       emis_new, z0_new, qsurf)
42
43    USE mod_surf_para
44    USE mod_synchro_omp
45   
46USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, fco2_land_inst, fco2_lu_inst
47
48!   
49! Cette routine sert d'interface entre le modele atmospherique et le
50! modele de sol continental. Appel a sechiba
51!
52! L. Fairhead 02/2000
53!
54! input:
55!   itime        numero du pas de temps
56!   dtime        pas de temps de la physique (en s)
57!   nisurf       index de la surface a traiter (1 = sol continental)
58!   knon         nombre de points de la surface a traiter
59!   knindex      index des points de la surface a traiter
60!   rlon         longitudes de la grille entiere
61!   rlat         latitudes de la grille entiere
62!   pctsrf       tableau des fractions de surface de chaque maille
63!   debut        logical: 1er appel a la physique (lire les restart)
64!   lafin        logical: dernier appel a la physique (ecrire les restart)
65!                     (si false calcul simplifie des fluxs sur les continents)
66!   plev         hauteur de la premiere couche (Pa)     
67!   u1_lay       vitesse u 1ere couche
68!   v1_lay       vitesse v 1ere couche
69!   temp_air     temperature de l'air 1ere couche
70!   spechum      humidite specifique 1ere couche
71!   epot_air     temp pot de l'air
72!   ccanopy      concentration CO2 canopee, correspond au co2_send de
73!                carbon_cycle_mod ou valeur constant co2_ppm
74!   tq_cdrag     cdrag
75!   petAcoef     coeff. A de la resolution de la CL pour t
76!   peqAcoef     coeff. A de la resolution de la CL pour q
77!   petBcoef     coeff. B de la resolution de la CL pour t
78!   peqBcoef     coeff. B de la resolution de la CL pour q
79!   precip_rain  precipitation liquide
80!   precip_snow  precipitation solide
81!   lwdown       flux IR descendant a la surface
82!   swnet        flux solaire net
83!   swdown       flux solaire entrant a la surface
84!   ps           pression au sol
85!   radsol       rayonnement net aus sol (LW + SW)
86!   
87!
88! output:
89!   evap         evaporation totale
90!   fluxsens     flux de chaleur sensible
91!   fluxlat      flux de chaleur latente
92!   tsol_rad     
93!   tsurf_new    temperature au sol
94!   alb1_new     albedo in visible SW interval
95!   alb2_new     albedo in near IR interval
96!   emis_new     emissivite
97!   z0_new       surface roughness
98!   qsurf        air moisture at surface
99!
100    INCLUDE "indicesol.h"
101    INCLUDE "temps.h"
102    INCLUDE "YOMCST.h"
103    INCLUDE "iniprint.h"
104    INCLUDE "dimensions.h"
105 
106!
107! Parametres d'entree
108!****************************************************************************************
109    INTEGER, INTENT(IN)                       :: itime
110    REAL, INTENT(IN)                          :: dtime
111    REAL, INTENT(IN)                          :: date0
112    INTEGER, INTENT(IN)                       :: knon
113    INTEGER, DIMENSION(klon), INTENT(IN)      :: knindex
114    LOGICAL, INTENT(IN)                       :: debut, lafin
115    REAL, DIMENSION(klon,nbsrf), INTENT(IN)   :: pctsrf
116    REAL, DIMENSION(klon), INTENT(IN)         :: rlon, rlat
117    REAL, DIMENSION(klon), INTENT(IN)         :: plev
118    REAL, DIMENSION(klon), INTENT(IN)         :: u1_lay, v1_lay
119    REAL, DIMENSION(klon), INTENT(IN)         :: temp_air, spechum
120    REAL, DIMENSION(klon), INTENT(IN)         :: epot_air, ccanopy
121    REAL, DIMENSION(klon), INTENT(IN)         :: tq_cdrag
122    REAL, DIMENSION(klon), INTENT(IN)         :: petAcoef, peqAcoef
123    REAL, DIMENSION(klon), INTENT(IN)         :: petBcoef, peqBcoef
124    REAL, DIMENSION(klon), INTENT(IN)         :: precip_rain, precip_snow
125    REAL, DIMENSION(klon), INTENT(IN)         :: lwdown, swnet, swdown, ps
126    REAL, DIMENSION(klon), INTENT(IN)         :: q2m, t2m
127
128! Parametres de sortie
129!****************************************************************************************
130    REAL, DIMENSION(klon), INTENT(OUT)        :: evap, fluxsens, fluxlat, qsurf
131    REAL, DIMENSION(klon), INTENT(OUT)        :: tsol_rad, tsurf_new
132    REAL, DIMENSION(klon), INTENT(OUT)        :: alb1_new, alb2_new
133    REAL, DIMENSION(klon), INTENT(OUT)        :: emis_new, z0_new
134
135! Local
136!****************************************************************************************
137    INTEGER                                   :: ij, jj, igrid, ireal, index
138    INTEGER                                   :: error
139    REAL, DIMENSION(klon)                     :: swdown_vrai
140    REAL, DIMENSION(klon)                     :: fco2_land_comp  ! sur grille compresse
141    REAL, DIMENSION(klon)                     :: fco2_lu_comp    ! sur grille compresse
142    CHARACTER (len = 20)                      :: modname = 'surf_land_orchidee'
143    CHARACTER (len = 80)                      :: abort_message
144    LOGICAL,SAVE                              :: check = .FALSE.
145    !$OMP THREADPRIVATE(check)
146
147! type de couplage dans sechiba
148!  character (len=10)   :: coupling = 'implicit'
149! drapeaux controlant les appels dans SECHIBA
150!  type(control_type), save   :: control_in
151! Preserved albedo
152    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: albedo_keep, zlev
153    !$OMP THREADPRIVATE(albedo_keep,zlev)
154! coordonnees geographiques
155    REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: lalo
156    !$OMP THREADPRIVATE(lalo)
157! pts voisins
158    INTEGER,ALLOCATABLE, DIMENSION(:,:), SAVE :: neighbours
159    !$OMP THREADPRIVATE(neighbours)
160! fractions continents
161    REAL,ALLOCATABLE, DIMENSION(:), SAVE      :: contfrac
162    !$OMP THREADPRIVATE(contfrac)
163! resolution de la grille
164    REAL, ALLOCATABLE, DIMENSION (:,:), SAVE  :: resolution
165    !$OMP THREADPRIVATE(resolution)
166
167    REAL, ALLOCATABLE, DIMENSION (:,:), SAVE  :: lon_scat, lat_scat 
168    !$OMP THREADPRIVATE(lon_scat,lat_scat)
169
170    LOGICAL, SAVE                             :: lrestart_read = .TRUE.
171    !$OMP THREADPRIVATE(lrestart_read)
172    LOGICAL, SAVE                             :: lrestart_write = .FALSE.
173    !$OMP THREADPRIVATE(lrestart_write)
174
175    REAL, DIMENSION(knon,2)                   :: albedo_out
176
177! Pb de nomenclature
178    REAL, DIMENSION(klon)                     :: petA_orc, peqA_orc
179    REAL, DIMENSION(klon)                     :: petB_orc, peqB_orc
180! Pb de correspondances de grilles
181    INTEGER, DIMENSION(:), SAVE, ALLOCATABLE  :: ig, jg
182    !$OMP THREADPRIVATE(ig,jg)
183    INTEGER :: indi, indj
184    INTEGER, SAVE, ALLOCATABLE,DIMENSION(:)   :: ktindex
185    !$OMP THREADPRIVATE(ktindex)
186
187! Essai cdrag
188    REAL, DIMENSION(klon)                     :: cdrag
189    INTEGER,SAVE                              :: offset
190    !$OMP THREADPRIVATE(offset)
191
192    REAL, DIMENSION(klon_glo)                 :: rlon_g,rlat_g
193    INTEGER, SAVE                             :: orch_comm
194    !$OMP THREADPRIVATE(orch_comm)
195
196    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: coastalflow
197    !$OMP THREADPRIVATE(coastalflow)
198    REAL, ALLOCATABLE, DIMENSION(:), SAVE     :: riverflow
199    !$OMP THREADPRIVATE(riverflow)
200   
201    INTEGER :: orch_omp_rank
202    INTEGER :: orch_omp_size
203!
204! Fin definition
205!****************************************************************************************
206
207    IF (check) WRITE(lunout,*)'Entree ', modname
208 
209! Initialisation
210 
211    IF (debut) THEN
212! Test of coherence between variable ok_veget and cpp key CPP_VEGET
213#ifndef CPP_VEGET
214       abort_message='Pb de coherence: ok_veget = .true. mais CPP_VEGET = .false.'
215       CALL abort_gcm(modname,abort_message,1)
216#endif
217
218       CALL Init_surf_para(knon)
219       ALLOCATE(ktindex(knon))
220       IF ( .NOT. ALLOCATED(albedo_keep)) THEN
221!ym          ALLOCATE(albedo_keep(klon))
222!ym bizarre que non alloué en knon precedement
223          ALLOCATE(albedo_keep(knon))
224          ALLOCATE(zlev(knon))
225       ENDIF
226! Pb de correspondances de grilles
227       ALLOCATE(ig(klon))
228       ALLOCATE(jg(klon))
229       ig(1) = 1
230       jg(1) = 1
231       indi = 0
232       indj = 2
233       DO igrid = 2, klon - 1
234          indi = indi + 1
235          IF ( indi > iim) THEN
236             indi = 1
237             indj = indj + 1
238          ENDIF
239          ig(igrid) = indi
240          jg(igrid) = indj
241       ENDDO
242       ig(klon) = 1
243       jg(klon) = jjm + 1
244
245       IF ((.NOT. ALLOCATED(lalo))) THEN
246          ALLOCATE(lalo(knon,2), stat = error)
247          IF (error /= 0) THEN
248             abort_message='Pb allocation lalo'
249             CALL abort_gcm(modname,abort_message,1)
250          ENDIF
251       ENDIF
252       IF ((.NOT. ALLOCATED(lon_scat))) THEN
253          ALLOCATE(lon_scat(iim,jjm+1), stat = error)
254          IF (error /= 0) THEN
255             abort_message='Pb allocation lon_scat'
256             CALL abort_gcm(modname,abort_message,1)
257          ENDIF
258       ENDIF
259       IF ((.NOT. ALLOCATED(lat_scat))) THEN
260          ALLOCATE(lat_scat(iim,jjm+1), stat = error)
261          IF (error /= 0) THEN
262             abort_message='Pb allocation lat_scat'
263             CALL abort_gcm(modname,abort_message,1)
264          ENDIF
265       ENDIF
266       lon_scat = 0.
267       lat_scat = 0.
268       DO igrid = 1, knon
269          index = knindex(igrid)
270          lalo(igrid,2) = rlon(index)
271          lalo(igrid,1) = rlat(index)
272       ENDDO
273
274       
275       
276       CALL Gather(rlon,rlon_g)
277       CALL Gather(rlat,rlat_g)
278
279       IF (is_mpi_root) THEN
280          index = 1
281          DO jj = 2, jjm
282             DO ij = 1, iim
283                index = index + 1
284                lon_scat(ij,jj) = rlon_g(index)
285                lat_scat(ij,jj) = rlat_g(index)
286             ENDDO
287          ENDDO
288          lon_scat(:,1) = lon_scat(:,2)
289          lat_scat(:,1) = rlat_g(1)
290          lon_scat(:,jjm+1) = lon_scat(:,2)
291          lat_scat(:,jjm+1) = rlat_g(klon_glo)
292       ENDIF
293   
294       CALL bcast(lon_scat)
295       CALL bcast(lat_scat)
296!
297! Allouer et initialiser le tableau des voisins et des fraction de continents
298!
299       IF ( (.NOT.ALLOCATED(neighbours))) THEN
300          ALLOCATE(neighbours(knon,8), stat = error)
301          IF (error /= 0) THEN
302             abort_message='Pb allocation neighbours'
303             CALL abort_gcm(modname,abort_message,1)
304          ENDIF
305       ENDIF
306       neighbours = -1.
307       IF (( .NOT. ALLOCATED(contfrac))) THEN
308          ALLOCATE(contfrac(knon), stat = error)
309          IF (error /= 0) THEN
310             abort_message='Pb allocation contfrac'
311             CALL abort_gcm(modname,abort_message,1)
312          ENDIF
313       ENDIF
314
315       DO igrid = 1, knon
316          ireal = knindex(igrid)
317          contfrac(igrid) = pctsrf(ireal,is_ter)
318       ENDDO
319
320
321       CALL Init_neighbours(knon,neighbours,knindex,pctsrf(:,is_ter))
322
323!
324!  Allocation et calcul resolutions
325       IF ( (.NOT.ALLOCATED(resolution))) THEN
326          ALLOCATE(resolution(knon,2), stat = error)
327          IF (error /= 0) THEN
328             abort_message='Pb allocation resolution'
329             CALL abort_gcm(modname,abort_message,1)
330          ENDIF
331       ENDIF
332       DO igrid = 1, knon
333          ij = knindex(igrid)
334          resolution(igrid,1) = cuphy(ij)
335          resolution(igrid,2) = cvphy(ij)
336       ENDDO
337     
338       ALLOCATE(coastalflow(klon), stat = error)
339       IF (error /= 0) THEN
340          abort_message='Pb allocation coastalflow'
341          CALL abort_gcm(modname,abort_message,1)
342       ENDIF
343       
344       ALLOCATE(riverflow(klon), stat = error)
345       IF (error /= 0) THEN
346          abort_message='Pb allocation riverflow'
347          CALL abort_gcm(modname,abort_message,1)
348       ENDIF
349!
350! Allocate variables needed for carbon_cycle_mod
351!
352       IF (carbon_cycle_cpl) THEN
353          IF (.NOT. ALLOCATED(fco2_land_inst)) THEN
354             ALLOCATE(fco2_land_inst(klon),stat=error)
355             IF (error /= 0)  CALL abort_gcm(modname,'Pb in allocation fco2_land_inst',1)
356             
357             ALLOCATE(fco2_lu_inst(klon),stat=error)
358             IF(error /=0) CALL abort_gcm(modname,'Pb in allocation fco2_lu_inst',1)
359          END IF
360       END IF
361       
362    ENDIF                          ! (fin debut)
363 
364
365!
366! Appel a la routine sols continentaux
367!
368    IF (lafin) lrestart_write = .TRUE.
369    IF (check) WRITE(lunout,*)'lafin ',lafin,lrestart_write
370     
371    petA_orc(1:knon) = petBcoef(1:knon) * dtime
372    petB_orc(1:knon) = petAcoef(1:knon)
373    peqA_orc(1:knon) = peqBcoef(1:knon) * dtime
374    peqB_orc(1:knon) = peqAcoef(1:knon)
375
376    cdrag = 0.
377    cdrag(1:knon) = tq_cdrag(1:knon)
378
379! zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/287.05*temp_air(1:knon))*9.80665)
380    zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/RD*temp_air(1:knon))*RG)
381
382
383! PF et PASB
384!   where(cdrag > 0.01)
385!     cdrag = 0.01
386!   endwhere
387!  write(*,*)'Cdrag = ',minval(cdrag),maxval(cdrag)
388
389 
390    IF (debut) THEN
391       CALL Init_orchidee_index(knon,knindex,offset,ktindex)
392       CALL Get_orchidee_communicator(orch_comm,orch_omp_size,orch_omp_rank)
393       CALL Init_synchro_omp
394       
395       IF (knon > 0) THEN
396#ifdef CPP_VEGET
397         CALL Init_intersurf(nbp_lon,nbp_lat,knon,ktindex,offset,orch_omp_size,orch_omp_rank,orch_comm)
398#endif
399       ENDIF
400
401       
402       IF (knon > 0) THEN
403
404#ifdef CPP_VEGET
405          CALL intersurf_main (itime+itau_phy-1, iim, jjm+1, knon, ktindex, dtime, &
406               lrestart_read, lrestart_write, lalo, &
407               contfrac, neighbours, resolution, date0, &
408               zlev,  u1_lay, v1_lay, spechum, temp_air, epot_air, ccanopy, &
409               cdrag, petA_orc, peqA_orc, petB_orc, peqB_orc, &
410               precip_rain, precip_snow, lwdown, swnet, swdown, ps, &
411               evap, fluxsens, fluxlat, coastalflow, riverflow, &
412               tsol_rad, tsurf_new, qsurf, albedo_out, emis_new, z0_new, &
413               lon_scat, lat_scat, q2m, t2m)
414#endif         
415       ENDIF
416
417       CALL Synchro_omp
418
419       albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
420
421    ENDIF
422
423   
424!  swdown_vrai(1:knon) = swnet(1:knon)/(1. - albedo_keep(1:knon))
425    swdown_vrai(1:knon) = swdown(1:knon)
426
427    IF (knon > 0) THEN
428#ifdef CPP_VEGET   
429       CALL intersurf_main (itime+itau_phy, iim, jjm+1, knon, ktindex, dtime,  &
430            lrestart_read, lrestart_write, lalo, &
431            contfrac, neighbours, resolution, date0, &
432            zlev,  u1_lay(1:knon), v1_lay(1:knon), spechum(1:knon), temp_air(1:knon), epot_air(1:knon), ccanopy(1:knon), &
433            cdrag(1:knon), petA_orc(1:knon), peqA_orc(1:knon), petB_orc(1:knon), peqB_orc(1:knon), &
434            precip_rain(1:knon), precip_snow(1:knon), lwdown(1:knon), swnet(1:knon), swdown_vrai(1:knon), ps(1:knon), &
435            evap(1:knon), fluxsens(1:knon), fluxlat(1:knon), coastalflow(1:knon), riverflow(1:knon), &
436            tsol_rad(1:knon), tsurf_new(1:knon), qsurf(1:knon), albedo_out(1:knon,:), emis_new(1:knon), z0_new(1:knon), &
437            lon_scat, lat_scat, q2m, t2m)
438#endif       
439    ENDIF
440
441    CALL Synchro_omp
442   
443    albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
444
445!* Send to coupler
446!
447    IF (type_ocean=='couple') THEN
448       CALL cpl_send_land_fields(itime, knon, knindex, &
449            riverflow, coastalflow)
450    ENDIF
451
452    alb1_new(1:knon) = albedo_out(1:knon,1)
453    alb2_new(1:knon) = albedo_out(1:knon,2)
454
455! Convention orchidee: positif vers le haut
456    fluxsens(1:knon) = -1. * fluxsens(1:knon)
457    fluxlat(1:knon)  = -1. * fluxlat(1:knon)
458   
459!  evap     = -1. * evap
460
461    IF (debut) lrestart_read = .FALSE.
462   
463    IF (debut) CALL Finalize_surf_para
464
465   
466! JG : TEMPORAIRE!!!! Les variables fco2_land_comp et fco2_lu_comp seront plus tard en sortie d'ORCHIDEE
467!      ici mis a valeur quelquonque pour test. Ces variables sont sur la grille compresse avec uniquement des points de terres
468
469    fco2_land_comp(:) = 1.
470    fco2_lu_comp(:)   = 10.
471
472! Decompress variables for the module carbon_cycle_mod
473    IF (carbon_cycle_cpl) THEN
474       fco2_land_inst(:)=0.
475       fco2_lu_inst(:)=0.
476       
477       DO igrid = 1, knon
478          ireal = knindex(igrid)
479          fco2_land_inst(ireal) = fco2_land_comp(igrid)
480          fco2_lu_inst(ireal)   = fco2_lu_comp(igrid)
481       END DO
482    END IF
483
484  END SUBROUTINE surf_land_orchidee
485!
486!****************************************************************************************
487!
488  SUBROUTINE Init_orchidee_index(knon,knindex,offset,ktindex)
489  USE mod_surf_para
490  USE mod_grid_phy_lmdz
491 
492    INTEGER,INTENT(IN)    :: knon
493    INTEGER,INTENT(IN)    :: knindex(klon)   
494    INTEGER,INTENT(OUT)   :: offset
495    INTEGER,INTENT(OUT)   :: ktindex(klon)
496   
497    INTEGER               :: ktindex_glo(knon_glo)
498    INTEGER               :: offset_para(0:omp_size*mpi_size-1)
499    INTEGER               :: LastPoint
500    INTEGER               :: task
501   
502    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
503   
504    CALL gather_surf(ktindex(1:knon),ktindex_glo)
505   
506    IF (is_mpi_root .AND. is_omp_root) THEN
507      LastPoint=0
508      DO Task=0,mpi_size*omp_size-1
509        IF (knon_glo_para(Task)>0) THEN
510           offset_para(task)= LastPoint-MOD(LastPoint,nbp_lon)
511           LastPoint=ktindex_glo(knon_glo_end_para(task))
512        ENDIF
513      ENDDO
514    ENDIF
515   
516    CALL bcast(offset_para)
517   
518    offset=offset_para(omp_size*mpi_rank+omp_rank)
519   
520    ktindex(1:knon)=ktindex(1:knon)-offset
521
522  END SUBROUTINE Init_orchidee_index
523
524!
525!************************* ***************************************************************
526!
527
528  SUBROUTINE Get_orchidee_communicator(orch_comm,orch_omp_size,orch_omp_rank)
529  USE  mod_surf_para
530     
531#ifdef CPP_MPI
532    INCLUDE 'mpif.h'
533#endif   
534
535    INTEGER,INTENT(OUT) :: orch_comm
536    INTEGER,INTENT(OUT) :: orch_omp_size
537    INTEGER,INTENT(OUT) :: orch_omp_rank
538    INTEGER             :: color
539    INTEGER             :: i,ierr
540!
541! End definition
542!****************************************************************************************
543   
544   
545    IF (is_omp_root) THEN         
546     
547      IF (knon_mpi==0) THEN
548         color = 0
549      ELSE
550         color = 1
551      ENDIF
552   
553#ifdef CPP_MPI   
554      CALL MPI_COMM_SPLIT(COMM_LMDZ_PHY,color,mpi_rank,orch_comm,ierr)
555#endif
556   
557    ENDIF
558    CALL bcast_omp(orch_comm)
559   
560    IF (knon_mpi /= 0) THEN
561      orch_omp_size=0
562      DO i=0,omp_size-1
563        IF (knon_omp_para(i) /=0) THEN
564          orch_omp_size=orch_omp_size+1
565          IF (i==omp_rank) orch_omp_rank=orch_omp_size-1
566        ENDIF
567      ENDDO
568    ENDIF
569   
570   
571  END SUBROUTINE Get_orchidee_communicator
572!
573!****************************************************************************************
574
575
576  SUBROUTINE Init_neighbours(knon,neighbours,knindex,pctsrf)
577    USE mod_grid_phy_lmdz
578    USE mod_surf_para   
579    INCLUDE "indicesol.h"
580
581#ifdef CPP_MPI
582    INCLUDE 'mpif.h'
583#endif   
584
585! Input arguments
586!****************************************************************************************
587    INTEGER, INTENT(IN)                     :: knon
588    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
589    REAL, DIMENSION(klon), INTENT(IN)       :: pctsrf
590   
591! Output arguments
592!****************************************************************************************
593    INTEGER, DIMENSION(knon,8), INTENT(OUT) :: neighbours
594
595! Local variables
596!****************************************************************************************
597    INTEGER                              :: i, igrid, jj, ij, iglob
598    INTEGER                              :: ierr, ireal, index
599    INTEGER, DIMENSION(8,3)              :: off_ini
600    INTEGER, DIMENSION(8)                :: offset 
601    INTEGER, DIMENSION(nbp_lon,nbp_lat)  :: correspond
602    INTEGER, DIMENSION(knon_glo)         :: ktindex_glo
603    INTEGER, DIMENSION(knon_glo,8)       :: neighbours_glo
604    REAL, DIMENSION(klon_glo)            :: pctsrf_glo
605    INTEGER                              :: ktindex(klon)
606!
607! End definition
608!****************************************************************************************
609
610    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
611   
612    CALL gather_surf(ktindex(1:knon),ktindex_glo)
613    CALL gather(pctsrf,pctsrf_glo)
614   
615    IF (is_mpi_root .AND. is_omp_root) THEN
616      neighbours_glo(:,:)=-1
617!  Initialisation des offset   
618!
619! offset bord ouest
620       off_ini(1,1) = - nbp_lon   ; off_ini(2,1) = - nbp_lon + 1     ; off_ini(3,1) = 1
621       off_ini(4,1) = nbp_lon + 1 ; off_ini(5,1) = nbp_lon           ; off_ini(6,1) = 2 * nbp_lon - 1
622       off_ini(7,1) = nbp_lon -1  ; off_ini(8,1) = - 1
623! offset point normal
624       off_ini(1,2) = - nbp_lon   ; off_ini(2,2) = - nbp_lon + 1     ; off_ini(3,2) = 1
625       off_ini(4,2) = nbp_lon + 1 ; off_ini(5,2) = nbp_lon           ; off_ini(6,2) = nbp_lon - 1
626       off_ini(7,2) = -1          ; off_ini(8,2) = - nbp_lon - 1
627! offset bord   est
628       off_ini(1,3) = - nbp_lon   ; off_ini(2,3) = - 2 * nbp_lon + 1 ; off_ini(3,3) = - nbp_lon + 1
629       off_ini(4,3) =  1          ; off_ini(5,3) = nbp_lon           ; off_ini(6,3) = nbp_lon - 1
630       off_ini(7,3) = -1          ; off_ini(8,3) = - nbp_lon - 1
631!
632!
633! Attention aux poles
634!
635       DO igrid = 1, knon_glo
636          index = ktindex_glo(igrid)
637          jj = INT((index - 1)/nbp_lon) + 1
638          ij = index - (jj - 1) * nbp_lon
639          correspond(ij,jj) = igrid
640       ENDDO
641       
642       DO igrid = 1, knon_glo
643          iglob = ktindex_glo(igrid)
644         
645          IF (MOD(iglob, nbp_lon) == 1) THEN
646             offset = off_ini(:,1)
647          ELSE IF(MOD(iglob, nbp_lon) == 0) THEN
648             offset = off_ini(:,3)
649          ELSE
650             offset = off_ini(:,2)
651          ENDIF
652         
653          DO i = 1, 8
654             index = iglob + offset(i)
655             ireal = (MIN(MAX(1, index - nbp_lon + 1), klon_glo))
656             IF (pctsrf_glo(ireal) > EPSFRA) THEN
657                jj = INT((index - 1)/nbp_lon) + 1
658                ij = index - (jj - 1) * nbp_lon
659                neighbours_glo(igrid, i) = correspond(ij, jj)
660             ENDIF
661          ENDDO
662       ENDDO
663
664    ENDIF
665   
666    DO i = 1, 8
667      CALL scatter_surf(neighbours_glo(:,i),neighbours(1:knon,i))
668    ENDDO
669  END SUBROUTINE Init_neighbours
670
671!
672!****************************************************************************************
673!
674#endif
675END MODULE surf_land_orchidee_mod
Note: See TracBrowser for help on using the repository browser.