source: LMDZ4/trunk/libf/phylmd/surf_land_orchidee_noopenmp_mod.F90 @ 1467

Last change on this file since 1467 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

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