source: LMDZ5/trunk/libf/phylmd/surf_land_orchidee_nofrein_mod.F90 @ 5314

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