Context Navigation

source: tags/ORCHIDEE_4_1/ORCHIDEE/src_sechiba/intersurf.f90 @ 7761

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

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

Download in other formats: