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

Last change on this file since 4799 was 2954, checked in by Laurent Fairhead, 7 years ago

Bug correction for new copses drag

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