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
RevLine 
[781]1!
[2571]2MODULE surf_land_orchidee_noz0h_mod
3#ifdef ORCHIDEE_NOZ0H
[781]4!
[2410]5! This module controles the interface towards the model ORCHIDEE.
[781]6!
[2410]7! Compatibility with ORCHIDIEE :
[2571]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.
[2410]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!
[781]15! Subroutines in this module : surf_land_orchidee
16!                              Init_orchidee_index
17!                              Get_orchidee_communicator
18!                              Init_neighbours
19
20  USE dimphy
[1067]21#ifdef CPP_VEGET
[781]22  USE intersurf     ! module d'ORCHIDEE
[1067]23#endif
[781]24  USE cpl_mod,      ONLY : cpl_send_land_fields
[996]25  USE surface_data, ONLY : type_ocean
[2351]26  USE geometry_mod, ONLY : dx, dy
[781]27  USE mod_grid_phy_lmdz
[2351]28  USE mod_phys_lmdz_para, mpi_root_rank=>mpi_master
[781]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, &
[2410]40       knindex, rlon, rlat, yrmu0, pctsrf, &
[781]41       debut, lafin, &
[2240]42       plev,  u1_lay, v1_lay, gustiness, temp_air, spechum, epot_air, ccanopy, &
[781]43       tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
44       precip_rain, precip_snow, lwdown, swnet, swdown, &
[1146]45       ps, q2m, t2m, &
[781]46       evap, fluxsens, fluxlat, &             
[888]47       tsol_rad, tsurf_new, alb1_new, alb2_new, &
[2952]48       emis_new, z0_new, z0h_new, qsurf, &
49       veget, lai, height)
[1279]50
51    USE mod_surf_para
52    USE mod_synchro_omp
[1454]53    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl
[1785]54    USE indice_sol_mod
[2311]55    USE print_control_mod, ONLY: lunout
[2344]56    USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat
[2348]57#ifdef CPP_VEGET
58    USE time_phylmdz_mod, ONLY: itau_phy
59#endif
[781]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
[1279]84!   ccanopy      concentration CO2 canopee, correspond au co2_send de
85!                carbon_cycle_mod ou valeur constant co2_ppm
[781]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
[888]106!   alb1_new     albedo in visible SW interval
107!   alb2_new     albedo in near IR interval
[781]108!   emis_new     emissivite
109!   z0_new       surface roughness
[2571]110!   z0h_new      surface roughness, it is a copy from z0_new
[781]111!   qsurf        air moisture at surface
112!
[793]113    INCLUDE "YOMCST.h"
[2952]114    INCLUDE "dimpft.h" 
[781]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
[2410]126    REAL, DIMENSION(klon), INTENT(IN)         :: yrmu0 ! cosine of solar zenith angle
[781]127    REAL, DIMENSION(klon), INTENT(IN)         :: plev
[2240]128    REAL, DIMENSION(klon), INTENT(IN)         :: u1_lay, v1_lay, gustiness
[781]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
[1146]136    REAL, DIMENSION(klon), INTENT(IN)         :: q2m, t2m
[781]137
138! Parametres de sortie
139!****************************************************************************************
140    REAL, DIMENSION(klon), INTENT(OUT)        :: evap, fluxsens, fluxlat, qsurf
[888]141    REAL, DIMENSION(klon), INTENT(OUT)        :: tsol_rad, tsurf_new
142    REAL, DIMENSION(klon), INTENT(OUT)        :: alb1_new, alb2_new
[2571]143    REAL, DIMENSION(klon), INTENT(OUT)        :: emis_new, z0_new, z0h_new
[2952]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
[781]147
148! Local
149!****************************************************************************************
150    INTEGER                                   :: ij, jj, igrid, ireal, index
151    INTEGER                                   :: error
[1146]152    REAL, DIMENSION(klon)                     :: swdown_vrai
[781]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)
[987]211   
212    INTEGER :: orch_omp_rank
213    INTEGER :: orch_omp_size
[781]214!
215! Fin definition
216!****************************************************************************************
217
218    IF (check) WRITE(lunout,*)'Entree ', modname
[2952]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
[781]224 
225! Initialisation
226 
227    IF (debut) THEN
[1067]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.'
[2311]231       CALL abort_physic(modname,abort_message,1)
[1067]232#endif
233
[987]234       CALL Init_surf_para(knon)
[781]235       ALLOCATE(ktindex(knon))
236       IF ( .NOT. ALLOCATED(albedo_keep)) THEN
[987]237!ym          ALLOCATE(albedo_keep(klon))
238!ym bizarre que non alloué en knon precedement
239          ALLOCATE(albedo_keep(knon))
[781]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
[2344]251          IF ( indi > nbp_lon) THEN
[781]252             indi = 1
253             indj = indj + 1
254          ENDIF
255          ig(igrid) = indi
256          jg(igrid) = indj
257       ENDDO
258       ig(klon) = 1
[2344]259       jg(klon) = nbp_lat
[781]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'
[2311]265             CALL abort_physic(modname,abort_message,1)
[781]266          ENDIF
267       ENDIF
268       IF ((.NOT. ALLOCATED(lon_scat))) THEN
[2344]269          ALLOCATE(lon_scat(nbp_lon,nbp_lat), stat = error)
[781]270          IF (error /= 0) THEN
271             abort_message='Pb allocation lon_scat'
[2311]272             CALL abort_physic(modname,abort_message,1)
[781]273          ENDIF
274       ENDIF
275       IF ((.NOT. ALLOCATED(lat_scat))) THEN
[2344]276          ALLOCATE(lat_scat(nbp_lon,nbp_lat), stat = error)
[781]277          IF (error /= 0) THEN
278             abort_message='Pb allocation lat_scat'
[2311]279             CALL abort_physic(modname,abort_message,1)
[781]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
[2344]297          DO jj = 2, nbp_lat-1
298             DO ij = 1, nbp_lon
[781]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)
[2344]306          lon_scat(:,nbp_lat) = lon_scat(:,2)
307          lat_scat(:,nbp_lat) = rlat_g(klon_glo)
[781]308       ENDIF
309   
[1023]310       CALL bcast(lon_scat)
311       CALL bcast(lat_scat)
[781]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'
[2311]319             CALL abort_physic(modname,abort_message,1)
[781]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'
[2311]327             CALL abort_physic(modname,abort_message,1)
[781]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'
[2311]345             CALL abort_physic(modname,abort_message,1)
[781]346          ENDIF
347       ENDIF
348       DO igrid = 1, knon
349          ij = knindex(igrid)
[2351]350          resolution(igrid,1) = dx(ij)
351          resolution(igrid,2) = dy(ij)
[781]352       ENDDO
353     
354       ALLOCATE(coastalflow(klon), stat = error)
355       IF (error /= 0) THEN
356          abort_message='Pb allocation coastalflow'
[2311]357          CALL abort_physic(modname,abort_message,1)
[781]358       ENDIF
359       
360       ALLOCATE(riverflow(klon), stat = error)
361       IF (error /= 0) THEN
362          abort_message='Pb allocation riverflow'
[2311]363          CALL abort_physic(modname,abort_message,1)
[781]364       ENDIF
[1279]365!
[1454]366! carbon_cycle_cpl not possible with this interface and version of ORHCHIDEE
[1279]367!
368       IF (carbon_cycle_cpl) THEN
[1454]369          abort_message='carbon_cycle_cpl not yet possible with this interface of ORCHIDEE'
[2311]370          CALL abort_physic(modname,abort_message,1)
[1279]371       END IF
372       
[781]373    ENDIF                          ! (fin debut)
[1279]374 
[781]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
[987]381     
[781]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)
[2011]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)
[781]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
[987]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       
[1023]407       IF (knon > 0) THEN
[1067]408#ifdef CPP_VEGET
[987]409         CALL Init_intersurf(nbp_lon,nbp_lat,knon,ktindex,offset,orch_omp_size,orch_omp_rank,orch_comm)
[1067]410#endif
[987]411       ENDIF
[802]412
[987]413       
414       IF (knon > 0) THEN
415
[1067]416#ifdef CPP_VEGET
[2344]417          CALL intersurf_main (itime+itau_phy-1, nbp_lon, nbp_lat, knon, ktindex, dtime, &
[802]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, &
[2410]425               lon_scat, lat_scat, q2m, t2m, coszang=yrmu0)
[1067]426#endif         
[781]427       ENDIF
428
[987]429       CALL Synchro_omp
430
[781]431       albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
432
433    ENDIF
434
[987]435   
[781]436!  swdown_vrai(1:knon) = swnet(1:knon)/(1. - albedo_keep(1:knon))
437    swdown_vrai(1:knon) = swdown(1:knon)
438
[987]439    IF (knon > 0) THEN
[1067]440#ifdef CPP_VEGET   
[2344]441       CALL intersurf_main (itime+itau_phy, nbp_lon, nbp_lat, knon, ktindex, dtime,  &
[802]442            lrestart_read, lrestart_write, lalo, &
443            contfrac, neighbours, resolution, date0, &
[781]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), &
[2410]449            lon_scat, lat_scat, q2m, t2m, coszang=yrmu0(1:knon))
[1067]450#endif       
[781]451    ENDIF
[2571]452    z0h_new(:)=z0_new(:)
[781]453
[987]454    CALL Synchro_omp
455   
[781]456    albedo_keep(1:knon) = (albedo_out(1:knon,1)+albedo_out(1:knon,2))/2.
457
458!* Send to coupler
459!
[996]460    IF (type_ocean=='couple') THEN
[781]461       CALL cpl_send_land_fields(itime, knon, knindex, &
462            riverflow, coastalflow)
463    ENDIF
464
[888]465    alb1_new(1:knon) = albedo_out(1:knon,1)
466    alb2_new(1:knon) = albedo_out(1:knon,2)
[781]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   
[987]476    IF (debut) CALL Finalize_surf_para
[1279]477
[987]478   
[781]479  END SUBROUTINE surf_land_orchidee
480!
481!****************************************************************************************
482!
[987]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)
[781]491   
[987]492    INTEGER               :: ktindex_glo(knon_glo)
493    INTEGER               :: offset_para(0:omp_size*mpi_size-1)
494    INTEGER               :: LastPoint
495    INTEGER               :: task
[781]496   
[987]497    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
[781]498   
[987]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
[781]509    ENDIF
510   
[987]511    CALL bcast(offset_para)
[781]512   
[987]513    offset=offset_para(omp_size*mpi_rank+omp_rank)
514   
515    ktindex(1:knon)=ktindex(1:knon)-offset
[781]516
[987]517  END SUBROUTINE Init_orchidee_index
518
[781]519!
[987]520!************************* ***************************************************************
[781]521!
[987]522
523  SUBROUTINE Get_orchidee_communicator(orch_comm,orch_omp_size,orch_omp_rank)
524  USE  mod_surf_para
525     
[1001]526#ifdef CPP_MPI
[781]527    INCLUDE 'mpif.h'
528#endif   
529
530    INTEGER,INTENT(OUT) :: orch_comm
[987]531    INTEGER,INTENT(OUT) :: orch_omp_size
532    INTEGER,INTENT(OUT) :: orch_omp_rank
[781]533    INTEGER             :: color
[987]534    INTEGER             :: i,ierr
[802]535!
536! End definition
537!****************************************************************************************
[781]538   
[987]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   
[1001]548#ifdef CPP_MPI   
[987]549      CALL MPI_COMM_SPLIT(COMM_LMDZ_PHY,color,mpi_rank,orch_comm,ierr)
[781]550#endif
[1023]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
[987]564   
[781]565   
566  END SUBROUTINE Get_orchidee_communicator
567!
568!****************************************************************************************
569
[987]570
571  SUBROUTINE Init_neighbours(knon,neighbours,knindex,pctsrf)
572    USE mod_grid_phy_lmdz
573    USE mod_surf_para   
[1785]574    USE indice_sol_mod
[987]575
[1001]576#ifdef CPP_MPI
[781]577    INCLUDE 'mpif.h'
578#endif   
579
[802]580! Input arguments
581!****************************************************************************************
[781]582    INTEGER, INTENT(IN)                     :: knon
[987]583    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
[781]584    REAL, DIMENSION(klon), INTENT(IN)       :: pctsrf
585   
[802]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 
[987]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)
[802]601!
602! End definition
603!****************************************************************************************
604
[987]605    ktindex(1:knon)=knindex(1:knon)+(klon_mpi_begin-1)+(klon_omp_begin-1)+nbp_lon-1
[781]606   
[987]607    CALL gather_surf(ktindex(1:knon),ktindex_glo)
608    CALL gather(pctsrf,pctsrf_glo)
[781]609   
[987]610    IF (is_mpi_root .AND. is_omp_root) THEN
611      neighbours_glo(:,:)=-1
[781]612!  Initialisation des offset   
613!
614! offset bord ouest
[987]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
[781]618! offset point normal
[987]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
[781]622! offset bord   est
[987]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
[781]626!
627!
628! Attention aux poles
629!
[987]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
[781]634          correspond(ij,jj) = igrid
635       ENDDO
[2126]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
[781]644       
[987]645       DO igrid = 1, knon_glo
646          iglob = ktindex_glo(igrid)
647         
648          IF (MOD(iglob, nbp_lon) == 1) THEN
[781]649             offset = off_ini(:,1)
[987]650          ELSE IF(MOD(iglob, nbp_lon) == 0) THEN
[781]651             offset = off_ini(:,3)
652          ELSE
653             offset = off_ini(:,2)
654          ENDIF
[987]655         
[781]656          DO i = 1, 8
657             index = iglob + offset(i)
[987]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)
[781]663             ENDIF
664          ENDDO
665       ENDDO
[2126]666       ENDIF !fin knon_glo == 1
[781]667
668    ENDIF
669   
[987]670    DO i = 1, 8
671      CALL scatter_surf(neighbours_glo(:,i),neighbours(1:knon,i))
[781]672    ENDDO
673  END SUBROUTINE Init_neighbours
[987]674
[781]675!
676!****************************************************************************************
677!
[1146]678#endif
[2571]679END MODULE surf_land_orchidee_noz0h_mod
Note: See TracBrowser for help on using the repository browser.