Context Navigation

source: tags/ORCHIDEE_2_0/ORCHIDEE/src_sechiba/intersurf.f90 @ 5391

Last change on this file since 5391 was 5391, checked in by josefine.ghattas, 6 years ago
Integrated changeset [5389] done on the trunk in the tag. This is needed for ESM configurations with transported carbon. P Cadule
Property svn:keywords set to `Revision Date HeadURL Date Author Revision`
File size: 88.3 KB

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

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

Download in other formats: