source: LMDZ5/branches/LMDZ5V1.0-dev/libf/phylmd/surf_land_orchidee_mod.F90 @ 4855

Last change on this file since 4855 was 1444, checked in by jghattas, 14 years ago

Modifications for carbon tracers :

  • add possibility to read source at different frequency
  • add dynamic varaiable to send fluxes in interface between ORCHIDEE and LMDZ : this is still under developpment under cpp key ORCH_NEW. No compatible official ORCHIDEE version does yet exist.
  • add variable co2_send in restart file.
  • clean up and bug fixes.


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