Context Navigation

source: branches/publications/ORCHIDEE_DFv1.0_site/src_sechiba/intersurf.f90 @ 6715

Last change on this file since 6715 was 5396, checked in by josefine.ghattas, 6 years ago

Added for ESM carbon coupling.

Cadule

Property svn:keywords set to Revision Date HeadURL Date Author Revision

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

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

Download in other formats: