New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

domain.F90 in trunk/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

Context Navigation

source: trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90 @ 6981

Last change on this file since 6981 was 6981, checked in by nicolasmartin, 8 years ago
Bugfixes for compilation issues : last commit for trcstp and gfortran doesn't allow starting lines with semicolon
Property svn:keywords set to `Id`
File size: 25.1 KB

Rev	Line
[3]	1	MODULE domain
	2	!!==============================================================================
	3	!! * MODULE domain *
	4	!! Ocean initialization : domain initialization
	5	!!==============================================================================
[1438]	6	!! History : OPA ! 1990-10 (C. Levy - G. Madec) Original code
	7	!! ! 1992-01 (M. Imbard) insert time step initialization
	8	!! ! 1996-06 (G. Madec) generalized vertical coordinate
	9	!! ! 1997-02 (G. Madec) creation of domwri.F
	10	!! ! 2001-05 (E.Durand - G. Madec) insert closed sea
	11	!! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module
	12	!! 2.0 ! 2005-11 (V. Garnier) Surface pressure gradient organization
[2528]	13	!! 3.3 ! 2010-11 (G. Madec) initialisation in C1D configuration
[4152]	14	!! 3.6 ! 2013 ( J. Simeon, C. Calone, G. Madec, C. Ethe ) Online coarsening of outputs
[6140]	15	!! 3.7 ! 2015-11 (G. Madec, A. Coward) time varying zgr by default
[3]	16	!!----------------------------------------------------------------------
[1438]	17
	18	!!----------------------------------------------------------------------
[3]	19	!! dom_init : initialize the space and time domain
	20	!! dom_nam : read and contral domain namelists
	21	!! dom_ctl : control print for the ocean domain
[5836]	22	!! dom_stiff : diagnose maximum grid stiffness/hydrostatic consistency (s-coordinate)
[3]	23	!!----------------------------------------------------------------------
[2528]	24	USE oce ! ocean variables
	25	USE dom_oce ! domain: ocean
[888]	26	USE sbc_oce ! surface boundary condition: ocean
[719]	27	USE phycst ! physical constants
[1601]	28	USE closea ! closed seas
[3]	29	USE domhgr ! domain: set the horizontal mesh
	30	USE domzgr ! domain: set the vertical mesh
	31	USE domstp ! domain: set the time-step
	32	USE dommsk ! domain: set the mask system
	33	USE domwri ! domain: write the meshmask file
[592]	34	USE domvvl ! variable volume
[2528]	35	USE c1d ! 1D vertical configuration
	36	USE dyncor_c1d ! Coriolis term (c1d case) (cor_c1d routine)
[5836]	37	!
	38	USE in_out_manager ! I/O manager
[6140]	39	USE wrk_nemo ! Memory Allocation
[5836]	40	USE lib_mpp ! distributed memory computing library
	41	USE lbclnk ! ocean lateral boundary condition (or mpp link)
[3294]	42	USE timing ! Timing
[3]	43
	44	IMPLICIT NONE
	45	PRIVATE
	46
[1438]	47	PUBLIC dom_init ! called by opa.F90
[3]	48
[1438]	49	!!-------------------------------------------------------------------------
[2528]	50	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[888]	51	!! $Id$
[2528]	52	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1438]	53	!!-------------------------------------------------------------------------
[3]	54	CONTAINS
	55
	56	SUBROUTINE dom_init
	57	!!----------------------------------------------------------------------
	58	!! * ROUTINE dom_init *
	59	!!
	60	!! ** Purpose : Domain initialization. Call the routines that are
[1601]	61	!! required to create the arrays which define the space
	62	!! and time domain of the ocean model.
[3]	63	!!
[1601]	64	!! ** Method : - dom_msk: compute the masks from the bathymetry file
	65	!! - dom_hgr: compute or read the horizontal grid-point position
	66	!! and scale factors, and the coriolis factor
	67	!! - dom_zgr: define the vertical coordinate and the bathymetry
	68	!! - dom_stp: defined the model time step
	69	!! - dom_wri: create the meshmask file if nmsh=1
[2528]	70	!! - 1D configuration, move Coriolis, u and v at T-point
[3]	71	!!----------------------------------------------------------------------
[6140]	72	INTEGER :: jk ! dummy loop indices
[3764]	73	INTEGER :: iconf = 0 ! local integers
[6140]	74	REAL(wp), POINTER, DIMENSION(:,:) :: z1_hu_0, z1_hv_0
[3]	75	!!----------------------------------------------------------------------
[1601]	76	!
[3764]	77	IF( nn_timing == 1 ) CALL timing_start('dom_init')
[3294]	78	!
[3]	79	IF(lwp) THEN
	80	WRITE(numout,*)
	81	WRITE(numout,*) 'dom_init : domain initialization'
	82	WRITE(numout,*) '~~~~~~~~'
	83	ENDIF
[1601]	84	!
[6140]	85	! !== Reference coordinate system ==!
[4490]	86	!
[6140]	87	CALL dom_nam ! read namelist ( namrun, namdom )
	88	CALL dom_clo ! Closed seas and lake
	89	CALL dom_hgr ! Horizontal mesh
	90	CALL dom_zgr ! Vertical mesh and bathymetry
	91	CALL dom_msk ! Masks
	92	IF( ln_sco ) CALL dom_stiff ! Maximum stiffness ratio/hydrostatic consistency
	93	!
	94	ht_0(:,:) = e3t_0(:,:,1) * tmask(:,:,1) ! Reference ocean thickness
	95	hu_0(:,:) = e3u_0(:,:,1) * umask(:,:,1)
	96	hv_0(:,:) = e3v_0(:,:,1) * vmask(:,:,1)
	97	DO jk = 2, jpk
[4490]	98	ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
	99	hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
	100	hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
	101	END DO
	102	!
[6140]	103	! !== time varying part of coordinate system ==!
[1601]	104	!
[6140]	105	IF( ln_linssh ) THEN ! Fix in time : set to the reference one for all
	106	! before ! now ! after !
[6981]	107	gdept_b = gdept_0 ; gdept_n = gdept_0 ! --- ! depth of grid-points
	108	gdepw_b = gdepw_0 ; gdepw_n = gdepw_0 ! --- !
	109	gde3w_n = gde3w_0 ! --- !
[6140]	110	!
[6981]	111	e3t_b = e3t_0 ; e3t_n = e3t_0 ; e3t_a = e3t_0 ! scale factors
	112	e3u_b = e3u_0 ; e3u_n = e3u_0 ; e3u_a = e3u_0 !
	113	e3v_b = e3v_0 ; e3v_n = e3v_0 ; e3v_a = e3v_0 !
	114	e3f_n = e3f_0 ! --- !
	115	e3w_b = e3w_0 ; e3w_n = e3w_0 ! --- !
	116	e3uw_b = e3uw_0 ; e3uw_n = e3uw_0 ! --- !
	117	e3vw_b = e3vw_0 ; e3vw_n = e3vw_0 ! --- !
[6140]	118	!
	119	CALL wrk_alloc( jpi,jpj, z1_hu_0, z1_hv_0 )
	120	!
	121	z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) ) ! _i mask due to ISF
	122	z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) )
	123	!
	124	! before ! now ! after !
[6981]	125	ht_n = ht_0 ! ! water column thickness
	126	hu_b = hu_0 ; hu_n = hu_0 ; hu_a = hu_0 !
	127	hv_b = hv_0 ; hv_n = hv_0 ; hv_a = hv_0 !
	128	r1_hu_b = z1_hu_0 ; r1_hu_n = z1_hu_0 ; r1_hu_a = z1_hu_0 ! inverse of water column thickness
	129	r1_hv_b = z1_hv_0 ; r1_hv_n = z1_hv_0 ; r1_hv_a = z1_hv_0 !
[6140]	130	!
	131	CALL wrk_dealloc( jpi,jpj, z1_hu_0, z1_hv_0 )
	132	!
	133	ELSE ! time varying : initialize before/now/after variables
	134	!
	135	CALL dom_vvl_init
	136	!
	137	ENDIF
[2528]	138	!
[6140]	139	IF( lk_c1d ) CALL cor_c1d ! 1D configuration: Coriolis set at T-point
[4370]	140	!
[6140]	141	CALL dom_stp ! time step
	142	IF( nmsh /= 0 .AND. .NOT. ln_iscpl ) CALL dom_wri ! Create a domain file
	143	IF( nmsh /= 0 .AND. ln_iscpl .AND. .NOT. ln_rstart ) CALL dom_wri ! Create a domain file
	144	IF( .NOT.ln_rstart ) CALL dom_ctl ! Domain control
[1438]	145	!
[3764]	146	IF( nn_timing == 1 ) CALL timing_stop('dom_init')
[3294]	147	!
[3]	148	END SUBROUTINE dom_init
	149
	150
	151	SUBROUTINE dom_nam
	152	!!----------------------------------------------------------------------
	153	!! * ROUTINE dom_nam *
	154	!!
	155	!! ** Purpose : read domaine namelists and print the variables.
	156	!!
	157	!! ** input : - namrun namelist
	158	!! - namdom namelist
[2528]	159	!! - namnc4 namelist ! "key_netcdf4" only
[3]	160	!!----------------------------------------------------------------------
	161	USE ioipsl
[6140]	162	NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list, &
	163	nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl , &
	164	& nn_it000, nn_itend , nn_date0 , nn_time0 , nn_leapy , nn_istate , &
	165	& nn_stock, nn_write , ln_mskland , ln_clobber , nn_chunksz, nn_euler , &
	166	& ln_cfmeta, ln_iscpl
	167	NAMELIST/namdom/ nn_bathy, rn_bathy , rn_e3zps_min, rn_e3zps_rat, nn_msh, rn_hmin, rn_isfhmin, &
	168	& rn_atfp , rn_rdt , nn_closea , ln_crs , jphgr_msh , &
	169	& ppglam0, ppgphi0, ppe1_deg, ppe2_deg, ppe1_m, ppe2_m, &
	170	& ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, &
[4147]	171	& ppa2, ppkth2, ppacr2
[2528]	172	#if defined key_netcdf4
	173	NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
	174	#endif
[4147]	175	INTEGER :: ios ! Local integer output status for namelist read
[3]	176	!!----------------------------------------------------------------------
	177
[4147]	178	REWIND( numnam_ref ) ! Namelist namrun in reference namelist : Parameters of the run
	179	READ ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
[5836]	180	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
[4147]	181
	182	REWIND( numnam_cfg ) ! Namelist namrun in configuration namelist : Parameters of the run
	183	READ ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
[5836]	184	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
[4624]	185	IF(lwm) WRITE ( numond, namrun )
[1601]	186	!
	187	IF(lwp) THEN ! control print
[3]	188	WRITE(numout,*)
	189	WRITE(numout,*) 'dom_nam : domain initialization through namelist read'
	190	WRITE(numout,*) '~~~~~~~ '
[1601]	191	WRITE(numout,*) ' Namelist namrun'
	192	WRITE(numout,*) ' job number nn_no = ', nn_no
	193	WRITE(numout,*) ' experiment name for output cn_exp = ', cn_exp
[4147]	194	WRITE(numout,*) ' file prefix restart input cn_ocerst_in= ', cn_ocerst_in
[5341]	195	WRITE(numout,*) ' restart input directory cn_ocerst_indir= ', cn_ocerst_indir
[4147]	196	WRITE(numout,*) ' file prefix restart output cn_ocerst_out= ', cn_ocerst_out
[5341]	197	WRITE(numout,*) ' restart output directory cn_ocerst_outdir= ', cn_ocerst_outdir
[1601]	198	WRITE(numout,*) ' restart logical ln_rstart = ', ln_rstart
[4370]	199	WRITE(numout,*) ' start with forward time step nn_euler = ', nn_euler
[1604]	200	WRITE(numout,*) ' control of time step nn_rstctl = ', nn_rstctl
[1601]	201	WRITE(numout,*) ' number of the first time step nn_it000 = ', nn_it000
	202	WRITE(numout,*) ' number of the last time step nn_itend = ', nn_itend
	203	WRITE(numout,*) ' initial calendar date aammjj nn_date0 = ', nn_date0
[6140]	204	WRITE(numout,*) ' initial time of day in hhmm nn_time0 = ', nn_time0
[1601]	205	WRITE(numout,*) ' leap year calendar (0/1) nn_leapy = ', nn_leapy
	206	WRITE(numout,*) ' initial state output nn_istate = ', nn_istate
[5341]	207	IF( ln_rst_list ) THEN
	208	WRITE(numout,*) ' list of restart dump times nn_stocklist =', nn_stocklist
	209	ELSE
	210	WRITE(numout,*) ' frequency of restart file nn_stock = ', nn_stock
	211	ENDIF
[1601]	212	WRITE(numout,*) ' frequency of output file nn_write = ', nn_write
	213	WRITE(numout,*) ' mask land points ln_mskland = ', ln_mskland
[5363]	214	WRITE(numout,*) ' additional CF standard metadata ln_cfmeta = ', ln_cfmeta
[1601]	215	WRITE(numout,*) ' overwrite an existing file ln_clobber = ', ln_clobber
	216	WRITE(numout,*) ' NetCDF chunksize (bytes) nn_chunksz = ', nn_chunksz
[6140]	217	WRITE(numout,*) ' IS coupling at the restart step ln_iscpl = ', ln_iscpl
[3]	218	ENDIF
	219
[1601]	220	no = nn_no ! conversion DOCTOR names into model names (this should disappear soon)
	221	cexper = cn_exp
	222	nrstdt = nn_rstctl
	223	nit000 = nn_it000
	224	nitend = nn_itend
	225	ndate0 = nn_date0
	226	nleapy = nn_leapy
	227	ninist = nn_istate
	228	nstock = nn_stock
[5341]	229	nstocklist = nn_stocklist
[1601]	230	nwrite = nn_write
[4370]	231	neuler = nn_euler
[5341]	232	IF ( neuler == 1 .AND. .NOT. ln_rstart ) THEN
[4370]	233	WRITE(ctmp1,*) 'ln_rstart =.FALSE., nn_euler is forced to 0 '
	234	CALL ctl_warn( ctmp1 )
	235	neuler = 0
	236	ENDIF
[3]	237
[1601]	238	! ! control of output frequency
[1335]	239	IF ( nstock == 0 .OR. nstock > nitend ) THEN
[1601]	240	WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
[783]	241	CALL ctl_warn( ctmp1 )
[1335]	242	nstock = nitend
[3]	243	ENDIF
	244	IF ( nwrite == 0 ) THEN
[1601]	245	WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
[783]	246	CALL ctl_warn( ctmp1 )
	247	nwrite = nitend
[3]	248	ENDIF
	249
[2528]	250	#if defined key_agrif
[1601]	251	IF( Agrif_Root() ) THEN
[2528]	252	#endif
	253	SELECT CASE ( nleapy ) ! Choose calendar for IOIPSL
	254	CASE ( 1 )
	255	CALL ioconf_calendar('gregorian')
	256	IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "gregorian", i.e. leap year'
	257	CASE ( 0 )
	258	CALL ioconf_calendar('noleap')
	259	IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "noleap", i.e. no leap year'
	260	CASE ( 30 )
	261	CALL ioconf_calendar('360d')
	262	IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "360d", i.e. 360 days in a year'
	263	END SELECT
	264	#if defined key_agrif
[1601]	265	ENDIF
[2528]	266	#endif
[3]	267
[4147]	268	REWIND( numnam_ref ) ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
	269	READ ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
	270	903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
[4152]	271
	272	!
[4147]	273	REWIND( numnam_cfg ) ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
	274	READ ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
	275	904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
[4624]	276	IF(lwm) WRITE ( numond, namdom )
[5836]	277	!
[3]	278	IF(lwp) THEN
[72]	279	WRITE(numout,*)
[1601]	280	WRITE(numout,*) ' Namelist namdom : space & time domain'
	281	WRITE(numout,*) ' flag read/compute bathymetry nn_bathy = ', nn_bathy
[4245]	282	WRITE(numout,*) ' Depth (if =0 bathy=jpkm1) rn_bathy = ', rn_bathy
[2528]	283	WRITE(numout,*) ' min depth of the ocean (>0) or rn_hmin = ', rn_hmin
	284	WRITE(numout,*) ' min number of ocean level (<0) '
[6140]	285	WRITE(numout,*) ' treshold to open the isf cavity rn_isfhmin = ', rn_isfhmin, ' (m)'
[1601]	286	WRITE(numout,*) ' minimum thickness of partial rn_e3zps_min = ', rn_e3zps_min, ' (m)'
	287	WRITE(numout,*) ' step level rn_e3zps_rat = ', rn_e3zps_rat
	288	WRITE(numout,*) ' create mesh/mask file(s) nn_msh = ', nn_msh
[2528]	289	WRITE(numout,*) ' = 0 no file created '
	290	WRITE(numout,*) ' = 1 mesh_mask '
	291	WRITE(numout,*) ' = 2 mesh and mask '
	292	WRITE(numout,*) ' = 3 mesh_hgr, msh_zgr and mask'
[4152]	293	WRITE(numout,*) ' ocean time step rn_rdt = ', rn_rdt
	294	WRITE(numout,*) ' asselin time filter parameter rn_atfp = ', rn_atfp
	295	WRITE(numout,*) ' suppression of closed seas (=0) nn_closea = ', nn_closea
	296	WRITE(numout,*) ' online coarsening of dynamical fields ln_crs = ', ln_crs
[4147]	297	WRITE(numout,*) ' type of horizontal mesh jphgr_msh = ', jphgr_msh
	298	WRITE(numout,*) ' longitude of first raw and column T-point ppglam0 = ', ppglam0
	299	WRITE(numout,*) ' latitude of first raw and column T-point ppgphi0 = ', ppgphi0
	300	WRITE(numout,*) ' zonal grid-spacing (degrees) ppe1_deg = ', ppe1_deg
	301	WRITE(numout,*) ' meridional grid-spacing (degrees) ppe2_deg = ', ppe2_deg
	302	WRITE(numout,*) ' zonal grid-spacing (degrees) ppe1_m = ', ppe1_m
	303	WRITE(numout,*) ' meridional grid-spacing (degrees) ppe2_m = ', ppe2_m
	304	WRITE(numout,*) ' ORCA r4, r2 and r05 coefficients ppsur = ', ppsur
	305	WRITE(numout,*) ' ppa0 = ', ppa0
	306	WRITE(numout,*) ' ppa1 = ', ppa1
	307	WRITE(numout,*) ' ppkth = ', ppkth
	308	WRITE(numout,*) ' ppacr = ', ppacr
	309	WRITE(numout,*) ' Minimum vertical spacing ppdzmin = ', ppdzmin
	310	WRITE(numout,*) ' Maximum depth pphmax = ', pphmax
	311	WRITE(numout,*) ' Use double tanf function for vertical coordinates ldbletanh = ', ldbletanh
	312	WRITE(numout,*) ' Double tanh function parameters ppa2 = ', ppa2
	313	WRITE(numout,*) ' ppkth2 = ', ppkth2
	314	WRITE(numout,*) ' ppacr2 = ', ppacr2
[223]	315	ENDIF
[5836]	316	!
[1601]	317	ntopo = nn_bathy ! conversion DOCTOR names into model names (this should disappear soon)
	318	e3zps_min = rn_e3zps_min
	319	e3zps_rat = rn_e3zps_rat
	320	nmsh = nn_msh
	321	atfp = rn_atfp
	322	rdt = rn_rdt
	323
[2528]	324	#if defined key_netcdf4
	325	! ! NetCDF 4 case ("key_netcdf4" defined)
[4147]	326	REWIND( numnam_ref ) ! Namelist namnc4 in reference namelist : NETCDF
	327	READ ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
	328	907 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
	329
	330	REWIND( numnam_cfg ) ! Namelist namnc4 in configuration namelist : NETCDF
	331	READ ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
	332	908 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
[4624]	333	IF(lwm) WRITE( numond, namnc4 )
[4147]	334
[2528]	335	IF(lwp) THEN ! control print
	336	WRITE(numout,*)
	337	WRITE(numout,*) ' Namelist namnc4 - Netcdf4 chunking parameters'
	338	WRITE(numout,*) ' number of chunks in i-dimension nn_nchunks_i = ', nn_nchunks_i
	339	WRITE(numout,*) ' number of chunks in j-dimension nn_nchunks_j = ', nn_nchunks_j
	340	WRITE(numout,*) ' number of chunks in k-dimension nn_nchunks_k = ', nn_nchunks_k
	341	WRITE(numout,*) ' apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
	342	ENDIF
[1601]	343
[2528]	344	! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
	345	! Note the chunk size in the unlimited (time) dimension will be fixed at 1
	346	snc4set%ni = nn_nchunks_i
	347	snc4set%nj = nn_nchunks_j
	348	snc4set%nk = nn_nchunks_k
	349	snc4set%luse = ln_nc4zip
	350	#else
	351	snc4set%luse = .FALSE. ! No NetCDF 4 case
	352	#endif
[1438]	353	!
[3]	354	END SUBROUTINE dom_nam
	355
	356
	357	SUBROUTINE dom_ctl
	358	!!----------------------------------------------------------------------
	359	!! * ROUTINE dom_ctl *
	360	!!
	361	!! ** Purpose : Domain control.
	362	!!
	363	!! ** Method : compute and print extrema of masked scale factors
	364	!!----------------------------------------------------------------------
	365	INTEGER :: iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
[1601]	366	INTEGER, DIMENSION(2) :: iloc !
[3]	367	REAL(wp) :: ze1min, ze1max, ze2min, ze2max
	368	!!----------------------------------------------------------------------
[1601]	369	!
	370	IF(lk_mpp) THEN
[4990]	371	CALL mpp_minloc( e1t(:,:), tmask_i(:,:), ze1min, iimi1,ijmi1 )
	372	CALL mpp_minloc( e2t(:,:), tmask_i(:,:), ze2min, iimi2,ijmi2 )
	373	CALL mpp_maxloc( e1t(:,:), tmask_i(:,:), ze1max, iima1,ijma1 )
	374	CALL mpp_maxloc( e2t(:,:), tmask_i(:,:), ze2max, iima2,ijma2 )
[181]	375	ELSE
[4990]	376	ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
	377	ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
	378	ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
	379	ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
[32]	380
[4990]	381	iloc = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]	382	iimi1 = iloc(1) + nimpp - 1
	383	ijmi1 = iloc(2) + njmpp - 1
[4990]	384	iloc = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]	385	iimi2 = iloc(1) + nimpp - 1
	386	ijmi2 = iloc(2) + njmpp - 1
[4990]	387	iloc = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]	388	iima1 = iloc(1) + nimpp - 1
	389	ijma1 = iloc(2) + njmpp - 1
[4990]	390	iloc = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]	391	iima2 = iloc(1) + nimpp - 1
	392	ijma2 = iloc(2) + njmpp - 1
[32]	393	ENDIF
[3]	394	IF(lwp) THEN
[1601]	395	WRITE(numout,*)
	396	WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
	397	WRITE(numout,*) '~~~~~~~'
[181]	398	WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
	399	WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
	400	WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
	401	WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
[3]	402	ENDIF
[1438]	403	!
[3]	404	END SUBROUTINE dom_ctl
	405
[5836]	406
[3680]	407	SUBROUTINE dom_stiff
	408	!!----------------------------------------------------------------------
	409	!! * ROUTINE dom_stiff *
	410	!!
	411	!! ** Purpose : Diagnose maximum grid stiffness/hydrostatic consistency
	412	!!
	413	!! ** Method : Compute Haney (1991) hydrostatic condition ratio
	414	!! Save the maximum in the vertical direction
	415	!! (this number is only relevant in s-coordinates)
	416	!!
	417	!! Haney, R. L., 1991: On the pressure gradient force
	418	!! over steep topography in sigma coordinate ocean models.
	419	!! J. Phys. Oceanogr., 21, 610???619.
	420	!!----------------------------------------------------------------------
	421	INTEGER :: ji, jj, jk
	422	REAL(wp) :: zrxmax
[6140]	423	REAL(wp), DIMENSION(4) :: zr1
[3680]	424	!!----------------------------------------------------------------------
[5836]	425	rx1(:,:) = 0._wp
	426	zrxmax = 0._wp
	427	zr1(:) = 0._wp
	428	!
[3680]	429	DO ji = 2, jpim1
	430	DO jj = 2, jpjm1
	431	DO jk = 1, jpkm1
[6140]	432	zr1(1) = ABS( ( gdepw_0(ji ,jj,jk )-gdepw_0(ji-1,jj,jk ) &
	433	& +gdepw_0(ji ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) ) &
	434	& / ( gdepw_0(ji ,jj,jk )+gdepw_0(ji-1,jj,jk ) &
	435	& -gdepw_0(ji ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) + rsmall ) ) * umask(ji-1,jj,jk)
	436	zr1(2) = ABS( ( gdepw_0(ji+1,jj,jk )-gdepw_0(ji ,jj,jk ) &
	437	& +gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji ,jj,jk+1) ) &
	438	& / ( gdepw_0(ji+1,jj,jk )+gdepw_0(ji ,jj,jk ) &
	439	& -gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji ,jj,jk+1) + rsmall ) ) * umask(ji ,jj,jk)
	440	zr1(3) = ABS( ( gdepw_0(ji,jj+1,jk )-gdepw_0(ji,jj ,jk ) &
	441	& +gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj ,jk+1) ) &
	442	& / ( gdepw_0(ji,jj+1,jk )+gdepw_0(ji,jj ,jk ) &
	443	& -gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj ,jk+1) + rsmall ) ) * vmask(ji,jj ,jk)
	444	zr1(4) = ABS( ( gdepw_0(ji,jj ,jk )-gdepw_0(ji,jj-1,jk ) &
	445	& +gdepw_0(ji,jj ,jk+1)-gdepw_0(ji,jj-1,jk+1) ) &
	446	& / ( gdepw_0(ji,jj ,jk )+gdepw_0(ji,jj-1,jk ) &
	447	& -gdepw_0(ji,jj ,jk+1)-gdepw_0(ji,jj-1,jk+1) + rsmall ) ) * vmask(ji,jj-1,jk)
	448	zrxmax = MAXVAL( zr1(1:4) )
	449	rx1(ji,jj) = MAX( rx1(ji,jj) , zrxmax )
[3680]	450	END DO
	451	END DO
	452	END DO
	453	CALL lbc_lnk( rx1, 'T', 1. )
[5836]	454	!
	455	zrxmax = MAXVAL( rx1 )
	456	!
[3680]	457	IF( lk_mpp ) CALL mpp_max( zrxmax ) ! max over the global domain
[5836]	458	!
[3680]	459	IF(lwp) THEN
	460	WRITE(numout,*)
	461	WRITE(numout,*) 'dom_stiff : maximum grid stiffness ratio: ', zrxmax
	462	WRITE(numout,*) '~~~~~~~~~'
	463	ENDIF
[5836]	464	!
[3680]	465	END SUBROUTINE dom_stiff
	466
[3]	467	!!======================================================================
	468	END MODULE domain

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

Download in other formats: