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.

traqsr.F90 in branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

Context Navigation

source: branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90 @ 4147

Last change on this file since 4147 was 4147, checked in by cetlod, 10 years ago
merge in dev_LOCEAN_2013, the 1st development branch dev_r3853_CNRS9_Confsetting, from its starting point ( r3853 ) on the trunk: see ticket #1169
Property svn:keywords set to `Id`
File size: 27.0 KB

Line
1	MODULE traqsr
2	!!======================================================================
3	!! * MODULE traqsr *
4	!! Ocean physics: solar radiation penetration in the top ocean levels
5	!!======================================================================
6	!! History : OPA ! 1990-10 (B. Blanke) Original code
7	!! 7.0 ! 1991-11 (G. Madec)
8	!! ! 1996-01 (G. Madec) s-coordinates
9	!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
10	!! - ! 2005-11 (G. Madec) zco, zps, sco coordinate
11	!! 3.2 ! 2009-04 (G. Madec & NEMO team)
12	!!----------------------------------------------------------------------
13
14	!!----------------------------------------------------------------------
15	!! tra_qsr : trend due to the solar radiation penetration
16	!! tra_qsr_init : solar radiation penetration initialization
17	!!----------------------------------------------------------------------
18	USE oce ! ocean dynamics and active tracers
19	USE dom_oce ! ocean space and time domain
20	USE sbc_oce ! surface boundary condition: ocean
21	USE trc_oce ! share SMS/Ocean variables
22	USE trdmod_oce ! ocean variables trends
23	USE trdtra ! ocean active tracers trends
24	USE in_out_manager ! I/O manager
25	USE phycst ! physical constants
26	USE prtctl ! Print control
27	USE iom ! I/O manager
28	USE fldread ! read input fields
29	USE lib_mpp ! MPP library
30	USE wrk_nemo ! Memory Allocation
31	USE timing ! Timing
32
33
34	IMPLICIT NONE
35	PRIVATE
36
37	PUBLIC tra_qsr ! routine called by step.F90 (ln_traqsr=T)
38	PUBLIC tra_qsr_init ! routine called by opa.F90
39
40	! !!* Namelist namtra_qsr: penetrative solar radiation
41	LOGICAL , PUBLIC :: ln_traqsr !: light absorption (qsr) flag
42	LOGICAL , PUBLIC :: ln_qsr_rgb !: Red-Green-Blue light absorption flag
43	LOGICAL , PUBLIC :: ln_qsr_2bd !: 2 band light absorption flag
44	LOGICAL , PUBLIC :: ln_qsr_bio !: bio-model light absorption flag
45	INTEGER , PUBLIC :: nn_chldta !: use Chlorophyll data (=1) or not (=0)
46	REAL(wp), PUBLIC :: rn_abs !: fraction absorbed in the very near surface (RGB & 2 bands)
47	REAL(wp), PUBLIC :: rn_si0 !: very near surface depth of extinction (RGB & 2 bands)
48	REAL(wp), PUBLIC :: rn_si1 !: deepest depth of extinction (water type I) (2 bands)
49
50	! Module variables
51	REAL(wp) :: xsi0r !: inverse of rn_si0
52	REAL(wp) :: xsi1r !: inverse of rn_si1
53	TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_chl ! structure of input Chl (file informations, fields read)
54	INTEGER, PUBLIC :: nksr ! levels below which the light cannot penetrate ( depth larger than 391 m)
55	REAL(wp), DIMENSION(3,61) :: rkrgb !: tabulated attenuation coefficients for RGB absorption
56
57	!! * Substitutions
58	# include "domzgr_substitute.h90"
59	# include "vectopt_loop_substitute.h90"
60	!!----------------------------------------------------------------------
61	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
62	!! $Id$
63	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
64	!!----------------------------------------------------------------------
65	CONTAINS
66
67	SUBROUTINE tra_qsr( kt )
68	!!----------------------------------------------------------------------
69	!! * ROUTINE tra_qsr *
70	!!
71	!! ** Purpose : Compute the temperature trend due to the solar radiation
72	!! penetration and add it to the general temperature trend.
73	!!
74	!! ** Method : The profile of the solar radiation within the ocean is defined
75	!! through 2 wavebands (rn_si0,rn_si1) or 3 wavebands (RGB) and a ratio rn_abs
76	!! Considering the 2 wavebands case:
77	!! I(k) = Qsr( rn_absEXP(z(k)/rn_si0) + (1.-rn_abs)*EXP(z(k)/rn_si1) )
78	!! The temperature trend associated with the solar radiation penetration
79	!! is given by : zta = 1/e3t dk[ I ] / (rau0*Cp)
80	!! At the bottom, boudary condition for the radiation is no flux :
81	!! all heat which has not been absorbed in the above levels is put
82	!! in the last ocean level.
83	!! In z-coordinate case, the computation is only done down to the
84	!! level where I(k) < 1.e-15 W/m2. In addition, the coefficients
85	!! used for the computation are calculated one for once as they
86	!! depends on k only.
87	!!
88	!! ** Action : - update ta with the penetrative solar radiation trend
89	!! - save the trend in ttrd ('key_trdtra')
90	!!
91	!! Reference : Jerlov, N. G., 1968 Optical Oceanography, Elsevier, 194pp.
92	!! Lengaigne et al. 2007, Clim. Dyn., V28, 5, 503-516.
93	!!----------------------------------------------------------------------
94	!
95	INTEGER, INTENT(in) :: kt ! ocean time-step
96	!
97	INTEGER :: ji, jj, jk ! dummy loop indices
98	INTEGER :: irgb ! local integers
99	REAL(wp) :: zchl, zcoef, zfact ! local scalars
100	REAL(wp) :: zc0, zc1, zc2, zc3 ! - -
101	REAL(wp) :: zz0, zz1, z1_e3t ! - -
102	REAL(wp), POINTER, DIMENSION(:,: ) :: zekb, zekg, zekr
103	REAL(wp), POINTER, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea, ztrdt
104	!!----------------------------------------------------------------------
105	!
106	IF( nn_timing == 1 ) CALL timing_start('tra_qsr')
107	!
108	CALL wrk_alloc( jpi, jpj, zekb, zekg, zekr )
109	CALL wrk_alloc( jpi, jpj, jpk, ze0, ze1, ze2, ze3, zea )
110	!
111	IF( kt == nit000 ) THEN
112	IF(lwp) WRITE(numout,*)
113	IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation'
114	IF(lwp) WRITE(numout,*) '~~~~~~~'
115	IF( .NOT.ln_traqsr ) RETURN
116	ENDIF
117
118	IF( l_trdtra ) THEN ! Save ta and sa trends
119	CALL wrk_alloc( jpi, jpj, jpk, ztrdt )
120	ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
121	ENDIF
122
123	! Set before qsr tracer content field
124	! ***********************************
125	IF( kt == nit000 ) THEN ! Set the forcing field at nit000 - 1
126	! ! -----------------------------------
127	IF( ln_rstart .AND. & ! Restart: read in restart file
128	& iom_varid( numror, 'qsr_hc_b', ldstop = .FALSE. ) > 0 ) THEN
129	IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field red in the restart file'
130	zfact = 0.5e0
131	CALL iom_get( numror, jpdom_autoglo, 'qsr_hc_b', qsr_hc_b ) ! before heat content trend due to Qsr flux
132	ELSE ! No restart or restart not found: Euler forward time stepping
133	zfact = 1.e0
134	qsr_hc_b(:,:,:) = 0.e0
135	ENDIF
136	ELSE ! Swap of forcing field
137	! ! ---------------------
138	zfact = 0.5e0
139	qsr_hc_b(:,:,:) = qsr_hc(:,:,:)
140	ENDIF
141	! Compute now qsr tracer content field
142	! ************************************
143
144	! ! ============================================== !
145	IF( lk_qsr_bio .AND. ln_qsr_bio ) THEN ! bio-model fluxes : all vertical coordinates !
146	! ! ============================================== !
147	DO jk = 1, jpkm1
148	qsr_hc(:,:,jk) = r1_rau0_rcp * ( etot3(:,:,jk) - etot3(:,:,jk+1) )
149	END DO
150	! Add to the general trend
151	DO jk = 1, jpkm1
152	DO jj = 2, jpjm1
153	DO ji = fs_2, fs_jpim1 ! vector opt.
154	z1_e3t = zfact / fse3t(ji,jj,jk)
155	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) * z1_e3t
156	END DO
157	END DO
158	END DO
159	CALL iom_put( 'qsr3d', etot3 ) ! Shortwave Radiation 3D distribution
160	! ! ============================================== !
161	ELSE ! Ocean alone :
162	! ! ============================================== !
163	!
164	! ! ------------------------- !
165	IF( ln_qsr_rgb) THEN ! R-G-B light penetration !
166	! ! ------------------------- !
167	! Set chlorophyl concentration
168	IF( nn_chldta == 1 .OR. lk_vvl ) THEN !* Variable Chlorophyll or ocean volume
169	!
170	IF( nn_chldta == 1 ) THEN !* Variable Chlorophyll
171	!
172	CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step
173	!
174	!CDIR COLLAPSE
175	!CDIR NOVERRCHK
176	DO jj = 1, jpj ! Separation in R-G-B depending of the surface Chl
177	!CDIR NOVERRCHK
178	DO ji = 1, jpi
179	zchl = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) )
180	irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 )
181	zekb(ji,jj) = rkrgb(1,irgb)
182	zekg(ji,jj) = rkrgb(2,irgb)
183	zekr(ji,jj) = rkrgb(3,irgb)
184	END DO
185	END DO
186	ELSE ! Variable ocean volume but constant chrlorophyll
187	zchl = 0.05 ! constant chlorophyll
188	irgb = NINT( 41 + 20.*LOG10( zchl ) + 1.e-15 )
189	zekb(:,:) = rkrgb(1,irgb) ! Separation in R-G-B depending of the chlorophyll
190	zekg(:,:) = rkrgb(2,irgb)
191	zekr(:,:) = rkrgb(3,irgb)
192	ENDIF
193	!
194	zcoef = ( 1. - rn_abs ) / 3.e0 ! equi-partition in R-G-B
195	ze0(:,:,1) = rn_abs * qsr(:,:)
196	ze1(:,:,1) = zcoef * qsr(:,:)
197	ze2(:,:,1) = zcoef * qsr(:,:)
198	ze3(:,:,1) = zcoef * qsr(:,:)
199	zea(:,:,1) = qsr(:,:)
200	!
201	DO jk = 2, nksr+1
202	!CDIR NOVERRCHK
203	DO jj = 1, jpj
204	!CDIR NOVERRCHK
205	DO ji = 1, jpi
206	zc0 = ze0(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * xsi0r )
207	zc1 = ze1(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * zekb(ji,jj) )
208	zc2 = ze2(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * zekg(ji,jj) )
209	zc3 = ze3(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * zekr(ji,jj) )
210	ze0(ji,jj,jk) = zc0
211	ze1(ji,jj,jk) = zc1
212	ze2(ji,jj,jk) = zc2
213	ze3(ji,jj,jk) = zc3
214	zea(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk)
215	END DO
216	END DO
217	END DO
218	!
219	DO jk = 1, nksr ! compute and add qsr trend to ta
220	qsr_hc(:,:,jk) = r1_rau0_rcp * ( zea(:,:,jk) - zea(:,:,jk+1) )
221	END DO
222	zea(:,:,nksr+1:jpk) = 0.e0 ! below 400m set to zero
223	CALL iom_put( 'qsr3d', zea ) ! Shortwave Radiation 3D distribution
224	!
225	ELSE !* Constant Chlorophyll
226	DO jk = 1, nksr
227	qsr_hc(:,:,jk) = etot3(:,:,jk) * qsr(:,:)
228	END DO
229	ENDIF
230
231	ENDIF
232	! ! ------------------------- !
233	IF( ln_qsr_2bd ) THEN ! 2 band light penetration !
234	! ! ------------------------- !
235	!
236	IF( lk_vvl ) THEN !* variable volume
237	zz0 = rn_abs * r1_rau0_rcp
238	zz1 = ( 1. - rn_abs ) * r1_rau0_rcp
239	DO jk = 1, nksr ! solar heat absorbed at T-point in the top 400m
240	DO jj = 1, jpj
241	DO ji = 1, jpi
242	zc0 = zz0 * EXP( -fsdepw(ji,jj,jk )xsi0r ) + zz1 EXP( -fsdepw(ji,jj,jk )*xsi1r )
243	zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)xsi0r ) + zz1 EXP( -fsdepw(ji,jj,jk+1)*xsi1r )
244	qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0tmask(ji,jj,jk) - zc1tmask(ji,jj,jk+1) )
245	END DO
246	END DO
247	END DO
248	ELSE !* constant volume: coef. computed one for all
249	DO jk = 1, nksr
250	DO jj = 2, jpjm1
251	DO ji = fs_2, fs_jpim1 ! vector opt.
252	qsr_hc(ji,jj,jk) = etot3(ji,jj,jk) * qsr(ji,jj)
253	END DO
254	END DO
255	END DO
256	!
257	ENDIF
258	!
259	ENDIF
260	!
261	! Add to the general trend
262	DO jk = 1, nksr
263	DO jj = 2, jpjm1
264	DO ji = fs_2, fs_jpim1 ! vector opt.
265	z1_e3t = zfact / fse3t(ji,jj,jk)
266	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) * z1_e3t
267	END DO
268	END DO
269	END DO
270	!
271	ENDIF
272	!
273	IF( lrst_oce ) THEN ! Write in the ocean restart file
274	! *******************************
275	IF(lwp) WRITE(numout,*)
276	IF(lwp) WRITE(numout,*) 'qsr tracer content forcing field written in ocean restart file ', &
277	& 'at it= ', kt,' date= ', ndastp
278	IF(lwp) WRITE(numout,*) '~~~~'
279	CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b', qsr_hc )
280	!
281	ENDIF
282
283	IF( l_trdtra ) THEN ! qsr tracers trends saved for diagnostics
284	ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
285	CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_qsr, ztrdt )
286	CALL wrk_dealloc( jpi, jpj, jpk, ztrdt )
287	ENDIF
288	! ! print mean trends (used for debugging)
289	IF(ln_ctl) CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' qsr - Ta: ', mask1=tmask, clinfo3='tra-ta' )
290	!
291	CALL wrk_dealloc( jpi, jpj, zekb, zekg, zekr )
292	CALL wrk_dealloc( jpi, jpj, jpk, ze0, ze1, ze2, ze3, zea )
293	!
294	IF( nn_timing == 1 ) CALL timing_stop('tra_qsr')
295	!
296	END SUBROUTINE tra_qsr
297
298
299	SUBROUTINE tra_qsr_init
300	!!----------------------------------------------------------------------
301	!! * ROUTINE tra_qsr_init *
302	!!
303	!! ** Purpose : Initialization for the penetrative solar radiation
304	!!
305	!! ** Method : The profile of solar radiation within the ocean is set
306	!! from two length scale of penetration (rn_si0,rn_si1) and a ratio
307	!! (rn_abs). These parameters are read in the namtra_qsr namelist. The
308	!! default values correspond to clear water (type I in Jerlov'
309	!! (1968) classification.
310	!! called by tra_qsr at the first timestep (nit000)
311	!!
312	!! ** Action : - initialize rn_si0, rn_si1 and rn_abs
313	!!
314	!! Reference : Jerlov, N. G., 1968 Optical Oceanography, Elsevier, 194pp.
315	!!----------------------------------------------------------------------
316	!
317	INTEGER :: ji, jj, jk ! dummy loop indices
318	INTEGER :: irgb, ierror, ioptio, nqsr ! local integer
319	INTEGER :: ios ! Local integer output status for namelist read
320	REAL(wp) :: zz0, zc0 , zc1, zcoef ! local scalars
321	REAL(wp) :: zz1, zc2 , zc3, zchl ! - -
322	REAL(wp), POINTER, DIMENSION(:,: ) :: zekb, zekg, zekr
323	REAL(wp), POINTER, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea
324	!
325	CHARACTER(len=100) :: cn_dir ! Root directory for location of ssr files
326	TYPE(FLD_N) :: sn_chl ! informations about the chlorofyl field to be read
327	!!
328	NAMELIST/namtra_qsr/ sn_chl, cn_dir, ln_traqsr, ln_qsr_rgb, ln_qsr_2bd, ln_qsr_bio, &
329	& nn_chldta, rn_abs, rn_si0, rn_si1
330	!!----------------------------------------------------------------------
331
332	!
333	IF( nn_timing == 1 ) CALL timing_start('tra_qsr_init')
334	!
335	CALL wrk_alloc( jpi, jpj, zekb, zekg, zekr )
336	CALL wrk_alloc( jpi, jpj, jpk, ze0, ze1, ze2, ze3, zea )
337	!
338
339	REWIND( numnam_ref ) ! Namelist namtra_qsr in reference namelist : Ratio and length of penetration
340	READ ( numnam_ref, namtra_qsr, IOSTAT = ios, ERR = 901)
341	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtra_qsr in reference namelist', lwp )
342
343	REWIND( numnam_cfg ) ! Namelist namtra_qsr in configuration namelist : Ratio and length of penetration
344	READ ( numnam_cfg, namtra_qsr, IOSTAT = ios, ERR = 902 )
345	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtra_qsr in configuration namelist', lwp )
346	WRITE ( numond, namtra_qsr )
347	!
348	IF(lwp) THEN ! control print
349	WRITE(numout,*)
350	WRITE(numout,*) 'tra_qsr_init : penetration of the surface solar radiation'
351	WRITE(numout,*) '~~~~~~~~~~~~'
352	WRITE(numout,*) ' Namelist namtra_qsr : set the parameter of penetration'
353	WRITE(numout,*) ' Light penetration (T) or not (F) ln_traqsr = ', ln_traqsr
354	WRITE(numout,*) ' RGB (Red-Green-Blue) light penetration ln_qsr_rgb = ', ln_qsr_rgb
355	WRITE(numout,*) ' 2 band light penetration ln_qsr_2bd = ', ln_qsr_2bd
356	WRITE(numout,*) ' bio-model light penetration ln_qsr_bio = ', ln_qsr_bio
357	WRITE(numout,*) ' RGB : Chl data (=1) or cst value (=0) nn_chldta = ', nn_chldta
358	WRITE(numout,*) ' RGB & 2 bands: fraction of light (rn_si1) rn_abs = ', rn_abs
359	WRITE(numout,*) ' RGB & 2 bands: shortess depth of extinction rn_si0 = ', rn_si0
360	WRITE(numout,*) ' 2 bands: longest depth of extinction rn_si1 = ', rn_si1
361	ENDIF
362
363	IF( ln_traqsr ) THEN ! control consistency
364	!
365	IF( .NOT.lk_qsr_bio .AND. ln_qsr_bio ) THEN
366	CALL ctl_warn( 'No bio model : force ln_qsr_bio = FALSE ' )
367	ln_qsr_bio = .FALSE.
368	ENDIF
369	!
370	ioptio = 0 ! Parameter control
371	IF( ln_qsr_rgb ) ioptio = ioptio + 1
372	IF( ln_qsr_2bd ) ioptio = ioptio + 1
373	IF( ln_qsr_bio ) ioptio = ioptio + 1
374	!
375	IF( ioptio /= 1 ) &
376	CALL ctl_stop( ' Choose ONE type of light penetration in namelist namtra_qsr', &
377	& ' 2 bands, 3 RGB bands or bio-model light penetration' )
378	!
379	IF( ln_qsr_rgb .AND. nn_chldta == 0 ) nqsr = 1
380	IF( ln_qsr_rgb .AND. nn_chldta == 1 ) nqsr = 2
381	IF( ln_qsr_2bd ) nqsr = 3
382	IF( ln_qsr_bio ) nqsr = 4
383	!
384	IF(lwp) THEN ! Print the choice
385	WRITE(numout,*)
386	IF( nqsr == 1 ) WRITE(numout,*) ' R-G-B light penetration - Constant Chlorophyll'
387	IF( nqsr == 2 ) WRITE(numout,*) ' R-G-B light penetration - Chl data '
388	IF( nqsr == 3 ) WRITE(numout,*) ' 2 bands light penetration'
389	IF( nqsr == 4 ) WRITE(numout,*) ' bio-model light penetration'
390	ENDIF
391	!
392	ENDIF
393	! ! ===================================== !
394	IF( ln_traqsr ) THEN ! Initialisation of Light Penetration !
395	! ! ===================================== !
396	!
397	xsi0r = 1.e0 / rn_si0
398	xsi1r = 1.e0 / rn_si1
399	! ! ---------------------------------- !
400	IF( ln_qsr_rgb ) THEN ! Red-Green-Blue light penetration !
401	! ! ---------------------------------- !
402	!
403	CALL trc_oce_rgb( rkrgb ) !* tabulated attenuation coef.
404	!
405	! !* level of light extinction
406	IF( ln_sco ) THEN ; nksr = jpkm1
407	ELSE ; nksr = trc_oce_ext_lev( r_si2, 0.33e2 )
408	ENDIF
409
410	IF(lwp) WRITE(numout,*) ' level of light extinction = ', nksr, ' ref depth = ', gdepw_0(nksr+1), ' m'
411	!
412	IF( nn_chldta == 1 ) THEN !* Chl data : set sf_chl structure
413	IF(lwp) WRITE(numout,*)
414	IF(lwp) WRITE(numout,*) ' Chlorophyll read in a file'
415	ALLOCATE( sf_chl(1), STAT=ierror )
416	IF( ierror > 0 ) THEN
417	CALL ctl_stop( 'tra_qsr_init: unable to allocate sf_chl structure' ) ; RETURN
418	ENDIF
419	ALLOCATE( sf_chl(1)%fnow(jpi,jpj,1) )
420	IF( sn_chl%ln_tint )ALLOCATE( sf_chl(1)%fdta(jpi,jpj,1,2) )
421	! ! fill sf_chl with sn_chl and control print
422	CALL fld_fill( sf_chl, (/ sn_chl /), cn_dir, 'tra_qsr_init', &
423	& 'Solar penetration function of read chlorophyll', 'namtra_qsr' )
424	!
425	ELSE !* constant Chl : compute once for all the distribution of light (etot3)
426	IF(lwp) WRITE(numout,*)
427	IF(lwp) WRITE(numout,*) ' Constant Chlorophyll concentration = 0.05'
428	IF( lk_vvl ) THEN ! variable volume
429	IF(lwp) WRITE(numout,*) ' key_vvl: light distribution will be computed at each time step'
430	ELSE ! constant volume: computes one for all
431	IF(lwp) WRITE(numout,*) ' fixed volume: light distribution computed one for all'
432	!
433	zchl = 0.05 ! constant chlorophyll
434	irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 )
435	zekb(:,:) = rkrgb(1,irgb) ! Separation in R-G-B depending of the chlorophyll
436	zekg(:,:) = rkrgb(2,irgb)
437	zekr(:,:) = rkrgb(3,irgb)
438	!
439	zcoef = ( 1. - rn_abs ) / 3.e0 ! equi-partition in R-G-B
440	ze0(:,:,1) = rn_abs
441	ze1(:,:,1) = zcoef
442	ze2(:,:,1) = zcoef
443	ze3(:,:,1) = zcoef
444	zea(:,:,1) = tmask(:,:,1) ! = ( ze0+ze1+z2+ze3 ) * tmask
445
446	DO jk = 2, nksr+1
447	!CDIR NOVERRCHK
448	DO jj = 1, jpj
449	!CDIR NOVERRCHK
450	DO ji = 1, jpi
451	zc0 = ze0(ji,jj,jk-1) * EXP( - fse3t_0(ji,jj,jk-1) * xsi0r )
452	zc1 = ze1(ji,jj,jk-1) * EXP( - fse3t_0(ji,jj,jk-1) * zekb(ji,jj) )
453	zc2 = ze2(ji,jj,jk-1) * EXP( - fse3t_0(ji,jj,jk-1) * zekg(ji,jj) )
454	zc3 = ze3(ji,jj,jk-1) * EXP( - fse3t_0(ji,jj,jk-1) * zekr(ji,jj) )
455	ze0(ji,jj,jk) = zc0
456	ze1(ji,jj,jk) = zc1
457	ze2(ji,jj,jk) = zc2
458	ze3(ji,jj,jk) = zc3
459	zea(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk)
460	END DO
461	END DO
462	END DO
463	!
464	DO jk = 1, nksr
465	etot3(:,:,jk) = r1_rau0_rcp * ( zea(:,:,jk) - zea(:,:,jk+1) )
466	END DO
467	etot3(:,:,nksr+1:jpk) = 0.e0 ! below 400m set to zero
468	ENDIF
469	ENDIF
470	!
471	ENDIF
472	! ! ---------------------------------- !
473	IF( ln_qsr_2bd ) THEN ! 2 bands light penetration !
474	! ! ---------------------------------- !
475	!
476	! ! level of light extinction
477	nksr = trc_oce_ext_lev( rn_si1, 1.e2 )
478	IF(lwp) THEN
479	WRITE(numout,*)
480	IF(lwp) WRITE(numout,*) ' level of light extinction = ', nksr, ' ref depth = ', gdepw_0(nksr+1), ' m'
481	ENDIF
482	!
483	IF( lk_vvl ) THEN ! variable volume
484	IF(lwp) WRITE(numout,*) ' key_vvl: light distribution will be computed at each time step'
485	ELSE ! constant volume: computes one for all
486	zz0 = rn_abs * r1_rau0_rcp
487	zz1 = ( 1. - rn_abs ) * r1_rau0_rcp
488	DO jk = 1, nksr !* solar heat absorbed at T-point computed once for all
489	DO jj = 1, jpj ! top 400 meters
490	DO ji = 1, jpi
491	zc0 = zz0 * EXP( -fsdepw(ji,jj,jk )xsi0r ) + zz1 EXP( -fsdepw(ji,jj,jk )*xsi1r )
492	zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)xsi0r ) + zz1 EXP( -fsdepw(ji,jj,jk+1)*xsi1r )
493	etot3(ji,jj,jk) = ( zc0 * tmask(ji,jj,jk) - zc1 * tmask(ji,jj,jk+1) )
494	END DO
495	END DO
496	END DO
497	etot3(:,:,nksr+1:jpk) = 0.e0 ! below 400m set to zero
498	!
499	ENDIF
500	ENDIF
501	! ! ===================================== !
502	ELSE ! No light penetration !
503	! ! ===================================== !
504	IF(lwp) THEN
505	WRITE(numout,*)
506	WRITE(numout,*) 'tra_qsr_init : NO solar flux penetration'
507	WRITE(numout,*) '~~~~~~~~~~~~'
508	ENDIF
509	ENDIF
510	!
511	CALL wrk_dealloc( jpi, jpj, zekb, zekg, zekr )
512	CALL wrk_dealloc( jpi, jpj, jpk, ze0, ze1, ze2, ze3, zea )
513	!
514	IF( nn_timing == 1 ) CALL timing_stop('tra_qsr_init')
515	!
516	END SUBROUTINE tra_qsr_init
517
518	!!======================================================================
519	END MODULE traqsr

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

Download in other formats: