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 branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

Context Navigation

source: branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90 @ 5802

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

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

Download in other formats: