Context Navigation

limit_read_mod.F90 @ 5006

Last change on this file since 5006 was 3406, checked in by jghattas, 6 years ago
Added all modifications in the model code that were used for the simulations with DYANMICO during the Grand Challeng 2018. Modifications done by Y. Meurdesoif, L. Fairhead and A.K. Traore
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 17.8 KB

Rev	Line
[996]	1	!
[3312]	2	! $Id: limit_read_mod.F90 3406 2018-10-23 09:16:41Z abarral $
[996]	3	!
	4	MODULE limit_read_mod
	5	!
	6	! This module reads the fichier "limit.nc" containing fields for surface forcing.
	7	!
	8	! Module subroutines :
	9	! limit_read_frac : call limit_read_tot and return the fractions
	10	! limit_read_rug_alb : return rugosity and albedo, if coupled ocean call limit_read_tot first
	11	! limit_read_sst : return sea ice temperature
	12	! limit_read_tot : read limit.nc and store the fields in local modules variables
	13	!
	14	IMPLICIT NONE
	15
	16	REAL, ALLOCATABLE, DIMENSION(:,:), SAVE, PRIVATE :: pctsrf
[1001]	17	!$OMP THREADPRIVATE(pctsrf)
[996]	18	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: rugos
[1001]	19	!$OMP THREADPRIVATE(rugos)
[996]	20	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: albedo
[1001]	21	!$OMP THREADPRIVATE(albedo)
[996]	22	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: sst
[1001]	23	!$OMP THREADPRIVATE(sst)
	24	LOGICAL,SAVE :: read_continents=.FALSE.
	25	!$OMP THREADPRIVATE(read_continents)
[996]	26
	27	CONTAINS
	28	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	29	!!
	30	!! Public subroutines :
	31	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	32
[3312]	33
	34	SUBROUTINE init_limit_read(first_day)
	35	USE mod_grid_phy_lmdz
	36	USE surface_data
	37	USE mod_phys_lmdz_para
	38	#ifdef CPP_XIOS
	39	USE XIOS
	40	#endif
	41	IMPLICIT NONE
	42	INTEGER, INTENT(IN) :: first_day
	43
	44
	45	IF ( type_ocean /= 'couple') THEN
	46	IF (grid_type==unstructured) THEN
	47	#ifdef CPP_XIOS
	48	IF (is_omp_master) CALL xios_set_file_attr("limit_read",enabled=.TRUE.,record_offset=first_day)
	49	#endif
	50	ENDIF
	51	ENDIF
	52
	53	END SUBROUTINE init_limit_read
	54
[996]	55	SUBROUTINE limit_read_frac(itime, dtime, jour, pctsrf_new, is_modified)
	56	!
	57	! This subroutine is called from "change_srf_frac" for case of
	58	! ocean=force or from ocean_slab_frac for ocean=slab.
	59	! The fraction for all sub-surfaces at actual time step is returned.
	60
	61	USE dimphy
[1785]	62	USE indice_sol_mod
[996]	63
	64	! Input arguments
	65	!****************************************************************************************
	66	INTEGER, INTENT(IN) :: itime ! time step
	67	INTEGER, INTENT(IN) :: jour ! current day
	68	REAL , INTENT(IN) :: dtime ! length of time step
	69
	70	! Output arguments
	71	!****************************************************************************************
	72	REAL, DIMENSION(klon,nbsrf), INTENT(OUT) :: pctsrf_new ! sub surface fractions
	73	LOGICAL, INTENT(OUT) :: is_modified ! true if pctsrf is modified at this time step
	74
	75	! End declaration
	76	!****************************************************************************************
	77
	78	! 1) Read file limit.nc
	79	CALL limit_read_tot(itime, dtime, jour, is_modified)
	80
	81	! 2) Return the fraction read in limit_read_tot
	82	pctsrf_new(:,:) = pctsrf(:,:)
	83
	84	END SUBROUTINE limit_read_frac
	85
	86	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	87
	88	SUBROUTINE limit_read_rug_alb(itime, dtime, jour, &
	89	knon, knindex, &
	90	rugos_out, alb_out)
	91	!
	92	! This subroutine is called from surf_land_bucket.
	93	! The flag "ok_veget" must can not be true. If coupled run, "ocean=couple"
	94	! then this routine will call limit_read_tot.
	95	!
	96	USE dimphy
	97	USE surface_data
	98
	99	! Input arguments
	100	!****************************************************************************************
	101	INTEGER, INTENT(IN) :: itime ! numero du pas de temps courant
	102	INTEGER, INTENT(IN) :: jour ! jour a lire dans l'annee
	103	REAL , INTENT(IN) :: dtime ! pas de temps de la physique (en s)
	104	INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid
	105	INTEGER, DIMENSION(klon), INTENT(IN) :: knindex ! grid point number for compressed grid
	106	! Output arguments
	107	!****************************************************************************************
	108	REAL, DIMENSION(klon), INTENT(OUT) :: rugos_out
	109	REAL, DIMENSION(klon), INTENT(OUT) :: alb_out
	110
	111	! Local variables
	112	!****************************************************************************************
	113	INTEGER :: i
	114	LOGICAL :: is_modified
	115	!****************************************************************************************
	116
[2209]	117	IF (type_ocean == 'couple'.OR. &
	118	(type_ocean == 'slab' .AND. version_ocean == 'sicINT')) THEN
[996]	119	! limit.nc has not yet been read. Do it now!
	120	CALL limit_read_tot(itime, dtime, jour, is_modified)
	121	END IF
	122
	123	DO i=1,knon
	124	rugos_out(i) = rugos(knindex(i))
	125	alb_out(i) = albedo(knindex(i))
	126	END DO
	127
	128	END SUBROUTINE limit_read_rug_alb
	129
	130	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	131
	132	SUBROUTINE limit_read_sst(knon, knindex, sst_out)
	133	!
	134	! This subroutine returns the sea surface temperature already read from limit.nc.
	135	!
	136	USE dimphy, ONLY : klon
	137
	138	INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid
	139	INTEGER, DIMENSION(klon), INTENT(IN) :: knindex ! grid point number for compressed grid
	140	REAL, DIMENSION(klon), INTENT(OUT) :: sst_out
	141
	142	INTEGER :: i
	143
	144	DO i = 1, knon
	145	sst_out(i) = sst(knindex(i))
	146	END DO
	147
	148	END SUBROUTINE limit_read_sst
	149
	150	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	151	!!
	152	!! Private subroutine :
	153	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	154
	155	SUBROUTINE limit_read_tot(itime, dtime, jour, is_modified)
	156	!
	157	! Read everything needed from limit.nc
	158	!
	159	! 0) Initialize
	160	! 1) Open the file limit.nc, if it is time
	161	! 2) Read fraction, if not type_ocean=couple
	162	! 3) Read sea surface temperature, if not type_ocean=couple
	163	! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
	164	! 5) Close file and distribuate variables to all processus
	165
	166	USE dimphy
	167	USE mod_grid_phy_lmdz
	168	USE mod_phys_lmdz_para
	169	USE surface_data, ONLY : type_ocean, ok_veget
	170	USE netcdf
[1785]	171	USE indice_sol_mod
[2665]	172	USE phys_cal_mod, ONLY : calend, year_len
	173	USE print_control_mod, ONLY: lunout, prt_level
[3312]	174	#ifdef CPP_XIOS
	175	USE XIOS, ONLY: xios_recv_field
	176	#endif
	177
[996]	178	IMPLICIT NONE
	179
	180	! In- and ouput arguments
	181	!****************************************************************************************
	182	INTEGER, INTENT(IN) :: itime ! numero du pas de temps courant
	183	INTEGER, INTENT(IN) :: jour ! jour a lire dans l'annee
	184	REAL , INTENT(IN) :: dtime ! pas de temps de la physique (en s)
	185
	186	LOGICAL, INTENT(OUT) :: is_modified ! true if pctsrf is modified at this time step
	187
	188	! Locals variables with attribute SAVE
	189	!****************************************************************************************
	190	! frequence de lecture des conditions limites (en pas de physique)
[1001]	191	INTEGER,SAVE :: lmt_pas
	192	!$OMP THREADPRIVATE(lmt_pas)
[996]	193	LOGICAL, SAVE :: first_call=.TRUE.
[1577]	194	!$OMP THREADPRIVATE(first_call)
	195	INTEGER, SAVE :: jour_lu = -1
	196	!$OMP THREADPRIVATE(jour_lu)
[996]	197	! Locals variables
	198	!****************************************************************************************
[2665]	199	INTEGER :: nid, nvarid, ndimid, nn
[996]	200	INTEGER :: ii, ierr
	201	INTEGER, DIMENSION(2) :: start, epais
	202	REAL, DIMENSION(klon_glo,nbsrf) :: pct_glo ! fraction at global grid
	203	REAL, DIMENSION(klon_glo) :: sst_glo ! sea-surface temperature at global grid
	204	REAL, DIMENSION(klon_glo) :: rug_glo ! rugosity at global grid
	205	REAL, DIMENSION(klon_glo) :: alb_glo ! albedo at global grid
[3312]	206
	207	REAL, DIMENSION(klon_mpi,nbsrf) :: pct_mpi ! fraction at global grid
	208	REAL, DIMENSION(klon_mpi) :: sst_mpi ! sea-surface temperature at global grid
	209	REAL, DIMENSION(klon_mpi) :: rug_mpi ! rugosity at global grid
	210	REAL, DIMENSION(klon_mpi) :: alb_mpi ! albedo at global grid
	211
[996]	212	CHARACTER(len=20) :: modname='limit_read_mod'
[2665]	213	CHARACTER(LEN=99) :: abort_message, calendar, str
[996]	214
	215	! End declaration
	216	!****************************************************************************************
	217
	218	!****************************************************************************************
	219	! 0) Initialization
	220	!
	221	!****************************************************************************************
	222	IF (first_call) THEN
[2770]	223	first_call=.FALSE.
[996]	224	! calculate number of time steps for one day
	225	lmt_pas = NINT(86400./dtime * 1.0)
	226
	227	! Allocate module save variables
	228	IF ( type_ocean /= 'couple' ) THEN
	229	ALLOCATE(pctsrf(klon,nbsrf), sst(klon), stat=ierr)
[2311]	230	IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating pctsrf and sst',1)
[996]	231	END IF
	232
	233	IF ( .NOT. ok_veget ) THEN
	234	ALLOCATE(rugos(klon), albedo(klon), stat=ierr)
[2311]	235	IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating rugos and albedo',1)
[996]	236	END IF
	237
	238	!$OMP MASTER ! Only master thread
	239	IF (is_mpi_root) THEN ! Only master processus
	240	ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
[2311]	241	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
[996]	242	'Pb d''ouverture du fichier de conditions aux limites',1)
[2665]	243
	244	!--- WARNING IF CALENDAR IS KNOWN AND DOES NOT MATCH THE ONE OF LMDZ
	245	ierr=NF90_INQ_VARID(nid, 'TEMPS', nvarid)
	246	ierr=NF90_GET_ATT(nid, nvarid, 'calendar', calendar)
	247	IF(ierr==NF90_NOERR.AND.calendar/=calend.AND.prt_level>=1) THEN
	248	WRITE(lunout,*)'BEWARE: gcm and limit.nc calendars differ: '
	249	WRITE(lunout,*)' '//TRIM(calend)//' for gcm'
	250	WRITE(lunout,*)' '//TRIM(calendar)//' for limit.nc file'
	251	END IF
	252
[3312]	253	!--- ERROR IF FILE RECORDS NUMBER IS NOT EQUAL TO EXPECTED NUMBER OF DAYS
	254	IF (grid_type==unstructured) THEN
	255	ierr=NF90_INQ_DIMID(nid,"time_year",ndimid)
	256	ELSE
	257	ierr=NF90_INQUIRE(nid, UnlimitedDimID=ndimid)
	258	ENDIF
[2665]	259	ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
[3406]	260	WRITE(abort_message,'(a,2(i0,a))')'limit.nc records number (',nn,') does no'//&
[2768]	261	't match year length (',year_len,')'
[2665]	262	IF(nn/=year_len) CALL abort_physic(modname,abort_message,1)
	263
	264	!--- ERROR IF FILES AND LMDZ HORIZONTAL RESOLUTIONS DO NOT MATCH
[3312]	265	IF (grid_type==unstructured) THEN
	266	ierr=NF90_INQ_DIMID(nid, 'cell', ndimid)
	267	ELSE
	268	ierr=NF90_INQ_DIMID(nid, 'points_physiques', ndimid)
	269	ENDIF
[2665]	270	ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
[2768]	271	WRITE(abort_message,'(a,2(i0,a))')'limit.nc horizontal number of cells (',nn, &
	272	') does not match LMDZ klon_glo (',klon_glo,')'
[2665]	273	IF(nn/=klon_glo) CALL abort_physic(modname,abort_message,1)
	274
[2768]	275	ierr = NF90_CLOSE(nid)
	276	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
	277	END IF ! is_mpi_root
	278	!$OMP END MASTER
	279	!$OMP BARRIER
	280	END IF
	281
	282	!****************************************************************************************
	283	! 1) Open the file limit.nc if it is the right moment to read, once a day.
	284	! The file is read only by the master thread of the master mpi process(is_mpi_root)
	285	! Check by the way if the number of records is correct.
	286	!
	287	!****************************************************************************************
	288
	289	is_modified = .FALSE.
[3406]	290	!ym IF (MOD(itime-1, lmt_pas) == 0 .OR. jour_lu /= jour ) THEN ! time to read
	291	! not REALLY PERIODIC
	292	IF (MOD(itime-1, lmt_pas) == 0) THEN ! time to read
[2768]	293	jour_lu = jour
	294	is_modified = .TRUE.
[3312]	295
	296	IF (grid_type==unstructured) THEN
	297
	298	#ifdef CPP_XIOS
	299	IF ( type_ocean /= 'couple') THEN
	300
	301	IF (is_omp_master) CALL xios_recv_field("foce_limin",pct_mpi(:,is_oce))
	302	IF (is_omp_master) CALL xios_recv_field("fsic_limin",pct_mpi(:,is_sic))
	303	! IF (read_continents .OR. itime == 1) THEN
	304	IF (is_omp_master) CALL xios_recv_field("fter_limin",pct_mpi(:,is_ter))
	305	IF (is_omp_master) CALL xios_recv_field("flic_limin",pct_mpi(:,is_lic))
	306	! ENDIF
	307	ENDIF! type_ocean /= couple
	308
	309	IF ( type_ocean /= 'couple') THEN
	310	IF (is_omp_master) CALL xios_recv_field("sst_limin",sst_mpi)
	311	ENDIF
	312
	313	IF (.NOT. ok_veget) THEN
	314	IF (is_omp_master) CALL xios_recv_field("alb_limin",alb_mpi)
	315	IF (is_omp_master) CALL xios_recv_field("rug_limin",rug_mpi)
	316	ENDIF
	317
	318	IF ( type_ocean /= 'couple') THEN
	319	CALL Scatter_omp(sst_mpi,sst)
	320	CALL Scatter_omp(pct_mpi(:,is_oce),pctsrf(:,is_oce))
	321	CALL Scatter_omp(pct_mpi(:,is_sic),pctsrf(:,is_sic))
	322	! IF (read_continents .OR. itime == 1) THEN
	323	CALL Scatter_omp(pct_mpi(:,is_ter),pctsrf(:,is_ter))
	324	CALL Scatter_omp(pct_mpi(:,is_lic),pctsrf(:,is_lic))
	325	! END IF
	326	END IF
	327
	328	IF (.NOT. ok_veget) THEN
	329	CALL Scatter_omp(alb_mpi, albedo)
	330	CALL Scatter_omp(rug_mpi, rugos)
	331	END IF
	332	#endif
	333
	334
	335	ELSE ! grid_type==regular
	336
[2768]	337	!$OMP MASTER ! Only master thread
[3406]	338	IF (is_mpi_root) THEN ! Only master processus!
[2768]	339
	340	ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
	341	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
	342	'Pb d''ouverture du fichier de conditions aux limites',1)
	343
[996]	344	! La tranche de donnees a lire:
	345	start(1) = 1
	346	start(2) = jour
	347	epais(1) = klon_glo
	348	epais(2) = 1
	349
	350
	351	!****************************************************************************************
	352	! 2) Read fraction if not type_ocean=couple
	353	!
	354	!****************************************************************************************
	355
	356	IF ( type_ocean /= 'couple') THEN
	357	!
	358	! Ocean fraction
	359	ierr = NF90_INQ_VARID(nid, 'FOCE', nvarid)
[2311]	360	IF (ierr /= NF90_NOERR) CALL abort_physic(modname, 'Le champ <FOCE> est absent',1)
[996]	361
	362	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_oce),start,epais)
[2311]	363	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FOCE>' ,1)
[996]	364	!
	365	! Sea-ice fraction
	366	ierr = NF90_INQ_VARID(nid, 'FSIC', nvarid)
[2311]	367	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FSIC> est absent',1)
[996]	368
	369	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_sic),start,epais)
[2311]	370	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FSIC>' ,1)
[996]	371
	372
	373	! Read land and continentals fraction only if asked for
	374	IF (read_continents .OR. itime == 1) THEN
	375	!
	376	! Land fraction
	377	ierr = NF90_INQ_VARID(nid, 'FTER', nvarid)
[2311]	378	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FTER> est absent',1)
[996]	379
	380	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_ter),start,epais)
[2311]	381	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FTER>',1)
[996]	382	!
	383	! Continentale ice fraction
	384	ierr = NF90_INQ_VARID(nid, 'FLIC', nvarid)
[2311]	385	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FLIC> est absent',1)
[996]	386
	387	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_lic),start,epais)
[2311]	388	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FLIC>',1)
[996]	389	END IF
	390
	391	END IF ! type_ocean /= couple
	392
	393	!****************************************************************************************
	394	! 3) Read sea-surface temperature, if not coupled ocean
	395	!
	396	!****************************************************************************************
	397	IF ( type_ocean /= 'couple') THEN
	398
	399	ierr = NF90_INQ_VARID(nid, 'SST', nvarid)
[2311]	400	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <SST> est absent',1)
[996]	401
	402	ierr = NF90_GET_VAR(nid,nvarid,sst_glo,start,epais)
[2311]	403	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <SST>',1)
[996]	404
	405	END IF
	406
	407	!****************************************************************************************
	408	! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
	409	!
	410	!****************************************************************************************
	411
	412	IF (.NOT. ok_veget) THEN
	413	!
	414	! Read albedo
	415	ierr = NF90_INQ_VARID(nid, 'ALB', nvarid)
[2311]	416	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <ALB> est absent',1)
[996]	417
	418	ierr = NF90_GET_VAR(nid,nvarid,alb_glo,start,epais)
[2311]	419	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <ALB>',1)
[996]	420	!
	421	! Read rugosity
	422	ierr = NF90_INQ_VARID(nid, 'RUG', nvarid)
[2311]	423	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <RUG> est absent',1)
[996]	424
	425	ierr = NF90_GET_VAR(nid,nvarid,rug_glo,start,epais)
[2311]	426	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <RUG>',1)
[996]	427
	428	END IF
	429
	430	!****************************************************************************************
	431	! 5) Close file and distribuate variables to all processus
	432	!
	433	!****************************************************************************************
	434	ierr = NF90_CLOSE(nid)
[2311]	435	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
[996]	436	ENDIF ! is_mpi_root
	437
	438	!$OMP END MASTER
	439	!$OMP BARRIER
	440
	441	IF ( type_ocean /= 'couple') THEN
	442	CALL Scatter(sst_glo,sst)
	443	CALL Scatter(pct_glo(:,is_oce),pctsrf(:,is_oce))
	444	CALL Scatter(pct_glo(:,is_sic),pctsrf(:,is_sic))
	445	IF (read_continents .OR. itime == 1) THEN
	446	CALL Scatter(pct_glo(:,is_ter),pctsrf(:,is_ter))
	447	CALL Scatter(pct_glo(:,is_lic),pctsrf(:,is_lic))
	448	END IF
	449	END IF
	450
	451	IF (.NOT. ok_veget) THEN
	452	CALL Scatter(alb_glo, albedo)
	453	CALL Scatter(rug_glo, rugos)
	454	END IF
	455
[3312]	456	ENDIF ! Grid type
	457
[996]	458	ENDIF ! time to read
	459
	460	END SUBROUTINE limit_read_tot
	461
	462	END MODULE limit_read_mod

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

Download in other formats:

Original Format