source: trunk/WRF.COMMON/INTERFACES/dynphy_wrf_titan_lmd/update_inputs_physiq_mod.F

Last change on this file was 2743, checked in by aslmd, 3 years ago

made the vertical discretization for Titan LES similar to that of Mars (simple linear law for znw and refinement close to the surface). no more need to prepare a file levels, just put ztop and number of levels in namelist. a file levels is saved at the end of the znw computation so that it can be read in update_inputs_physiq_mod in the case of Titan.

added MESOSCALE precompiling flags in call_profilgases to set constant methane abundance. therefore no need for profile.def file and the number of vertical levels is set by namelist -- no need to set the same number in both levels and profile.def

outputing ustar from vdifc

ending hardcoding of ptop in update_inputs_physiq_mod for Titan and generic LES/mesoscale. streamlining grod%ptop through module_lmd_driver as ptopwrf within variables_mod module.

adding flux outputs to check surface and atmosphere energy budgets
File size: 21.8 KB
Line 
1MODULE update_inputs_physiq_mod
2
3CONTAINS
4
5!SUBROUTINE update_inputs_physiq_time
6!SUBROUTINE update_inputs_physiq_tracers
7!SUBROUTINE update_inputs_physiq_constants
8!SUBROUTINE update_inputs_physiq_geom
9!SUBROUTINE update_inputs_physiq_surf
10!SUBROUTINE update_inputs_physiq_soil
11!SUBROUTINE update_inputs_physiq_turb
12!SUBROUTINE update_inputs_physiq_rad
13!SUBROUTINE update_inputs_physiq_slope
14
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
17SUBROUTINE update_inputs_physiq_time(&
18            JULYR,JULDAY,GMT,&
19            elaps,&
20            lct_input,lon_input,ls_input,&
21            MY)
22
23  USE variables_mod, only: JD_cur,JH_cur_split,phour_ini
24  !use callkeys_mod, only : tlocked
25  !! JD_cur <> pday ! Julian day
26  !! JH_cur_split <> ptime ! Julian hour (fraction of day)
27
28  INTEGER, INTENT(IN) :: JULDAY, JULYR
29  REAL, INTENT(IN) :: GMT,elaps,lon_input,ls_input,lct_input
30  REAL,INTENT(OUT) :: MY
31
32  IF (JULYR .ne. 9999) THEN
33    !
34    ! specified
35    !
36    JH_cur_split = (GMT + elaps/57420.) !! universal time (0<JH_cur_split<1): JH_cur_split=0.5 at 12:00 UT
37    JH_cur_split = MODULO(JH_cur_split,24.)   !! the two arguments of MODULO must be of the same type
38    JH_cur_split = JH_cur_split / 24.
39    JD_cur = (JULDAY - 1 + INT((57420.*GMT+elaps)/1378080.))
40    JD_cur = MODULO(int(JD_cur),2)
41    MY = (JULYR-2000) + (1378080*(JULDAY - 1)+57420.*GMT+elaps)/2756160.
42    MY = INT(MY)
43  ELSE
44    !
45    ! idealized
46    !
47    JH_cur_split = lct_input - lon_input / 15. + elaps/57420.0
48    JH_cur_split = MODULO(JH_cur_split,24.)
49    JH_cur_split = JH_cur_split / 24.
50    JD_cur = INT((57420.0*(lct_input - lon_input / 15.) + elaps)/1378080) !! ls2sol???
51    JD_cur = MODULO(int(pday),365) !! 365 to be changed
52    MY = 2024
53    !day_ini = floor(ls2sol(ls_input)) !! JD_cur at firstcall is day_ini
54  ENDIF
55  print *,'** Titan ** TIME IS', JD_cur, JH_cur_split*24.
56
57END SUBROUTINE update_inputs_physiq_time
58
59!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
60!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
61SUBROUTINE update_inputs_physiq_tracers(nq,MARS_MODE)
62
63  use tracer_h, only: noms,nqtot_p
64  !use phys_state_var_mod , only : qsurf
65  !USE dimphy, only : klon,klev
66  INTEGER, INTENT(IN) :: nq,MARS_MODE
67  INTEGER :: i,k
68!  CHARACTER(len=20), INTENT(INOUT) :: tname(nq) ! tracer names
69  logical :: end_of_file
70
71  !! TRACERS
72  IF (.not.ALLOCATED(noms)) ALLOCATE(noms(nq)) !! est fait dans initracer normalement
73     
74                 !! tableau dans tracer_mod.F90
75  print*,'nq',nq
76  print*,'MARS_MODE',MARS_MODE
77  nqtot=nq
78  noms(:)="zolbxs"
79  !noms(:)=tname(:)
80  print*,'noms',noms
81  !!---------------------!!
82  !! OUTPUT FOR CHECKING !!
83  !!---------------------!!
84
85END SUBROUTINE update_inputs_physiq_tracers
86
87!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
88!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
89SUBROUTINE update_inputs_physiq_constants
90
91   !USE module_model_constants
92   !use comcstfi_h, only: omeg,mugaz
93   !use planete_h,  only: year_day,periheli,aphelie, &
94   !                       peri_day,obliquit,emin_turb, &
95   !                       lmixmin
96   !use planete_mod, only: year_day, periastr, apoastr, peri_day,&
97   !                    obliquit, z0, lmixmin, emin_turb
98   !use surfdat_h,  only: emissiv,iceradius, &
99   !                      emisice,dtemisice
100   !                       z0_default
101   !use comsoil_h, only: volcapa
102   !use comcstfi_mod, only: rad, omeg, g, mugaz, rcp, cpp, r
103   !! comcstfi_h
104   !use phys_state_var_mod, only :cloudfrac,totcloudfrac,hice,rnat,pctsrf_sic,tsea_ice
105   !use iniorbit
106   !use time_phylmdz_mod, only: dtphys, daysec,day_ini
107 
108 
109     !open(17,file='controle.txt',form='formatted',status='old')
110     !rewind(17)
111     !read(17,*)
112     !read(17,*)
113     !read(17,*) day_ini !(tab0+3)
114     !read(17,*)
115     !read(17,*) !tab0+5)
116     !read(17,*) omeg !(tab0+6)
117     !read(17,*) !(tab0+7)
118     !read(17,*)  !(tab0+8)
119     !read(17,*)  !(tab0+9)
120     !read(17,*) daysec
121     !read(17,*) dtphys !tab0+11)
122     !read(17,*)
123     !read(17,*)
124     !read(17,*) year_day !(tab0+14)
125     !read(17,*) periastr !tab0+15)
126     !read(17,*) apoastr !tab0+16)
127     !read(17,*) peri_day !tab0+17)
128     !read(17,*) obliquit !tab0+18)
129     !read(17,*) z0
130     !read(17,*)
131     !read(17,*)
132     !read(17,*)
133     !read(17,*)
134     !read(17,*) emisice(1)
135     !read(17,*) emisice(2)
136     !read(17,*) emissiv
137     !read(17,*)
138     !read(17,*)
139     !read(17,*)
140     !read(17,*)
141     !read(17,*) iceradius(1)
142     !read(17,*) iceradius(2)
143     !read(17,*) dtemisice(1)
144     !read(17,*) dtemisice(2)
145     !close(17)
146     !cpp=(8.314511/(mugaz/1000.0))/rcp
147     !print*,'cpp',cpp
148     !print*,'g',g
149
150     !emissiv(:)=EMIS
151     !cloudfrac(:,:)=0.5
152     !totcloudfrac(:)=0.5
153     !hice(:)=0.
154     !rnat(:)=0.
155     !pctsrf_sic(:)=0.
156     !tsea_ice(:)=0.
157     !qsurf(:,:) = 0.
158     !print*,'iceradius',iceradius,'dtemisice',dtemisice
159     !print*,'apoastr,periastr,year_day,peri_day,obliq',apoastr,periastr,year_day,peri_day,obliquit
160     !print*,'emissiv',emissiv
161 
162     ! SUBROUTINE iniorbit(apoastr,periastr,year_day,peri_day,obliq)
163
164END SUBROUTINE update_inputs_physiq_constants
165
166!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
167!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
168SUBROUTINE update_inputs_physiq_geom( &
169            ims,ime,jms,jme,&
170            ips,ipe,jps,jpe,&
171            JULYR,ngrid,nlayer,&
172            DX,DY,MSFT,&
173            lat_input, lon_input,&
174            XLAT,XLONG)
175
176   ! in WRF (share)
177   USE module_model_constants, only: DEGRAD,p0
178   ! in LMD (phymars)
179   !use comgeomfi_h, only: ini_fillgeom
180   ! in LMD (phy_common)
181   USE mod_grid_phy_lmdz, ONLY: init_grid_phy_lmdz
182   USE geometry_mod, ONLY: latitude,latitude_deg,&
183                           longitude,longitude_deg,&
184                           cell_area
185   use comdiurn_h, only: sinlat, coslat, sinlon, coslon
186   !use planetwide_mod, only: planetwide_sumval
187   use comgeomfi_h, only: totarea, totarea_planet
188   use planete_mod, only: ini_planete_mod
189   ! in LMD-WRF interface
190   USE variables_mod, only: ptopwrf
191
192   INTEGER, INTENT(IN) :: ims,ime,jms,jme
193   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe,JULYR,ngrid,nlayer
194   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: &
195     MSFT,XLAT,XLONG
196   REAL, INTENT(IN) :: dx,dy
197   REAL, INTENT(IN) :: lat_input, lon_input
198   INTEGER :: i,j,subs,ig,k
199   REAL,DIMENSION(ngrid) :: plon, plat, parea
200   REAL, DIMENSION(nlayer+1) :: znw
201   REAL*8, DIMENSION(nlayer+1) :: apdyn,bpdyn
202   REAL*8 :: PP0
203   DO j = jps,jpe
204   DO i = ips,ipe
205
206    !-----------------------------------!
207    ! 1D subscript for physics "cursor" !
208    !-----------------------------------!
209    subs = (j-jps)*(ipe-ips+1)+(i-ips+1)
210
211    !----------------------------------------!
212    ! Surface of each part of the grid (m^2) ! 
213    !----------------------------------------!
214    !parea(subs) = dx*dy                           !! 1. idealized cases - computational grid
215    parea(subs) = (dx/msft(i,j))*(dy/msft(i,j))    !! 2. WRF map scale factors - assume that msfx=msfy (msf=covariance)
216    !parea(subs)=dx*dy/msfu(i,j)                   !! 3. special for Mercator GCM-like simulations
217
218    !---------------------------------------------!
219    ! Mass-point latitude and longitude (radians) !
220    !---------------------------------------------!
221    IF (JULYR .ne. 9999) THEN
222     plat(subs) = XLAT(i,j)*DEGRAD
223     plon(subs) = XLONG(i,j)*DEGRAD
224    ELSE
225     !!! IDEALIZED CASE
226     IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION lat: ',lat_input
227     IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION lon: ',lon_input
228     plat(subs) = lat_input*DEGRAD
229     plon(subs) = lon_input*DEGRAD
230    ENDIF
231
232   ENDDO
233   ENDDO
234
235   !! FILL GEOMETRICAL ARRAYS !!
236   !call ini_fillgeom(ngrid,plat,plon,parea)
237
238   !!! ----------------------------------------------------------
239   !!! --- initializing geometry in phy_common
240   !!! --- (this is quite planet-independent)
241   !!! ----------------------------------------------------------
242   ! initialize mod_grid_phy_lmdz
243   CALL init_grid_phy_lmdz(1,1,ipe-ips+1,jpe-jps+1,nlayer)
244   ! fill in geometry_mod variables
245   ! ... copy over local grid longitudes and latitudes
246   ! ... partly what is done in init_geometry
247   IF(.not.ALLOCATED(longitude)) ALLOCATE(longitude(ngrid))
248   IF(.not.ALLOCATED(longitude_deg)) ALLOCATE(longitude_deg(ngrid))
249   IF(.not.ALLOCATED(latitude)) ALLOCATE(latitude(ngrid))
250   IF(.not.ALLOCATED(latitude_deg)) ALLOCATE(latitude_deg(ngrid))
251   IF(.not.ALLOCATED(cell_area)) ALLOCATE(cell_area(ngrid))
252   longitude(:) = plon(:)
253   latitude(:) = plat(:)
254   longitude_deg(:) = plon(:)/DEGRAD
255   latitude_deg(:) = plat(:)/DEGRAD
256   cell_area(:) = parea(:)
257   !call planetwide_sumval(parea,totarea_planet)
258   !print*,'parea',parea(1)
259   !totarea=SSUM(ngrid,parea,1)
260   totarea=ngrid*parea(1)
261   !totarea_planet=SSUM(ngrid,parea,1)
262   totarea_planet=ngrid*parea(1)
263
264   IF (.not.ALLOCATED(sinlat)) ALLOCATE(sinlat(ngrid))
265   IF (.not.ALLOCATED(coslat)) ALLOCATE(coslat(ngrid))
266   IF (.not.ALLOCATED(sinlon)) ALLOCATE(sinlon(ngrid))
267   IF (.not.ALLOCATED(coslon)) ALLOCATE(coslon(ngrid))
268   DO ig=1,ngrid
269     sinlat(ig)=sin(plat(ig))
270     coslat(ig)=cos(plat(ig))
271     sinlon(ig)=sin(plon(ig))
272     coslon(ig)=cos(plon(ig))
273   ENDDO
274   
275   open(unit=12,file='levels',form='formatted',status='old')
276   rewind(12)
277   DO k=1, nlayer
278   read(12,*) znw(k)
279   !write(6,*) 'read level ', k,grid%znw(k)
280   ENDDO
281   close(12)
282   
283   apdyn=ptopwrf*(1-znw)
284   bpdyn=znw
285   PP0=p0
286   CALL ini_planete_mod(nlayer,PP0,apdyn,bpdyn)
287   !!! ----------------------------------------------------------
288
289END SUBROUTINE update_inputs_physiq_geom
290
291!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
292!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
293SUBROUTINE update_inputs_physiq_surf( &
294            ims,ime,jms,jme,&
295            ips,ipe,jps,jpe,&
296            JULYR,MARS_MODE,&
297            M_ALBEDO,CST_AL,&
298            M_TSURF,M_EMISS,M_CO2ICE,&
299            M_GW,M_Z0,CST_Z0,&
300            M_H2OICE,&
301            phisfi_val)
302
303   use surfdat_h, only: phisfi, albedodat,  &
304                       zmea, zstd, zsig, zgam, zthe
305   use planete_mod, only: z0
306   use phys_state_var_mod, only : tsurf, emis, qsurf
307
308   INTEGER, INTENT(IN) :: ims,ime,jms,jme
309   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe,JULYR,MARS_MODE
310   INTEGER :: i,j,subs,nlast,iq
311   REAL, INTENT(IN  ) :: CST_AL, phisfi_val, CST_Z0
312   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: &
313     M_ALBEDO,M_TSURF,M_EMISS,M_CO2ICE,M_H2OICE,M_Z0
314   REAL, DIMENSION( ims:ime, 5, jms:jme ), INTENT(IN   )  :: M_GW 
315
316   !! TITAN: constant surface roughness
317   !! see below commented part if varying is needed
318   z0 = CST_Z0
319   IF (z0 == 0.) THEN
320       PRINT *, 'WELL, z0 is 0, this is no good. Setting to default value from Tokano 2006'
321       z0 = 0.005 ! (m) value for Huygens landing site from Tokano 2006
322   ENDIF
323
324   !print*,'ALLOCATED(phisfi)',ALLOCATED(phisfi)
325   !print*,'size phisfi',size(phisfi)
326   DO j = jps,jpe
327   DO i = ips,ipe
328
329     !-----------------------------------!
330     ! 1D subscript for physics "cursor" !
331     !-----------------------------------!
332     subs = (j-jps)*(ipe-ips+1)+(i-ips+1)
333
334     !---------------------!
335     ! Ground geopotential !
336     !---------------------!
337     phisfi(subs) = phisfi_val
338!     print*,'size phisfi',size(phisfi)
339     !print*,'phisfi',phisfi(subs)
340     !---------------!
341     ! Ground albedo !
342     !---------------!
343     IF (JULYR .ne. 9999) THEN
344      IF (CST_AL == 0) THEN
345       albedodat(subs)=M_ALBEDO(i,j)
346      ELSE
347       albedodat(subs)=CST_AL
348       IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** SET CONSTANT ALBEDO ', albedodat
349      ENDIF
350     ELSE
351      IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION albedo: ', CST_AL
352      albedodat(subs)=CST_AL
353     ENDIF
354
355     !-----------------------------------------!
356     ! Gravity wave parametrization            !
357     ! NB: usually 0 in mesoscale applications !
358     !-----------------------------------------!
359     zmea(subs)=0.
360     zstd(subs)=0.
361     zsig(subs)=0.
362     zgam(subs)=0.
363     zthe(subs)=0.
364
365     !----------------------------!
366     ! Variable surface roughness !
367     !----------------------------!
368     !IF (JULYR .ne. 9999) THEN
369     ! IF (CST_Z0 == 0) THEN
370     !  z0(subs) = M_Z0(i,j)
371     ! ELSE
372     !  z0(subs) = CST_Z0
373     !  IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** SET CONSTANT SURF ROUGHNESS (m) ',CST_Z0
374     ! ENDIF
375     !ELSE
376     ! IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION z0 (m) ', CST_Z0
377     ! z0(subs)=CST_Z0
378     !ENDIF
379     !!!!! ADDITIONAL SECURITY. THIS MIGHT HAPPEN WITH OLD INIT FILES.
380     !IF (z0(subs) == 0.) THEN
381     !  IF ( (i == ips) .AND. (j == jps) ) PRINT *, 'WELL, z0 is 0, this is no good. Setting to old defaults value 0.01 m'
382     !  z0(subs) = 0.01
383     !ENDIF
384     !!!!! ADDITIONAL SECURITY. INTERP+SMOOTH IN GEOGRID MIGHT YIELD NEGATIVE Z0 !!!
385     !IF (z0(subs) < 0.) THEN
386     !  PRINT *, 'WELL, z0 is NEGATIVE, this is impossible. better stop here.'
387     !  PRINT *, 'advice --> correct interpolation / smoothing of z0 in WPS'
388     !  PRINT *, '           -- or check the constant value set in namelist.input'
389     !  STOP
390     !ENDIF
391
392     !-----------------------------------------------!
393     ! Ground temperature, emissivity, CO2 ice cover !
394     !-----------------------------------------------!
395     tsurf(subs) = M_TSURF(i,j)
396     emis(subs) = M_EMISS(i,j)     
397     !enddo
398     !-------------------!
399     ! Tracer at surface !
400     !-------------------!
401     qsurf(subs,:)=0. ! default case
402     SELECT CASE (MARS_MODE)
403       CASE(1)
404       qsurf(subs,2)=M_H2OICE(i,j)    !! logique avec noms(2) = 'h2o_ice' defini ci-dessus
405                                      !! ----- retrocompatible ancienne physique
406                                      !! ----- [H2O ice is last tracer in qsurf in LMD physics]
407       CASE(2)   
408       qsurf(subs,1)=0.                !! not coupled with lifting for the moment [non remobilise]
409       !CASE(3)
410       !qsurf(subs,1)=q_prof(1,1)                !!! temporaire, a definir       
411       !qsurf(subs,2)=q_prof(1,2)     
412       CASE(11)
413       qsurf(subs,2)=M_H2OICE(i,j)    !! logique avec noms(2) = 'h2o_ice' defini ci-dessus
414       qsurf(subs,3)=0.                !! not coupled with lifting for the moment [non remobilise]
415       CASE(12)
416       qsurf(subs,2)=M_H2OICE(i,j)    !! logique avec noms(2) = 'h2o_ice' defini ci-dessus
417       qsurf(subs,3)=0.                !! not coupled with lifting for the moment [non remobilise]
418     END SELECT
419
420   ENDDO
421   ENDDO
422 
423   !!---------------------!!
424   !! OUTPUT FOR CHECKING !!
425   !!---------------------!!
426   nlast = (ipe-ips+1)*(jpe-jps+1)
427   print*,"check: phisfi",phisfi(1),phisfi(nlast)
428   print*,"check: albedodat",albedodat(1),albedodat(nlast)
429   print*,"check: zmea",zmea(1),zmea(nlast)
430   print*,"check: zstd",zstd(1),zstd(nlast)
431   print*,"check: zsig",zsig(1),zsig(nlast)
432   print*,"check: zgam",zgam(1),zgam(nlast)
433   print*,"check: zthe",zthe(1),zthe(nlast)
434   print*,"check: z0",z0
435   print*,"check: tsurf",tsurf(1),tsurf(nlast)
436   print*,"check: emis",emis(1),emis(nlast)
437   print*,"check: qsurf",qsurf(1,:),qsurf(nlast,:)
438
439END SUBROUTINE update_inputs_physiq_surf
440
441!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
442!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
443SUBROUTINE update_inputs_physiq_soil( &
444            ims,ime,jms,jme,&
445            ips,ipe,jps,jpe,&
446            JULYR,nsoil,&
447            M_TI,CST_TI,&
448            M_ISOIL,M_DSOIL,&
449            M_TSOIL,M_TSURF)
450
451   use comsoil_h, only: inertiedat,mlayer,layer,volcapa
452   use phys_state_var_mod, only: tsoil
453
454   INTEGER, INTENT(IN) :: ims,ime,jms,jme
455   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe,JULYR,nsoil
456   INTEGER :: i,j,subs,nlast
457   REAL, INTENT(IN  ) :: CST_TI
458   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: &
459     M_TI, M_TSURF
460   REAL, DIMENSION( ims:ime, nsoil, jms:jme ), INTENT(IN)  :: &
461     M_TSOIL, M_ISOIL, M_DSOIL
462   REAL :: inertiedat_val
463   REAL :: lay1,alpha
464
465   DO j = jps,jpe
466   DO i = ips,ipe
467
468     !-----------------------------------!
469     ! 1D subscript for physics "cursor" !
470     !-----------------------------------!
471     subs = (j-jps)*(ipe-ips+1)+(i-ips+1)
472
473     !-----------------!
474     ! Thermal Inertia !
475     !-----------------!
476     IF (JULYR .ne. 9999) THEN
477      IF (CST_TI == 0) THEN
478       inertiedat_val=M_TI(i,j)
479      ELSE
480       inertiedat_val=CST_TI
481       IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** SET CONSTANT THERMAL INERTIA ', inertiedat_val
482      ENDIF
483     ELSE
484      IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION inertia: ',CST_TI
485      inertiedat_val=CST_TI
486     ENDIF
487     !inertiedat(subs) = inertiedat_val
488     !--pb de dimensions???!!???
489     IF (JULYR .ne. 9999) THEN
490       inertiedat(subs,:)=M_ISOIL(i,:,j) !! verifier que cest bien hires TI en surface
491       mlayer(0:nsoil-1)=M_DSOIL(i,:,j)
492     ELSE
493        IF ( nsoil .lt. 18 ) THEN
494            PRINT *,'** Mars ** WRONG NUMBER OF SOIL LAYERS. SHOULD BE 18 AND IT IS ',nsoil
495            STOP
496        ENDIF
497        IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION isoil and dsoil standard'
498        do k=1,nsoil
499         inertiedat(subs,k) = inertiedat_val
500         !mlayer(k-1) = sqrt(887.75/3.14)*((2.**(k-0.5))-1.) * inertiedat_val / wvolcapa    ! old setting
501         mlayer(k-1) = 2.E-4 * (2.**(k-0.5-1.))                                            ! new gcm settings
502        enddo
503     ENDIF
504     IF ( (i == ips) .AND. (j == jps) ) &
505         PRINT *,'** Mars ** TI and depth profiles are',inertiedat(subs,:)!,mlayer(0:nsoil-1)
506
507     !!!!!!!!!!!!!!!!! DONE in soil_setting.F
508     ! 1.5 Build layer(); following the same law as mlayer()
509     ! Assuming layer distribution follows mid-layer law:
510     ! layer(k)=lay1*alpha**(k-1)
511     lay1=sqrt(mlayer(0)*mlayer(1))
512     alpha=mlayer(1)/mlayer(0)
513     do k=1,nsoil
514       layer(k)=lay1*(alpha**(k-1))
515     enddo
516
517     !------------------------!
518     ! Deep soil temperatures !
519     !------------------------!
520     IF (M_TSOIL(i,1,j) .gt. 0. .and. JULYR .ne. 9999) THEN
521       tsoil(subs,:)=M_TSOIL(i,:,j)
522     ELSE
523       IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** no tsoil. set it to tsurf.'
524       do k=1,nsoil
525        !print*,'M_TSURF(i,j)',M_TSURF(1,:)
526        !print*,'size M_TSURF',size(M_TSURF)
527        !print*,'size tsoil',size(tsoil)
528        tsoil(subs,k) = M_TSURF(i,j)
529        !print*,'tsoil(subs,k)',tsoil(subs,k)
530       enddo
531     ENDIF
532
533   ENDDO
534   ENDDO
535
536   volcapa=1.e6   
537
538   print*,'zolbxs'
539   !!---------------------!!
540   !! OUTPUT FOR CHECKING !!
541   !!---------------------!!
542   nlast = (ipe-ips+1)*(jpe-jps+1)
543   print*,"check: inertiedat",inertiedat(1,:),inertiedat(nlast,:)
544   print*,"check: mlayer",mlayer(:)
545   print*,"check: layer",layer(:)
546   print*,"check: tsoil",tsoil(1,:),tsoil(nlast,:)
547
548END SUBROUTINE update_inputs_physiq_soil
549
550!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
551!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
552SUBROUTINE update_inputs_physiq_turb( &
553            ims,ime,jms,jme,kms,kme,&
554            ips,ipe,jps,jpe,&
555            RESTART,isles,&
556            M_Q2,M_WSTAR)
557
558   use turb_mod, only: q2,wstar,turb_resolved
559   !use phys_state_var_mod, only : q2,
560
561   INTEGER, INTENT(IN) :: ims,ime,jms,jme,kms,kme
562   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe
563   INTEGER :: i,j,subs,nlast     
564   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: M_WSTAR
565   REAL, DIMENSION( ims:ime, kms:kme+1, jms:jme ), INTENT(IN) :: M_Q2
566   LOGICAL, INTENT(IN ) :: RESTART,isles
567
568   !! to know if this is turbulence-resolving run or not
569
570   turb_resolved = isles
571   print*,'isles',isles
572
573   DO j = jps,jpe
574   DO i = ips,ipe
575
576     !-----------------------------------!
577     ! 1D subscript for physics "cursor" !
578     !-----------------------------------!
579     subs = (j-jps)*(ipe-ips+1)+(i-ips+1)
580
581     !PBL wind variance
582     IF (.not. restart) THEN
583      q2(subs,:) = 1.E-6     
584      wstar(subs)=0.
585     ELSE
586      q2(subs,:)=M_Q2(i,:,j)!
587      !q2(subs,:) = 1.e-3
588      wstar(subs)=M_WSTAR(i,j)
589     ENDIF
590
591   ENDDO
592   ENDDO
593 
594   !!---------------------!!
595   !! OUTPUT FOR CHECKING !!
596   !!---------------------!!
597   nlast = (ipe-ips+1)*(jpe-jps+1)
598   print*,"check: q2",q2(1,1)!,q2(nlast,:)
599   print*,"check: wstar",wstar(1),wstar(nlast)
600
601END SUBROUTINE update_inputs_physiq_turb
602
603!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
604!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
605SUBROUTINE update_inputs_physiq_rad( &
606            ims,ime,jms,jme,&
607            ips,ipe,jps,jpe,&
608            RESTART,&
609            M_FLUXRAD)
610
611   !use dimradmars_mod, only: fluxrad
612   use phys_state_var_mod, only : fluxrad
613
614   INTEGER, INTENT(IN) :: ims,ime,jms,jme
615   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe
616   INTEGER :: i,j,subs,nlast     
617   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: M_FLUXRAD
618   LOGICAL, INTENT(IN ) :: RESTART
619
620   DO j = jps,jpe
621   DO i = ips,ipe
622
623     !-----------------------------------!
624     ! 1D subscript for physics "cursor" !
625     !-----------------------------------!
626     subs = (j-jps)*(ipe-ips+1)+(i-ips+1)
627
628     ! fluxrad_save
629     IF (.not. restart) THEN
630      fluxrad(subs)=0.
631     ELSE
632      fluxrad(subs)=M_FLUXRAD(i,j)
633     ENDIF
634     !! et fluxrad_sky ???!???
635
636   ENDDO
637   ENDDO
638 
639   !!---------------------!!
640   !! OUTPUT FOR CHECKING !!
641   !!---------------------!!
642   nlast = (ipe-ips+1)*(jpe-jps+1)
643   print*,"check: fluxrad",fluxrad(1),fluxrad(nlast)
644
645END SUBROUTINE update_inputs_physiq_rad
646
647!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
648!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
649SUBROUTINE update_inputs_physiq_slope( &
650            ims,ime,jms,jme,&
651            ips,ipe,jps,jpe,&
652            JULYR,&
653            SLPX,SLPY)
654
655   !USE module_model_constants, only: DEGRAD
656   !USE slope_mod, ONLY: theta_sl, psi_sl
657
658   INTEGER, INTENT(IN) :: ims,ime,jms,jme
659   INTEGER, INTENT(IN) :: ips,ipe,jps,jpe,JULYR
660   REAL, DIMENSION( ims:ime, jms:jme ), INTENT(IN)  :: SLPX,SLPY
661   INTEGER :: i,j,subs,nlast
662
663END SUBROUTINE update_inputs_physiq_slope
664
665END MODULE update_inputs_physiq_mod
Note: See TracBrowser for help on using the repository browser.