Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

intersurf.f90 @ 7329

Last change on this file since 7329 was 7191, checked in by haicheng.zhang, 3 years ago
Implementing latral transfers of sediment and POC from land to ocean through river in ORCHILEAK
Property svn:keywords set to `Revision Date HeadURL Date Author Revision`
File size: 92.2 KB

Line
1	! ================================================================================================================================
2	! MODULE : intersurf
3	!
4	! CONTACT : orchidee-help _at_ listes.ipsl.fr
5	!
6	! LICENCE : IPSL (2006)
7	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8	!
9	!>\BRIEF Subroutines for interfacing between the driver dim2_driver and sechiba and between LMDZ and sechiba
10	!!
11	!!\n DESCRIPTION: This module contains subroutines for the interfacing between ORCHIDEE and LMDZ and between the driver
12	!! dim2_driver and the reste of ORCHIDEE. Follwoing subroutines exists:
13	!!
14	!! - intersurf_initialize_2d : Intitialization and call to sechiba_initialize, called from dim2_driver
15	!! before the first call to intersurf_main_2d.
16	!! - intersurf_main_2d : Main subroutine will call sechiba_main, called from dim2_driver at each time-step
17	!! - init_intersurf : Initialize grid information when coupling with LMDZ. This subroutine is called
18	!! from LMDZ before call to intersurf_initialize_gathered
19	!! - intersurf_initialize_gathered : Intitialization and call to sechiba_initialize, called from LMDZ
20	!! before the first call to intersurf_main_2d.
21	!! - intersurf_main_gathered : Main subroutine will call sechiba_main, called from LMDZ at each time-step
22	!! - intsurf_time : Initialize and update time information, called in the initialization phase
23	!! and at each time step from the different intersurf subroutines.
24	!!
25	!! RECENT CHANGE(S):
26	!!
27	!! REFERENCE(S) : None
28	!!
29	!! SVN :
30	!! $HeadURL$
31	!! $Date$
32	!! $Revision$
33	!! \n
34	!_ ================================================================================================================================
35	MODULE intersurf
36
37	USE IOIPSL
38	USE xios_orchidee
39	USE ioipsl_para
40	USE defprec
41	USE sechiba
42	USE constantes
43	USE constantes_soil
44	USE control
45	USE pft_parameters
46	USE mod_orchidee_para
47	USE solar
48	USE grid
49	USE time, ONLY : one_year, one_day, dt_sechiba, time_initialize, time_nextstep, julian_diff
50	USE time, ONLY : year_end, month_end, day_end, sec_end
51	USE thermosoilc, ONLY : thermosoilc_levels
52	USE ioipslctrl, ONLY : ioipslctrl_history, ioipslctrl_restini, ok_histsync, max_hist_level, dw
53
54	IMPLICIT NONE
55
56	PRIVATE
57	PUBLIC :: Init_intersurf, intersurf_main, intersurf_main_2d, intersurf_main_gathered
58	PUBLIC :: intersurf_initialize_2d, intersurf_initialize_gathered, intersurf_clear
59
60
61	! Interface called from LMDZ in older verisons
62	INTERFACE intersurf_main
63	MODULE PROCEDURE intersurf_main_gathered
64	END INTERFACE
65
66	LOGICAL, SAVE :: lstep_init_intersurf=.TRUE.!! Initialisation has to be done one time
67	!$OMP THREADPRIVATE(lstep_init_intersurf)
68	INTEGER(i_std), SAVE :: printlev_loc !! Write level to this module
69	!$OMP THREADPRIVATE(printlev_loc)
70	INTEGER(i_std), SAVE :: hist_id, rest_id !! IDs for history and restart files
71	!$OMP THREADPRIVATE(hist_id, rest_id)
72	INTEGER(i_std), SAVE :: hist2_id !! ID for the second history files (Hi-frequency ?)
73	!$OMP THREADPRIVATE(hist2_id)
74	INTEGER(i_std), SAVE :: hist_id_stom, hist_id_stom_IPCC, rest_id_stom !! Dito for STOMATE
75	!$OMP THREADPRIVATE(hist_id_stom, hist_id_stom_IPCC, rest_id_stom)
76	INTEGER(i_std), SAVE :: itau_offset !! This offset is used to phase the
77	!! calendar of the GCM or the driver.
78	!$OMP THREADPRIVATE(itau_offset)
79	REAL(r_std), SAVE :: date0_shifted
80	!$OMP THREADPRIVATE(date0_shifted)
81	REAL(r_std), SAVE :: julian0 !! first day of this year
82	!$OMP THREADPRIVATE(julian0)
83	LOGICAL, SAVE :: ok_q2m_t2m=.TRUE. !! Flag ok if the variables are present in the call in gcm.
84	!$OMP THREADPRIVATE(ok_q2m_t2m)
85	LOGICAL, SAVE :: fatmco2 !! Flag to force the value of atmospheric CO2 for vegetation.
86	!$OMP THREADPRIVATE(fatmco2)
87	REAL(r_std), SAVE :: atmco2 !! atmospheric CO2
88	!$OMP THREADPRIVATE(atmco2)
89	REAL(r_std), SAVE :: coeff_rel !! Coefficient for time filter on riverflow and coastalflow
90	!$OMP THREADPRIVATE(coeff_rel)
91	REAL(r_std), ALLOCATABLE, DIMENSION(:,:), SAVE :: riverflow_cpl0 !! Value from previous time step for riverflow
92	!$OMP THREADPRIVATE(riverflow_cpl0)
93	REAL(r_std), ALLOCATABLE, DIMENSION(:,:), SAVE :: coastalflow_cpl0 !! Value from previous time step for coastalflow
94	!$OMP THREADPRIVATE(coastalflow_cpl0)
95	INTEGER(i_std), SAVE :: nb_fields_in=-1 !! Number of fields to give to the GCM
96	!$OMP THREADPRIVATE(nb_fields_in)
97	INTEGER(i_std), SAVE :: nb_fields_out=-1 !! Number of fields to get from the GCM
98	!$OMP THREADPRIVATE(nb_fields_out)
99
100	PUBLIC lstep_init_intersurf
101
102	CONTAINS
103
104	!! =============================================================================================================================
105	!! SUBROUTINE: intersurf_initialize_2d
106	!!
107	!>\BRIEF Initialization and call to sechiba_initialize
108	!!
109	!! DESCRIPTION: Initialization of module variables, read options from parameter file, initialize output files and call to
110	!! sechiba_initialize.
111	!!
112	!! This subroutine is called from dim2_driver before the first call to intersurf_main_2d.
113	!!
114	!! \n
115	!_ ==============================================================================================================================
116
117	SUBROUTINE intersurf_initialize_2d (kjit, iim, jjm, kjpindex, kindex, xrdt, &
118	lrestart_read, lrestart_write, lon, lat, zcontfrac, zresolution, date0, &
119	zlev, u, v, qair, temp_air, epot_air, ccanopy, &
120	cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
121	precip_rain, precip_snow, lwdown, swnet, swdown, pb, &
122	vevapp, fluxsens, fluxlat, coastalflow0, riverflow0, &
123	tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m)
124
125	IMPLICIT NONE
126
127	!! 0. Variable and parameter declaration
128	!! 0.1 Input variables
129	INTEGER(i_std),INTENT (in) :: kjit !! Time step number
130	INTEGER(i_std),INTENT (in) :: iim, jjm !! Dimension of input fields
131	INTEGER(i_std),INTENT (in) :: kjpindex !! Number of continental points
132	REAL(r_std),INTENT (in) :: xrdt !! Time step in seconds
133	LOGICAL, INTENT (in) :: lrestart_read !! Logical for _restart_ file to read
134	LOGICAL, INTENT (in) :: lrestart_write !! Logical for _restart_ file to write'
135	REAL(r_std), INTENT (in) :: date0 !! Date at which kjit = 0
136	INTEGER(i_std),DIMENSION (kjpindex), INTENT (in) :: kindex !! Index for continental points
137	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: u !! Lowest level wind speed
138	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: v !! Lowest level wind speed
139	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: zlev !! Height of first layer
140	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: qair !! Lowest level specific humidity
141	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: precip_rain !! Rain precipitation
142	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: precip_snow !! Snow precipitation
143	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: lwdown !! Down-welling long-wave flux
144	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: swnet !! Net surface short-wave flux
145	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: swdown !! Downwelling surface short-wave flux
146	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: temp_air !! Air temperature in Kelvin
147	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: epot_air !! Air potential energy
148	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: ccanopy !! CO2 concentration in the canopy
149	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: petAcoef !! Coeficients A from the PBL resolution
150	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: peqAcoef !! One for T and another for q
151	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: petBcoef !! Coeficients B from the PBL resolution
152	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: peqBcoef !! One for T and another for q
153	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: pb !! Surface pressure
154	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: lon, lat !! Geographical coordinates
155	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: zcontfrac !! Fraction of continent in the grid
156	REAL(r_std),DIMENSION (iim,jjm,2), INTENT(in) :: zresolution !! resolution in x and y dimensions
157
158	!! 0.2 Output variables
159	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: z0m !! Surface roughness
160	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: coastalflow0 !! Diffuse flow of water into the ocean (m^3/s)
161	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: riverflow0 !! Largest rivers flowing into the ocean (m^3/s)
162	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: tsol_rad !! Radiative surface temperature
163	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: vevapp !! Total of evaporation
164	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: temp_sol_new !! New soil temperature
165	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: qsurf !! Surface specific humidity
166	REAL(r_std),DIMENSION (iim,jjm,2), INTENT(out) :: albedo !! Albedo
167	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: fluxsens !! Sensible chaleur flux
168	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: fluxlat !! Latent chaleur flux
169	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: emis !! Emissivity
170
171	!! 0.3 Modified variables
172	REAL(r_std),DIMENSION (iim,jjm), INTENT(inout) :: cdrag !! Cdrag
173
174	!! 0.4 Local variables
175	REAL(r_std),DIMENSION (kjpindex) :: zu !! Work array to keep u
176	REAL(r_std),DIMENSION (kjpindex) :: zv !! Work array to keep v
177	REAL(r_std),DIMENSION (kjpindex) :: zzlev !! Work array to keep zlev
178	REAL(r_std),DIMENSION (kjpindex) :: zqair !! Work array to keep qair
179	REAL(r_std),DIMENSION (kjpindex) :: zprecip_rain !! Work array to keep precip_rain
180	REAL(r_std),DIMENSION (kjpindex) :: zprecip_snow !! Work array to keep precip_snow
181	REAL(r_std),DIMENSION (kjpindex) :: zlwdown !! Work array to keep lwdown
182	REAL(r_std),DIMENSION (kjpindex) :: zswnet !! Work array to keep swnet
183	REAL(r_std),DIMENSION (kjpindex) :: zswdown !! Work array to keep swdown
184	REAL(r_std),DIMENSION (kjpindex) :: ztemp_air !! Work array to keep temp_air
185	REAL(r_std),DIMENSION (kjpindex) :: zepot_air !! Work array to keep epot_air
186	REAL(r_std),DIMENSION (kjpindex) :: zccanopy !! Work array to keep ccanopy
187	REAL(r_std),DIMENSION (kjpindex) :: zpetAcoef !! Work array to keep petAcoef
188	REAL(r_std),DIMENSION (kjpindex) :: zpeqAcoef !! Work array to keep peqAcoef
189	REAL(r_std),DIMENSION (kjpindex) :: zpetBcoef !! Work array to keep petBcoef
190	REAL(r_std),DIMENSION (kjpindex) :: zpeqBcoef !! Work array to keep peqVcoef
191	REAL(r_std),DIMENSION (kjpindex) :: zcdrag !! Work array to keep cdrag
192	REAL(r_std),DIMENSION (kjpindex) :: zpb !! Work array to keep pb
193	REAL(r_std),DIMENSION (kjpindex) :: zz0m, zz0h !! Work array to keep zz0m, zz0h
194	REAL(r_std),DIMENSION (kjpindex,nflow) :: zcoastal !! Work array to keep coastalflow (m^3/dt)
195	REAL(r_std),DIMENSION (kjpindex,nflow) :: zriver !! Work array to keep riverflow (m^3/dt)
196	REAL(r_std),DIMENSION (kjpindex) :: znetco2 !! Work array to keep netco2flux
197	REAL(r_std),DIMENSION (kjpindex) :: zcarblu !! Work array to keep fco2_land_use
198	REAL(r_std),DIMENSION (kjpindex) :: ztsol_rad !! Work array to keep tsol_rad
199	REAL(r_std),DIMENSION (kjpindex) :: zvevapp !! Work array to keep vevapp
200	REAL(r_std),DIMENSION (kjpindex) :: ztemp_sol_new !! Work array to keep temp_sol_new
201	REAL(r_std),DIMENSION (kjpindex) :: zqsurf !! Work array to keep qsurf
202	REAL(r_std),DIMENSION (kjpindex,2) :: zalbedo !! Work array to keep albedo
203	REAL(r_std),DIMENSION (kjpindex) :: zfluxsens !! Work array to keep fluxsens
204	REAL(r_std),DIMENSION (kjpindex) :: zfluxlat !! Work array to keep fluxlat
205	REAL(r_std),DIMENSION (kjpindex) :: zemis !! Work array to keep emis
206	REAL(r_std),DIMENSION (iim,jjm,nflow) :: coastalflow !! Diffuse flow of water into the ocean (m^3/s)
207	REAL(r_std),DIMENSION (iim,jjm,nflow) :: riverflow !! Largest rivers flowing into the ocean (m^3/s)
208	REAL(r_std),ALLOCATABLE, DIMENSION (:) :: soilth_lev !! Vertical soil axis for thermal scheme (m)
209	INTEGER(i_std) :: i, j, ik
210	INTEGER(i_std) :: ier
211	INTEGER(i_std) :: itau_sechiba
212	INTEGER :: old_fileout !! old Logical Int for std IO output
213
214
215	IF (printlev >= 2) WRITE(numout,*) 'Start intersurf_initialize_2d'
216	CALL ipslnlf_p(new_number=numout,old_number=old_fileout)
217
218	! Initialize variables in module time
219	CALL time_initialize(kjit, date0, xrdt, "END")
220
221	! Configuration of SSL specific parameters
222	CALL control_initialize
223
224	! Initialize specific write level
225	printlev_loc=get_printlev('instersurf')
226
227	OFF_LINE_MODE = .TRUE.
228
229	DO ik=1,kjpindex
230
231	j = ((kindex(ik)-1)/iim) + 1
232	i = (kindex(ik) - (j-1)*iim)
233
234	!- Create the internal coordinate table
235	!-
236	lalo(ik,1) = lat(i,j)
237	lalo(ik,2) = lon(i,j)
238	!
239	!- Store the fraction of the continents only once so that the user
240	!- does not change them afterwards.
241	!-
242	contfrac(ik) = zcontfrac(i,j)
243	ENDDO
244	CALL gather(contfrac,contfrac_g)
245	CALL gather(lalo,lalo_g)
246	CALL gather2D_mpi(lon,lon_g)
247	CALL gather2D_mpi(lat,lat_g)
248
249	CALL ioipslctrl_restini(kjit, date0, xrdt, rest_id, rest_id_stom, itau_offset, date0_shifted)
250	itau_sechiba = kjit + itau_offset
251
252	!!- Initialize module for output with XIOS
253	!
254	! Get the vertical soil levels for the thermal scheme, to be used in xios_orchidee_init
255	ALLOCATE(soilth_lev(ngrnd), stat=ier)
256	IF (ier /= 0) THEN
257	CALL ipslerr_p(3,'intersurf_main_2d', 'Error in allocation of soilth_lev','','')
258	END IF
259	IF (hydrol_cwrr) THEN
260	soilth_lev(1:ngrnd) = znt(:)
261	ELSE
262	soilth_lev(1:ngrnd) = thermosoilc_levels()
263	END IF
264
265	CALL xios_orchidee_init( MPI_COMM_ORCH, &
266	date0, year_end, month_end, day_end, julian_diff, &
267	lon, lat, soilth_lev )
268
269	!- Initialize IOIPSL sechiba output files
270	CALL ioipslctrl_history(iim, jjm, lon, lat, kindex, kjpindex, itau_sechiba, date0_shifted, xrdt, hist_id, &
271	hist2_id, hist_id_stom, hist_id_stom_IPCC)
272
273	!
274	! Shift the time step to phase the two models
275	!
276	itau_sechiba = kjit + itau_offset
277
278	! Update the calendar in xios by sending the new time step
279	! Special case : the model is only in initialization phase and the current itau_sechiba is not a real time step.
280	! Therefor give itau_sechiba+1 to xios to have a correct time axis in output files.
281	CALL xios_orchidee_update_calendar(itau_sechiba+1)
282
283	!
284	! 1. gather input fields from kindex array
285	! Warning : I'm not sure this interface with one dimension array is the good one
286	!
287	!WRITE(numout,*) 'SUMinterfini0',iim,jjm,kjpindex
288	DO ik=1, kjpindex
289
290	j = ((kindex(ik)-1)/iim) + 1
291	i = (kindex(ik) - (j-1)*iim)
292
293	zu(ik) = u(i,j)
294	zv(ik) = v(i,j)
295	zzlev(ik) = zlev(i,j)
296	zqair(ik) = qair(i,j)
297	zprecip_rain(ik) = precip_rain(i,j)*xrdt
298	zprecip_snow(ik) = precip_snow(i,j)*xrdt
299	zlwdown(ik) = lwdown(i,j)
300	zswnet(ik) = swnet(i,j)
301	zswdown(ik) = swdown(i,j)
302	ztemp_air(ik) = temp_air(i,j)
303	zepot_air(ik) = epot_air(i,j)
304	zccanopy(ik) = ccanopy(i,j)
305	zpetAcoef(ik) = petAcoef(i,j)
306	zpeqAcoef(ik) = peqAcoef(i,j)
307	zpetBcoef(ik) = petBcoef(i,j)
308	zpeqBcoef(ik) = peqBcoef(i,j)
309	zcdrag(ik) = cdrag(i,j)
310	zpb(ik) = pb(i,j)
311	!WRITE(numout,*) 'SUMinterfini0',i,j,ik,zlev(i,j),zzlev(ik)
312	ENDDO
313	!
314	! 2. save the grid
315	!
316	CALL histwrite_p(hist_id, 'LandPoints', itau_sechiba+1, (/ ( REAL(ik), ik=1,kjpindex ) /), kjpindex, kindex)
317	CALL histwrite_p(hist_id, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
318	IF ( ok_stomate ) THEN
319	CALL histwrite_p(hist_id_stom, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
320	CALL histwrite_p(hist_id_stom_IPCC, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
321	ENDIF
322
323	CALL histwrite_p(hist_id, 'Contfrac', itau_sechiba+1, contfrac, kjpindex, kindex)
324	IF ( is_omp_root .AND. hist_id > 0 ) THEN
325	! Always syncronize output after initialization
326	CALL histsync(hist_id)
327	END IF
328
329	CALL histwrite_p(hist2_id, 'LandPoints', itau_sechiba+1, (/ ( REAL(ik), ik=1,kjpindex ) /), kjpindex, kindex)
330	CALL histwrite_p(hist2_id, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
331	CALL histwrite_p(hist2_id, 'Contfrac', itau_sechiba+1, contfrac, kjpindex, kindex)
332	IF ( is_omp_root .AND. hist2_id > 0 ) THEN
333	! Always syncronize output after initialization
334	CALL histsync(hist2_id)
335	ENDIF
336
337	!
338	! 3. call sechiba for continental points only
339	!
340	IF (printlev_loc >= 3) WRITE(numout,*) 'Before call to sechiba_initialize'
341	!WRITE(numout,*) 'Intsurfini1_tem_soil',MINVAL(zcoastal(:,1)),MAXVAL(zcoastal(:,1)),MINVAL(zriver(:,1)),MAXVAL(zriver(:,1))
342	CALL sechiba_initialize( &
343	itau_sechiba, iim*jjm, kjpindex, kindex, &
344	lalo, contfrac, neighbours, resolution, zzlev, &
345	zu, zv, zqair, ztemp_air, ztemp_air, &
346	zpetAcoef, zpeqAcoef, zpetBcoef, zpeqBcoef, &
347	zprecip_rain, zprecip_snow, zlwdown, zswnet, zswdown, &
348	zpb, rest_id, hist_id, hist2_id, &
349	rest_id_stom, hist_id_stom, hist_id_stom_IPCC, &
350	zcoastal, zriver, ztsol_rad, zvevapp, zqsurf, &
351	zz0m, zz0h, zalbedo, zfluxsens, zfluxlat, &
352	zemis, znetco2, zcarblu, ztemp_sol_new, zcdrag)
353	!WRITE(numout,*) 'Intsurfini2_tem_soil',MINVAL(zcoastal(:,1)),MAXVAL(zcoastal(:,1)),MINVAL(zriver(:,1)),MAXVAL(zriver(:,1))
354	IF (printlev_loc >= 3) WRITE(numout,*) 'After call to sechiba_initialize'
355	!
356	! 5. scatter output fields
357	!
358	z0m(:,:) = undef_sechiba
359	coastalflow(:,:,:)= undef_sechiba
360	riverflow(:,:,:) = undef_sechiba
361	tsol_rad(:,:) = undef_sechiba
362	vevapp(:,:) = undef_sechiba
363	temp_sol_new(:,:) = undef_sechiba
364	qsurf(:,:) = undef_sechiba
365	albedo(:,:,:) = undef_sechiba
366	fluxsens(:,:) = undef_sechiba
367	fluxlat(:,:) = undef_sechiba
368	emis(:,:) = undef_sechiba
369	cdrag(:,:) = undef_sechiba
370
371	DO ik=1, kjpindex
372	j = ((kindex(ik)-1)/iim) + 1
373	i = (kindex(ik) - (j-1)*iim)
374
375	z0m(i,j) = zz0m(ik)
376	coastalflow(i,j,:)= zcoastal(ik,:)
377	riverflow(i,j,:) = zriver(ik,:)
378	tsol_rad(i,j) = ztsol_rad(ik)
379	vevapp(i,j) = zvevapp(ik)
380	temp_sol_new(i,j) = ztemp_sol_new(ik)
381	qsurf(i,j) = zqsurf(ik)
382	albedo(i,j,1) = zalbedo(ik,1)
383	albedo(i,j,2) = zalbedo(ik,2)
384	fluxsens(i,j) = zfluxsens(ik)
385	fluxlat(i,j) = zfluxlat(ik)
386	emis(i,j) = zemis(ik)
387	cdrag(i,j) = zcdrag(ik)
388	!WRITE(numout,*) 'SUMinterfini3',i,j,ik,z0m(i,j),zz0m(ik),zlev(i,j)
389	ENDDO
390
391	!
392	! 6. Transform the water fluxes into Kg/m^2s and m^3/s
393	!
394	DO ik=1, kjpindex
395
396	j = ((kindex(ik)-1)/iim) + 1
397	i = (kindex(ik) - (j-1)*iim)
398
399	vevapp(i,j) = vevapp(i,j)/xrdt
400	coastalflow(i,j,:) = coastalflow(i,j,:)/xrdt
401	riverflow(i,j,:) = riverflow(i,j,:)/xrdt
402
403	ENDDO
404	coastalflow0(:,:) = coastalflow(:,:,1)
405	riverflow0(:,:) = riverflow(:,:,1)
406
407	IF (is_root_prc) CALL getin_dump
408
409	lstep_init_intersurf = .FALSE.
410	CALL ipslnlf_p(new_number=old_fileout)
411	IF (printlev_loc >= 1) WRITE (numout,*) 'End intersurf_initialize_2d'
412
413	DO ik=1, kjpindex
414	j = ((kindex(ik)-1)/iim) + 1
415	i = (kindex(ik) - (j-1)*iim)
416	!WRITE(numout,*) 'SUMinterfini4',i,j,ik,z0m(i,j),zz0m(ik),zlev(i,j)
417	ENDDO
418	END SUBROUTINE intersurf_initialize_2d
419
420
421	!! =============================================================================================================================
422	!! SUBROUTINE: intersurf_main_2d
423	!!
424	!>\BRIEF Main subroutine to call ORCHIDEE from dim2_driver using variables on a 2d grid.
425	!!
426	!! DESCRIPTION: This subroutine is the main interface for ORCHIDEE when it is called from the offline driver dim2_driver.
427	!! The variables are all on the 2D grid including ocean points. intersurf_initialize_2d should be called before
428	!! this subroutine is called. This subroutine is called at each time step.
429	!!
430	!! \n
431	!_ ==============================================================================================================================
432
433	SUBROUTINE intersurf_main_2d (kjit, iim, jjm, kjpindex, kindex, xrdt, &
434	lrestart_read, lrestart_write, lon, lat, zcontfrac, zresolution, date0, &
435	zlev, u, v, qair, temp_air, epot_air, ccanopy, &
436	cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
437	precip_rain, precip_snow, lwdown, swnet, swdown, pb, &
438	vevapp, fluxsens, fluxlat, coastalflow0, riverflow0, &
439	tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m, &
440	coszang)
441
442	IMPLICIT NONE
443
444	!! 0. Variable and parameter declaration
445	!! 0.1 Input variables
446	INTEGER(i_std),INTENT (in) :: kjit !! Time step number
447	INTEGER(i_std),INTENT (in) :: iim, jjm !! Dimension of input fields
448	INTEGER(i_std),INTENT (in) :: kjpindex !! Number of continental points
449	REAL(r_std),INTENT (in) :: xrdt !! Time step in seconds
450	LOGICAL, INTENT (in) :: lrestart_read !! Logical for _restart_ file to read
451	LOGICAL, INTENT (in) :: lrestart_write !! Logical for _restart_ file to write'
452	REAL(r_std), INTENT (in) :: date0 !! Date at which kjit = 0
453	INTEGER(i_std),DIMENSION (kjpindex), INTENT (in) :: kindex !! Index for continental points
454	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: u !! Lowest level wind speed
455	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: v !! Lowest level wind speed
456	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: zlev !! Height of first layer
457	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: qair !! Lowest level specific humidity
458	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: precip_rain !! Rain precipitation
459	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: precip_snow !! Snow precipitation
460	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: lwdown !! Down-welling long-wave flux
461	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: swnet !! Net surface short-wave flux
462	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: swdown !! Downwelling surface short-wave flux
463	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: coszang !! Cosine of the solar zenith angle (unitless)
464	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: temp_air !! Air temperature in Kelvin
465	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: epot_air !! Air potential energy
466	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: ccanopy !! CO2 concentration in the canopy
467	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: petAcoef !! Coeficients A from the PBL resolution
468	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: peqAcoef !! One for T and another for q
469	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: petBcoef !! Coeficients B from the PBL resolution
470	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: peqBcoef !! One for T and another for q
471	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: pb !! Surface pressure
472	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: lon, lat !! Geographical coordinates
473	REAL(r_std),DIMENSION (iim,jjm), INTENT(in) :: zcontfrac !! Fraction of continent in the grid
474	REAL(r_std),DIMENSION (iim,jjm,2), INTENT(in) :: zresolution !! resolution in x and y dimensions
475
476	!! 0.2 Output variables
477	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: z0m !! Surface roughness for momemtum
478	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: coastalflow0 !! Diffuse flow of water into the ocean (m^3/s)
479	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: riverflow0 !! Largest rivers flowing into the ocean (m^3/s)
480	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: tsol_rad !! Radiative surface temperature
481	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: vevapp !! Total of evaporation
482	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: temp_sol_new !! New soil temperature
483	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: qsurf !! Surface specific humidity
484	REAL(r_std),DIMENSION (iim,jjm,2), INTENT(out) :: albedo !! Albedo
485	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: fluxsens !! Sensible chaleur flux
486	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: fluxlat !! Latent chaleur flux
487	REAL(r_std),DIMENSION (iim,jjm), INTENT(out) :: emis !! Emissivity
488
489	!! 0.3 Modified variables
490	REAL(r_std),DIMENSION (iim,jjm), INTENT(inout) :: cdrag !! Cdrag
491
492	!! 0.4 Local variables
493	REAL(r_std),DIMENSION (kjpindex) :: zu !! Work array to keep u
494	REAL(r_std),DIMENSION (kjpindex) :: zv !! Work array to keep v
495	REAL(r_std),DIMENSION (kjpindex) :: zzlev !! Work array to keep zlev
496	REAL(r_std),DIMENSION (kjpindex) :: zqair !! Work array to keep qair
497	REAL(r_std),DIMENSION (kjpindex) :: zprecip_rain !! Work array to keep precip_rain
498	REAL(r_std),DIMENSION (kjpindex) :: zprecip_snow !! Work array to keep precip_snow
499	REAL(r_std),DIMENSION (kjpindex) :: zlwdown !! Work array to keep lwdown
500	REAL(r_std),DIMENSION (kjpindex) :: zswnet !! Work array to keep swnet
501	REAL(r_std),DIMENSION (kjpindex) :: zswdown !! Work array to keep swdown
502	REAL(r_std),DIMENSION (kjpindex) :: zcoszang !! Work array to keep coszang
503	REAL(r_std),DIMENSION (kjpindex) :: ztemp_air !! Work array to keep temp_air
504	REAL(r_std),DIMENSION (kjpindex) :: zepot_air !! Work array to keep epot_air
505	REAL(r_std),DIMENSION (kjpindex) :: zccanopy !! Work array to keep ccanopy
506	REAL(r_std),DIMENSION (kjpindex) :: zpetAcoef !! Work array to keep petAcoef
507	REAL(r_std),DIMENSION (kjpindex) :: zpeqAcoef !! Work array to keep peqAcoef
508	REAL(r_std),DIMENSION (kjpindex) :: zpetBcoef !! Work array to keep petBcoef
509	REAL(r_std),DIMENSION (kjpindex) :: zpeqBcoef !! Work array to keep peqVcoef
510	REAL(r_std),DIMENSION (kjpindex) :: zcdrag !! Work array to keep cdrag
511	REAL(r_std),DIMENSION (kjpindex) :: zpb !! Work array to keep pb
512	REAL(r_std),DIMENSION (kjpindex) :: zz0m, zz0h !! Work array to keep zz0m, zz0h
513	REAL(r_std),DIMENSION (kjpindex,nvm) :: zveget !! Work array to keep veget
514	REAL(r_std),DIMENSION (kjpindex,nvm) :: zlai !! Work array to keep lai
515	REAL(r_std),DIMENSION (kjpindex,nvm) :: zheight !! Work array to keep height
516	REAL(r_std),DIMENSION (kjpindex,nflow) :: zcoastal !! Work array to keep coastalflow (m^3/dt)
517	REAL(r_std),DIMENSION (kjpindex,nflow) :: zriver !! Work array to keep riverflow (m^3/dt)
518	REAL(r_std),DIMENSION (kjpindex) :: znetco2 !! Work array to keep netco2flux
519	REAL(r_std),DIMENSION (kjpindex) :: zcarblu !! Work array to keep fco2_land_use
520	REAL(r_std),DIMENSION (kjpindex) :: ztsol_rad !! Work array to keep tsol_rad
521	REAL(r_std),DIMENSION (kjpindex) :: zvevapp !! Work array to keep vevapp
522	REAL(r_std),DIMENSION (kjpindex) :: ztemp_sol_new !! Work array to keep temp_sol_new
523	REAL(r_std),DIMENSION (kjpindex) :: zqsurf !! Work array to keep qsurf
524	REAL(r_std),DIMENSION (kjpindex,2) :: zalbedo !! Work array to keep albedo
525	REAL(r_std),DIMENSION (kjpindex) :: zfluxsens !! Work array to keep fluxsens
526	REAL(r_std),DIMENSION (kjpindex) :: zfluxlat !! Work array to keep fluxlat
527	REAL(r_std),DIMENSION (kjpindex) :: zemis !! Work array to keep emis
528	REAL(r_std),DIMENSION (iim,jjm,nflow) :: coastalflow !! Diffuse flow of water into the ocean (m^3/s)
529	REAL(r_std),DIMENSION (iim,jjm,nflow) :: riverflow !! Largest rivers flowing into the ocean (m^3/s)
530	REAL(r_std),ALLOCATABLE, DIMENSION (:) :: soilth_lev !! Vertical soil axis for thermal scheme (m)
531	INTEGER(i_std) :: i, j, ik
532	INTEGER(i_std) :: ier
533	INTEGER(i_std) :: itau_sechiba
534	INTEGER :: old_fileout !! old Logical Int for std IO output
535
536	IF (printlev_loc >= 3) WRITE(numout,*) 'Start intersurf_main_2d'
537	CALL ipslnlf_p(new_number=numout,old_number=old_fileout)
538
539	!
540	! Shift the time step to phase the two models
541	!
542	itau_sechiba = kjit + itau_offset
543	!
544
545	! Update the calendar in xios by sending the new time step
546	CALL xios_orchidee_update_calendar(itau_sechiba)
547
548	! Update the calendar and all time variables in module time
549	CALL time_nextstep(itau_sechiba)
550	!
551	! 1. gather input fields from kindex array
552	! Warning : I'm not sure this interface with one dimension array is the good one
553	!
554	DO ik=1, kjpindex
555
556	j = ((kindex(ik)-1)/iim) + 1
557	i = (kindex(ik) - (j-1)*iim)
558
559	zu(ik) = u(i,j)
560	zv(ik) = v(i,j)
561	zzlev(ik) = zlev(i,j)
562	zqair(ik) = qair(i,j)
563	zprecip_rain(ik) = precip_rain(i,j)*xrdt
564	zprecip_snow(ik) = precip_snow(i,j)*xrdt
565	zlwdown(ik) = lwdown(i,j)
566	zswnet(ik) = swnet(i,j)
567	zswdown(ik) = swdown(i,j)
568	zcoszang(ik) = coszang(i,j)
569	ztemp_air(ik) = temp_air(i,j)
570	zepot_air(ik) = epot_air(i,j)
571	zccanopy(ik) = ccanopy(i,j)
572	zpetAcoef(ik) = petAcoef(i,j)
573	zpeqAcoef(ik) = peqAcoef(i,j)
574	zpetBcoef(ik) = petBcoef(i,j)
575	zpeqBcoef(ik) = peqBcoef(i,j)
576	zcdrag(ik) = cdrag(i,j)
577	zpb(ik) = pb(i,j)
578	!WRITE(numout,*) 'SUMinterf',ik,i,j
579	!WRITE(numout,*) 'SUMvu',u(i,j),v(i,j)
580	!WRITE(numout,*) 'SUMtemp',temp_air(i,j),qair(i,j)
581	ENDDO
582
583	!
584	! 3. call sechiba for continental points only
585	!
586	IF (printlev_loc >= 3) WRITE(numout,*) 'Before call to sechiba_main from intersurf_main_2d'
587	!WRITE(numout,*) 'ZHC_Intersurf1_tem_soil',MINVAL(ztemp_sol_new(:)),MAXVAL(ztemp_sol_new(:))
588	!WRITE(numout,*) 'UV1',MINVAL(zu(:)),MAXVAL(zu(:)),MINVAL(zv(:)),MAXVAL(zv(:))
589	!WRITE(numout,*) 'UV11',MINVAL(u(:,:)),MAXVAL(u(:,:)),MINVAL(v(:,:)),MAXVAL(v(:,:))
590	!WRITE(numout,*) 'ZHC_Intersurf1_flow',MINVAL(zcoastal(:,:)),MAXVAL(zcoastal(:,:)),MINVAL(zriver(:,:)),MAXVAL(zriver(:,:))
591	CALL sechiba_main (itau_sechiba, iim*jjm, kjpindex, kindex, &
592	lrestart_read, lrestart_write, &
593	lalo, contfrac, neighbours, resolution, &
594	zzlev, zu, zv, zqair, zqair, ztemp_air, ztemp_air, zepot_air, zccanopy, &
595	zcdrag, zpetAcoef, zpeqAcoef, zpetBcoef, zpeqBcoef, &
596	zprecip_rain ,zprecip_snow, zlwdown, zswnet, zswdown, zcoszang, zpb, &
597	zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, &
598	ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0m, zz0h,&
599	zveget, zlai, zheight, &
600	rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC )
601	!WRITE(numout,*) 'ZHC_Intersurf2_tem_soil',MINVAL(ztemp_sol_new(:)),MAXVAL(ztemp_sol_new(:))
602	!WRITE(numout,*) 'ZHC_Intersurf2_flow',MINVAL(zcoastal(:,:)),MAXVAL(zcoastal(:,:)),MINVAL(zriver(:,:)),MAXVAL(zriver(:,:))
603	IF (printlev_loc >= 3) WRITE(numout,*) 'After call to sechiba_main'
604
605	!
606	! 5. scatter output fields
607	!
608	z0m(:,:) = undef_sechiba
609	coastalflow(:,:,:)= undef_sechiba
610	riverflow(:,:,:) = undef_sechiba
611	tsol_rad(:,:) = undef_sechiba
612	vevapp(:,:) = undef_sechiba
613	temp_sol_new(:,:) = undef_sechiba
614	qsurf(:,:) = undef_sechiba
615	albedo(:,:,:) = undef_sechiba
616	fluxsens(:,:) = undef_sechiba
617	fluxlat(:,:) = undef_sechiba
618	emis(:,:) = undef_sechiba
619	cdrag(:,:) = undef_sechiba
620	!
621	DO ik=1, kjpindex
622
623
624	j = ((kindex(ik)-1)/iim) + 1
625	i = (kindex(ik) - (j-1)*iim)
626
627	z0m(i,j) = zz0m(ik)
628	coastalflow(i,j,:) = zcoastal(ik,:)
629	riverflow(i,j,:) = zriver(ik,:)
630	tsol_rad(i,j) = ztsol_rad(ik)
631	vevapp(i,j) = zvevapp(ik)
632	temp_sol_new(i,j) = ztemp_sol_new(ik)
633	qsurf(i,j) = zqsurf(ik)
634	albedo(i,j,1) = zalbedo(ik,1)
635	albedo(i,j,2) = zalbedo(ik,2)
636	fluxsens(i,j) = zfluxsens(ik)
637	fluxlat(i,j) = zfluxlat(ik)
638	emis(i,j) = zemis(ik)
639	cdrag(i,j) = zcdrag(ik)
640
641	ENDDO
642
643	CALL xios_orchidee_send_field("LandPoints" ,(/ ( REAL(ik), ik=1,kjpindex ) /))
644	CALL xios_orchidee_send_field("areas", area)
645	CALL xios_orchidee_send_field("contfrac",contfrac)
646	CALL xios_orchidee_send_field("temp_air",ztemp_air)
647	CALL xios_orchidee_send_field("qair",zqair)
648	CALL xios_orchidee_send_field("swnet",zswnet)
649	CALL xios_orchidee_send_field("swdown",zswdown)
650	CALL xios_orchidee_send_field("pb",zpb)
651
652	IF ( .NOT. almaoutput ) THEN
653	!
654	! scattered during the writing
655	!
656	CALL histwrite_p (hist_id, 'evap', itau_sechiba, zvevapp, kjpindex, kindex)
657	CALL histwrite_p (hist_id, 'coastalflow',itau_sechiba, zcoastal(:,ih2o), kjpindex, kindex)
658	CALL histwrite_p (hist_id, 'riverflow',itau_sechiba, zriver(:,ih2o), kjpindex, kindex)
659	CALL histwrite_p (hist_id, 'temp_sol', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
660	CALL histwrite_p (hist_id, 'tsol_max', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
661	CALL histwrite_p (hist_id, 'tsol_min', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
662	CALL histwrite_p (hist_id, 'fluxsens', itau_sechiba, zfluxsens, kjpindex, kindex)
663	CALL histwrite_p (hist_id, 'fluxlat', itau_sechiba, zfluxlat, kjpindex, kindex)
664	CALL histwrite_p (hist_id, 'swnet', itau_sechiba, zswnet, kjpindex, kindex)
665	CALL histwrite_p (hist_id, 'swdown', itau_sechiba, zswdown, kjpindex, kindex)
666	CALL histwrite_p (hist_id, 'alb_vis', itau_sechiba, zalbedo(:,1), kjpindex, kindex)
667	CALL histwrite_p (hist_id, 'alb_nir', itau_sechiba, zalbedo(:,2), kjpindex, kindex)
668	CALL histwrite_p (hist_id, 'tair', itau_sechiba, ztemp_air, kjpindex, kindex)
669	CALL histwrite_p (hist_id, 'qair', itau_sechiba, zqair, kjpindex, kindex)
670	CALL histwrite_p (hist_id, 'q2m', itau_sechiba, zqair, kjpindex, kindex)
671	CALL histwrite_p (hist_id, 't2m', itau_sechiba, ztemp_air, kjpindex, kindex)
672	CALL histwrite_p (hist2_id, 'evap', itau_sechiba, zvevapp, kjpindex, kindex)
673	CALL histwrite_p (hist2_id, 'coastalflow',itau_sechiba, zcoastal(:,ih2o), kjpindex, kindex)
674	CALL histwrite_p (hist2_id, 'riverflow',itau_sechiba, zriver(:,ih2o), kjpindex, kindex)
675	CALL histwrite_p (hist2_id, 'temp_sol', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
676	CALL histwrite_p (hist2_id, 'tsol_max', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
677	CALL histwrite_p (hist2_id, 'tsol_min', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
678	CALL histwrite_p (hist2_id, 'fluxsens', itau_sechiba, zfluxsens, kjpindex, kindex)
679	CALL histwrite_p (hist2_id, 'fluxlat', itau_sechiba, zfluxlat, kjpindex, kindex)
680	CALL histwrite_p (hist2_id, 'swnet', itau_sechiba, zswnet, kjpindex, kindex)
681	CALL histwrite_p (hist2_id, 'swdown', itau_sechiba, zswdown, kjpindex, kindex)
682	CALL histwrite_p (hist2_id, 'alb_vis', itau_sechiba, zalbedo(:,1), kjpindex, kindex)
683	CALL histwrite_p (hist2_id, 'alb_nir', itau_sechiba, zalbedo(:,2), kjpindex, kindex)
684	CALL histwrite_p (hist2_id, 'tair', itau_sechiba, ztemp_air, kjpindex, kindex)
685	CALL histwrite_p (hist2_id, 'qair', itau_sechiba, zqair, kjpindex, kindex)
686	CALL histwrite_p (hist2_id, 'q2m', itau_sechiba, zqair, kjpindex, kindex)
687	CALL histwrite_p (hist2_id, 't2m', itau_sechiba, ztemp_air, kjpindex, kindex)
688	ELSE
689	CALL histwrite_p (hist_id, 'Evap', itau_sechiba, zvevapp, kjpindex, kindex)
690	CALL histwrite_p (hist_id, 'SWnet', itau_sechiba, zswnet, kjpindex, kindex)
691	CALL histwrite_p (hist_id, 'Qh', itau_sechiba, zfluxsens, kjpindex, kindex)
692	CALL histwrite_p (hist_id, 'Qle', itau_sechiba, zfluxlat, kjpindex, kindex)
693	CALL histwrite_p (hist_id, 'AvgSurfT', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
694	CALL histwrite_p (hist_id, 'RadT', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
695	CALL histwrite_p (hist_id, 'Tair', itau_sechiba, ztemp_air, kjpindex, kindex)
696	CALL histwrite_p (hist_id, 'Qair', itau_sechiba, zqair, kjpindex, kindex)
697	CALL histwrite_p (hist2_id, 'Evap', itau_sechiba, zvevapp, kjpindex, kindex)
698	CALL histwrite_p (hist2_id, 'SWnet', itau_sechiba, zswnet, kjpindex, kindex)
699	CALL histwrite_p (hist2_id, 'Qh', itau_sechiba, zfluxsens, kjpindex, kindex)
700	CALL histwrite_p (hist2_id, 'Qle', itau_sechiba, zfluxlat, kjpindex, kindex)
701	CALL histwrite_p (hist2_id, 'AvgSurfT', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
702	CALL histwrite_p (hist2_id, 'RadT', itau_sechiba, ztemp_sol_NEW, kjpindex, kindex)
703	ENDIF
704	!
705	IF ( is_omp_root ) THEN
706	IF ( (dw .EQ. xrdt) .AND. hist_id > 0 ) THEN
707	! Syncronize output but only if flag ok_histsync is set to true
708	IF (ok_histsync) CALL histsync(hist_id)
709	ENDIF
710	END IF
711
712	!
713	! 6. Transform the water fluxes into Kg/m^2s and m^3/s
714	!
715	DO ik=1, kjpindex
716
717	j = ((kindex(ik)-1)/iim) + 1
718	i = (kindex(ik) - (j-1)*iim)
719
720	vevapp(i,j) = vevapp(i,j)/xrdt
721	coastalflow(i,j,:) = coastalflow(i,j,:)/xrdt
722	riverflow(i,j,:) = riverflow(i,j,:)/xrdt
723
724	ENDDO
725	coastalflow0(:,:) = coastalflow(:,:,1)
726	riverflow0(:,:) = riverflow(:,:,1)
727
728	CALL ipslnlf_p(new_number=old_fileout)
729	IF (printlev_loc >= 3) WRITE (numout,*) 'End intersurf_main_2d'
730
731	END SUBROUTINE intersurf_main_2d
732
733
734	!! =============================================================================================================================
735	!! SUBROUTINE: init_intersurf
736	!!
737	!>\BRIEF Initialize grid information
738	!!
739	!! DESCRIPTION: This subroutine is called from LMDZ before first call to intersurf_main_gathered or
740	!! intersurf_initialize_gathered
741	!!
742	!! \n
743	!_ ==============================================================================================================================
744
745	SUBROUTINE init_intersurf(nbp_l_lon,nbp_l_lat,kjpindex,kindex,orch_offset,orch_omp_size,orch_omp_rank,COMM)
746
747	USE mod_orchidee_para
748	USE timer
749	IMPLICIT NONE
750
751	INTEGER,INTENT(IN) :: nbp_l_lon
752	INTEGER,INTENT(IN) :: nbp_l_lat
753	INTEGER,INTENT(IN) :: kjpindex
754	INTEGER,INTENT(IN) :: kindex(:)
755	INTEGER,INTENT(IN) :: orch_offset
756	INTEGER,INTENT(IN) :: COMM
757	INTEGER,INTENT(IN) :: orch_omp_size
758	INTEGER,INTENT(IN) :: orch_omp_rank
759
760	INTEGER,DIMENSION(kjpindex) :: kindex_offset
761
762	IF (printlev >= 1) WRITE(,) 'Start ORCHIDEE'
763
764	IF (orch_omp_rank==0) THEN
765	CALL Init_timer
766	CALL start_timer(timer_mpi)
767	CALL grid_set_glo(nbp_l_lon,nbp_l_lat)
768	ENDIF
769	CALL barrier2_omp()
770	CALL init_orchidee_data_para(kjpindex,kindex,orch_offset,orch_omp_size,orch_omp_rank,COMM)
771	CALL Set_stdout_file('out_orchidee')
772
773	IF (printlev >= 1) WRITE(numout,*) 'Start ORCHIDEE intitalization phase'
774
775	IF (is_omp_root) CALL grid_allocate_glo(4)
776	CALL barrier2_omp()
777	CALL init_ioipsl_para
778
779	kindex_offset(:)=kindex(:)+offset
780	CALL gather(kindex_offset,index_g)
781	CALL bcast(index_g)
782
783	IF (printlev_loc >= 2) THEN
784	WRITE(numout,*) "kjpindex = ",kjpindex
785	WRITE(numout,*) "offset for OMP = ",offset_omp
786	WRITE(numout,*) "Index num local for continental points = ",kindex
787	WRITE(numout,*) "Index num global for continental points = ",kindex_offset
788	IF (is_omp_root) THEN
789	WRITE(numout,*) "ROOT OMP, Index global MPI : ",kindex_mpi(:)
790	ENDIF
791	IF (is_root_prc) THEN
792	WRITE(numout,*) "ROOT global, Index global : ",index_g(:)
793	ENDIF
794	END IF
795
796	! Allocation of grid variables
797	CALL grid_init ( kjpindex, 4, "RegLonLat", "2DGrid" )
798
799
800	END SUBROUTINE init_intersurf
801
802	!! =============================================================================================================================
803	!! SUBROUTINE: intersurf_initialize_gathered
804	!!
805	!>\BRIEF Initialization and call to sechiba_initialize
806	!!
807	!! DESCRIPTION: Initialization of module variables, read options from parameter file, initialize output files and call to
808	!! sechiba_initialize.
809	!!
810	!! This subroutine can be called directly from GCM(LMDZ). If it is not called before the first call to
811	!! intersurf_main_gathered, then it will be done from there. This possibility is done to keep backward
812	!! compatibility with LMDZ.
813	!!
814	!! \n
815	!_ ==============================================================================================================================
816
817	SUBROUTINE intersurf_initialize_gathered (kjit, iim_glo, jjm_glo, kjpindex, kindex, xrdt, &
818	lrestart_read, lrestart_write, latlon, zcontfrac, zneighbours, zresolution, date0, &
819	zlev, u, v, qair, temp_air, epot_air, &
820	cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
821	precip_rain, precip_snow, lwdown, swnet, swdown, pb, &
822	vevapp, fluxsens, fluxlat, coastalflow_cpl, riverflow_cpl, &
823	tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m, lon_scat_g, lat_scat_g, &
824	q2m, t2m, z0h, nvm_out, &
825	field_out_names, fields_out, field_in_names, fields_in)
826
827	USE mod_orchidee_para
828	IMPLICIT NONE
829
830	!! 0. Variable and parameter declaration
831	!! 0.1 Input
832	INTEGER(i_std),INTENT (in) :: kjit !! Time step number
833	INTEGER(i_std),INTENT (in) :: iim_glo !! Dimension of global fields
834	INTEGER(i_std),INTENT (in) :: jjm_glo !! Dimension of global fields
835	INTEGER(i_std),INTENT (in) :: kjpindex !! Number of continental points
836	REAL(r_std),INTENT (in) :: xrdt !! Time step in seconds
837	LOGICAL, INTENT (in) :: lrestart_read !! Logical for _restart_ file to read
838	LOGICAL, INTENT (in) :: lrestart_write !! Logical for _restart_ file to write'
839	REAL(r_std), INTENT (in) :: date0 !! Date at which kjit = 0
840	INTEGER(i_std),DIMENSION (kjpindex), INTENT (in) :: kindex !! Index for continental points
841	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: u !! Lowest level wind speed
842	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: v !! Lowest level wind speed
843	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: zlev !! Height of first layer
844	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: qair !! Lowest level specific humidity
845	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: precip_rain !! Rain precipitation
846	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: precip_snow !! Snow precipitation
847	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: lwdown !! Down-welling long-wave flux
848	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: swnet !! Net surface short-wave flux
849	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: swdown !! Downwelling surface short-wave flux
850	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: temp_air !! Air temperature in Kelvin
851	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: epot_air !! Air potential energy
852	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: petAcoef !! Coeficients A from the PBL resolution
853	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: peqAcoef !! One for T and another for q
854	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: petBcoef !! Coeficients B from the PBL resolution
855	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: peqBcoef !! One for T and another for q
856	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: pb !! Surface pressure
857	REAL(r_std),DIMENSION (kjpindex,2), INTENT(in) :: latlon !! Geographical coordinates
858	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: zcontfrac !! Fraction of continent
859	INTEGER(i_std),DIMENSION (kjpindex,NbNeighb), INTENT(in):: zneighbours !! neighbours
860	REAL(r_std),DIMENSION (kjpindex,2), INTENT(in) :: zresolution !! size of the grid box
861	REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN) :: lon_scat_g !! The scattered values for longitude
862	REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN) :: lat_scat_g !! The scattered values for latitudes
863
864	!! 0.2 Output variables
865	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: z0m !! Surface roughness (momentum)
866	REAL(r_std),DIMENSION (kjpindex,nflow), INTENT(out) :: coastalflow_cpl!! Diffuse flow of water into the ocean (m^3/s)
867	REAL(r_std),DIMENSION (kjpindex,nflow), INTENT(out) :: riverflow_cpl !! Largest rivers flowing into the ocean (m^3/s)
868	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: tsol_rad !! Radiative surface temperature
869	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: vevapp !! Total of evaporation
870	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: temp_sol_new !! New soil temperature
871	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: qsurf !! Surface specific humidity
872	REAL(r_std),DIMENSION (kjpindex,2), INTENT(out) :: albedo !! Albedo
873	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: fluxsens !! Sensible chaleur flux
874	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: fluxlat !! Latent chaleur flux
875	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: emis !! Emissivity
876
877	!! 0.3 Modified variables
878	REAL(r_std),DIMENSION (kjpindex), INTENT(inout) :: cdrag !! Cdrag
879
880	!! 0.4 Optional input and output variables
881	REAL(r_std),DIMENSION (kjpindex), INTENT(in), OPTIONAL :: q2m !! Surface specific humidity
882	REAL(r_std),DIMENSION (kjpindex), INTENT(in), OPTIONAL :: t2m !! Surface air temperature
883	REAL(r_std),DIMENSION (kjpindex), INTENT(out),OPTIONAL :: z0h !! Surface roughness (heat)
884	INTEGER(i_std), INTENT(out), OPTIONAL :: nvm_out !! Number of vegetation types, PFTs
885	CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_in_names !! Names for deposit variables to be transported
886	!! from chemistry model by GCM to ORCHIDEE
887	REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(IN) :: fields_in !! Fields for deposit variables to be transported
888	!! from chemistry model by GCM to ORCHIDEE
889	CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_out_names!! Names for emission variables to be transported
890	!! to chemistry model by GCM from ORCHIDEE
891	REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: fields_out !! Fields for emission variables to be transported
892	!! to chemistry model by GCM from ORCHIDEE
893
894	!! 0.5 Local variables
895	REAL(r_std),DIMENSION (kjpindex) :: zprecip_rain !! Work array to keep precip_rain
896	REAL(r_std),DIMENSION (kjpindex) :: zprecip_snow !! Work array to keep precip_snow
897	REAL(r_std),DIMENSION (kjpindex) :: zz0m !! Work array to keep zz0m
898	REAL(r_std),DIMENSION (kjpindex) :: zz0h !! Work array to keep zz0h
899	REAL(r_std),DIMENSION (kjpindex) :: zcdrag !! Work array for surface drag
900	REAL(r_std),DIMENSION (kjpindex,nflow) :: zcoastal !! Work array to keep coastal flow
901	REAL(r_std),DIMENSION (kjpindex,nflow) :: zriver !! Work array to keep river out flow
902	REAL(r_std),DIMENSION (kjpindex) :: znetco2 !! Work array to keep netco2flux
903	REAL(r_std),DIMENSION (kjpindex) :: zcarblu !! Work array to keep fco2_land_use
904	REAL(r_std),DIMENSION (kjpindex) :: ztsol_rad !! Work array to keep tsol_rad
905	REAL(r_std),DIMENSION (kjpindex) :: zvevapp !! Work array to keep vevapp
906	REAL(r_std),DIMENSION (kjpindex) :: ztemp_sol_new !! Work array to keep temp_sol_new
907	REAL(r_std),DIMENSION (kjpindex) :: zqsurf !! Work array to keep qsurf
908	REAL(r_std),DIMENSION (kjpindex,2) :: zalbedo !! Work array to keep albedo
909	REAL(r_std),DIMENSION (kjpindex) :: zfluxsens !! Work array to keep fluxsens
910	REAL(r_std),DIMENSION (kjpindex) :: zfluxlat !! Work array to keep fluxlat
911	REAL(r_std),DIMENSION (kjpindex) :: zemis !! Work array to keep emis
912	REAL(r_std),ALLOCATABLE, DIMENSION (:) :: soilth_lev !! Vertical soil axis for thermal scheme (m)
913	REAL(r_std),DIMENSION (kjpindex) :: q2m_loc !! Work array for q2m or qair
914	REAL(r_std),DIMENSION (kjpindex) :: t2m_loc !! Work array for t2m or temp_air
915	REAL(r_std),DIMENSION (:,:),ALLOCATABLE :: lon_scat !! The scattered values for longitude
916	REAL(r_std),DIMENSION (:,:),ALLOCATABLE :: lat_scat !! The scattered values for latitude
917	INTEGER(i_std) :: i, j, ik
918	INTEGER(i_std) :: ier
919	INTEGER(i_std) :: itau_sechiba
920	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: tmp_lon, tmp_lat, tmp_lev
921	INTEGER :: old_fileout !! old Logical Int for std IO output
922	REAL,ALLOCATABLE,DIMENSION(:,:) :: lalo_mpi
923	REAL(r_std),DIMENSION (kjpindex) :: landpoints !! Land point vector
924	REAL(r_std) :: tau_outflow !! Number of days for temoral filter on river- and coastalflow
925
926	! Initialize specific write level
927	printlev_loc=get_printlev('instersurf')
928
929	IF (printlev_loc >= 1) WRITE(numout,*) 'Entering intersurf_initialize_gathered'
930	IF (printlev_loc >= 2) WRITE(numout,*) 'using printlev_loc for intersurf:', printlev_loc
931
932	OFF_LINE_MODE = .FALSE.
933	CALL ipslnlf_p(new_number=numout, old_number=old_fileout)
934
935	! Initialize variables in module time
936	CALL time_initialize( kjit, date0, xrdt, "END" )
937
938	! Configuration of SSL specific parameters
939	CALL control_initialize
940
941	! Set the variable ok_q2m_t2m=true if q2m and t2m are present in the call from the gcm.
942	! If one of the variables are not present in this subroutine, set ok_q2m_t2m=.FALSE.
943	! Otherwise do not change the current value. Note that the current value for ok_q2m_t2m comes
944	! either from default initialization (true) or from intersurf_main_gathered.
945	IF (.NOT. PRESENT(q2m) .OR. .NOT. PRESENT(t2m)) THEN
946	ok_q2m_t2m=.FALSE.
947	END IF
948
949	IF (ok_q2m_t2m) THEN
950	t2m_loc=t2m
951	q2m_loc=q2m
952	ELSE
953	t2m_loc=temp_air
954	q2m_loc=qair
955	END IF
956
957
958	IF (is_omp_root) THEN
959	ALLOCATE(lon_scat(iim_g,jj_nb))
960	ALLOCATE(lat_scat(iim_g,jj_nb))
961	ELSE
962	ALLOCATE(lon_scat(1,1))
963	ALLOCATE(lat_scat(1,1))
964	ENDIF
965
966	! Create the internal coordinate table
967	!
968	lalo(:,:) = latlon(:,:)
969	CALL gather(lalo,lalo_g)
970	!
971	!-
972	!- Store variable to help describe the grid
973	!- once the points are gathered.
974	!-
975	neighbours(:,:) = zneighbours(:,:)
976	CALL gather(neighbours,neighbours_g)
977	!
978	resolution(:,:) = zresolution(:,:)
979	CALL gather(resolution,resolution_g)
980	!
981	area(:) = resolution(:,1)*resolution(:,2)
982	CALL gather(area,area_g)
983	!
984	!- Store the fraction of the continents only once so that the user
985	!- does not change them afterwards.
986	!
987	contfrac(:) = zcontfrac(:)
988	CALL gather(contfrac,contfrac_g)
989	!
990	!
991	! Create the internal coordinate table
992	!
993	IF ( (.NOT.ALLOCATED(tmp_lon))) THEN
994	ALLOCATE(tmp_lon(iim_g,jj_nb))
995	ENDIF
996	IF ( (.NOT.ALLOCATED(tmp_lat))) THEN
997	ALLOCATE(tmp_lat(iim_g,jj_nb))
998	ENDIF
999	IF ( (.NOT.ALLOCATED(tmp_lev))) THEN
1000	ALLOCATE(tmp_lev(iim_g,jj_nb))
1001	ENDIF
1002	!
1003	! Either we have the scattered coordinates as arguments or
1004	! we have to do the work here.
1005	!
1006	IF (is_omp_root) THEN
1007	lon_scat(:,:)=zero
1008	lat_scat(:,:)=zero
1009	CALL scatter2D_mpi(lon_scat_g,lon_scat)
1010	CALL scatter2D_mpi(lat_scat_g,lat_scat)
1011	lon_scat(:,1)=lon_scat(:,2)
1012	lon_scat(:,jj_nb)=lon_scat(:,2)
1013	lat_scat(:,1)=lat_scat(iim_g,1)
1014	lat_scat(:,jj_nb)=lat_scat(1,jj_nb)
1015
1016	tmp_lon(:,:) = lon_scat(:,:)
1017	tmp_lat(:,:) = lat_scat(:,:)
1018
1019	IF (is_mpi_root) THEN
1020	lon_g(:,:) = lon_scat_g(:,:)
1021	lat_g(:,:) = lat_scat_g(:,:)
1022	ENDIF
1023	ENDIF
1024
1025	!Config Key = FORCE_CO2_VEG
1026	!Config Desc = Flag to force the value of atmospheric CO2 for vegetation.
1027	!Config If = Only in coupled mode
1028	!Config Def = FALSE
1029	!Config Help = If this flag is set to true, the ATM_CO2 parameter is used
1030	!Config to prescribe the atmospheric CO2.
1031	!Config This Flag is only use in couple mode.
1032	!Config Units = [FLAG]
1033	fatmco2=.FALSE.
1034	CALL getin_p('FORCE_CO2_VEG',fatmco2)
1035	!
1036	! Next flag is only use in couple mode with a gcm in intersurf.
1037	! In forced mode, it has already been read and set in driver.
1038	IF ( fatmco2 ) THEN
1039	atmco2=350.
1040	CALL getin_p('ATM_CO2',atmco2)
1041	WRITE(numout,*) 'atmco2 ',atmco2
1042	ENDIF
1043
1044	CALL ioipslctrl_restini(kjit, date0, xrdt, rest_id, rest_id_stom, itau_offset, date0_shifted)
1045	itau_sechiba = kjit + itau_offset
1046
1047	!!- Initialize module for output with XIOS
1048	!
1049	! Get the vertical soil levels for the thermal scheme, to be used in xios_orchidee_init
1050	ALLOCATE(soilth_lev(ngrnd), stat=ier)
1051	IF (ier /= 0) THEN
1052	CALL ipslerr_p(3,'intersurf_main_gathered', 'Error in allocation of soilth_lev','','')
1053	END IF
1054	IF (hydrol_cwrr) THEN
1055	soilth_lev(1:ngrnd) = znt(:)
1056	ELSE
1057	soilth_lev(1:ngrnd) = thermosoilc_levels()
1058	END IF
1059
1060	CALL xios_orchidee_init( MPI_COMM_ORCH, &
1061	date0, year_end, month_end, day_end, julian_diff, &
1062	tmp_lon, tmp_lat, soilth_lev )
1063
1064	!- Initialize IOIPSL sechiba output files
1065	CALL ioipslctrl_history(iim_g, jj_nb, tmp_lon, tmp_lat, kindex, kjpindex, itau_sechiba, &
1066	date0_shifted, xrdt, hist_id, hist2_id, hist_id_stom, hist_id_stom_IPCC)
1067
1068	CALL bcast_omp(hist_id)
1069	CALL bcast_omp(hist2_id)
1070	CALL bcast_omp(hist_id_stom)
1071	CALL bcast_omp(hist_id_stom_IPCC)
1072
1073	! Count number of extra output fields to the GCM if it is not already done.
1074	IF (nb_fields_out == -1) THEN
1075	! nb_fields_out is not yet calculated. Do it now.
1076	! This means that the call is done directly from GCM.
1077	IF (PRESENT(field_out_names)) THEN
1078	nb_fields_out=SIZE(field_out_names)
1079	ELSE
1080	nb_fields_out=0
1081	ENDIF
1082	END IF
1083
1084	! Count number of extra input fields to the GCM if it is not already done.
1085	IF (nb_fields_in == -1) THEN
1086	! nb_fields_in is not yet calculated. Do it now.
1087	! This means that the call is done directly from GCM.
1088	IF (PRESENT(field_in_names)) THEN
1089	nb_fields_in=SIZE(field_in_names)
1090	ELSE
1091	nb_fields_in=0
1092	ENDIF
1093	END IF
1094
1095
1096	!
1097	!! Change to be in the orchidee context for XIOS
1098	!
1099	CALL xios_orchidee_change_context("orchidee")
1100
1101	!
1102	! Shift the time step to phase the two models
1103	!
1104	itau_sechiba = kjit + itau_offset
1105
1106	! Update the calendar in xios by sending the new time step
1107	! Special case : the model is only in initialization phase and the current itau_sechiba is not a real time step.
1108	! Therefor give itau_sechiba+1 to xios to have a correct time axis in output files.
1109	CALL xios_orchidee_update_calendar(itau_sechiba+1)
1110
1111	!
1112	! 1. Just change the units of some input fields
1113	!
1114	DO ik=1, kjpindex
1115
1116	zprecip_rain(ik) = precip_rain(ik)*xrdt
1117	zprecip_snow(ik) = precip_snow(ik)*xrdt
1118	zcdrag(ik) = cdrag(ik)
1119
1120	ENDDO
1121
1122	! Fields for deposit variables : to be transport from chemistry model by GCM to ORCHIDEE.
1123	! There are currently no fields to be transported into ORCHIDEE in this way
1124	DO i = 1, nb_fields_in
1125	WRITE(numout,*) i," Champ = ",TRIM(field_in_names(i))
1126	SELECT CASE(TRIM(field_in_names(i)))
1127	CASE DEFAULT
1128	CALL ipslerr_p (3,'intsurf_gathered', &
1129	'You ask in GCM an unknown field '//TRIM(field_in_names(i))//&
1130	' to give to ORCHIDEE for this specific version.',&
1131	'This model won''t be able to continue.', &
1132	'(check your tracer parameters in GCM)')
1133	END SELECT
1134	ENDDO
1135
1136	!
1137	! 3. save the grid
1138	!
1139	landpoints(:)=(/ ( REAL(ik), ik=1,kjpindex ) /)
1140	CALL histwrite_p(hist_id, 'LandPoints', itau_sechiba+1, landpoints, kjpindex, kindex)
1141	CALL histwrite_p(hist_id, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
1142	IF ( ok_stomate ) THEN
1143	CALL histwrite_p(hist_id_stom, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
1144	IF ( hist_id_stom_ipcc > 0 ) &
1145	CALL histwrite_p(hist_id_stom_IPCC, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
1146	ENDIF
1147	CALL histwrite_p(hist_id, 'Contfrac', itau_sechiba+1, contfrac, kjpindex, kindex)
1148
1149	! Syncronize output but only if flag ok_histsync is set to true
1150	IF (ok_histsync) THEN
1151	IF (is_omp_root .AND. hist_id > 0) THEN
1152	CALL histsync(hist_id)
1153	END IF
1154	END IF
1155
1156	IF ( hist2_id > 0 ) THEN
1157	CALL histwrite_p(hist2_id, 'LandPoints', itau_sechiba+1, landpoints, kjpindex, kindex)
1158	CALL histwrite_p(hist2_id, 'Areas', itau_sechiba+1, area, kjpindex, kindex)
1159	CALL histwrite_p(hist2_id, 'Contfrac', itau_sechiba+1, contfrac, kjpindex, kindex)
1160
1161	! Syncronize output but only if flag ok_histsync is set to true
1162	IF (ok_histsync .AND. is_omp_root) THEN
1163	CALL histsync(hist2_id)
1164	ENDIF
1165	ENDIF
1166	!
1167
1168	!
1169	! 4. call sechiba for continental points only
1170	!
1171	IF ( printlev_loc>=3 ) WRITE(numout,*) 'Before call to sechiba_initialize'
1172	!WRITE(numout,*) 'SUMinterfini4',SUM(z0m(:)),SUM(zz0m(:)),SUM(zlev(:))
1173	CALL sechiba_initialize( &
1174	itau_sechiba, iim_g*jj_nb, kjpindex, kindex, &
1175	lalo, contfrac, neighbours, resolution, zlev, &
1176	u, v, qair, t2m_loc, temp_air, &
1177	petAcoef, peqAcoef, petBcoef, peqBcoef, &
1178	zprecip_rain, zprecip_snow, lwdown, swnet, swdown, &
1179	pb, rest_id, hist_id, hist2_id, &
1180	rest_id_stom, hist_id_stom, hist_id_stom_IPCC, &
1181	zcoastal, zriver, ztsol_rad, zvevapp, zqsurf, &
1182	zz0m, zz0h, zalbedo, zfluxsens, zfluxlat, zemis, &
1183	znetco2, zcarblu, ztemp_sol_new, zcdrag)
1184
1185	IF ( printlev_loc>=3 ) WRITE(numout,*) 'After call to sechiba_initialize'
1186	!WRITE(numout,*) 'SUMinterfini5',SUM(z0m(:)),SUM(zz0m(:)),SUM(zlev(:))
1187	!
1188	! 6. scatter output fields
1189	!
1190	DO ik=1, kjpindex
1191	!WRITE(numout,*) 'SUMinterfini6',ik,z0m(ik),zz0m(ik)
1192	z0m(ik) = zz0m(ik)
1193	tsol_rad(ik) = ztsol_rad(ik)
1194	vevapp(ik) = zvevapp(ik)/xrdt ! Transform into kg/m^2/s
1195	temp_sol_new(ik) = ztemp_sol_new(ik)
1196	qsurf(ik) = zqsurf(ik)
1197	albedo(ik,1) = zalbedo(ik,1)
1198	albedo(ik,2) = zalbedo(ik,2)
1199	fluxsens(ik) = zfluxsens(ik)
1200	fluxlat(ik) = zfluxlat(ik)
1201	emis(ik) = zemis(ik)
1202	cdrag(ik) = zcdrag(ik)
1203	!WRITE(numout,*) 'SUMinterfini7',ik,z0m(ik),zz0m(ik)
1204	ENDDO
1205	! z0h is a optional output variable. Check first if it is present in the call from LMDZ.
1206	IF ( PRESENT(z0h) ) z0h(:) = zz0h(:)
1207
1208
1209	!Config Key = TAU_OUTFLOW
1210	!Config Desc = Number of days over which the coastal- and riverflow will be distributed
1211	!Config If = Only in coupled mode
1212	!Config Def = 0
1213	!Config Help = The default value 0 makes the distribution instanteneous
1214	!Config Units = [days]
1215	tau_outflow = 0
1216	CALL getin_p('TAU_OUTFLOW',tau_outflow)
1217	IF (tau_outflow <=xrdt/one_day) THEN
1218	coeff_rel = 1.0
1219	ELSE
1220	coeff_rel = (1.0 - exp(-xrdt/(tau_outflow*one_day)))
1221	END IF
1222	IF (printlev_loc >=2) WRITE(numout,*) 'tau_outflow, coeff_rel = ', tau_outflow, coeff_rel
1223
1224	! Allocate and read riverflow_cpl0 from restart file. Initialize to 0 if the variable was not found in the restart file.
1225	ALLOCATE(riverflow_cpl0(kjpindex,nflow), stat=ier)
1226	IF (ier /= 0) CALL ipslerr_p(3,'intersurf_initialize_gathered', 'Error in allocation of riverflow_cpl0','','')
1227	CALL restget_p (rest_id, 'riverflow_cpl0', nbp_glo, nflow, 1, kjit, .TRUE., riverflow_cpl0, "gather", nbp_glo, index_g)
1228	IF ( ALL(riverflow_cpl0(:,:) == val_exp) ) riverflow_cpl0(:,:)=0
1229
1230	! Allocate and read coastalflow_cpl0 from restart file. Initialize to 0 if the variable was not found in the restart file.
1231	ALLOCATE(coastalflow_cpl0(kjpindex,nflow), stat=ier)
1232	IF (ier /= 0) CALL ipslerr_p(3,'intersurf_initialize_gathered', 'Error in allocation of coastalflow_cpl0','','')
1233	CALL restget_p (rest_id, 'coastalflow_cpl0', nbp_glo, nflow, 1, kjit, .TRUE., coastalflow_cpl0, "gather", nbp_glo, index_g)
1234	IF ( ALL(coastalflow_cpl0(:,:) == val_exp) ) coastalflow_cpl0(:,:)=0
1235
1236	! Do not applay the filter now in initialization phase.
1237	! These variables will never be used anyway in the initialization phase.
1238	! Transform into m^3/s
1239	riverflow_cpl = zriver/xrdt
1240	coastalflow_cpl = zcoastal/xrdt
1241
1242	! Fields for emission variables : to be transport by GCM to chemistry model.
1243	DO i = 1, nb_fields_out
1244	SELECT CASE(TRIM(field_out_names(i)))
1245	CASE("fCO2_land")
1246	fields_out(:,i)=znetco2(:)
1247	CASE("fCO2_land_use")
1248	fields_out(:,i)=zcarblu(:)
1249	CASE DEFAULT
1250	CALL ipslerr_p (3,'intsurf_gathered', &
1251	'You ask from GCM an unknown field '//TRIM(field_out_names(i))//&
1252	' to ORCHIDEE for this specific version.',&
1253	'This model won''t be able to continue.', &
1254	'(check your tracer parameters in GCM)')
1255	END SELECT
1256	END DO
1257
1258	! Copy the number of vegetation types to local output variable
1259	IF (PRESENT(nvm_out)) nvm_out=nvm
1260
1261	IF(is_root_prc) CALL getin_dump
1262	lstep_init_intersurf = .FALSE.
1263
1264	CALL ipslnlf_p(new_number=old_fileout)
1265	!
1266	!! Change back to be in the LMDZ context for XIOS
1267	!
1268	CALL xios_orchidee_change_context("LMDZ")
1269
1270	IF (printlev_loc >= 2) WRITE (numout,*) 'End intersurf_initialize_gathered'
1271	IF (printlev_loc >= 1) WRITE (numout,*) 'Initialization phase for ORCHIDEE is finished.'
1272
1273	END SUBROUTINE intersurf_initialize_gathered
1274
1275
1276	!! =============================================================================================================================
1277	!! SUBROUTINE: intersurf_main_gathered
1278	!!
1279	!>\BRIEF Main subroutine to call ORCHIDEE from the gcm (LMDZ) using variables on a 1D grid with only land points.
1280	!!
1281	!! DESCRIPTION: This subroutine is the main interface for ORCHIDEE when it is called from the gcm (LMDZ).
1282	!! The variables are all gathered before entering this subroutine on the 1D grid with only landpoints.
1283	!!
1284	!! \n
1285	!_ ==============================================================================================================================
1286
1287	SUBROUTINE intersurf_main_gathered (kjit, iim_glo, jjm_glo, kjpindex, kindex, xrdt, &
1288	lrestart_read, lrestart_write, latlon, zcontfrac, zneighbours, zresolution, date0, &
1289	zlev, u, v, qair, temp_air, epot_air, ccanopy, &
1290	cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
1291	precip_rain, precip_snow, lwdown, swnet, swdown, pb, &
1292	vevapp, fluxsens, fluxlat, coastalflow_cpl, riverflow_cpl, &
1293	tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m,lon_scat_g, lat_scat_g, q2m, t2m, z0h, &
1294	veget, lai, height, &
1295	field_out_names, fields_out, field_in_names, fields_in, &
1296	coszang)
1297
1298	USE mod_orchidee_para
1299	IMPLICIT NONE
1300
1301	!! 0. Variable and parameter declaration
1302	!! 0.1 Input variables
1303	INTEGER(i_std),INTENT (in) :: kjit !! Time step number
1304	INTEGER(i_std),INTENT (in) :: iim_glo !! Dimension of global fields
1305	INTEGER(i_std),INTENT (in) :: jjm_glo !! Dimension of global fields
1306	INTEGER(i_std),INTENT (in) :: kjpindex !! Number of continental points
1307	REAL(r_std),INTENT (in) :: xrdt !! Time step in seconds
1308	LOGICAL, INTENT (in) :: lrestart_read !! Logical for _restart_ file to read
1309	LOGICAL, INTENT (in) :: lrestart_write !! Logical for _restart_ file to write'
1310	REAL(r_std), INTENT (in) :: date0 !! Date at which kjit = 0
1311	INTEGER(i_std),DIMENSION (kjpindex), INTENT (in) :: kindex !! Index for continental points
1312	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: u !! Lowest level wind speed
1313	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: v !! Lowest level wind speed
1314	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: zlev !! Height of first layer
1315	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: qair !! Lowest level specific humidity
1316	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: precip_rain !! Rain precipitation
1317	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: precip_snow !! Snow precipitation
1318	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: lwdown !! Down-welling long-wave flux
1319	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: swnet !! Net surface short-wave flux
1320	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: swdown !! Downwelling surface short-wave flux
1321	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: temp_air !! Air temperature in Kelvin
1322	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: epot_air !! Air potential energy
1323	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: ccanopy !! CO2 concentration in the canopy
1324	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: petAcoef !! Coeficients A from the PBL resolution
1325	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: peqAcoef !! One for T and another for q
1326	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: petBcoef !! Coeficients B from the PBL resolution
1327	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: peqBcoef !! One for T and another for q
1328	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: pb !! Surface pressure
1329	REAL(r_std),DIMENSION (kjpindex,2), INTENT(in) :: latlon !! Geographical coordinates
1330	REAL(r_std),DIMENSION (kjpindex), INTENT(in) :: zcontfrac !! Fraction of continent
1331	INTEGER(i_std),DIMENSION (kjpindex,NbNeighb), INTENT(in):: zneighbours !! neighbours
1332	REAL(r_std),DIMENSION (kjpindex,2), INTENT(in) :: zresolution !! size of the grid box
1333	REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN) :: lon_scat_g !! The scattered values for longitude
1334	REAL(r_std),DIMENSION (iim_glo,jjm_glo), INTENT(IN) :: lat_scat_g !! The scattered values for latitude
1335
1336	!! 0.2 Output variables
1337	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: z0m !! Surface roughness for momentum
1338	REAL(r_std),DIMENSION (kjpindex,nflow), INTENT(out) :: coastalflow_cpl !! Diffuse flow of water into the ocean, time filtered (m^3/s)
1339	REAL(r_std),DIMENSION (kjpindex,nflow), INTENT(out) :: riverflow_cpl !! Largest rivers flowing into the ocean, time filtered (m^3/s)
1340	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: tsol_rad !! Radiative surface temperature
1341	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: vevapp !! Total of evaporation
1342	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: temp_sol_new !! New soil temperature
1343	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: qsurf !! Surface specific humidity
1344	REAL(r_std),DIMENSION (kjpindex,2), INTENT(out) :: albedo !! Albedo
1345	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: fluxsens !! Sensible chaleur flux
1346	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: fluxlat !! Latent chaleur flux
1347	REAL(r_std),DIMENSION (kjpindex), INTENT(out) :: emis !! Emissivity
1348
1349	!! 0.3 Modified variables
1350	REAL(r_std),DIMENSION (kjpindex), INTENT(inout) :: cdrag !! Cdrag
1351
1352	!! 0.4 Optional input and output variables
1353	REAL(r_std),DIMENSION (kjpindex), INTENT(in), OPTIONAL:: q2m !! Surface specific humidity
1354	REAL(r_std),DIMENSION (kjpindex), INTENT(in), OPTIONAL:: t2m !! Surface air temperature
1355	REAL(r_std),DIMENSION (kjpindex), INTENT(out),OPTIONAL:: z0h !! Surface roughness for heat
1356	REAL(r_std),DIMENSION (kjpindex,nvm), INTENT(out), OPTIONAL :: veget !! Fraction of vegetation type (unitless, 0-1)
1357	REAL(r_std),DIMENSION (kjpindex,nvm), INTENT(out), OPTIONAL :: lai !! Leaf area index (m^2 m^{-2}
1358	REAL(r_std),DIMENSION (kjpindex,nvm), INTENT(out), OPTIONAL :: height !! Vegetation Height (m)
1359	REAL(r_std), DIMENSION(kjpindex), OPTIONAL, INTENT(in):: coszang !! Cosine of the solar zenith angle (unitless)
1360	CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_in_names !! Names for deposit variables to be transported
1361	!! from chemistry model by GCM to ORCHIDEE
1362	REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(IN) :: fields_in !! Fields for deposit variables to be transported
1363	!! from chemistry model by GCM to ORCHIDEE
1364	CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_out_names !! Names for emission variables to be transported
1365	!! to chemistry model by GCM from ORCHIDEE
1366	REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: fields_out !! Fields for emission variables to be transported
1367	!! to chemistry model by GCM from ORCHIDEE
1368
1369	!! 0.5 Local variables
1370	REAL(r_std),DIMENSION (kjpindex) :: zccanopy !! Work array to keep ccanopy
1371	REAL(r_std),DIMENSION (kjpindex) :: zprecip_rain !! Work array to keep precip_rain
1372	REAL(r_std),DIMENSION (kjpindex) :: zprecip_snow !! Work array to keep precip_snow
1373	REAL(r_std),DIMENSION (kjpindex) :: zz0m, zz0h !! Work array to keep zz0m, zz0h
1374	REAL(r_std),DIMENSION (kjpindex) :: zcdrag !! Work array for surface drag
1375	REAL(r_std),DIMENSION (kjpindex,nflow) :: zcoastal !! Work array to keep coastal flow
1376	REAL(r_std),DIMENSION (kjpindex,nflow) :: zriver !! Work array to keep river out flow
1377	REAL(r_std),DIMENSION (kjpindex) :: znetco2 !! Work array to keep netco2flux
1378	REAL(r_std),DIMENSION (kjpindex) :: zcarblu !! Work array to keep fco2_land_use
1379	REAL(r_std),DIMENSION (kjpindex) :: ztsol_rad !! Work array to keep tsol_rad
1380	REAL(r_std),DIMENSION (kjpindex) :: zvevapp !! Work array to keep vevapp
1381	REAL(r_std),DIMENSION (kjpindex) :: ztemp_sol_new !! Work array to keep temp_sol_new
1382	REAL(r_std),DIMENSION (kjpindex) :: zqsurf !! Work array to keep qsurf
1383	REAL(r_std),DIMENSION (kjpindex,2) :: zalbedo !! Work array to keep albedo
1384	REAL(r_std),DIMENSION (kjpindex) :: zfluxsens !! Work array to keep fluxsens
1385	REAL(r_std),DIMENSION (kjpindex) :: zfluxlat !! Work array to keep fluxlat
1386	REAL(r_std),DIMENSION (kjpindex) :: zemis !! Work array to keep emis
1387	REAL(r_std),DIMENSION (kjpindex) :: zcoszang !! Work array to keep coszang
1388	REAL(r_std),ALLOCATABLE, DIMENSION (:) :: soilth_lev !! Vertical soil axis for thermal scheme (m)
1389	REAL(r_std),DIMENSION (kjpindex) :: q2m_loc !! Work array for q2m or qair
1390	REAL(r_std),DIMENSION (kjpindex) :: t2m_loc !! Work array for t2m or temp_air
1391	REAL(r_std),DIMENSION (kjpindex,nvm) :: zveget !! Work array to keep veget
1392	REAL(r_std),DIMENSION (kjpindex,nvm) :: zlai !! Work array to keep lai
1393	REAL(r_std),DIMENSION (kjpindex,nvm) :: zheight !! Work array to keep height
1394	REAL(r_std),DIMENSION (:,:),ALLOCATABLE :: lon_scat !! The scattered values for longitude
1395	REAL(r_std),DIMENSION (:,:),ALLOCATABLE :: lat_scat !! The scattered values for latitude
1396	INTEGER(i_std) :: i, j, ik
1397	INTEGER(i_std) :: ier
1398	INTEGER(i_std) :: itau_sechiba
1399	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: tmp_lon, tmp_lat, tmp_lev
1400	INTEGER :: old_fileout !! old Logical Int for std IO output
1401	REAL,ALLOCATABLE,DIMENSION(:,:) :: lalo_mpi
1402	REAL(r_std),DIMENSION (kjpindex) :: landpoints !! Local landpoints vector
1403
1404
1405	IF (printlev_loc >= 3) WRITE(numout,*) 'Start intersurf_main_gathered'
1406	CALL ipslnlf_p(new_number=numout, old_number=old_fileout)
1407
1408	IF (lstep_init_intersurf) THEN
1409	! Test if q2m and t2m are present
1410	IF (PRESENT(q2m) .AND. PRESENT(t2m)) THEN
1411	ok_q2m_t2m=.TRUE.
1412	ELSE
1413	ok_q2m_t2m=.FALSE.
1414	ENDIF
1415
1416	! Test if field_out_names and field_in_names are present and if so, count
1417	! the number of extra fields to exchange.
1418	IF (PRESENT(field_out_names)) THEN
1419	nb_fields_out=SIZE(field_out_names)
1420	ELSE
1421	nb_fields_out=0
1422	ENDIF
1423
1424	IF (PRESENT(field_in_names)) THEN
1425	nb_fields_in=SIZE(field_in_names)
1426	ELSE
1427	nb_fields_in=0
1428	ENDIF
1429
1430	CALL intersurf_initialize_gathered (kjit, iim_glo, jjm_glo, kjpindex, kindex, xrdt, &
1431	lrestart_read, lrestart_write, latlon, zcontfrac, zneighbours, zresolution, date0, &
1432	zlev, u, v, qair, temp_air, epot_air, &
1433	cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
1434	precip_rain, precip_snow, lwdown, swnet, swdown, pb, &
1435	vevapp, fluxsens, fluxlat, coastalflow_cpl, riverflow_cpl, &
1436	tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m,lon_scat_g, lat_scat_g, &
1437	q2m, t2m, zz0h, &
1438	field_out_names=field_out_names, fields_out=fields_out, &
1439	field_in_names=field_in_names, fields_in=fields_in)
1440
1441	! z0h is a optional output variable. Check first if it is present in the call from LMDZ.
1442	IF ( PRESENT(z0h) ) z0h(:) = zz0h(:)
1443	! Return from subroutine intersurf_main_gathered
1444	RETURN
1445	END IF
1446
1447	!
1448	!! Change to be in the orchidee context for XIOS
1449	!
1450	CALL xios_orchidee_change_context("orchidee")
1451
1452	!
1453	! Shift the time step to phase the two models
1454	!
1455	itau_sechiba = kjit + itau_offset
1456
1457	! Update the calendar in xios by sending the new time step
1458	CALL xios_orchidee_update_calendar(itau_sechiba)
1459
1460	! Update the calendar and all time variables in module time
1461	CALL time_nextstep(itau_sechiba)
1462
1463	!
1464	! 1. Just change the units of some input fields
1465	!
1466	DO ik=1, kjpindex
1467
1468	zprecip_rain(ik) = precip_rain(ik)*xrdt
1469	zprecip_snow(ik) = precip_snow(ik)*xrdt
1470	zcdrag(ik) = cdrag(ik)
1471
1472	ENDDO
1473
1474	!>> VOC in coupled mode
1475	IF ( PRESENT(coszang) ) THEN
1476	zcoszang(:) = coszang(:)
1477	ELSE
1478	zcoszang(:) = zero
1479	ENDIF
1480
1481	! Fields for deposit variables : to be transport from chemistry model by GCM to ORCHIDEE.
1482	DO i = 1, nb_fields_in
1483	WRITE(numout,*) i," Champ = ",TRIM(field_in_names(i))
1484	SELECT CASE(TRIM(field_in_names(i)))
1485	CASE DEFAULT
1486	CALL ipslerr_p (3,'intsurf_gathered', &
1487	'You ask in GCM an unknown field '//TRIM(field_in_names(i))//&
1488	' to give to ORCHIDEE for this specific version.',&
1489	'This model won''t be able to continue.', &
1490	'(check your tracer parameters in GCM)')
1491	END SELECT
1492	ENDDO
1493
1494	!
1495	! 2. modification of co2
1496	!
1497	IF ( fatmco2 ) THEN
1498	zccanopy(:) = atmco2
1499	WRITE (numout,*) 'Modification of the ccanopy value. CO2 = ',atmco2
1500	ELSE
1501	zccanopy(:) = ccanopy(:)
1502	ENDIF
1503
1504	!
1505	! 4. call sechiba for continental points only
1506	!
1507	IF ( printlev_loc >= 3 ) WRITE(numout,*) 'Before call to sechiba_main from intersurf_main_gathered'
1508
1509
1510	IF (ok_q2m_t2m) THEN
1511	t2m_loc=t2m
1512	q2m_loc=q2m
1513	ELSE
1514	t2m_loc=temp_air
1515	q2m_loc=qair
1516	END IF
1517
1518	CALL sechiba_main (itau_sechiba, iim_g*jj_nb, kjpindex, kindex, &
1519	lrestart_read, lrestart_write, &
1520	lalo, contfrac, neighbours, resolution, &
1521	zlev, u, v, qair, q2m_loc, t2m_loc, temp_air, epot_air, zccanopy, &
1522	zcdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
1523	zprecip_rain ,zprecip_snow, lwdown, swnet, swdown, zcoszang, pb, &
1524	zvevapp, zfluxsens, zfluxlat, zcoastal, zriver, znetco2, zcarblu, &
1525	ztsol_rad, ztemp_sol_new, zqsurf, zalbedo, zemis, zz0m, zz0h, &
1526	zveget, zlai, zheight, &
1527	rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC )
1528
1529	IF ( printlev_loc>=3 ) WRITE(numout,*) 'After call to sechiba_main'
1530
1531	!
1532	! 6. scatter output fields
1533	!
1534	DO ik=1, kjpindex
1535	z0m(ik) = zz0m(ik)
1536	tsol_rad(ik) = ztsol_rad(ik)
1537	temp_sol_new(ik) = ztemp_sol_new(ik)
1538	qsurf(ik) = zqsurf(ik)
1539	albedo(ik,1) = zalbedo(ik,1)
1540	albedo(ik,2) = zalbedo(ik,2)
1541	fluxsens(ik) = zfluxsens(ik)
1542	fluxlat(ik) = zfluxlat(ik)
1543	emis(ik) = zemis(ik)
1544	cdrag(ik) = zcdrag(ik)
1545	! Transform the water fluxes into Kg/m^2/s
1546	vevapp(ik) = zvevapp(ik)/xrdt
1547	ENDDO
1548
1549	! Copy variables only if the optional variables are present in the call to the subroutine
1550	IF (PRESENT(veget)) veget(:,:) = zveget(:,:)
1551	IF (PRESENT(lai)) lai(:,:) = zlai(:,:)
1552	IF (PRESENT(height)) height(:,:) = zheight(:,:)
1553
1554	! Applay time filter to distribut the river- and coastalflow over a longer time period.
1555	! When coeff_rel=1(default case when tau_outflow=0), the distribution is instanteneous.
1556	! Use TAU_OUTFLOW in run.def to set the number of days of distribution.
1557	! The water fluxes zriver and zcoastal coming from sechiba are at the same time transfromed
1558	! from m^3/dt_sechiba into m^3/s by dividing with the sechiba time step (xrdt).
1559	riverflow_cpl = coeff_relzriver/xrdt + (1.-coeff_rel)riverflow_cpl0
1560	riverflow_cpl0 = riverflow_cpl
1561
1562	coastalflow_cpl = coeff_relzcoastal/xrdt + (1.-coeff_rel)coastalflow_cpl0
1563	coastalflow_cpl0 = coastalflow_cpl
1564
1565
1566	! z0h is a optional output variable. Check first if it is present in the call from LMDZ.
1567	IF ( PRESENT(z0h) ) z0h(:) = zz0h(:)
1568
1569	CALL xios_orchidee_send_field("LandPoints" ,(/ ( REAL(ik), ik=1,kjpindex ) /))
1570	CALL xios_orchidee_send_field("areas", area)
1571	CALL xios_orchidee_send_field("contfrac",contfrac)
1572	CALL xios_orchidee_send_field("temp_air",temp_air)
1573	CALL xios_orchidee_send_field("qair",qair)
1574	CALL xios_orchidee_send_field("swnet",swnet)
1575	CALL xios_orchidee_send_field("swdown",swdown)
1576	CALL xios_orchidee_send_field("pb",pb)
1577	CALL xios_orchidee_send_field("riverflow_cpl",riverflow_cpl(:,ih2o))
1578	CALL xios_orchidee_send_field("coastalflow_cpl",coastalflow_cpl(:,ih2o))
1579
1580
1581	IF (ok_q2m_t2m) THEN
1582	CALL xios_orchidee_send_field("t2m",t2m)
1583	CALL xios_orchidee_send_field("q2m",q2m)
1584	ELSE
1585	CALL xios_orchidee_send_field("t2m",temp_air)
1586	CALL xios_orchidee_send_field("q2m",qair)
1587	ENDIF
1588
1589	IF ( .NOT. almaoutput ) THEN
1590	!
1591	! scattered during the writing
1592	!
1593	CALL histwrite_p (hist_id, 'evap', itau_sechiba, zvevapp, kjpindex, kindex)
1594	CALL histwrite_p (hist_id, 'coastalflow',itau_sechiba, zcoastal(:,ih2o), kjpindex, kindex)
1595	CALL histwrite_p (hist_id, 'riverflow',itau_sechiba, zriver(:,ih2o), kjpindex, kindex)
1596	CALL histwrite_p (hist_id, 'temp_sol', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1597	CALL histwrite_p (hist_id, 'tsol_max', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1598	CALL histwrite_p (hist_id, 'tsol_min', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1599	CALL histwrite_p (hist_id, 'fluxsens', itau_sechiba, zfluxsens, kjpindex, kindex)
1600	CALL histwrite_p (hist_id, 'fluxlat', itau_sechiba, zfluxlat, kjpindex, kindex)
1601	CALL histwrite_p (hist_id, 'swnet', itau_sechiba, swnet, kjpindex, kindex)
1602	CALL histwrite_p (hist_id, 'swdown', itau_sechiba, swdown, kjpindex, kindex)
1603	CALL histwrite_p (hist_id, 'alb_vis', itau_sechiba, zalbedo(:,1), kjpindex, kindex)
1604	CALL histwrite_p (hist_id, 'alb_nir', itau_sechiba, zalbedo(:,2), kjpindex, kindex)
1605	CALL histwrite_p (hist_id, 'tair', itau_sechiba, temp_air, kjpindex, kindex)
1606	CALL histwrite_p (hist_id, 'qair', itau_sechiba, qair, kjpindex, kindex)
1607	IF (ok_q2m_t2m) THEN
1608	CALL histwrite_p (hist_id, 't2m', itau_sechiba, t2m, kjpindex, kindex)
1609	CALL histwrite_p (hist_id, 'q2m', itau_sechiba, q2m, kjpindex, kindex)
1610	ELSE
1611	CALL histwrite_p (hist_id, 't2m', itau_sechiba, temp_air, kjpindex, kindex)
1612	CALL histwrite_p (hist_id, 'q2m', itau_sechiba, qair, kjpindex, kindex)
1613	ENDIF
1614	IF ( hist2_id > 0 ) THEN
1615	CALL histwrite_p (hist2_id, 'evap', itau_sechiba, zvevapp, kjpindex, kindex)
1616	CALL histwrite_p (hist2_id, 'coastalflow',itau_sechiba, zcoastal(:,ih2o), kjpindex, kindex)
1617	CALL histwrite_p (hist2_id, 'riverflow',itau_sechiba, zriver(:,ih2o), kjpindex, kindex)
1618	CALL histwrite_p (hist2_id, 'temp_sol', itau_sechiba, temp_sol_new, kjpindex, kindex)
1619	CALL histwrite_p (hist2_id, 'tsol_max', itau_sechiba, temp_sol_new, kjpindex, kindex)
1620	CALL histwrite_p (hist2_id, 'tsol_min', itau_sechiba, temp_sol_new, kjpindex, kindex)
1621	CALL histwrite_p (hist2_id, 'fluxsens', itau_sechiba, zfluxsens, kjpindex, kindex)
1622	CALL histwrite_p (hist2_id, 'fluxlat', itau_sechiba, zfluxlat, kjpindex, kindex)
1623	CALL histwrite_p (hist2_id, 'swnet', itau_sechiba, swnet, kjpindex, kindex)
1624	CALL histwrite_p (hist2_id, 'swdown', itau_sechiba, swdown, kjpindex, kindex)
1625	CALL histwrite_p (hist2_id, 'alb_vis', itau_sechiba, zalbedo(:,1), kjpindex, kindex)
1626	CALL histwrite_p (hist2_id, 'alb_nir', itau_sechiba, zalbedo(:,2), kjpindex, kindex)
1627	CALL histwrite_p (hist2_id, 'tair', itau_sechiba, temp_air, kjpindex, kindex)
1628	CALL histwrite_p (hist2_id, 'qair', itau_sechiba, qair, kjpindex, kindex)
1629	IF (ok_q2m_t2m) THEN
1630	CALL histwrite_p (hist2_id, 't2m', itau_sechiba, t2m, kjpindex, kindex)
1631	CALL histwrite_p (hist2_id, 'q2m', itau_sechiba, q2m, kjpindex, kindex)
1632	ELSE
1633	CALL histwrite_p (hist2_id, 't2m', itau_sechiba, temp_air, kjpindex, kindex)
1634	CALL histwrite_p (hist2_id, 'q2m', itau_sechiba, qair, kjpindex, kindex)
1635	ENDIF
1636	ENDIF
1637	ELSE
1638	CALL histwrite_p (hist_id, 'Evap', itau_sechiba, zvevapp, kjpindex, kindex)
1639	CALL histwrite_p (hist_id, 'SWnet', itau_sechiba, swnet, kjpindex, kindex)
1640	CALL histwrite_p (hist_id, 'Qh', itau_sechiba, zfluxsens, kjpindex, kindex)
1641	CALL histwrite_p (hist_id, 'Qle', itau_sechiba, zfluxlat, kjpindex, kindex)
1642	CALL histwrite_p (hist_id, 'AvgSurfT', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1643	CALL histwrite_p (hist_id, 'RadT', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1644	IF ( hist2_id > 0 ) THEN
1645	CALL histwrite_p (hist2_id, 'Evap', itau_sechiba, zvevapp, kjpindex, kindex)
1646	CALL histwrite_p (hist2_id, 'SWnet', itau_sechiba, swnet, kjpindex, kindex)
1647	CALL histwrite_p (hist2_id, 'Qh', itau_sechiba, zfluxsens, kjpindex, kindex)
1648	CALL histwrite_p (hist2_id, 'Qle', itau_sechiba, zfluxlat, kjpindex, kindex)
1649	CALL histwrite_p (hist2_id, 'AvgSurfT', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1650	CALL histwrite_p (hist2_id, 'RadT', itau_sechiba, ztemp_sol_new, kjpindex, kindex)
1651	ENDIF
1652	ENDIF
1653
1654	! Syncronize output but only if flag ok_histsync is set to true
1655	IF (ok_histsync .AND. is_omp_root) THEN
1656	IF ( (dw .EQ. xrdt) .AND. hist_id > 0 ) THEN
1657	CALL histsync(hist_id)
1658	ENDIF
1659	ENDIF
1660
1661
1662	! Fields for emission variables : to be transport by GCM to chemistry model.
1663	DO i = 1, nb_fields_out
1664	SELECT CASE(TRIM(field_out_names(i)))
1665	CASE("fCO2_land")
1666	fields_out(:,i)=znetco2(:)
1667	CASE("fCO2_land_use")
1668	fields_out(:,i)=zcarblu(:)
1669	CASE DEFAULT
1670	CALL ipslerr_p (3,'intsurf_gathered', &
1671	& 'You ask from GCM an unknown field '//TRIM(field_out_names(i))//&
1672	& ' to ORCHIDEE for this specific version.',&
1673	& 'This model won''t be able to continue.', &
1674	& '(check your tracer parameters in GCM)')
1675	END SELECT
1676	ENDDO
1677
1678
1679	! Write variables to restart file
1680	IF (lrestart_write) THEN
1681	CALL restput_p (rest_id, 'riverflow_cpl0', nbp_glo, nflow, 1, kjit, riverflow_cpl0, 'scatter', nbp_glo, index_g)
1682	CALL restput_p (rest_id, 'coastalflow_cpl0', nbp_glo, nflow, 1, kjit, coastalflow_cpl0, 'scatter', nbp_glo, index_g)
1683	END IF
1684
1685	!
1686	CALL ipslnlf_p(new_number=old_fileout)
1687	!
1688
1689	!
1690	!! Finalize the XIOS orchidee context if it is the last call
1691	!
1692	IF (lrestart_write) THEN
1693	CALL xios_orchidee_context_finalize
1694	END IF
1695	!
1696	!! Change back to be in the LMDZ context for XIOS
1697	!
1698	CALL xios_orchidee_change_context("LMDZ")
1699
1700	IF (printlev_loc >= 3) WRITE (numout,*) 'End intersurf_main_gathered'
1701
1702	END SUBROUTINE intersurf_main_gathered
1703
1704	!! ================================================================================================================================
1705	!! SUBROUTINE : intersurf_clear
1706	!!
1707	!>\BRIEF Clear intersurf module and underlaying modules
1708	!!
1709	!! DESCRIPTION : Deallocate memory and reset initialization variables to there original values.
1710	!! Call the clearing for sechiba module.
1711	!!
1712	!_ ================================================================================================================================
1713	SUBROUTINE intersurf_clear
1714	CALL sechiba_clear
1715	END SUBROUTINE intersurf_clear
1716
1717	END MODULE intersurf

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

Context Navigation

source: branches/publications/ORCHIDEE-Clateral/src_sechiba/intersurf.f90 @ 7329

Download in other formats: