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

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

Parametrization of drag by copses
Need version 4465 of ORCHIDEE at least

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