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.

ldfdyn_c2d.h90 in trunk/NEMO/OPA_SRC/LDF – NEMO

Context Navigation

source: trunk/NEMO/OPA_SRC/LDF/ldfdyn_c2d.h90 @ 1152

Last change on this file since 1152 was 1152, checked in by rblod, 16 years ago
Convert cvs header to svn Id, step II
Property svn:eol-style set to `native` Property svn:executable set to ``* Property svn:keywords set to `Id`
File size: 11.1 KB

Line
1	!!----------------------------------------------------------------------
2	!! * ldfdyn_c2d.h90 *
3	!!----------------------------------------------------------------------
4	!! ldf_dyn_c2d : set the lateral viscosity coefficients
5	!! ldf_dyn_c2d_orca : specific case for orca r2 and r4
6	!!----------------------------------------------------------------------
7
8	!!----------------------------------------------------------------------
9	!! OPA 9.0 , LOCEAN-IPSL (2005)
10	!! $Id$
11	!! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
12	!!----------------------------------------------------------------------
13
14	SUBROUTINE ldf_dyn_c2d( ld_print )
15	!!----------------------------------------------------------------------
16	!! * ROUTINE ldf_dyn_c2d *
17	!!
18	!! ** Purpose : initializations of the horizontal ocean physics
19	!!
20	!! ** Method :
21	!! 2D eddy viscosity coefficients ( longitude, latitude )
22	!!
23	!! harmonic operator : ahm1 is defined at t-point
24	!! ahm2 is defined at f-point
25	!! + isopycnal : ahm3 is defined at u-point
26	!! or geopotential ahm4 is defined at v-point
27	!! iso-model level : ahm3, ahm4 not used
28	!!
29	!! biharmonic operator : ahm1 is defined at u-point
30	!! ahm2 is defined at v-point
31	!! : ahm3, ahm4 not used
32	!!
33	!!----------------------------------------------------------------------
34	!! * Arguments
35	LOGICAL, INTENT (in) :: ld_print ! If true, output arrays on numout
36
37	!! * Local variables
38	INTEGER :: ji, jj
39	REAL(wp) :: za00, zd_max, zetmax, zeumax, zefmax, zevmax
40	!!----------------------------------------------------------------------
41
42	IF(lwp) WRITE(numout,*)
43	IF(lwp) WRITE(numout,*) 'ldf_dyn_c2d : 2d lateral eddy viscosity coefficient'
44	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
45	IF(lwp) WRITE(numout,*)
46
47	! harmonic operator (ahm1, ahm2) : ( T- and F- points) (used for laplacian operators
48	! =============================== whatever its orientation is)
49	IF( ln_dynldf_lap ) THEN
50	! define ahm1 and ahm2 at the right grid point position
51	! (USER: modify ahm1 and ahm2 following your desiderata)
52
53	zd_max = MAX( MAXVAL( e1t(:,:) ), MAXVAL( e2t(:,:) ) )
54	IF( lk_mpp ) CALL mpp_max( zd_max ) ! max over the global domain
55
56	IF(lwp) WRITE(numout,*) ' laplacian operator: ahm proportional to e1'
57	IF(lwp) WRITE(numout,*) ' maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0
58
59	za00 = ahm0 / zd_max
60	DO jj = 1, jpj
61	DO ji = 1, jpi
62	zetmax = MAX( e1t(ji,jj), e2t(ji,jj) )
63	zefmax = MAX( e1f(ji,jj), e2f(ji,jj) )
64	ahm1(ji,jj) = za00 * zetmax
65	ahm2(ji,jj) = za00 * zefmax
66	END DO
67	END DO
68
69	IF( ln_dynldf_iso ) THEN
70	IF(lwp) WRITE(numout,*) ' Caution, as implemented now, the isopycnal part of momentum'
71	IF(lwp) WRITE(numout,*) ' mixing use aht0 as eddy viscosity coefficient. Thus, it is'
72	IF(lwp) WRITE(numout,*) ' uniform and you must be sure that your ahm is greater than'
73	IF(lwp) WRITE(numout,*) ' aht0 everywhere in the model domain.'
74	ENDIF
75
76	! Special case for ORCA R2 and R4 configurations (overwrite the value of ahm1 ahm2)
77	! ==============================================
78	IF( cp_cfg == "orca" .AND. ( jp_cfg == 2 .OR. jp_cfg == 4 ) ) CALL ldf_dyn_c2d_orca( ld_print )
79
80	! Control print
81	IF( lwp .AND. ld_print ) THEN
82	WRITE(numout,*)
83	WRITE(numout,*) 'inildf: 2D ahm1 array'
84	CALL prihre(ahm1,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
85	WRITE(numout,*)
86	WRITE(numout,*) 'inildf: 2D ahm2 array'
87	CALL prihre(ahm2,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
88	ENDIF
89	ENDIF
90
91	! biharmonic operator (ahm3, ahm4) : at U- and V-points (used for bilaplacian operator
92	! ================================= whatever its orientation is)
93	IF( ln_dynldf_bilap ) THEN
94	! (USER: modify ahm3 and ahm4 following your desiderata)
95	! Here: ahm is proportional to the cube of the maximum of the gridspacing
96	! in the to horizontal direction
97
98	zd_max = MAX( MAXVAL( e1u(:,:) ), MAXVAL( e2u(:,:) ) )
99	IF( lk_mpp ) CALL mpp_max( zd_max ) ! max over the global domain
100
101	IF(lwp) WRITE(numout,) ' bi-laplacian operator: ahm proportional to e1*3 '
102	IF(lwp) WRITE(numout,*) ' maximum grid-spacing = ', zd_max, ' maximum value for ahm = ', ahm0
103
104	za00 = ahm0 / ( zd_max * zd_max * zd_max )
105	DO jj = 1, jpj
106	DO ji = 1, jpi
107	zeumax = MAX( e1u(ji,jj), e2u(ji,jj) )
108	zevmax = MAX( e1v(ji,jj), e2v(ji,jj) )
109	ahm3(ji,jj) = za00 * zeumax * zeumax * zeumax
110	ahm4(ji,jj) = za00 * zevmax * zevmax * zevmax
111	END DO
112	END DO
113
114	! Control print
115	IF( lwp .AND. ld_print ) THEN
116	WRITE(numout,*)
117	WRITE(numout,*) 'inildf: ahm3 array'
118	CALL prihre(ahm3,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
119	WRITE(numout,*)
120	WRITE(numout,*) 'inildf: ahm4 array'
121	CALL prihre(ahm4,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
122	ENDIF
123	ENDIF
124
125
126	END SUBROUTINE ldf_dyn_c2d
127
128
129	SUBROUTINE ldf_dyn_c2d_orca( ld_print )
130	!!----------------------------------------------------------------------
131	!! * ROUTINE ldf_dyn_c2d *
132	!!
133	!! ** W A R N I N G **
134	!!
135	!! ORCA R2 and R4 configurations
136	!!
137	!! ** W A R N I N G **
138	!!
139	!! ** Purpose : initializations of the lateral viscosity for orca R2
140	!!
141	!! ** Method : blah blah blah...
142	!!
143	!!----------------------------------------------------------------------
144	!! * Modules used
145	USE ldftra_oce, ONLY : aht0
146
147	!! * Arguments
148	LOGICAL, INTENT (in) :: ld_print ! If true, output arrays on numout
149
150	!! * Local variables
151	INTEGER :: ji, jj, jn ! dummy loop indices
152	INTEGER :: inum ! temporary logical unit
153	INTEGER :: iim, ijm
154	INTEGER :: ifreq, il1, il2, ij, ii
155	INTEGER, DIMENSION(jpidta,jpidta) :: idata
156	INTEGER, DIMENSION(jpi ,jpj ) :: icof
157
158	REAL(wp) :: zahmeq, zcoft, zcoff, zmsk
159
160	CHARACTER (len=15) :: clexp
161	!!----------------------------------------------------------------------
162
163	IF(lwp) WRITE(numout,*)
164	IF(lwp) WRITE(numout,*) 'inildf: 2d eddy viscosity coefficient'
165	IF(lwp) WRITE(numout,*) '~~~~~~ --'
166	IF(lwp) WRITE(numout,*)
167	IF(lwp) WRITE(numout,*) ' orca ocean model'
168	IF(lwp) WRITE(numout,*)
169
170	#if defined key_antarctic
171	# include "ldfdyn_antarctic.h90"
172	#elif defined key_arctic
173	# include "ldfdyn_arctic.h90"
174	#else
175	! Read 2d integer array to specify western boundary increase in the
176	! ===================== equatorial strip (20N-20S) defined at t-points
177
178	CALL ctlopn( inum, 'ahmcoef', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', &
179	& 1, numout, lwp, 1 )
180	REWIND inum
181	READ(inum,9101) clexp, iim, ijm
182	READ(inum,'(/)')
183	ifreq = 40
184	il1 = 1
185	DO jn = 1, jpidta/ifreq+1
186	READ(inum,'(/)')
187	il2 = MIN( jpidta, il1+ifreq-1 )
188	READ(inum,9201) ( ii, ji = il1, il2, 5 )
189	READ(inum,'(/)')
190	DO jj = jpjdta, 1, -1
191	READ(inum,9202) ij, ( idata(ji,jj), ji = il1, il2 )
192	END DO
193	il1 = il1 + ifreq
194	END DO
195
196	DO jj = 1, nlcj
197	DO ji = 1, nlci
198	icof(ji,jj) = idata( mig(ji), mjg(jj) )
199	END DO
200	END DO
201	DO jj = nlcj+1, jpj
202	DO ji = 1, nlci
203	icof(ji,jj) = icof(ji,nlcj)
204	END DO
205	END DO
206	DO jj = 1, jpj
207	DO ji = nlci+1, jpi
208	icof(ji,jj) = icof(nlci,jj)
209	END DO
210	END DO
211
212	9101 FORMAT(1x,a15,2i8)
213	9201 FORMAT(3x,13(i3,12x))
214	9202 FORMAT(i3,41i3)
215
216
217	! Set ahm1 and ahm2 ( T- and F- points) (used for laplacian operator)
218	! =================
219	! define ahm1 and ahm2 at the right grid point position
220	! (USER: modify ahm1 and ahm2 following your desiderata)
221
222
223	! Decrease ahm to zahmeq m2/s in the tropics
224	! (from 90 to 20 degre: ahm = constant
225	! from 20 to 2.5 degre: ahm = decrease in (1-cos)/2
226	! from 2.5 to 0 degre: ahm = constant
227	! symmetric in the south hemisphere)
228
229	zahmeq = aht0
230
231	DO jj = 1, jpj
232	DO ji = 1, jpi
233	IF( ABS( gphif(ji,jj) ) >= 20. ) THEN
234	ahm2(ji,jj) = ahm0
235	ELSEIF( ABS( gphif(ji,jj) ) <= 2.5 ) THEN
236	ahm2(ji,jj) = zahmeq
237	ELSE
238	ahm2(ji,jj) = zahmeq + (ahm0-zahmeq)/2. &
239	* ( 1. - COS( rad * ( ABS(gphif(ji,jj))-2.5 ) * 180. / 17.5 ) )
240	ENDIF
241	IF( ABS( gphit(ji,jj) ) >= 20. ) THEN
242	ahm1(ji,jj) = ahm0
243	ELSEIF( ABS( gphit(ji,jj) ) <= 2.5 ) THEN
244	ahm1(ji,jj) = zahmeq
245	ELSE
246	ahm1(ji,jj) = zahmeq + (ahm0-zahmeq)/2. &
247	* ( 1. - COS( rad * ( ABS(gphit(ji,jj))-2.5 ) * 180. / 17.5 ) )
248	ENDIF
249	END DO
250	END DO
251
252	! increase along western boundaries of equatorial strip
253	! t-point
254	DO jj = 1, jpjm1
255	DO ji = 1, jpim1
256	zcoft = FLOAT( icof(ji,jj) ) / 100.
257	ahm1(ji,jj) = zcoft * ahm0 + (1.-zcoft) * ahm1(ji,jj)
258	END DO
259	END DO
260	! f-point
261	icof(:,:) = icof(:,:) * tmask(:,:,1)
262	DO jj = 1, jpjm1
263	DO ji = 1, jpim1
264	zmsk = tmask(ji,jj+1,1) + tmask(ji+1,jj+1,1) + tmask(ji,jj,1) + tmask(ji,jj+1,1)
265	IF( zmsk == 0. ) THEN
266	zcoff = 1.
267	ELSE
268	zcoff = FLOAT( icof(ji,jj+1) + icof(ji+1,jj+1) + icof(ji,jj) + icof(ji,jj+1) ) &
269	/ (zmsk * 100.)
270	ENDIF
271	ahm2(ji,jj) = zcoff * ahm0 + (1.-zcoff) * ahm2(ji,jj)
272	END DO
273	END DO
274	#endif
275
276	! Lateral boundary conditions on ( ahm1, ahm2 )
277	! ==============
278	CALL lbc_lnk( ahm1, 'T', 1. ) ! T-point, unchanged sign
279	CALL lbc_lnk( ahm2, 'F', 1. ) ! F-point, unchanged sign
280
281	! Control print
282	IF( lwp .AND. ld_print ) THEN
283	WRITE(numout,*)
284	WRITE(numout,*) 'inildf: 2D ahm1 array'
285	CALL prihre(ahm1,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
286	WRITE(numout,*)
287	WRITE(numout,*) 'inildf: 2D ahm2 array'
288	CALL prihre(ahm2,jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)
289	ENDIF
290
291	END SUBROUTINE ldf_dyn_c2d_orca

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

Download in other formats:

Original Format