source: LMDZ4/trunk/libf/phytherm/ocean_slab_mod.F90 @ 965

Last change on this file since 965 was 814, checked in by Laurent Fairhead, 17 years ago

Rajout de la physique utilisant les thermiques FH
LF
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 23.3 KB
RevLine 
[814]1!
2! $Header$
3!
4MODULE ocean_slab_mod
5!
6! This module is used for both surface ocean and sea-ice when using the slab ocean,
7! "ocean=slab".
8!
9  USE surface_data,     ONLY : tau_gl, calice, calsno
10  USE fonte_neige_mod,  ONLY : fonte_neige
11  USE calcul_fluxs_mod, ONLY : calcul_fluxs
12  USE dimphy
13 
14  IMPLICIT NONE
15
16  INTEGER, PRIVATE, SAVE                           :: lmt_pas, julien, idayvrai
17  !$OMP THREADPRIVATE(lmt_pas,julien,idayvrai)
18  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_seaice
19  !$OMP THREADPRIVATE(tmp_seaice)
20  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_tslab_loc
21  !$OMP THREADPRIVATE(tmp_tslab_loc)
22  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: slab_bils
23  !$OMP THREADPRIVATE(slab_bils)
24  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE , SAVE  :: lmt_bils
25  !$OMP THREADPRIVATE(lmt_bils)
26  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: tmp_pctsrf_slab
27  !$OMP THREADPRIVATE(tmp_pctsrf_slab)
28  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_tslab
29  !$OMP THREADPRIVATE(tmp_tslab)
30  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_radsol
31  !$OMP THREADPRIVATE(tmp_radsol)
32  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_flux_o, tmp_flux_g
33  !$OMP THREADPRIVATE(tmp_flux_o,tmp_flux_g)
34  LOGICAL, PRIVATE, SAVE                           :: check = .FALSE.
35  !$OMP THREADPRIVATE(check)
36
37CONTAINS
38!
39!****************************************************************************************
40!
41  SUBROUTINE ocean_slab_init(dtime, tslab_rst, seaice_rst, pctsrf_rst)
42
43    INCLUDE "indicesol.h"
44    INCLUDE "iniprint.h"
45
46! Input variables
47!****************************************************************************************
48    REAL, INTENT(IN)                         :: dtime
49! Variables read from restart file
50    REAL, DIMENSION(klon), INTENT(IN)        :: tslab_rst         
51    REAL, DIMENSION(klon), INTENT(IN)        :: seaice_rst
52    REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: pctsrf_rst
53
54
55! Local variables
56!****************************************************************************************
57    INTEGER                :: error
58    CHARACTER (len = 80)   :: abort_message
59    CHARACTER (len = 20)   :: modname = 'ocean_slab_intit'
60
61
62    WRITE(lunout,*) '************************'
63    WRITE(lunout,*) 'SLAB OCEAN est actif, prenez precautions !'
64    WRITE(lunout,*) '************************'   
65
66! Allocate variables initialize from restart fields
67    ALLOCATE(tmp_tslab(klon), stat = error)
68    IF (error /= 0) THEN
69       abort_message='Pb allocation tmp_tslab'
70       CALL abort_gcm(modname,abort_message,1)
71    ENDIF
72    tmp_tslab(:) = tslab_rst(:)
73
74    ALLOCATE(tmp_tslab_loc(klon), stat = error)
75    IF (error /= 0) THEN
76       abort_message='Pb allocation tmp_tslab_loc'
77       CALL abort_gcm(modname,abort_message,1)
78    ENDIF
79    tmp_tslab_loc(:) = tslab_rst(:)
80
81    ALLOCATE(tmp_seaice(klon), stat = error)
82    IF (error /= 0) THEN
83       abort_message='Pb allocation tmp_seaice'
84       CALL abort_gcm(modname,abort_message,1)
85    ENDIF
86    tmp_seaice(:) = seaice_rst(:)
87
88    ALLOCATE(tmp_pctsrf_slab(klon,nbsrf), stat = error)
89    IF (error /= 0) THEN
90       abort_message='Pb allocation tmp_pctsrf_slab'
91       CALL abort_gcm(modname,abort_message,1)
92    ENDIF
93    tmp_pctsrf_slab(:,:) = pctsrf_rst(:,:)
94   
95! Allocate some other variables internal in module mod_oceanslab
96    ALLOCATE(tmp_radsol(klon), stat = error)
97    IF (error /= 0) THEN
98       abort_message='Pb allocation tmp_radsol'
99       CALL abort_gcm(modname,abort_message,1)
100    ENDIF
101
102    ALLOCATE(tmp_flux_o(klon), stat = error)
103    IF (error /= 0) THEN
104       abort_message='Pb allocation tmp_flux_o'
105       CALL abort_gcm(modname,abort_message,1)
106    ENDIF
107   
108    ALLOCATE(tmp_flux_g(klon), stat = error)
109    IF (error /= 0) THEN
110       abort_message='Pb allocation tmp_flux_g'
111       CALL abort_gcm(modname,abort_message,1)
112    ENDIF
113
114! a mettre un slab_bils aussi en force !!!
115    ALLOCATE(slab_bils(klon), stat = error)
116    IF (error /= 0) THEN
117       abort_message='Pb allocation slab_bils'
118       CALL abort_gcm(modname,abort_message,1)
119    ENDIF
120    slab_bils(:) = 0.0   
121
122    ALLOCATE(lmt_bils(klon), stat = error)
123    IF (error /= 0) THEN
124       abort_message='Pb allocation lmt_bils'
125       CALL abort_gcm(modname,abort_message,1)
126    ENDIF
127
128
129! pour une lecture une fois par jour   
130    lmt_pas = NINT(86400./dtime * 1.0)
131
132  END SUBROUTINE ocean_slab_init
133!
134!****************************************************************************************
135!
136  SUBROUTINE ocean_slab_noice( &
137       dtime, knon, knindex, &
138       p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
139       petAcoef, peqAcoef, petBcoef, peqBcoef, &
140       ps, u1_lay, v1_lay, &
141       radsol, snow, qsurf, agesno, &
142       evap, fluxsens, fluxlat, &
143       tsurf_new, &
144       dflux_s, dflux_l, pctsrf_oce)
145
146    INCLUDE "indicesol.h"
147    INCLUDE "iniprint.h"
148
149! Input arguments
150!****************************************************************************************
151    INTEGER, INTENT(IN)                  :: knon
152    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
153    REAL, INTENT(IN)                     :: dtime
154    REAL, DIMENSION(klon), INTENT(IN)    :: p1lay
155    REAL, DIMENSION(klon), INTENT(IN)    :: tq_cdrag
156    REAL, DIMENSION(klon), INTENT(IN)    :: precip_rain, precip_snow
157    REAL, DIMENSION(klon), INTENT(IN)    :: temp_air, spechum
158    REAL, DIMENSION(klon), INTENT(IN)    :: petAcoef, peqAcoef
159    REAL, DIMENSION(klon), INTENT(IN)    :: petBcoef, peqBcoef
160    REAL, DIMENSION(klon), INTENT(IN)    :: ps
161    REAL, DIMENSION(klon), INTENT(IN)    :: u1_lay, v1_lay
162
163! In/Output arguments
164!****************************************************************************************
165    REAL, DIMENSION(klon), INTENT(INOUT) :: radsol
166    REAL, DIMENSION(klon), INTENT(INOUT) :: snow
167    REAL, DIMENSION(klon), INTENT(INOUT) :: agesno
168   
169! Output arguments
170!****************************************************************************************
171    REAL, DIMENSION(klon), INTENT(OUT)   :: qsurf
172    REAL, DIMENSION(klon), INTENT(OUT)   :: evap, fluxsens, fluxlat
173    REAL, DIMENSION(klon), INTENT(OUT)   :: tsurf_new
174    REAL, DIMENSION(klon), INTENT(OUT)   :: dflux_s, dflux_l     
175    REAL, DIMENSION(klon), INTENT(OUT)   :: pctsrf_oce
176
177! Local variables
178!****************************************************************************************
179    INTEGER                :: i
180    REAL, DIMENSION(klon)  :: cal, beta, dif_grnd
181    REAL, DIMENSION(klon)  :: alb_neig, tsurf_temp
182
183!****************************************************************************************
184    IF (check) WRITE(*,*)' Entering ocean_slab_noice'   
185
186    tsurf_new(1:knon) = tmp_tslab(knindex(1:knon))
187    pctsrf_oce(:)   = tmp_pctsrf_slab(:,is_oce)
188   
189    tsurf_temp(:) = tsurf_new(:)
190    cal = 0.
191    beta = 1.
192    dif_grnd = 0.
193    alb_neig(:) = 0.
194    agesno(:) = 0.
195   
196    CALL calcul_fluxs(knon, is_oce, dtime, &
197         tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
198         precip_rain, precip_snow, snow, qsurf,  &
199         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
200         petAcoef, peqAcoef, petBcoef, peqBcoef, &
201         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
202   
203    tmp_flux_o(:) = 0.0
204    tmp_radsol(:) = 0.0
205
206    DO i=1, knon
207       tmp_radsol(knindex(i))=radsol(i)
208       
209       IF (pctsrf_oce(knindex(i)) .GT. epsfra) &
210            tmp_flux_o(knindex(i)) = fluxsens(i) + fluxlat(i)
211    ENDDO
212   
213  END SUBROUTINE ocean_slab_noice
214!
215!****************************************************************************************
216!
217  SUBROUTINE ocean_slab_ice(   &
218       itime, dtime, jour, knon, knindex, &
219       debut, &
220       tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
221       petAcoef, peqAcoef, petBcoef, peqBcoef, &
222       ps, u1_lay, v1_lay, &
223       radsol, snow, qsurf, qsol, agesno, &
224       tsoil, &
225       alblw, evap, fluxsens, fluxlat, &
226       tsurf_new, alb_new, dflux_s, dflux_l, pctsrf_sic)
227
228    INCLUDE "indicesol.h"
229    INCLUDE "dimsoil.h"
230    INCLUDE "YOMCST.h"
231    INCLUDE "iniprint.h"
232    INCLUDE "clesphys.h"
233
234! Input arguments 
235!****************************************************************************************
236    INTEGER, INTENT(IN)                  :: itime, jour, knon
237    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
238    REAL, INTENT(IN)                     :: dtime
239    REAL, DIMENSION(klon), INTENT(IN)    :: tsurf
240    REAL, DIMENSION(klon), INTENT(IN)    :: p1lay
241    REAL, DIMENSION(klon), INTENT(IN)    :: tq_cdrag
242    REAL, DIMENSION(klon), INTENT(IN)    :: precip_rain, precip_snow
243    REAL, DIMENSION(klon), INTENT(IN)    :: temp_air, spechum
244    REAL, DIMENSION(klon), INTENT(IN)    :: petAcoef, peqAcoef
245    REAL, DIMENSION(klon), INTENT(IN)    :: petBcoef, peqBcoef
246    REAL, DIMENSION(klon), INTENT(IN)    :: ps
247    REAL, DIMENSION(klon), INTENT(IN)    :: u1_lay, v1_lay
248    LOGICAL, INTENT(IN)                  :: debut
249
250!In/Output arguments
251!****************************************************************************************
252    REAL, DIMENSION(klon), INTENT(INOUT)          :: radsol
253    REAL, DIMENSION(klon), INTENT(INOUT)          :: snow, qsol
254    REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
255    REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
256
257! Output arguments
258!****************************************************************************************
259    REAL, DIMENSION(klon), INTENT(OUT)            :: qsurf
260    REAL, DIMENSION(klon), INTENT(OUT)            :: alblw
261    REAL, DIMENSION(klon), INTENT(OUT)            :: evap, fluxsens, fluxlat
262    REAL, DIMENSION(klon), INTENT(OUT)            :: tsurf_new, alb_new
263    REAL, DIMENSION(klon), INTENT(OUT)            :: dflux_s, dflux_l     
264    REAL, DIMENSION(klon), INTENT(OUT)            :: pctsrf_sic
265
266! Local variables
267!****************************************************************************************
268    INTEGER                              :: i
269    REAL, DIMENSION(klon)                :: cal, beta, dif_grnd, capsol
270    REAL, DIMENSION(klon)                :: alb_neig, tsurf_temp
271    REAL, DIMENSION(klon)                :: soilcap, soilflux
272    REAL, DIMENSION(klon)                :: zfra
273    REAL, PARAMETER                      :: t_grnd=271.35
274    REAL                                 :: amn, amx
275    REAL, DIMENSION(klon)                :: tslab
276    REAL, DIMENSION(klon)                :: seaice ! glace de mer (kg/m2)
277    REAL, DIMENSION(klon,nbsrf)          :: pctsrf_new
278
279!****************************************************************************************
280
281    IF (check) WRITE(*,*)'Entering surface_seaice, knon=',knon
282
283! Initialization of output variables
284    alb_new(:) = 0.0
285
286!****************************************************************************************
287!
288!
289!****************************************************************************************
290    IF ( ok_slab_sicOBS) THEN   
291       ! glace de mer observee, lecture conditions limites
292       CALL interfoce_lim(itime, dtime, jour, &
293            knon, knindex, &
294            debut, &
295            tsurf_new, pctsrf_new)
296
297       tmp_pctsrf_slab(:,:) = pctsrf_new(:,:)
298       WRITE(lunout,*) 'jour lecture pctsrf_new sic =',jour
299
300    ELSE
301       pctsrf_new=tmp_pctsrf_slab
302    ENDIF
303
304    DO i = 1, knon
305       tsurf_new(i) = tsurf(i)
306       IF (pctsrf_new(knindex(i),is_sic) < EPSFRA) THEN
307          snow(i) = 0.0
308          tsurf_new(i) = RTT - 1.8
309          IF (soil_model) tsoil(i,:) = RTT -1.8
310       ENDIF
311    ENDDO
312   
313    CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, capsol, dif_grnd)
314   
315    IF (soil_model) THEN
316       CALL soil(dtime, is_sic, knon, snow, tsurf_new, tsoil, soilcap, soilflux)
317       cal(1:knon) = RCPD / soilcap(1:knon)
318       radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)
319    ELSE
320       dif_grnd = 1.0 / tau_gl
321       cal = RCPD * calice
322       WHERE (snow > 0.0) cal = RCPD * calsno
323    ENDIF
324    tsurf_temp = tsurf_new
325    beta = 1.0
326
327!****************************************************************************************
328!
329!
330!****************************************************************************************
331    CALL calcul_fluxs(knon, is_sic, dtime, &
332         tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
333         precip_rain, precip_snow, snow, qsurf,  &
334         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
335         petAcoef, peqAcoef, petBcoef, peqBcoef, &
336         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
337   
338    CALL fonte_neige( knon, is_sic, knindex, dtime, &
339         tsurf_temp, precip_rain, precip_snow, &
340         snow, qsol, tsurf_new, evap)
341
342!****************************************************************************************
343!     calcul albedo
344!
345!****************************************************************************************
346    CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
347    WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
348    zfra(1:knon) = MAX(0.0,MIN(1.0,snow(1:knon)/(snow(1:knon)+10.0)))
349    alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
350         0.6 * (1.0-zfra(1:knon))
351    !!      alb_new(1 : knon) = 0.6
352
353!
354!IM: flux entre l'ocean et la glace de mer pour le "slab" ocean
355    tmp_flux_g(:) = 0.0
356    DO i = 1, knon
357!
358!IM: faire dependre le coefficient de conduction de la glace de mer
359!    de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
360!    actuel correspond a 3m de glace de mer, cf. L.Li
361!
362       IF ((cal(i).GT.1.0e-15) .AND. (pctsrf_new(knindex(i),is_sic) .GT. epsfra)) THEN
363          tmp_flux_g(knindex(i))=(tsurf_new(i)-t_grnd) &
364               * dif_grnd(i) *RCPD/cal(i)
365       ENDIF
366!
367!IM: Attention: ne pas initialiser le tmp_radsol puisque c'est deja fait sur is_oce;
368!IM:            tmp_radsol doit etre le flux solaire qui arrive sur l'ocean
369!IM:            et non pas celui qui arrive sur la glace de mer
370    ENDDO
371   
372!
373! 2eme appel a interfoce pour le cumul et le passage des flux a l'ocean
374!
375
376    IF (check) THEN
377       amn=MIN(tmp_tslab(1),1000.)
378       amx=MAX(tmp_tslab(1),-1000.)
379       DO i=2, klon
380          amn=MIN(tmp_tslab(i),amn)
381          amx=MAX(tmp_tslab(i),amx)
382       ENDDO
383       
384       WRITE(lunout,*) ' debut avant interfoce_slab min max tmp_tslab',amn,amx
385    ENDIF !(check) THEN
386
387    tslab = tmp_tslab   
388
389    CALL interfoce_slab(klon, debut, itime, dtime, jour, &
390         tmp_radsol, tmp_flux_o, tmp_flux_g, tmp_pctsrf_slab, &
391         tslab, seaice, pctsrf_new)
392   
393    tmp_pctsrf_slab=pctsrf_new
394
395    DO i=1, knon
396       tmp_tslab(knindex(i))=tslab(knindex(i))
397    ENDDO
398     
399!!$    z0_new = 0.002
400!!$    z0_new = SQRT(z0_new**2+rugoro**2)
401    alblw(1:knon) = alb_new(1:knon)
402       
403!****************************************************************************************
404! Return the fraction of sea-ice
405! NB! jg : Peut-etre un probleme, faut-il prend aussi tmp_pctsrf_slab(:,is_oce)???
406!****************************************************************************************
407    pctsrf_sic(:) =  tmp_pctsrf_slab(:,is_sic)
408
409
410  END SUBROUTINE ocean_slab_ice
411!
412!****************************************************************************************
413!
414  SUBROUTINE interfoce_slab(klon, debut, itap, dtime, ijour, &
415       radsol, fluxo, fluxg, pctsrf, &
416       tslab, seaice, pctsrf_slab)
417!
418! Cette routine calcule la temperature d'un slab ocean, la glace de mer
419! et les pourcentages de la maille couverte par l'ocean libre et/ou
420! la glace de mer pour un "slab" ocean de 50m
421!
422! Conception: Laurent Li
423! Re-ecriture + adaptation LMDZ4: I. Musat
424!
425! input:
426!   klon         nombre total de points de grille
427!   debut        logical: 1er appel a la physique
428!   itap         numero du pas de temps
429!   dtime        pas de temps de la physique (en s)
430!   ijour        jour dans l'annee en cours
431!   radsol       rayonnement net au sol (LW + SW)
432!   fluxo        flux turbulent (sensible + latent) sur les mailles oceaniques
433!   fluxg        flux de conduction entre la surface de la glace de mer et l'ocean
434!   pctsrf       tableau des pourcentages de surface de chaque maille
435! output:
436!   tslab        temperature de l'ocean libre
437!   seaice       glace de mer (kg/m2)
438!   pctsrf_slab  "pourcentages" (valeurs entre 0. et 1.) surfaces issus du slab
439
440    INCLUDE "indicesol.h"
441    INCLUDE "YOMCST.h"
442    INCLUDE "iniprint.h"
443    INCLUDE "clesphys.h"
444
445! Input arguments
446!****************************************************************************************
447    INTEGER, INTENT(IN)                       :: klon
448    LOGICAL, INTENT(IN)                       :: debut    ! not used
449    INTEGER, INTENT(IN)                       :: itap
450    REAL, INTENT(IN)                          :: dtime       ! not used
451    INTEGER, INTENT(IN)                       :: ijour
452    REAL, DIMENSION(klon), INTENT(IN)         :: radsol
453    REAL, DIMENSION(klon), INTENT(IN)         :: fluxo
454    REAL, DIMENSION(klon), INTENT(IN)         :: fluxg
455    REAL, DIMENSION(klon, nbsrf), INTENT(IN)  :: pctsrf
456
457! Output arguments
458!****************************************************************************************
459    REAL, DIMENSION(klon), INTENT(OUT)        :: tslab
460    REAL, DIMENSION(klon), INTENT(OUT)        :: seaice ! glace de mer (kg/m2)
461    REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: pctsrf_slab
462
463! Local variables
464!****************************************************************************************
465    REAL                    :: amn, amx
466    REAL, PARAMETER         :: unjour=86400.
467    REAL, PARAMETER         :: cyang=50.0 * 4.228e+06 ! capacite calorifique volumetrique de l'eau J/(m2 K)
468    REAL, PARAMETER         :: cbing=0.334e+05        ! J/kg
469    REAL, DIMENSION(klon)   :: siceh !hauteur de la glace de mer (m)
470    INTEGER                 :: i
471    REAL                    :: zz, za, zb
472!
473!****************************************************************************************
474!
475    julien = MOD(ijour,360)
476
477    IF (check ) THEN
478       amn=MIN(tmp_tslab_loc(1),1000.)
479       amx=MAX(tmp_tslab_loc(1),-1000.)
480       DO i=2, klon
481          amn=MIN(tmp_tslab_loc(i),amn)
482          amx=MAX(tmp_tslab_loc(i),amx)
483       ENDDO
484
485       WRITE(lunout,*) ' debut min max tslab',amn,amx
486       WRITE(lunout,*) ' itap,lmt_pas unjour',itap,lmt_pas,unjour
487    ENDIF
488
489    pctsrf_slab(1:klon,1:nbsrf) = pctsrf(1:klon,1:nbsrf)
490!
491! lecture du bilan au sol lmt_bils issu d'une simulation forcee en debut de journee
492!
493    IF (MOD(itap,lmt_pas) .EQ. 1) THEN
494       ! 1er pas de temps de la journee
495       idayvrai = ijour
496       CALL condsurf(julien,idayvrai, lmt_bils)
497    ENDIF
498
499    DO i = 1, klon
500       IF((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
501            (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
502!
503! fabriquer de la glace si congelation atteinte:
504!         
505          IF (tmp_tslab_loc(i).LT.(RTT-1.8)) THEN
506             zz =  (RTT-1.8)-tmp_tslab_loc(i)
507             tmp_seaice(i) = tmp_seaice(i) + cyang/cbing * zz
508             seaice(i) = tmp_seaice(i)
509             tmp_tslab_loc(i) = RTT-1.8
510          ENDIF
511!
512! faire fondre de la glace si temperature est superieure a 0:
513!
514          IF ((tmp_tslab_loc(i).GT.RTT) .AND. (tmp_seaice(i).GT.0.0)) THEN
515             zz = cyang/cbing * (tmp_tslab_loc(i)-RTT)
516             zz = MIN(zz,tmp_seaice(i))
517             tmp_seaice(i) = tmp_seaice(i) - zz
518             seaice(i) = tmp_seaice(i)
519             tmp_tslab_loc(i) = tmp_tslab_loc(i) - zz*cbing/cyang
520          ENDIF
521!
522! limiter la glace de mer a 10 metres (10000 kg/m2)
523!
524          IF(tmp_seaice(i).GT.45.) THEN
525             tmp_seaice(i) = MIN(tmp_seaice(i),10000.0)
526          ELSE
527             tmp_seaice(i) = 0.
528          ENDIF
529          seaice(i) = tmp_seaice(i)
530          siceh(i)=tmp_seaice(i)/1000. !en metres
531!
532! determiner la nature du sol (glace de mer ou ocean libre):
533!
534! on fait dependre la fraction de seaice "pctsrf(i,is_sic)"
535! de l'epaisseur de seaice :
536! pctsrf(i,is_sic)=1. si l'epaisseur de la glace de mer est >= a 20cm
537! et pctsrf(i,is_sic) croit lineairement avec seaice de 0. a 20cm d'epaisseur
538!
539
540          IF(.NOT.ok_slab_sicOBS) THEN
541             pctsrf_slab(i,is_sic)=MIN(siceh(i)/0.20, &
542                  1.-(pctsrf_slab(i,is_ter)+pctsrf_slab(i,is_lic)))
543             pctsrf_slab(i,is_oce)=1.0 - &
544                  (pctsrf_slab(i,is_ter)+pctsrf_slab(i,is_lic)+pctsrf_slab(i,is_sic))
545          ELSE
546             IF (i.EQ.1) WRITE(lunout,*) 'cas ok_slab_sicOBS TRUE : passe sur la modif.'
547          ENDIF !(.NOT.ok_slab_sicOBS) then
548       ENDIF !pctsrf
549    ENDDO
550!
551! Calculer le bilan du flux de chaleur au sol :
552!
553    DO i = 1, klon
554       za = radsol(i) + fluxo(i)
555       zb = fluxg(i)
556       IF((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
557            (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
558          slab_bils(i)=slab_bils(i)+(za*pctsrf_slab(i,is_oce) &
559               +zb*pctsrf_slab(i,is_sic))/ FLOAT(lmt_pas)
560       ENDIF
561    ENDDO !klon
562!
563! calcul tslab
564!
565    IF (MOD(itap,lmt_pas).EQ.0) THEN !fin de journee
566!
567! calcul tslab
568       amn=MIN(tmp_tslab_loc(1),1000.)
569       amx=MAX(tmp_tslab_loc(1),-1000.)
570       DO i = 1, klon
571          IF ((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
572               (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
573             tmp_tslab_loc(i) = tmp_tslab_loc(i) + &
574                  (slab_bils(i)-lmt_bils(i)) &
575                  /cyang*unjour
576
577! on remet l'accumulation a 0
578             slab_bils(i) = 0.
579          ENDIF !pctsrf
580!
581          IF (check) THEN
582             IF(i.EQ.1) THEN 
583                WRITE(lunout,*) 'i tmp_tslab_loc MOD(itap,lmt_pas).EQ.0 sicINTER',i,itap, &
584                     tmp_tslab_loc(i), tmp_seaice(i)
585             ENDIF
586             
587             amn=MIN(tmp_tslab_loc(i),amn)
588             amx=MAX(tmp_tslab_loc(i),amx)
589          ENDIF
590       ENDDO !klon
591    ENDIF !(MOD(itap,lmt_pas).EQ.0) THEN
592
593    IF ( check ) THEN
594       WRITE(lunout,*) 'fin min max tslab',amn,amx
595    ENDIF
596
597    tslab  = tmp_tslab_loc
598    seaice = tmp_seaice
599
600  END SUBROUTINE interfoce_slab
601!
602!**************************************************************************************** 
603!
604  SUBROUTINE ocean_slab_final(tslab_rst, seaice_rst)
605
606! This subroutine will send to phyredem the variables concerning the slab
607! ocean that should be written to restart file.
608
609!****************************************************************************************
610
611    REAL, DIMENSION(klon), INTENT(OUT) :: tslab_rst
612    REAL, DIMENSION(klon), INTENT(OUT) :: seaice_rst
613
614!****************************************************************************************
615! Set the output variables
616    tslab_rst(:)  = tmp_tslab(:)
617!    tslab_rst(:)  = tmp_tslab_loc(:)
618    seaice_rst(:) = tmp_seaice(:)
619
620! Deallocation of all variables in module
621    DEALLOCATE(tmp_tslab, tmp_tslab_loc, tmp_pctsrf_slab)
622    DEALLOCATE(tmp_seaice, tmp_radsol, tmp_flux_o, tmp_flux_g)
623    DEALLOCATE(slab_bils, lmt_bils)
624
625  END SUBROUTINE ocean_slab_final
626!
627!****************************************************************************************
628!
629  SUBROUTINE ocean_slab_get_vars(tslab_loc, seaice_loc, flux_o_loc, flux_g_loc)
630! "Get some variables from module ocean_slab_mod"
631! This subroutine prints variables to a external routine
632
633    REAL, DIMENSION(klon), INTENT(OUT) :: tslab_loc
634    REAL, DIMENSION(klon), INTENT(OUT) :: seaice_loc
635    REAL, DIMENSION(klon), INTENT(OUT) :: flux_o_loc
636    REAL, DIMENSION(klon), INTENT(OUT) :: flux_g_loc
637
638! Set the output variables
639    tslab_loc(:)  = tmp_tslab(:)
640    seaice_loc(:) = tmp_seaice(:)
641    flux_o_loc(:) = tmp_flux_o(:)
642    flux_g_loc(:) = tmp_flux_g(:)
643
644  END SUBROUTINE ocean_slab_get_vars
645!
646!****************************************************************************************
647!
648END MODULE ocean_slab_mod
Note: See TracBrowser for help on using the repository browser.