source: LMDZ5/trunk/libf/phylmd/surf_land_orchidee_noz0h_mod.F90 @ 2912

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