Context Navigation

surf_land_orchidee_mod.F90 @ 799

Last change on this file since 799 was 799, checked in by Laurent Fairhead, 17 years ago
Rajoute la possibilite d'utiliser l'ancienne interface d'ORCHIDEE selon que la cle CPP ORC_PREPAR est definie (interface pre orchidee_1_9) ou non (interface post orchide_1_9). LF
Property svn:eol-style set to `native` Property svn:executable set to ``* Property svn:keywords set to `Author Date Id Revision`
File size: 22.1 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format