Context Navigation

source: branches/publications/ORCHIDEE-ICE_SurfaceMassBalance/src_driver/orchideedriver.f90 @ 8398

Last change on this file since 8398 was 4287, checked in by josefine.ghattas, 7 years ago
Enhencement on clear functions. Added call to clear functions from all offline drivers. See ticket #232
File size: 38.8 KB

Line
1	! =================================================================================================================================
2	! PROGRAM : orchideedriver
3	!
4	! CONTACT : jan.polcher@lmd.jussieu.fr
5	!
6	! LICENCE : IPSL (2016)
7	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8	!
9	!>\BRIEF This is the main program for the new driver. This only organises the data and calls sechiba_main.
10	!! The main work is done in glogrid.f90 and forcing_tools.f90.
11	!!
12	!!\n DESCRIPTION: Call the various modules to get the forcing data and provide it to SECHIBA. The only complexity
13	!! is setting-up the domain decomposition and distributing the grid information.
14	!! The code is parallel from tip to toe using the domain decomposition inherited from LMDZ.
15	!!
16	!! RECENT CHANGE(S): None
17	!!
18	!! REFERENCE(S) :
19	!!
20	!! SVN :
21	!! $HeadURL: $
22	!! $Date: $
23	!! $Revision: $
24	!! \n
25	!_ ================================================================================================================================
26	!
27	PROGRAM orchidedriver
28	!---------------------------------------------------------------------
29	!-
30	!-
31	!---------------------------------------------------------------------
32	USE defprec
33	USE netcdf
34	!
35	!
36	USE ioipsl_para
37	USE mod_orchidee_para
38	!
39	USE grid
40	USE timer
41	!
42	USE forcing_tools
43	USE globgrd
44	!
45	USE sechiba
46	USE control
47	USE ioipslctrl
48	!
49	USE thermosoilc, ONLY : thermosoilc_levels
50	!
51	!-
52	IMPLICIT NONE
53	!-
54	CHARACTER(LEN=80) :: gridfilename
55	CHARACTER(LEN=80), DIMENSION(100) :: forfilename
56	INTEGER(i_std) :: nb_forcefile
57	CHARACTER(LEN=8) :: model_guess
58	INTEGER(i_std) :: iim_glo, jjm_glo, file_id
59	!-
60	INTEGER(i_std) :: nbseg
61	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: lon_glo, lat_glo, area_glo
62	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: mask_glo
63	REAL(r_std), ALLOCATABLE, DIMENSION(:,:,:,:) :: corners_glo
64	REAL(r_std), ALLOCATABLE, DIMENSION(:,:,:) :: corners_lon, corners_lat
65	INTEGER(i_std) :: nbindex_g, kjpindex
66	INTEGER(i_std), ALLOCATABLE, DIMENSION(:) :: kindex, kindex_g
67	REAL(r_std), DIMENSION(2) :: zoom_lon, zoom_lat
68	!
69	! Variables for the global grid available on all procs and used
70	! to fill the ORCHIDEE variable on the root_proc
71	!
72	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: lalo_glo
73	REAL(r_std), ALLOCATABLE, DIMENSION(:) :: contfrac_glo
74	CHARACTER(LEN=20) :: calendar
75	!-
76	!- Variables local to each processors.
77	!-
78	INTEGER(i_std) :: i, j, ik, nbdt, first_point
79	INTEGER(i_std) :: itau, itau_offset, itau_sechiba
80	REAL(r_std) :: date0, date0_shifted, dt, julian, julian0
81	REAL(r_std) :: date0_tmp, dt_tmp
82	INTEGER(i_std) :: nbdt_tmp
83	REAL(r_std) :: timestep_interval(2), timestep_int_next(2)
84	!
85	INTEGER(i_std) :: rest_id, rest_id_stom
86	INTEGER(i_std) :: hist_id, hist2_id, hist_id_stom, hist_id_stom_IPCC
87	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: lalo_loc
88	INTEGER(i_std) :: iim, jjm, ier
89	REAL(r_std), ALLOCATABLE, DIMENSION(:,:) :: lon, lat
90	REAL(r_std),ALLOCATABLE, DIMENSION (:) :: soilth_lev !! Vertical soil axis for thermal scheme (m)
91	!-
92	!- input fields
93	!-
94	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: u !! Lowest level wind speed
95	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: v !! Lowest level wind speed
96	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: zlev_uv !! Height of first layer
97	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: zlev_tq !! Height of first layer
98	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: qair !! Lowest level specific humidity
99	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: precip_rain !! Rain precipitation
100	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: precip_snow !! Snow precipitation
101	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: lwdown !! Down-welling long-wave flux
102	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: swdown !! Downwelling surface short-wave flux
103	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: sinang !! cosine of solar zenith angle
104	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: temp_air !! Air temperature in Kelvin
105	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: epot_air !! Air potential energy
106	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: ccanopy !! CO2 concentration in the canopy
107	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: petAcoef !! Coeficients A from the PBL resolution
108	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: peqAcoef !! One for T and another for q
109	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: petBcoef !! Coeficients B from the PBL resolution
110	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: peqBcoef !! One for T and another for q
111	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: cdrag !! Cdrag
112	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: pb !! Lowest level pressure
113	!-
114	!- output fields
115	!-
116	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: z0m !! Surface roughness for momentum (m)
117	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: z0h !! Surface roughness for heat (m)
118	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: coastalflow !! Diffuse flow of water into the ocean (m^3/dt)
119	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: riverflow !! Largest rivers flowing into the ocean (m^3/dt)
120	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: tsol_rad !! Radiative surface temperature
121	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: vevapp !! Total of evaporation
122	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: temp_sol_new !! New soil temperature
123	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: qsurf !! Surface specific humidity
124	REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: albedo !! Albedo
125	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: fluxsens !! Sensible chaleur flux
126	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: fluxlat !! Latent chaleur flux
127	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: emis !! Emissivity
128	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: netco2 !! netco2flux
129	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: carblu !! fco2_land_use
130	!-
131	!-
132	!-
133	REAL(r_std) :: atmco2
134	REAL(r_std), ALLOCATABLE, DIMENSION (:) :: u_tq, v_tq, swnet
135	LOGICAL :: lrestart_read = .TRUE. !! Logical for _restart_ file to read
136	LOGICAL :: lrestart_write = .FALSE. !! Logical for _restart_ file to write'
137	!
138	! Timer variables
139	!
140	LOGICAL, PARAMETER :: timemeasure=.TRUE.
141	REAL(r_std) :: waitput_cputime=0.0, waitget_cputime=0.0, orchidee_cputime=0.0
142	REAL(r_std) :: waitput_walltime=0.0, waitget_walltime=0.0, orchidee_walltime=0.0
143	!
144	!
145	! Print point
146	!
147	!! REAL(r_std), DIMENSION(2) :: testpt=(/44.8,-25.3/)
148	!! REAL(r_std), DIMENSION(2) :: testpt=(/44.8,-18.3/)
149	!! REAL(r_std), DIMENSION(2) :: testpt=(/-60.25,-5.25/)
150	!! REAL(r_std), DIMENSION(2) :: testpt=(/46.7,10.3/)
151	!! REAL(r_std), DIMENSION(2) :: testpt=(/0.25,49.25/)
152	! Case when no ouput is desired.
153	REAL(r_std), DIMENSION(2) :: testpt=(/9999.99,9999.99/)
154	INTEGER(i_std) :: ktest
155
156	OFF_LINE_MODE = .TRUE.
157
158	!-
159	!---------------------------------------------------------------------------------------
160	!-
161	!- Define MPI communicator
162	!-
163	!---------------------------------------------------------------------------------------
164	!-
165	!
166	! Set parallel processing in ORCHIDEE
167	!
168	CALL Init_orchidee_para()
169	!
170	!====================================================================================
171	!
172	! Start timer now that the paralelisation is initialized.
173	!
174	IF ( timemeasure ) THEN
175	CALL init_timer
176	CALL start_timer(timer_global)
177	CALL start_timer(timer_mpi)
178	ENDIF
179	!
180	!
181	!---------------------------------------------------------------------------------------
182	!-
183	!- Start the getconf processes
184	!-
185	!---------------------------------------------------------------------------------------
186	!-
187	!! CALL getin_name("run.def")
188	!-
189	!Config Key = GRID_FILE
190	!Config Desc = Name of file containing the forcing data
191	!Config If = [-]
192	!Config Def = grid_file.nc
193	!Config Help = This is the name of the file from which we will read
194	!Config or write into it the description of the grid from
195	!Config the forcing file.
196	!Config compliant.
197	!Config Units = [FILE]
198	!-
199	gridfilename='NONE'
200	CALL getin_p('GRID_FILE', gridfilename)
201	!-
202	forfilename(:)=" "
203	forfilename(1)='forcing_file.nc'
204	CALL getin_p('FORCING_FILE', forfilename)
205	!-
206	!- Define the zoom
207	!-
208	zoom_lon=(/-180,180/)
209	zoom_lat=(/-90,90/)
210	!
211	!Config Key = LIMIT_WEST
212	!Config Desc = Western limit of region
213	!Config If = [-]
214	!Config Def = -180.
215	!Config Help = Western limit of the region we are
216	!Config interested in. Between -180 and +180 degrees
217	!Config The model will use the smalest regions from
218	!Config region specified here and the one of the forcing file.
219	!Config Units = [Degrees]
220	!-
221	CALL getin_p('LIMIT_WEST',zoom_lon(1))
222	!-
223	!Config Key = LIMIT_EAST
224	!Config Desc = Eastern limit of region
225	!Config If = [-]
226	!Config Def = 180.
227	!Config Help = Eastern limit of the region we are
228	!Config interested in. Between -180 and +180 degrees
229	!Config The model will use the smalest regions from
230	!Config region specified here and the one of the forcing file.
231	!Config Units = [Degrees]
232	!-
233	CALL getin_p('LIMIT_EAST',zoom_lon(2))
234	!-
235	!Config Key = LIMIT_NORTH
236	!Config Desc = Northern limit of region
237	!Config If = [-]
238	!Config Def = 90.
239	!Config Help = Northern limit of the region we are
240	!Config interested in. Between +90 and -90 degrees
241	!Config The model will use the smalest regions from
242	!Config region specified here and the one of the forcing file.
243	!Config Units = [Degrees]
244	!-
245	CALL getin_p('LIMIT_NORTH',zoom_lat(2))
246	!-
247	!Config Key = LIMIT_SOUTH
248	!Config Desc = Southern limit of region
249	!Config If = [-]
250	!Config Def = -90.
251	!Config Help = Southern limit of the region we are
252	!Config interested in. Between 90 and -90 degrees
253	!Config The model will use the smalest regions from
254	!Config region specified here and the one of the forcing file.
255	!Config Units = [Degrees]
256	!-
257	CALL getin_p('LIMIT_SOUTH',zoom_lat(1))
258	IF ( (zoom_lon(1)+180 < EPSILON(zoom_lon(1))) .AND. (zoom_lon(2)-180 < EPSILON(zoom_lon(2))) .AND.&
259	&(zoom_lat(1)+90 < EPSILON(zoom_lat(1))) .AND. (zoom_lat(2)-90 < EPSILON(zoom_lat(2))) ) THEN
260	!
261	!Config Key = WEST_EAST
262	!Config Desc = Longitude interval to use from the forcing data
263	!Config If = [-]
264	!Config Def = -180, 180
265	!Config Help = This function allows to zoom into the forcing data
266	!Config Units = [degrees east]
267	!-
268	CALL getin_p('WEST_EAST', zoom_lon)
269	!
270	!Config Key = SOUTH_NORTH
271	!Config Desc = Latitude interval to use from the forcing data
272	!Config If = [-]
273	!Config Def = -90, 90
274	!Config Help = This function allows to zoom into the forcing data
275	!Config Units = [degrees north]
276	!-
277	CALL getin_p('SOUTH_NORTH', zoom_lat)
278	ENDIF
279	!-
280	!-
281	!- Get some basic variables from the run.def
282	!-
283	atmco2=350.
284	CALL getin_p('ATM_CO2',atmco2)
285	!
286	!====================================================================================
287	!-
288	!-
289	!- Get the grid on all processors.
290	!-
291	!---------------------------------------------------------------------------------------
292	!-
293	!- Read the grid, only on the root proc. from the forcing file using tools in the globgrd module.
294	!- The grid is then broadcast to all other broadcast.
295	!
296	nb_forcefile = 0
297	DO ik=1,100
298	IF ( INDEX(forfilename(ik), '.nc') > 0 ) nb_forcefile = nb_forcefile+1
299	ENDDO
300	!
301	IF ( is_root_prc) THEN
302	CALL globgrd_getdomsz(gridfilename, iim_glo, jjm_glo, nbindex_g, model_guess, file_id, forfilename, zoom_lon, zoom_lat)
303	nbseg = 4
304	ENDIF
305	!
306	CALL bcast(iim_glo)
307	CALL bcast(jjm_glo)
308	CALL bcast(nbindex_g)
309	CALL bcast(nbseg)
310	!-
311	!- Allocation of memory
312	!- variables over the entire grid (thus in x,y)
313	ALLOCATE(lon_glo(iim_glo, jjm_glo))
314	ALLOCATE(lat_glo(iim_glo, jjm_glo))
315	ALLOCATE(mask_glo(iim_glo, jjm_glo))
316	ALLOCATE(area_glo(iim_glo, jjm_glo))
317	ALLOCATE(corners_glo(iim_glo, jjm_glo, nbseg, 2))
318	!
319	! Gathered variables
320	ALLOCATE(kindex_g(nbindex_g))
321	ALLOCATE(contfrac_glo(nbindex_g))
322	!-
323	IF ( is_root_prc) THEN
324	CALL globgrd_getgrid(file_id, iim_glo, jjm_glo, nbindex_g, model_guess, &
325	& lon_glo, lat_glo, mask_glo, area_glo, corners_glo,&
326	& kindex_g, contfrac_glo, calendar)
327	ENDIF
328	!
329	CALL bcast(lon_glo)
330	CALL bcast(lat_glo)
331	CALL bcast(mask_glo)
332	CALL bcast(area_glo)
333	CALL bcast(corners_glo)
334	CALL bcast(kindex_g)
335	CALL bcast(contfrac_glo)
336	CALL bcast(calendar)
337	CALL bcast(model_guess)
338	!
339	ALLOCATE(lalo_glo(nbindex_g,2))
340	DO ik=1,nbindex_g
341	!
342	j = ((kindex_g(ik)-1)/iim_glo)+1
343	i = (kindex_g(ik)-(j-1)*iim_glo)
344	!
345	IF ( i > iim_glo .OR. j > jjm_glo ) THEN
346	WRITE(100+mpi_rank,*) "Error in the indexing (ik, kindex, i, j) : ", ik, kindex(ik), i, j
347	STOP "ERROR in orchideedriver"
348	ENDIF
349	!
350	lalo_glo(ik,1) = lat_glo(i,j)
351	lalo_glo(ik,2) = lon_glo(i,j)
352	!
353	ENDDO
354	!
355	WRITE(,) "Rank", mpi_rank, " Before parallel region All land points : ", nbindex_g
356	WRITE(,) "Rank", mpi_rank, " from ", iim_glo, " point in Lon. and ", jjm_glo, "in Lat."
357	!-
358	!---------------------------------------------------------------------------------------
359	!-
360	!- Now that the grid is distributed on all procs we can start
361	!- initialise the ORCHIDEE domain on each proc (longitude, latitude, indices)
362	!-
363	!---------------------------------------------------------------------------------------
364	!-
365	!- init_data_para also transfers kindex_g to index_g (the variable used in ORCHIDEE)
366	!-
367	CALL grid_set_glo(iim_glo, jjm_glo, nbindex_g)
368	CALL grid_allocate_glo(nbseg)
369	! Copy the list of indexes of land points into index_g used by ORCHIDEE and then broacast to all
370	! processors
371	CALL bcast(nbindex_g)
372	IF ( is_root_prc) index_g = kindex_g
373	CALL bcast(index_g)
374	!
375	WRITE(numout,*) "Rank", mpi_rank, "Into Init_orchidee_data_para_driver with ", nbindex_g
376	WRITE(numout,*) "Rank", mpi_rank, "Into ", index_g(1), index_g(nbindex_g)
377	!
378	CALL Init_orchidee_data_para_driver(nbindex_g,index_g)
379	CALL init_ioipsl_para
380	!
381	WRITE(numout,*) "Rank", mpi_rank, "After init_data_para global size : ", nbp_glo, SIZE(index_g), iim_g, iim_glo, jjm_g, jjm_glo
382	WRITE(numout,'("After init_data_para local : ij_nb, jj_nb",2I4)') iim_glo, jj_nb
383	!
384	! Allocate grid on the local processor
385	!
386	IF ( model_guess == "regular") THEN
387	CALL grid_init (nbp_loc, nbseg, "RegLonLat", "ForcingGrid")
388	ELSE IF ( model_guess == "WRF") THEN
389	CALL grid_init (nbp_loc, nbseg, "RegXY", "WRFGrid")
390	ELSE
391	CALL ipslerr(3, "orchidedriver", "The grid found in the GRID_FILE is not supported by ORCHIDEE", "", "")
392	ENDIF
393	!
394	!
395	! Transfer the global grid variables to the ORCHIDEE version on the root proc
396	! _glo -> _g
397	! Variables *_g were allocated with the CALL init_grid
398	!
399	IF ( is_root_prc) THEN
400	!
401	lalo_g(:,:) = lalo_glo(:,:)
402	lon_g(:,:) = lon_glo(:,:)
403	lat_g(:,:) = lat_glo(:,:)
404	!
405	ENDIF
406	!
407	!
408	! Set the local dimensions of the fields
409	!
410	iim = iim_glo
411	jjm = jj_nb
412	kjpindex = nbp_loc
413	!
414	WRITE(numout,*) mpi_rank, "DIMENSIONS of grid on processor : iim, jjm, kjpindex = ", iim, jjm, kjpindex
415	!
416	! Allocate the local arrays we need :
417	!
418	ALLOCATE(lon(iim,jjm), lat(iim,jjm))
419	ALLOCATE(kindex(kjpindex))
420	!
421	lon=lon_glo(:,jj_para_begin(mpi_rank):jj_para_end(mpi_rank))
422	lat=lat_glo(:,jj_para_begin(mpi_rank):jj_para_end(mpi_rank))
423	!
424	!
425	! Redistribute the indeces on all procs (apple distribution of land points)
426	!
427	CALL bcast(lon_g)
428	CALL bcast(lat_g)
429	CALL scatter(index_g, kindex)
430	!
431	!
432	! Apply the offset needed so that kindex refers to the index of the land point
433	! on the current region, i.e. the local lon lat domain.
434	!
435	kindex(1:kjpindex)=kindex(1:kjpindex)-(jj_begin-1)*iim_glo
436	!
437	! This routine transforms the global grid into a series of polygons for all land
438	! points identified by index_g.
439	!
440	CALL grid_stuff(nbindex_g, iim_g, jjm_g, lon_g, lat_g, index_g, contfrac_glo)
441	!
442	! Distribute the global lalo to the local processor level lalo
443	!
444	ALLOCATE(lalo_loc(kjpindex,2))
445	CALL scatter(lalo_glo, lalo_loc)
446	lalo(:,:) = lalo_loc(:,:)
447	!
448	!====================================================================================
449	!-
450	!- Prepare the time for the simulation
451	!-
452	!- Set the calendar and get some information
453	!-
454	CALL ioconf_calendar(calendar)
455	CALL ioget_calendar(one_year, one_day)
456	!-
457	!- get the time period for the run
458	!-
459	CALL forcing_integration_time(date0, dt, nbdt)
460	!
461	!
462	!
463	!====================================================================================
464	!-
465	!- Initialize the forcing files and prepare the time stepping through the data.
466	!-
467	!
468	CALL forcing_open(forfilename, iim_glo, jjm_glo, lon_glo, lat_glo, nbindex_g, zoom_lon, zoom_lat, &
469	& index_g, kjpindex, numout)
470	!
471	!
472	ALLOCATE(zlev_tq(kjpindex), zlev_uv(kjpindex))
473	ALLOCATE(u(kjpindex), v(kjpindex), pb(kjpindex))
474	ALLOCATE(temp_air(kjpindex))
475	ALLOCATE(qair(kjpindex))
476	ALLOCATE(petAcoef(kjpindex), peqAcoef(kjpindex), petBcoef(kjpindex), peqBcoef(kjpindex))
477	ALLOCATE(ccanopy(kjpindex))
478	ALLOCATE(cdrag(kjpindex))
479	ALLOCATE(precip_rain(kjpindex))
480	ALLOCATE(precip_snow(kjpindex))
481	ALLOCATE(swdown(kjpindex))
482	ALLOCATE(swnet(kjpindex))
483	ALLOCATE(lwdown(kjpindex))
484	ALLOCATE(sinang(kjpindex))
485	ALLOCATE(vevapp(kjpindex))
486	ALLOCATE(fluxsens(kjpindex))
487	ALLOCATE(fluxlat(kjpindex))
488	ALLOCATE(coastalflow(kjpindex))
489	ALLOCATE(riverflow(kjpindex))
490	ALLOCATE(netco2(kjpindex))
491	ALLOCATE(carblu(kjpindex))
492	ALLOCATE(tsol_rad(kjpindex))
493	ALLOCATE(temp_sol_new(kjpindex))
494	ALLOCATE(qsurf(kjpindex))
495	ALLOCATE(albedo(kjpindex,2))
496	ALLOCATE(emis(kjpindex))
497	ALLOCATE(epot_air(kjpindex))
498	ALLOCATE(u_tq(kjpindex), v_tq(kjpindex))
499	ALLOCATE(z0m(kjpindex))
500	ALLOCATE(z0h(kjpindex))
501	!-
502	!---------------------------------------------------------------------------------------
503	!-
504	!- Get a first set of forcing data
505	!-
506	!---------------------------------------------------------------------------------------
507	!-
508	!- Some default values so that the operations before the ORCHIDEE initialisation do not fail.
509	!-
510	z0m(:) = 0.1
511	albedo(:,:) = 0.13
512	!-
513	!====================================================================================
514	!-
515	!- Initialise the ORCHIDEE system in 4 steps :
516	!- 1 The control flags,
517	!- 2 the restart system of IOIPSL
518	!- 3 The history mechanism
519	!- 4 Finally the first call to SECHIBA will initialise all the internal variables
520	!
521	CALL control_initialize(dt)
522	!
523	itau = 0
524	!
525	CALL ioipslctrl_restini(itau, date0, dt, rest_id, rest_id_stom, itau_offset, date0_shifted)
526	WRITE(numout,*) "itau_offset : ", itau_offset, date0, date0_shifted
527	WRITE(numout,*) "itau_offset diff = ", date0_shifted, date0, date0_shifted-date0
528	!
529	! Get the vertical soil levels for the thermal scheme, to be used in xios_orchidee_init
530	ALLOCATE(soilth_lev(ngrnd), stat=ier)
531	IF (ier /= 0) CALL ipslerr_p(3,'orchideedriver', 'Error in allocation of soilth_lev','','')
532	IF (hydrol_cwrr) THEN
533	soilth_lev(1:ngrnd) = znt(:)
534	ELSE
535	soilth_lev(1:ngrnd) = thermosoilc_levels()
536	END IF
537	!
538	! To ensure that itau starts with 0 at date0 for the restart, we have to set an off-set to achieve this.
539	! itau_offset will get used to prduce itau_sechiba.
540	!
541	itau_offset=-itau_offset-1
542	!
543	! Get the date of the first time step
544	!
545	julian = date0 + 0.5*(dt/one_day)
546	CALL ju2ymds (julian, year, month, day, sec)
547	WRITE(,) "itau_offset : date0 : ", year, month, day, sec
548	!
549	CALL xios_orchidee_init( MPI_COMM_ORCH, &
550	date0, year, month, day, &
551	lon, lat, soilth_lev)
552	!
553	!- Initialize IOIPSL sechiba output files
554	itau_sechiba = itau+itau_offset
555	CALL ioipslctrl_history(iim, jjm, lon, lat, kindex, kjpindex, itau_sechiba, &
556	date0, dt, hist_id, hist2_id, hist_id_stom, hist_id_stom_IPCC)
557	WRITE(numout,*) "HISTORY : Defined for ", itau_sechiba, date0, dt
558	!
559	!-
560	!---------------------------------------------------------------------------------------
561	!-
562	!- Go into the time loop
563	!-
564	!---------------------------------------------------------------------------------------
565	!-
566	DO itau = 1,nbdt
567	!
568	timestep_interval(1) = date0 + (itau-1)*(dt/one_day)
569	timestep_interval(2) = date0 + itau*(dt/one_day)
570	julian = date0 + (itau-0.5)*(dt/one_day)
571	!
572	! Get the forcing data
573	!
574	CALL forcing_getvalues(timestep_interval, dt, zlev_tq, zlev_uv, temp_air, qair, &
575	& precip_rain, precip_snow, swdown, lwdown, sinang, u, v, pb)
576	!-
577	in_julian = itau2date(itau, date0, dt)
578	CALL ju2ymds (julian, year, month, day, sec)
579	CALL ymds2ju (year,1,1,zero, julian0)
580	julian_diff = in_julian-julian0
581	!
582	IF ( itau == nbdt ) lrestart_write = .TRUE.
583	!
584	! Adaptation of the forcing data to SECHIBA's needs
585	!
586	! Contrary to what the documentation says, ORCHIDEE expects surface pressure in hPa.
587	pb(:) = pb(:)/100.
588	epot_air(:) = cp_airtemp_air(:)+cte_gravzlev_tq(:)
589	ccanopy(:) = atmco2
590	cdrag(:) = 0.0
591	!
592	petBcoef(:) = epot_air(:)
593	peqBcoef(:) = qair(:)
594	petAcoef(:) = zero
595	peqAcoef(:) = zero
596	!
597	! Interpolate the wind (which is at hight zlev_uv) to the same height
598	! as the temperature and humidity (at zlev_tq).
599	!
600	u_tq(:) = u(:)*LOG(zlev_tq(:)/z0m(:))/LOG(zlev_uv(:)/z0m(:))
601	v_tq(:) = v(:)*LOG(zlev_tq(:)/z0m(:))/LOG(zlev_uv(:)/z0m(:))
602	!
603	!
604	swnet(:) =(1.-(albedo(:,1)+albedo(:,2))/2.)*swdown(:)
605	!
606	!
607	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, temp_air, "RECEIVED Air temperature")
608	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, qair, "RECEIVED Air humidity")
609	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, precip_rain*one_day, "RECEIVED Rainfall")
610	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, precip_snow*one_day, "RECEIVED Snowfall")
611	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, swnet, "RECEIVED net solar")
612	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, lwdown, "RECEIVED lwdown")
613	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, u, "RECEIVED East-ward wind")
614	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, v, "RECEIVED North-ward wind")
615	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, pb*100, "RECEIVED surface pressure")
616	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, zlev_uv, "RECEIVED UV height")
617	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, zlev_tq, "RECEIVED TQ height")
618	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, sinang, "RECEIVED sinang")
619	!
620	IF ( itau .NE. 1 ) THEN
621	IF ( timemeasure ) THEN
622	waitget_cputime = waitget_cputime + Get_cpu_Time(timer_global)
623	waitget_walltime = waitget_walltime + Get_real_Time(timer_global)
624	CALL stop_timer(timer_global)
625	CALL start_timer(timer_global)
626	ENDIF
627	ENDIF
628	!
629	!---------------------------------------------------------------------------------------
630	!-
631	!- IF first time step : Call to SECHIBA_initialize to set-up ORCHIDEE before doing an actual call
632	!- which will provide the first fluxes.
633	!-
634	!---------------------------------------------------------------------------------------
635	!
636	itau_sechiba = itau+itau_offset
637	!
638	! Update the calendar in xios by sending the new time step
639	CALL xios_orchidee_update_calendar(itau_sechiba)
640	!
641	IF ( itau == 1 ) THEN
642	!
643	IF ( timemeasure ) THEN
644	WRITE(numout,*) '------> CPU Time for start-up of main : ',Get_cpu_Time(timer_global)
645	WRITE(numout,*) '------> Real Time for start-up of main : ',Get_real_Time(timer_global)
646	CALL stop_timer(timer_global)
647	CALL start_timer(timer_global)
648	ENDIF
649	!
650	CALL sechiba_initialize( &
651	itau_sechiba, iim*jjm, kjpindex, kindex, &
652	lalo_loc, contfrac, neighbours, resolution, zlev_tq, &
653	u_tq, v_tq, qair, temp_air, temp_air, &
654	petAcoef, peqAcoef, petBcoef, peqBcoef, &
655	precip_rain, precip_snow, lwdown, swnet, swdown, &
656	pb, rest_id, hist_id, hist2_id, &
657	rest_id_stom, hist_id_stom, hist_id_stom_IPCC, &
658	coastalflow, riverflow, tsol_rad, vevapp, qsurf, &
659	z0m, z0h, albedo, fluxsens, fluxlat, emis, &
660	netco2, carblu, temp_sol_new, cdrag)
661	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, temp_sol_new, "Init temp_sol_new")
662	!
663	! Net solar and the wind at the right hight are recomputed with the correct values.
664	!
665	swnet(:) =(1.-(albedo(:,1)+albedo(:,2))/2.)*swdown(:)
666	u_tq(:) = u(:)*LOG(zlev_tq(:)/z0m(:))/LOG(zlev_uv(:)/z0m(:))
667	v_tq(:) = v(:)*LOG(zlev_tq(:)/z0m(:))/LOG(zlev_uv(:)/z0m(:))
668	!
669	lrestart_read = .FALSE.
670	!
671	CALL histwrite_p(hist_id, 'LandPoints', itau+1, (/ REAL(kindex) /), kjpindex, kindex)
672	CALL histwrite_p(hist_id, 'Areas', itau+1, area, kjpindex, kindex)
673	CALL histwrite_p(hist_id, 'Contfrac', itau+1, contfrac, kjpindex, kindex)
674	!
675	IF ( timemeasure ) THEN
676	WRITE(numout,*) '------> CPU Time for set-up of ORCHIDEE : ',Get_cpu_Time(timer_global)
677	WRITE(numout,*) '------> Real Time for set-up of ORCHIDEE : ',Get_real_Time(timer_global)
678	CALL stop_timer(timer_global)
679	CALL start_timer(timer_global)
680	ENDIF
681	!
682	ENDIF
683	!
684	!---------------------------------------------------------------------------------------
685	!-
686	!- Main call to SECHIBA
687	!-
688	!---------------------------------------------------------------------------------------
689	!
690	!
691	CALL sechiba_main (itau_sechiba, iim*jjm, kjpindex, kindex, &
692	& lrestart_read, lrestart_write, &
693	& lalo_loc, contfrac, neighbours, resolution, &
694	! First level conditions
695	& zlev_tq, u_tq, v_tq, qair, qair, temp_air, temp_air, epot_air, ccanopy, &
696	! Variables for the implicit coupling
697	& cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
698	! Rain, snow, radiation and surface pressure
699	& precip_rain ,precip_snow, lwdown, swnet, swdown, sinang, pb, &
700	! Output : Fluxes
701	& vevapp, fluxsens, fluxlat, coastalflow, riverflow, netco2, carblu, &
702	! Surface temperatures and surface properties
703	& tsol_rad, temp_sol_new, qsurf, albedo, emis, z0m, z0h, &
704	! File ids
705	& rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC)
706	!
707	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, temp_sol_new, "Produced temp_sol_new")
708	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, fluxsens, "Produced fluxsens")
709	CALL forcing_printpoint(julian, testpt(1), testpt(2), kjpindex, lalo_loc, fluxlat, "Produced fluxlat")
710	!
711	IF ( timemeasure ) THEN
712	orchidee_cputime = orchidee_cputime + Get_cpu_Time(timer_global)
713	orchidee_walltime = orchidee_walltime + Get_real_Time(timer_global)
714	CALL stop_timer(timer_global)
715	CALL start_timer(timer_global)
716	ENDIF
717	!
718	!---------------------------------------------------------------------------------------
719	!-
720	!- Write diagnostics
721	!-
722	!---------------------------------------------------------------------------------------
723	!
724	CALL xios_orchidee_send_field("LandPoints" ,(/ ( REAL(ik), ik=1,kjpindex ) /))
725	CALL xios_orchidee_send_field("areas", area)
726	CALL xios_orchidee_send_field("contfrac",contfrac)
727	CALL xios_orchidee_send_field("temp_air",temp_air)
728	CALL xios_orchidee_send_field("qair",qair)
729	CALL xios_orchidee_send_field("swnet",swnet)
730	CALL xios_orchidee_send_field("swdown",swdown)
731	! zpb in hPa, output in Pa
732	CALL xios_orchidee_send_field("pb",pb)
733	!
734	IF ( .NOT. almaoutput ) THEN
735	!
736	! ORCHIDEE INPUT variables
737	!
738	CALL histwrite_p (hist_id, 'swdown', itau_sechiba, swdown, kjpindex, kindex)
739	CALL histwrite_p (hist_id, 'tair', itau_sechiba, temp_air, kjpindex, kindex)
740	CALL histwrite_p (hist_id, 'qair', itau_sechiba, qair, kjpindex, kindex)
741	CALL histwrite_p (hist_id, 'evap', itau_sechiba, vevapp, kjpindex, kindex)
742	CALL histwrite_p (hist_id, 'coastalflow',itau_sechiba, coastalflow, kjpindex, kindex)
743	CALL histwrite_p (hist_id, 'riverflow',itau_sechiba, riverflow, kjpindex, kindex)
744	!
745	CALL histwrite_p (hist_id, 'temp_sol', itau_sechiba, temp_sol_new, kjpindex, kindex)
746	CALL histwrite_p (hist_id, 'tsol_max', itau_sechiba, temp_sol_new, kjpindex, kindex)
747	CALL histwrite_p (hist_id, 'tsol_min', itau_sechiba, temp_sol_new, kjpindex, kindex)
748	CALL histwrite_p (hist_id, 'fluxsens', itau_sechiba, fluxsens, kjpindex, kindex)
749	CALL histwrite_p (hist_id, 'fluxlat', itau_sechiba, fluxlat, kjpindex, kindex)
750	CALL histwrite_p (hist_id, 'swnet', itau_sechiba, swnet, kjpindex, kindex)
751	CALL histwrite_p (hist_id, 'alb_vis', itau_sechiba, albedo(:,1), kjpindex, kindex)
752	CALL histwrite_p (hist_id, 'alb_nir', itau_sechiba, albedo(:,2), kjpindex, kindex)
753	!
754	IF ( hist2_id > 0 ) THEN
755	CALL histwrite_p (hist2_id, 'swdown', itau_sechiba, swdown, kjpindex, kindex)
756	CALL histwrite_p (hist2_id, 'tair', itau_sechiba, temp_air, kjpindex, kindex)
757	CALL histwrite_p (hist2_id, 'qair', itau_sechiba, qair, kjpindex, kindex)
758	!
759	CALL histwrite_p (hist2_id, 'evap', itau_sechiba, vevapp, kjpindex, kindex)
760	CALL histwrite_p (hist2_id, 'coastalflow',itau_sechiba, coastalflow, kjpindex, kindex)
761	CALL histwrite_p (hist2_id, 'riverflow',itau_sechiba, riverflow, kjpindex, kindex)
762	!
763	CALL histwrite_p (hist2_id, 'temp_sol', itau_sechiba, temp_sol_new, kjpindex, kindex)
764	CALL histwrite_p (hist2_id, 'tsol_max', itau_sechiba, temp_sol_new, kjpindex, kindex)
765	CALL histwrite_p (hist2_id, 'tsol_min', itau_sechiba, temp_sol_new, kjpindex, kindex)
766	CALL histwrite_p (hist2_id, 'fluxsens', itau_sechiba, fluxsens, kjpindex, kindex)
767	CALL histwrite_p (hist2_id, 'fluxlat', itau_sechiba, fluxlat, kjpindex, kindex)
768	CALL histwrite_p (hist2_id, 'swnet', itau_sechiba, swnet, kjpindex, kindex)
769	!
770	CALL histwrite_p (hist2_id, 'alb_vis', itau_sechiba, albedo(:,1), kjpindex, kindex)
771	CALL histwrite_p (hist2_id, 'alb_nir', itau_sechiba, albedo(:,2), kjpindex, kindex)
772	ENDIF
773	ELSE
774	!
775	! Input variables
776	!
777	CALL histwrite_p (hist_id, 'SinAng', itau_sechiba, sinang, kjpindex, kindex)
778	CALL histwrite_p (hist_id, 'LWdown', itau_sechiba, lwdown, kjpindex, kindex)
779	CALL histwrite_p (hist_id, 'SWdown', itau_sechiba, swdown, kjpindex, kindex)
780	CALL histwrite_p (hist_id, 'Tair', itau_sechiba, temp_air, kjpindex, kindex)
781	CALL histwrite_p (hist_id, 'Qair', itau_sechiba, qair, kjpindex, kindex)
782	CALL histwrite_p (hist_id, 'SurfP', itau_sechiba, pb, kjpindex, kindex)
783	CALL histwrite_p (hist_id, 'Windu', itau_sechiba, u_tq, kjpindex, kindex)
784	CALL histwrite_p (hist_id, 'Windv', itau_sechiba, v_tq, kjpindex, kindex)
785	!
786	CALL histwrite_p (hist_id, 'Evap', itau_sechiba, vevapp, kjpindex, kindex)
787	CALL histwrite_p (hist_id, 'SWnet', itau_sechiba, swnet, kjpindex, kindex)
788	CALL histwrite_p (hist_id, 'Qh', itau_sechiba, fluxsens, kjpindex, kindex)
789	CALL histwrite_p (hist_id, 'Qle', itau_sechiba, fluxlat, kjpindex, kindex)
790	CALL histwrite_p (hist_id, 'AvgSurfT', itau_sechiba, temp_sol_new, kjpindex, kindex)
791	CALL histwrite_p (hist_id, 'RadT', itau_sechiba, temp_sol_new, kjpindex, kindex)
792	!
793	! There is a mess with the units passed to the coupler. To be checked with Marc
794	!
795	IF ( river_routing ) THEN
796	CALL histwrite_p (hist_id, 'CoastalFlow',itau_sechiba, coastalflow, kjpindex, kindex)
797	CALL histwrite_p (hist_id, 'RiverFlow',itau_sechiba, riverflow, kjpindex, kindex)
798	ENDIF
799	!
800	IF ( hist2_id > 0 ) THEN
801	CALL histwrite_p (hist2_id, 'Evap', itau_sechiba, vevapp, kjpindex, kindex)
802	CALL histwrite_p (hist2_id, 'SWnet', itau_sechiba, swnet, kjpindex, kindex)
803	CALL histwrite_p (hist2_id, 'Qh', itau_sechiba, fluxsens, kjpindex, kindex)
804	CALL histwrite_p (hist2_id, 'Qle', itau_sechiba, fluxlat, kjpindex, kindex)
805	CALL histwrite_p (hist2_id, 'AvgSurfT', itau_sechiba, temp_sol_new, kjpindex, kindex)
806	CALL histwrite_p (hist2_id, 'RadT', itau_sechiba, temp_sol_new, kjpindex, kindex)
807	ENDIF
808	ENDIF
809	!
810	!
811	ENDDO
812	!-
813	!-
814	!---------------------------------------------------------------------------------------
815	!-
816	!- Close eveything
817	!-
818	!--
819	!
820	CALL histclo
821	IF(is_root_prc) THEN
822	CALL restclo
823	CALL getin_dump
824	ENDIF
825
826	!-
827	!- Deallocate all variables and reset initialization variables
828	!-
829	CALL orchideedriver_clear()
830
831	WRITE(numout,*) "End at proc ", mpi_rank
832	!
833	!
834	!---------------------------------------------------------------------------------------
835	!-
836	!- Get time and close IOIPSL, OASIS and MPI
837	!-
838	!---------------------------------------------------------------------------------------
839	!-
840	IF ( timemeasure ) THEN
841	WRITE(numout,*) '------> Total CPU Time waiting to get forcing : ',waitget_cputime
842	WRITE(numout,*) '------> Total Real Time waiting to get forcing : ',waitget_walltime
843	WRITE(numout,*) '------> Total CPU Time for ORCHIDEE : ', orchidee_cputime
844	WRITE(numout,*) '------> Total Real Time for ORCHIDEE : ', orchidee_walltime
845	WRITE(numout,*) '------> Total CPU Time waiting to put fluxes : ',waitput_cputime
846	WRITE(numout,*) '------> Total Real Time waiting to put fluxes : ',waitput_walltime
847	WRITE(numout,*) '------> Total CPU Time for closing : ', Get_cpu_Time(timer_global)
848	WRITE(numout,*) '------> Total Real Time for closing : ', Get_real_Time(timer_global)
849	WRITE(numout,*) '------> Total without MPI : CPU Time : ', Get_cpu_Time(timer_mpi)
850	WRITE(numout,*) '------> Total without MPI : Real Time : ', Get_real_Time(timer_mpi)
851	CALL stop_timer(timer_global)
852	CALL stop_timer(timer_mpi)
853	ENDIF
854	!
855	CALL Finalize_mpi
856	CONTAINS
857
858	!! ================================================================================================================================
859	!! SUBROUTINE : orchideedriver_clear
860	!!
861	!>\BRIEF Clear orchideedriver
862	!!
863	!! DESCRIPTION : Deallocate memory and reset initialization variables to there original values
864	!! This subroutine calls forcing_tools_clear and sechiba_clear.
865	!!
866	!_ ================================================================================================================================
867	SUBROUTINE orchideedriver_clear
868
869	!- Deallocate all variables existing on all procs
870
871	IF ( ALLOCATED(lon_glo) ) DEALLOCATE(lon_glo)
872	IF ( ALLOCATED(lat_glo) ) DEALLOCATE(lat_glo)
873	IF ( ALLOCATED(mask_glo) ) DEALLOCATE(mask_glo)
874	IF ( ALLOCATED(area_glo) ) DEALLOCATE(area_glo)
875	IF ( ALLOCATED(corners_glo) ) DEALLOCATE(corners_glo)
876	IF ( ALLOCATED(kindex_g) ) DEALLOCATE(kindex_g)
877	IF ( ALLOCATED(contfrac_glo) ) DEALLOCATE(contfrac_glo)
878	IF ( ALLOCATED(lalo_glo) ) DEALLOCATE(lalo_glo)
879	IF ( ALLOCATED(lon) ) DEALLOCATE(lon)
880	IF ( ALLOCATED(lat) ) DEALLOCATE(lat)
881	IF ( ALLOCATED(kindex) ) DEALLOCATE(kindex)
882	IF ( ALLOCATED(lalo_loc) ) DEALLOCATE(lalo_loc)
883	IF ( ALLOCATED(zlev_tq) ) DEALLOCATE(zlev_tq)
884	IF ( ALLOCATED(zlev_uv) ) DEALLOCATE(zlev_uv)
885	IF ( ALLOCATED(u) ) DEALLOCATE(u)
886	IF ( ALLOCATED(v) ) DEALLOCATE(v)
887	IF ( ALLOCATED(pb) ) DEALLOCATE(pb)
888	IF ( ALLOCATED(temp_air) ) DEALLOCATE(temp_air)
889	IF ( ALLOCATED(qair) ) DEALLOCATE(qair)
890	IF ( ALLOCATED(precip_rain) ) DEALLOCATE(precip_rain)
891	IF ( ALLOCATED(precip_snow) ) DEALLOCATE(precip_snow)
892	IF ( ALLOCATED(swdown) ) DEALLOCATE(swdown)
893	IF ( ALLOCATED(swnet) ) DEALLOCATE(swnet)
894	IF ( ALLOCATED(lwdown) ) DEALLOCATE(lwdown)
895	IF ( ALLOCATED(sinang) ) DEALLOCATE(sinang)
896	IF ( ALLOCATED(epot_air) ) DEALLOCATE(epot_air)
897	IF ( ALLOCATED(ccanopy) ) DEALLOCATE(ccanopy)
898	IF ( ALLOCATED(cdrag) ) DEALLOCATE(cdrag)
899	IF ( ALLOCATED(swnet) ) DEALLOCATE(swnet)
900	IF ( ALLOCATED(petAcoef) ) DEALLOCATE(petAcoef)
901	IF ( ALLOCATED(peqAcoef) ) DEALLOCATE(peqAcoef)
902	IF ( ALLOCATED(petBcoef) ) DEALLOCATE(petBcoef)
903	IF ( ALLOCATED(peqBcoef) ) DEALLOCATE(peqBcoef)
904	IF ( ALLOCATED(u_tq) ) DEALLOCATE(u_tq)
905	IF ( ALLOCATED(v_tq) ) DEALLOCATE(v_tq)
906	IF ( ALLOCATED(vevapp) ) DEALLOCATE(vevapp)
907	IF ( ALLOCATED(fluxsens) ) DEALLOCATE(fluxsens)
908	IF ( ALLOCATED(fluxlat) ) DEALLOCATE(fluxlat)
909	IF ( ALLOCATED(coastalflow) ) DEALLOCATE(coastalflow)
910	IF ( ALLOCATED(riverflow) ) DEALLOCATE(riverflow)
911	IF ( ALLOCATED(netco2) ) DEALLOCATE(netco2)
912	IF ( ALLOCATED(carblu) ) DEALLOCATE(carblu)
913	IF ( ALLOCATED(tsol_rad) ) DEALLOCATE(tsol_rad)
914	IF ( ALLOCATED(temp_sol_new) ) DEALLOCATE(temp_sol_new)
915	IF ( ALLOCATED(qsurf) ) DEALLOCATE(qsurf)
916	IF ( ALLOCATED(albedo) ) DEALLOCATE(albedo)
917	IF ( ALLOCATED(emis) ) DEALLOCATE(emis)
918	IF ( ALLOCATED(z0m) ) DEALLOCATE(z0m)
919
920	! Deallocate and reset variables in forcing_tools module
921	CALL forcing_tools_clear()
922
923	! Deallocate and reset variables in sechiba module and underlaying modules
924	CALL sechiba_clear()
925
926	WRITE(numout,*) "Memory deallocated in orchideedriver_clear"
927
928	END SUBROUTINE orchideedriver_clear
929
930	END PROGRAM orchidedriver

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

Download in other formats: