Context Navigation

source: branches/publications/ORCHIDEE_gmd-2018-57/src_driver/orchideedriver.f90 @ 5143

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

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

Download in other formats: