source: LMDZ6/trunk/libf/phylmd/create_etat0_unstruct.F90 @ 3817

Last change on this file since 3817 was 3535, checked in by yann meurdesoif, 6 years ago

For unstructured grid (dynamico).
A file can be used to impose land fraction without using Relief.nc.
The land franction is available at 0.05°, 0.25° and 0.5° degree resolution.
If land_water_0.05.nc or land_water_0.25.nc or land_water_0.50.nc is present in the runing directory, it will be use to impose ter/oce fraction.

YM
File size: 7.8 KB
Line 
1MODULE create_etat0_unstruct_mod
2
3
4
5
6
7
8CONTAINS
9 
10  SUBROUTINE init_create_etat0_unstruct
11#ifdef CPP_XIOS
12  USE xios
13  USE netcdf
14  USE mod_phys_lmdz_para
15  IMPLICIT NONE
16  INTEGER :: file_id, iret
17 
18   ! for coupling activate ocean fraction reading from file "ocean_fraction.nc"
19    IF (is_omp_master) THEN
20
21      IF (NF90_OPEN("ocean_fraction.nc", NF90_NOWRITE, file_id)==NF90_NOERR) THEN
22        CALL xios_set_file_attr("frac_ocean",enabled=.TRUE.)
23        CALL xios_set_field_attr("mask",field_ref="frac_ocean_read")
24        iret=NF90_CLOSE(file_id)
25      ELSE IF (NF90_OPEN("land_water_0.05.nc", NF90_NOWRITE, file_id)==NF90_NOERR) THEN
26        CALL xios_set_file_attr("land_water",name="land_water_0.05",enabled=.TRUE.)
27        CALL xios_set_field_attr("mask",field_ref="land_water")
28        iret=NF90_CLOSE(file_id)
29      ELSE IF (NF90_OPEN("land_water_0.25.nc", NF90_NOWRITE, file_id)==NF90_NOERR) THEN
30        CALL xios_set_file_attr("land_water",name="land_water_0.25",enabled=.TRUE.)
31        CALL xios_set_field_attr("mask",field_ref="land_water")
32        iret=NF90_CLOSE(file_id)
33      ELSE IF (NF90_OPEN("land_water_0.50.nc", NF90_NOWRITE, file_id)==NF90_NOERR) THEN
34        CALL xios_set_file_attr("land_water",name="land_water_0.50",enabled=.TRUE.)
35        CALL xios_set_field_attr("mask",field_ref="land_water")
36        iret=NF90_CLOSE(file_id)
37      ENDIF
38
39    ENDIF
40
41#endif
42 
43  END SUBROUTINE init_create_etat0_unstruct
44 
45 
46  SUBROUTINE create_etat0_unstruct
47  USE dimphy
48#ifdef CPP_XIOS
49  USE xios
50  USE infotrac_phy
51  USE fonte_neige_mod
52  USE pbl_surface_mod
53  USE phys_state_var_mod
54  USE indice_sol_mod
55  USE mod_phys_lmdz_para
56  USE print_control_mod, ONLY: lunout
57  USE geometry_mod
58  USE ioipsl_getin_p_mod, ONLY: getin_p
59
60  IMPLICIT NONE
61  INCLUDE 'dimsoil.h'
62
63    LOGICAL :: no_ter_antartique   ! If true, no land points are allowed at Antartic
64    REAL,    DIMENSION(klon)                 :: tsol
65    REAL,    DIMENSION(klon)                 :: sn
66    REAL,    DIMENSION(klon)                 :: rugmer
67    REAL,    DIMENSION(klon)                 :: run_off_lic_0
68    REAL,    DIMENSION(klon)                 :: lic
69    REAL,    DIMENSION(klon)                 :: fder
70
71    REAL,    DIMENSION(klon,nbsrf)           :: qsolsrf, snsrf
72    REAL,    DIMENSION(klon,nsoilmx,nbsrf)   :: tsoil
73   
74    REAL,    DIMENSION(klon_mpi)             :: tsol_mpi, qsol_mpi, zmasq_mpi, lic_mpi
75    REAL,    DIMENSION(klon_mpi)             :: zmea_mpi, zstd_mpi, zsig_mpi, zgam_mpi, zthe_mpi
76    REAL,    DIMENSION(klon_mpi)             :: cell_area_mpi
77    REAL,    DIMENSION(klon_mpi,nbsrf)       :: pctsrf_mpi
78
79    INTEGER :: ji,j,i
80 
81    IF (is_omp_master) THEN
82      CALL xios_recv_field("ts",tsol_mpi)
83      CALL xios_recv_field("qs",qsol_mpi)
84      CALL xios_recv_field("mask",zmasq_mpi)
85      CALL xios_recv_field("landice",lic_mpi)
86      CALL xios_recv_field("zmea",zmea_mpi)
87      CALL xios_recv_field("zstd",zstd_mpi)
88      CALL xios_recv_field("zsig",zsig_mpi)
89      CALL xios_recv_field("zgam",zgam_mpi)
90      CALL xios_recv_field("zthe",zthe_mpi)
91    ENDIF
92    CALL scatter_omp(tsol_mpi,tsol)
93    CALL scatter_omp(qsol_mpi,qsol)
94    CALL scatter_omp(zmasq_mpi,zmasq)
95    CALL scatter_omp(lic_mpi,lic)
96    CALL scatter_omp(zmea_mpi,zmea)
97    CALL scatter_omp(zstd_mpi,zstd)
98    CALL scatter_omp(zsig_mpi,zsig)
99    CALL scatter_omp(zgam_mpi,zgam)
100    CALL scatter_omp(zthe_mpi,zthe)
101
102    radsol(:)   = 0.0
103    rugmer(:) = 0.001
104    sn(:)     = 0
105
106    WHERE(qsol(:)<0) qsol(:)=0
107       
108    WHERE(   zmasq(:)<EPSFRA) zmasq(:)=0.
109    WHERE(1.-zmasq(:)<EPSFRA) zmasq(:)=1.
110
111    pctsrf(:,:) = 0
112    pctsrf(:,is_lic)=lic
113    WHERE(pctsrf(:,is_lic)<EPSFRA) pctsrf(:,is_lic)=0.
114    WHERE(zmasq(:)<EPSFRA)         pctsrf(:,is_lic)=0.
115
116    pctsrf(:,is_ter)=zmasq(:)
117
118!--- Adequation with soil/sea mask
119    DO ji=1,klon
120      IF(zmasq(ji)>EPSFRA) THEN
121        IF(pctsrf(ji,is_lic)>=zmasq(ji)) THEN
122          pctsrf(ji,is_lic)=zmasq(ji)
123          pctsrf(ji,is_ter)=0.
124        ELSE
125          pctsrf(ji,is_ter)=zmasq(ji)-pctsrf(ji,is_lic)
126          IF(pctsrf(ji,is_ter)<EPSFRA) THEN
127            pctsrf(ji,is_ter)=0.
128            pctsrf(ji,is_lic)=zmasq(ji)
129          END IF
130        END IF
131      END IF
132    END DO
133
134
135  !--- Option no_ter_antartique removes all land fractions souther than 60S.
136  !--- Land ice is set instead of the land fractions on these latitudes.
137  !--- The ocean and sea-ice fractions are not changed.
138  no_ter_antartique=.FALSE.
139  CALL getin_p('no_ter_antartique',no_ter_antartique)
140  WRITE(lunout,*)"no_ter_antartique=",no_ter_antartique
141  IF (no_ter_antartique) THEN
142     ! Remove all land fractions souther than 60S and set land-ice instead
143     WRITE(lunout,*) "Remove land fractions souther than 60deg south by increasing"
144     WRITE(lunout,*) "the continental ice fractions. No land can now be found at Antartic."
145     DO ji=1, klon
146        IF (latitude_deg(ji)<-60.0) THEN
147           pctsrf(ji,is_lic) = pctsrf(ji,is_lic) + pctsrf(ji,is_ter)
148           pctsrf(ji,is_ter) = 0
149        END IF
150     END DO
151  END IF
152   
153! sub-surface ocean and sea ice (sea ice set to zero for start)
154!*******************************************************************************
155    pctsrf(:,is_oce)=(1.-zmasq(:))
156    WHERE(pctsrf(:,is_oce)<EPSFRA) pctsrf(:,is_oce)=0.
157   
158    zval(:)=max(0.,zmea-2*zstd(:))
159    zpic(:)=zmea+2*zstd(:)
160   
161!! WARNING    DON'T FORGET FOR LATER
162!!ym  IF(couple) pctsrf(:,is_oce)=ocemask_fi(:)
163!!
164   
165! Init: tsol, qsol, sn, evap, tsoil, rain_fall, snow_fall, solsw, sollw, frugs
166!*******************************************************************************
167    DO i=1,nbsrf
168     ftsol(:,i) = tsol
169    END DO
170 
171    DO i=1,nbsrf
172     snsrf(:,i) = sn
173    END DO
174!albedo SB >>>
175!ym error : the sub surface dimension is the third not second
176!    falb_dir(:,is_ter,:)=0.08; falb_dir(:,is_lic,:)=0.6
177!    falb_dir(:,is_oce,:)=0.5;  falb_dir(:,is_sic,:)=0.6
178    falb_dir(:,:,is_ter)=0.08; falb_dir(:,:,is_lic)=0.6
179    falb_dir(:,:,is_oce)=0.5;  falb_dir(:,:,is_sic)=0.6
180
181!ym falb_dif has been forgotten, initialize with defaukt value found in phyetat0 or 0 ?
182!ym probably the uninitialized value was 0 for standard (regular grid) case
183    falb_dif(:,:,:)=0
184
185!albedo SB <<<
186    fevap(:,:) = 0.
187    DO i=1,nbsrf
188     qsolsrf(:,i)=150.
189    END DO
190 
191    DO i=1,nbsrf
192      DO j=1,nsoilmx
193        tsoil(:,j,i) = tsol
194      END DO
195    END DO
196 
197    rain_fall = 0.; snow_fall = 0.
198    solsw = 165.;   sollw = -53.
199!ym warning missing init for sollwdown => set to 0
200  sollwdown  = 0.
201   
202   
203    t_ancien = 273.15
204    u_ancien=0
205    v_ancien=0
206    q_ancien = 0.
207    agesno = 0.
208
209    z0m(:,is_oce) = rugmer(:)
210
211   z0m(:,is_ter) = MAX(1.0e-05,zstd(:)*zsig(:)/2.0)
212   z0m(:,is_lic) = MAX(1.0e-05,zstd(:)*zsig(:)/2.0)
213
214   z0m(:,is_sic) = 0.001
215   z0h(:,:)=z0m(:,:)
216
217    fder = 0.0
218    clwcon = 0.0
219    rnebcon = 0.0
220    ratqs = 0.0
221    run_off_lic_0 = 0.0
222    rugoro = 0.0
223
224! Before phyredem calling, surface modules and values to be saved in startphy.nc
225! are initialized
226!*******************************************************************************
227    pbl_tke(:,:,:) = 1.e-8
228    zmax0(:) = 40.
229    f0(:) = 1.e-5
230    sig1(:,:) = 0.
231    w01(:,:) = 0.
232    wake_deltat(:,:) = 0.
233    wake_deltaq(:,:) = 0.
234    wake_s(:) = 0.
235    wake_cstar(:) = 0.
236    wake_fip(:) = 0.
237    wake_pe = 0.
238    fm_therm = 0.
239    entr_therm = 0.
240    detr_therm = 0.
241    ale_bl = 0.
242    ale_bl_trig =0.
243    alp_bl =0.
244    CALL fonte_neige_init(run_off_lic_0)
245    CALL pbl_surface_init( fder, snsrf, qsolsrf, tsoil )
246
247    CALL gather_omp(cell_area,cell_area_mpi)
248    CALL gather_omp(pctsrf,pctsrf_mpi)
249    IF (is_omp_master) THEN
250      CALL xios_send_field("area_ce0l",cell_area_mpi)
251      CALL xios_send_field("fract_oce_ce0l",pctsrf_mpi(:,is_oce))
252      CALL xios_send_field("fract_sic_ce0l",pctsrf_mpi(:,is_sic))
253    ENDIF
254   
255    CALL phyredem( "startphy.nc" )
256
257#endif
258  END SUBROUTINE create_etat0_unstruct
259
260
261END MODULE create_etat0_unstruct_mod
Note: See TracBrowser for help on using the repository browser.