40 USE modd_topodyn, ONLY : nncat, nmesht, xdmaxt, xdxt, xmpara, nnmc, xconn, nline,&
41 xslop, xdarea, xlambda
49 USE modi_write_file_vecmap
50 USE modi_write_file_map
52 USE yomhook
,ONLY : lhook, dr_hook
53 USE parkind1
,ONLY : jprb
59 REAL,
DIMENSION(:,:),
INTENT(IN) :: prw
60 REAL,
DIMENSION(:,:),
INTENT(OUT) :: pdef
61 REAL,
DIMENSION(:,:),
INTENT(OUT) :: pkappa
62 REAL,
DIMENSION(:),
INTENT(OUT) :: pkappac
63 LOGICAL,
DIMENSION(:),
INTENT(OUT) :: gtopd
71 REAL :: zkval, zkvalmin, zkvalmax
74 REAL :: zndmaxav,znkav
83 REAL,
DIMENSION(NMESHT) :: zdmax
84 REAL,
DIMENSION(NMESHT) :: zrw
85 REAL,
DIMENSION(NMESHT) :: zdini
86 REAL,
DIMENSION(NMESHT) :: zmask
87 REAL,
DIMENSION(NMESHT) :: zkappa_pack, zdmax_pack
97 REAL(KIND=JPRB) :: zhook_handle
99 IF (lhook) CALL dr_hook(
'TOPODYN_LAT',0,zhook_handle)
121 zdmax(:) = xdmaxt(jj,:)
131 IF ( zrw(j1)>0.0 .AND. zrw(j1)/=xundef)
THEN
139 CALL
flowdown(nnmc(jj),zmask,xconn(jj,:,:),nline(jj,:))
141 WHERE (zmask == 0.0) zmask = xundef
149 za = nnmc(jj) * zdx**2
154 IF (zmask(j1)/=xundef)
THEN
156 pkappa(jj,j1) = zrw(j1)
158 zdini(j1) = zdmax(j1) - zrw(j1)
160 IF ( zdini(j1) <0.0 )
THEN
162 ztmp = ztmp - zdini(j1)
167 zdav = zdav + zdini(j1)
171 pkappa(jj,j1) = xundef
179 WHERE ( zdini(:)>0. ) zdini(:) = zdini(:)-ztmp/(count(zdini(:)>0.))
182 IF (inpcon >= nnmc(jj)/1000)
THEN
190 CALL
flowdown(nnmc(jj),pkappa(jj,:),xconn(jj,:,:),nline(jj,:))
197 DO WHILE ( .NOT.gfound .AND. j2.LE.nnmc(jj) )
199 IF (zmask(j2)/=xundef)
THEN
202 zkval = pkappa(jj,j2) * exp(xlambda(jj,j2))
207 pkappa(jj,j2) = zkval
211 pkappa(jj,j2) = xundef
221 IF (zmask(j1)/=xundef)
THEN
223 zkval = pkappa(jj,j1) * exp(xlambda(jj,j1))
227 IF (zkval.GT.zkvalmax)
THEN
229 ELSEIF (zkval.LT.zkvalmin)
THEN
233 pkappa(jj,j1) = zkval
237 pkappa(jj,j1) = xundef
246 i_dim = count( zmask(1:nnmc(jj))/=xundef )
247 zkappa_pack(:) = xundef
248 zdmax_pack(:) = xundef
249 zkappa_pack(1:i_dim) = pack(pkappa(jj,1:nnmc(jj)),zmask(1:nnmc(jj))/=xundef)
250 zdmax_pack(1:i_dim) = pack(zdmax(1:nnmc(jj)),zmask(1:nnmc(jj))/=xundef)
252 inkappa = int((zkvalmax - zkvalmin) / xstepk)
256 zkval = zkvalmin + (xstepk * (j1-1))
264 IF ( zkappa_pack(j2).GE.zkval )
THEN
267 ELSEIF (zkappa_pack(j2).LE.( zkval-(zdmax_pack(j2)/zm)) )
THEN
270 zndmaxav = zndmaxav + zdmax_pack(j2)
272 znkav = znkav + zkappa_pack(j2)
280 zndmaxav = zndmaxav /
REAL(inad)
283 IF ( inpcon == inas .OR. inpcon == inad .OR. inpcon == (inad+inas))
THEN
286 znkav = znkav /
REAL(inpcon - inad - inas)
289 IF (inpcon /= 0)
THEN
290 znas =
REAL(INAS) /
REAL(inpcon)
291 znad =
REAL(INAD) /
REAL(inpcon)
294 zfunc = (1 - znas - znad) * ( zkval - znkav )
295 IF (zm /= 0.) zfunc = zfunc + (znad * (zndmaxav / zm))
297 zdif = abs( zfunc - zdav )
299 IF ( zdif.LT.zdifmin )
THEN
320 IF ( zas<1. .AND. zad<1. .AND. (zas + zad/=1.) )
THEN
322 zdav2 = (zdav - zdmaxav * zad) / (1 - zas - zad)
326 IF (zas>=1.)
WRITE(*,*)
'ALL THE AREA IS SATURATED'
327 IF (zad>=1.)
WRITE(*,*)
'ALL THE AREA HAS A MAXIMAL DEFICIT'
328 WRITE(*,*)
'ALL THE AREA',zas,zad
341 IF ( zmask(j1)/=xundef )
THEN
343 IF ( (pkappa(jj,j1).GT.(pkappac(jj) - zdmax(j1)/zm)) .AND. (pkappa(jj,j1).LT.pkappac(jj)) )
THEN
345 pdef(jj,j1) = zm * (zkav - pkappa(jj,j1)) + zdav2
347 IF (pdef(jj,j1) < 0.0) pdef(jj,j1) = 0.0
349 ELSEIF ( pkappa(jj,j1).GE.pkappac(jj) )
THEN
354 ELSEIF ( pkappa(jj,j1).LE.(pkappac(jj) - zdmax(j1)/zm) )
THEN
356 pdef(jj,j1) = zdmax(j1)
366 pdef(jj,j1) = zdmax(j1)
375 IF (pdef(jj,j1)<0.0)
THEN
376 WRITE(*,*)
'LAMBDA=',pkappa(jj,j1),
'LAMBDAC=',pkappac(jj)
388 pkappa(jj,:) = xundef
390 pdef(jj,:) = zdini(:)
397 IF (lhook) CALL dr_hook(
'TOPODYN_LAT',1,zhook_handle)
subroutine flowdown(KNMC, PVAR, PCONN, KLINE)
subroutine topodyn_lat(PRW, PDEF, PKAPPA, PKAPPAC, GTOPD)