The TELEMAC-MASCARET system  trunk
front2.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE front2
3 ! *****************
4 !
5  &(nfrliq,lihbor,liubor,
6  & x,y,nbor,kp1bor,dejavu,npoin,nptfr,klog,listin,numliq,maxfro)
7 !
8 !***********************************************************************
9 ! BIEF V7P1
10 !***********************************************************************
11 !
12 !brief IDENTIFIES AND NUMBERS THE LIQUID AND SOLID BOUNDARIES.
13 !
14 !note SOLID BOUNDARIES ARE INDICATED WITH LIHBOR(K) = KLOG
15 !+ FOR A BOUNDARY NODE NUMBER K.
16 !+ A SEGMENT CONNECTING A LIQUID AND A SOLID NODE IS
17 !+ CONSIDERED TO BE SOLID.
18 !
19 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
20 !+ 13/07/2010
21 !+ V6P0
22 !+ Translation of French comments within the FORTRAN sources into
23 !+ English comments
24 !
25 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
26 !+ 21/08/2010
27 !+ V6P0
28 !+ Creation of DOXYGEN tags for automated documentation and
29 !+ cross-referencing of the FORTRAN sources
30 !
31 !history J-M HERVOUET
32 !+ 15/09/2015
33 !+ V7P1
34 !+ Now checking the maximum number of boundaries before using the
35 !+ would-be undersized arrays.
36 !
37 !history R. ATA
38 !+ 19/04/2016
39 !+ v7p2
40 !+ Bug correction:LIQF and SOLF were not changed for cases where
41 !+ the 3 nodes (prevK, K,kp1bor(k)) are all liquid
42 !
43 !history S.E.BOURBAN (HRW))
44 !+ 02/02/2017
45 !+ v7p2
46 !+ Bug correction: Allowing for circular liquid boundaries
47 !
48 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
49 !| DEJAVU |<->| WORK ARRAY
50 !| FINLIQ |<--| END OF LIQUID BOUNDARIES
51 !| FINSOL |<--| END OF SOLID BOUNDARIES
52 !| KLOG |-->| LIHBOR(K)=KLOG : SOLID BOUNDARY
53 !| KP1BOR |-->| GIVES THE NEXT BOUNDARY POINT IN A CONTOUR
54 !| LIHBOR |-->| TYPE OF BOUNDARY CONDITIONS ON DEPTH
55 !| LISTIN |-->| IF YES, PRINTING ON LISTING
56 !| LIUBOR |-->| TYPE OF BOUNDARY CONDITIONS ON U
57 !| MAXFRO |-->| MAXIMUM NUMBER OF LIQUID OR SOLID BOUNDARIES
58 !| NBOR |-->| GLOBAL NUMBER OF BOUNDARY POINTS
59 !| NFRLIQ |<--| NUMBER OF LIQUID BOUNDARIES
60 !| NPOIN |-->| NUMBER OF POINTS
61 !| NPTFR |-->| NUMBER OF BOUNDARY POINTS
62 !| NUMLIQ |-->| BOUNDARY NUMBER OF BOUNDARY POINTS
63 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
64 !
66  IMPLICIT NONE
67 !
68 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
69 !
70  INTEGER, INTENT(IN) :: NPOIN,NPTFR,KLOG,MAXFRO
71  INTEGER, INTENT(OUT) :: NFRLIQ
72  INTEGER , INTENT(IN) :: LIHBOR(nptfr),LIUBOR(nptfr)
73  DOUBLE PRECISION, INTENT(IN) :: X(npoin) , Y(npoin)
74  INTEGER, INTENT(IN) :: NBOR(nptfr),KP1BOR(nptfr)
75  INTEGER, INTENT(OUT) :: DEJAVU(nptfr)
76  LOGICAL, INTENT(IN) :: LISTIN
77  INTEGER, INTENT(OUT) :: NUMLIQ(nptfr)
78 !
79 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
80 !
81  INTEGER K,KPREV,IDEP,SOL1,LIQ1,L1,L2,L3,NILE
82 !
83  LOGICAL SOLF,LIQF,SOLD,LIQD
84 !
85  DOUBLE PRECISION MINNS,MAXNS,EPS,YMIN,NS
86  DOUBLE PRECISION, PARAMETER :: EPSS=1.d-6
87  INTEGER :: MAXFROLIQ
88  INTEGER :: MAXFROSOL
89  INTEGER, ALLOCATABLE :: TMP(:)
90  INTEGER, ALLOCATABLE :: DEBLIQ(:),FINLIQ(:)
91  INTEGER, ALLOCATABLE :: DEBSOL(:),FINSOL(:)
92  INTEGER :: NFRSOL
93  INTEGER :: IERR
94 !
95  INTRINSIC abs
96 !
97 !-----------------------------------------------------------------------
98 !
99 ! INITIALISES
100 !
101 ! DEJAVU : MARKS WITH 1 THE POINTS THAT HAVE ALREADY BEEN TREATED
102 ! NILE : NUMBER OF ISLANDS
103 !
104  DO k=1,nptfr
105  dejavu(k) = 0
106  ENDDO ! K
107  maxfroliq = maxfro
108  maxfrosol = maxfro
109  ALLOCATE(debliq(maxfroliq),stat=ierr)
110  CALL check_allocate(ierr,'DEBLIQ')
111  ALLOCATE(finliq(maxfroliq),stat=ierr)
112  CALL check_allocate(ierr,'FINLIQ')
113  ALLOCATE(debsol(maxfrosol),stat=ierr)
114  CALL check_allocate(ierr,'DEBSOL')
115  ALLOCATE(finsol(maxfrosol),stat=ierr)
116  CALL check_allocate(ierr,'FINSOL')
117 !
118  nile = 0
119  idep = 1
120  nfrliq = 0
121  nfrsol = 0
122  debliq = 0
123  finliq = 0
124  debsol = 0
125  finsol = 0
126 !
127 !-----------------------------------------------------------------------
128 !
129 ! COMES BACK TO LABEL 20 IF THERE IS AT LEAST 1 ISLAND
130 !
131 20 CONTINUE
132 !
133 ! LOOKS FOR THE SOUTH-WESTERNMOST POINT (THERE CAN BE MORE THAN 1)
134 !
135  minns = x(nbor(idep)) + y(nbor(idep))
136  maxns = minns
137  ymin = y(nbor(idep))
138 !
139  DO k = 1 , nptfr
140  IF(dejavu(k).EQ.0) THEN
141  ns = x(nbor(k)) + y(nbor(k))
142  IF(ns.LT.minns) THEN
143  idep = k
144  minns = ns
145  ymin = y(nbor(k))
146  ENDIF
147  IF(ns.GT.maxns) maxns = ns
148  ENDIF
149  ENDDO ! K
150 !
151  eps = (maxns-minns) * 1.d-4
152 !
153 ! SELECTS THE SOUTHERNMOST POINT FROM THE SOUTH-WESTERNMOST CANDIDATES
154 !
155  DO k = 1 , nptfr
156  IF(dejavu(k).EQ.0) THEN
157  ns = x(nbor(k)) + y(nbor(k))
158  IF(abs(minns-ns).LT.eps) THEN
159  IF(y(nbor(k)).LT.ymin) THEN
160  idep = k
161  ymin = y(nbor(k))
162  ENDIF
163  ENDIF
164  ENDIF
165  ENDDO ! K
166 !
167 !-----------------------------------------------------------------------
168 !
169 ! NUMBERS AND LOCATES THE CONTOUR BOUNDARIES STARTING
170 ! AT POINT IDEP
171 !
172 ! SOLD = .TRUE. : THE BOUNDARY STARTING AT IDEP IS SOLID
173 ! LIQD = .TRUE. : THE BOUNDARY STARTING AT IDEP IS LIQUID
174 ! SOLF = .TRUE. : THE BOUNDARY ENDING AT IDEP IS SOLID
175 ! LIQF = .TRUE. : THE BOUNDARY ENDING AT IDEP IS LIQUID
176 ! LIQ1 : NUMBER OF THE 1ST LIQUID BOUNDARY OF THE CONTOUR
177 ! SOL1 : NUMBER OF THE 1ST SOLID BOUNDARY OF THE CONTOUR
178 !
179  k = idep
180 !
181  sol1 = 0
182  liq1 = 0
183  liqf = .false.
184  solf = .false.
185 !
186 ! TYPE OF THE 1ST SEGMENT
187 !
188 ! LAW OF PREDOMINANCE SOLID OVER LIQUID
189  IF(lihbor(k).EQ.klog.OR.lihbor(kp1bor(k)).EQ.klog) THEN
190 ! THE 1ST SEGMENT IS SOLID
191  nfrsol = nfrsol + 1
192  IF(nfrsol.GT.maxfrosol) THEN
193  ! REALLOCATE ARRAY
194  ALLOCATE(tmp(maxfrosol),stat=ierr)
195  CALL check_allocate(ierr,'TMP')
196 
197  tmp = debsol
198  DEALLOCATE(debsol)
199  ALLOCATE(debsol(maxfrosol*10),stat=ierr)
200  CALL check_allocate(ierr,'DEBSOL_UP')
201  debsol(1:maxfrosol) = tmp
202 
203  tmp = finsol
204  DEALLOCATE(finsol)
205  ALLOCATE(finsol(maxfrosol*10),stat=ierr)
206  CALL check_allocate(ierr,'FINSOL_UP')
207  finsol(1:maxfrosol) = tmp
208 
209  DEALLOCATE(tmp)
210  maxfrosol = maxfrosol*10
211  ENDIF
212 
213  sol1 = nfrsol
214  sold = .true.
215  liqd = .false.
216  debsol(nfrsol) = k
217  finsol(nfrsol) = k ! IN CASE OF CIRCULAR BOUNDARY
218  ELSE
219 ! THE 1ST SEGMENT IS LIQUID
220  nfrliq = nfrliq + 1
221  IF(nfrliq.GT.maxfroliq) THEN
222  ! REALLOCATE ARRAY
223  ALLOCATE(tmp(maxfroliq),stat=ierr)
224  CALL check_allocate(ierr,'TMP')
225 
226  tmp = debliq
227  DEALLOCATE(debliq)
228  ALLOCATE(debliq(maxfroliq*10),stat=ierr)
229  CALL check_allocate(ierr,'DEBLIQ_UP')
230  debliq(1:maxfroliq) = tmp
231 
232  tmp = finliq
233  DEALLOCATE(finliq)
234  ALLOCATE(finliq(maxfroliq*10),stat=ierr)
235  CALL check_allocate(ierr,'FINLIQ_UP')
236  finliq(1:maxfroliq) = tmp
237 
238  DEALLOCATE(tmp)
239  maxfroliq = maxfroliq*10
240  ENDIF
241  liq1 = nfrliq
242  liqd = .true.
243  sold = .false.
244  debliq(nfrliq) = k
245  finliq(nfrliq) = k ! IN CASE OF CIRCULAR BOUNDARY
246  ENDIF
247 !
248  dejavu(k) = 1
249  kprev = k
250  k = kp1bor(k)
251 !
252 50 CONTINUE
253 !
254 ! LOOKS FOR TRANSITION POINTS FROM THE POINT FOLLOWING IDEB
255 !
256 ! ALSO LOOKS FOR ISOLATED POINTS TO DETECT THE ERRORS IN
257 ! THE DATA
258 !
259  l1 = lihbor(kprev)
260  l2 = lihbor(k)
261  l3 = lihbor(kp1bor(k))
262 !
263  IF(l1.EQ.klog.AND.l2.NE.klog.AND.l3.NE.klog) THEN
264 ! SOLID-LIQUID TRANSITION AT POINT K
265  nfrliq = nfrliq + 1
266  IF(nfrliq.GT.maxfroliq) THEN
267  ! REALLOCATE ARRAY
268  ALLOCATE(tmp(maxfroliq),stat=ierr)
269  CALL check_allocate(ierr,'TMP')
270 
271  tmp = debliq
272  DEALLOCATE(debliq)
273  ALLOCATE(debliq(maxfroliq*10),stat=ierr)
274  CALL check_allocate(ierr,'DEBLIQ_UP')
275  debliq(1:maxfroliq) = tmp
276 
277  tmp = finliq
278  DEALLOCATE(finliq)
279  ALLOCATE(finliq(maxfroliq*10),stat=ierr)
280  CALL check_allocate(ierr,'FINLIQ_UP')
281  finliq(1:maxfroliq) = tmp
282 
283  DEALLOCATE(tmp)
284  maxfroliq = maxfroliq*10
285  ENDIF
286  finsol(nfrsol) = k
287  debliq(nfrliq) = k
288  liqf = .true.
289  solf = .false.
290  ELSEIF(l1.NE.klog.AND.l2.NE.klog.AND.l3.EQ.klog) THEN
291 ! LIQUID-SOLID TRANSITION AT POINT K
292  nfrsol = nfrsol + 1
293  IF(nfrsol.GT.maxfrosol) THEN
294  ! REALLOCATE ARRAY
295  ALLOCATE(tmp(maxfrosol),stat=ierr)
296  CALL check_allocate(ierr,'TMP')
297 
298  tmp = debsol
299  DEALLOCATE(debsol)
300  ALLOCATE(debsol(maxfrosol*10),stat=ierr)
301  CALL check_allocate(ierr,'DEBSOL_UP')
302  debsol(1:maxfrosol) = tmp
303 
304  tmp = finsol
305  DEALLOCATE(finsol)
306  ALLOCATE(finsol(maxfrosol*10),stat=ierr)
307  CALL check_allocate(ierr,'FINSOL_UP')
308  finsol(1:maxfrosol) = tmp
309 
310  DEALLOCATE(tmp)
311  maxfrosol = maxfrosol*10
312  ENDIF
313  finliq(nfrliq) = k
314  debsol(nfrsol) = k
315  liqf = .false.
316  solf = .true.
317  ELSEIF(l1.NE.klog.AND.l2.NE.klog.AND.l3.NE.klog) THEN
318  liqf = .true.
319  solf = .false.
320 ! LIQUID-LIQUID TRANSITIONS AT POINT K
321  IF(l2.NE.l3.OR.liubor(k).NE.liubor(kp1bor(k))) THEN
322  finliq(nfrliq) = k
323  nfrliq = nfrliq + 1
324  IF(nfrliq.GT.maxfroliq) THEN
325  ! REALLOCATE ARRAY
326  ALLOCATE(tmp(maxfroliq),stat=ierr)
327  CALL check_allocate(ierr,'TMP')
328 
329  tmp = debliq
330  DEALLOCATE(debliq)
331  ALLOCATE(debliq(maxfroliq*10),stat=ierr)
332  CALL check_allocate(ierr,'DEBLIQ_UP')
333  debliq(1:maxfroliq) = tmp
334 
335  tmp = finliq
336  DEALLOCATE(finliq)
337  ALLOCATE(finliq(maxfroliq*10),stat=ierr)
338  CALL check_allocate(ierr,'FINLIQ_UP')
339  finliq(1:maxfroliq) = tmp
340 
341  DEALLOCATE(tmp)
342  maxfroliq = maxfroliq*10
343  debliq(nfrliq) = kp1bor(k)
344  ENDIF
345  debliq(nfrliq) = kp1bor(k)
346  ENDIF
347  ELSEIF(l1.EQ.klog.AND.l2.NE.klog.AND.l3.EQ.klog) THEN
348 ! ERROR IN THE DATA
349  WRITE(lu,103) k, nbor(k)
350  CALL plante(1)
351  stop
352  ELSEIF(l1.NE.klog.AND.l2.EQ.klog.AND.l3.NE.klog) THEN
353 ! ERROR IN THE DATA
354  WRITE(lu,105) k, nbor(k)
355  CALL plante(1)
356  stop
357  ENDIF
358 !
359  dejavu(k) = 1
360  kprev = k
361  k = kp1bor(k)
362  IF(k.NE.idep) GO TO 50
363 !
364 ! CASE WHERE THE BOUNDARY TYPE CHANGES AT IDEP
365 !
366  IF(solf) THEN
367 ! THE LAST CONTOUR BOUNDARY WAS SOLID
368  IF(sold) THEN
369 ! THE FIRST CONTOUR BOUNDARY WAS SOLID
370  IF( sol1.NE.nfrsol )THEN
371  debsol(sol1) = debsol(nfrsol)
372  nfrsol = nfrsol - 1
373  ENDIF
374  ELSEIF(liqd) THEN
375 ! THE FIRST CONTOUR BOUNDARY WAS LIQUID
376  debliq(liq1) = idep
377  finsol(nfrsol) = idep
378  ENDIF
379 !
380  ELSEIF(liqf) THEN
381 ! THE LAST CONTOUR BOUNDARY WAS LIQUID
382  IF(liqd) THEN
383 ! THE FIRST CONTOUR BOUNDARY WAS LIQUID
384  IF( liq1.NE.nfrliq )THEN
385  debliq(liq1) = debliq(nfrliq)
386  nfrliq = nfrliq - 1
387  ENDIF
388  ELSEIF(sold) THEN
389 ! THE FIRST CONTOUR BOUNDARY WAS SOLID
390  debsol(sol1) = idep
391  finliq(nfrliq) = idep
392  ENDIF
393 !
394  ELSE
395 ! CASE WHERE THE WHOLE CONTOUR HAS THE SAME TYPE
396  IF(sol1.NE.0) THEN
397  debsol(sol1) = idep
398  finsol(sol1) = idep
399  ELSEIF(liq1.NE.0) THEN
400  debliq(liq1) = idep
401  finliq(liq1) = idep
402  ELSE
403  IF(listin) THEN
404  WRITE(lu,'(1X,A)') 'IMPOSSIBLE CASE IN FRONT2'
405  ENDIF
406  CALL plante(1)
407  stop
408  ENDIF
409  ENDIF
410 !
411 !-----------------------------------------------------------------------
412 !
413 ! CHECKS WHETHER THERE ARE OTHER CONTOURS LEFT:
414 !
415  DO k = 1 , nptfr
416  IF(dejavu(k).EQ.0) THEN
417  idep = k
418  nile = nile + 1
419  GO TO 20
420  ENDIF
421  ENDDO ! K
422 !
423 !-----------------------------------------------------------------------
424 !
425  DO k=1,nptfr
426  numliq(k)=0
427  ENDDO ! K
428 !
429 !-----------------------------------------------------------------------
430 !
431 ! PRINTS OUT THE RESULTS AND COMPUTES NUMLIQ
432 !
433  IF(nile.NE.0.AND.listin) WRITE(lu,169) nile
434 !
435  IF(nfrliq.NE.0) THEN
436  IF(listin) WRITE(lu,170) nfrliq
437  DO k = 1, nfrliq
438 !
439 ! MARKS THE NUMBERS OF THE LIQUID BOUNDARIES
440 !
441  l1=debliq(k)
442  numliq(l1)=k
443 707 l1=kp1bor(l1)
444  numliq(l1)=k
445  IF(l1.NE.finliq(k)) GO TO 707
446 !
447 ! END OF MARKING
448 !
449  IF(listin) WRITE(lu,190)
450  & k,debliq(k),nbor(debliq(k)),
451  & x(nbor(debliq(k))),y(nbor(debliq(k))),
452  & finliq(k),nbor(finliq(k)),
453  & x(nbor(finliq(k))),y(nbor(finliq(k)))
454  ENDDO ! K
455  ENDIF
456 !
457  IF(nfrsol.NE.0) THEN
458  IF(listin) WRITE(lu,101) nfrsol
459  DO k = 1, nfrsol
460  IF(listin) WRITE(lu,190)
461  & k,debsol(k),nbor(debsol(k)),
462  & x(nbor(debsol(k))),y(nbor(debsol(k))),
463  & finsol(k),nbor(finsol(k)),
464  & x(nbor(finsol(k))),y(nbor(finsol(k)))
465  ENDDO ! K
466  ENDIF
467 !
468 !-----------------------------------------------------------------------
469 !
470  GO TO 2000
471 !
472 !-----------------------------------------------------------------------
473 !
474 ! FORMATS
475 !
476 169 FORMAT(/,1x,'THERE IS ',1i5,' ISLAND(S) IN THE DOMAIN')
477 170 FORMAT(/,1x,'THERE IS ',1i5,' LIQUID BOUNDARIES:')
478 101 FORMAT(/,1x,'THERE IS ',1i5,' SOLID BOUNDARIES:')
479 103 FORMAT(/,1x,'FRONT2 : ERROR AT BOUNDARY POINT ',1i8,
480  & /,1x,' LOCAL NUMBER ',1i8,
481  & /,1x,' LIQUID POINT BETWEEN TWO SOLID POINTS')
482 105 FORMAT(/,1x,'FRONT2 : ERROR AT BOUNDARY POINT ',1i8,
483  & /,1x,' LOCAL NUMBER ',1i8,
484  & /,1x,' SOLID POINT BETWEEN TWO LIQUID POINTS')
485 190 FORMAT(/,1x,'BOUNDARY ',1i4,' : ',/,1x,
486  & ' BEGINS AT BOUNDARY POINT: ',1i8,
487  & ' , WITH GLOBAL NUMBER: ',1i9,/,1x,
488  & ' AND COORDINATES: ',g16.7,3x,g16.7,
489  & /,1x,' ENDS AT BOUNDARY POINT: ',1i8,
490  & ' , WITH GLOBAL NUMBER: ',1i9,/,1x,
491  & ' AND COORDINATES: ',g16.7,3x,g16.7)
492 !
493 !-----------------------------------------------------------------------
494 !
495 2000 CONTINUE
496  DEALLOCATE(debsol)
497  DEALLOCATE(finsol)
498  DEALLOCATE(debliq)
499  DEALLOCATE(finliq)
500 !
501 !-----------------------------------------------------------------------
502 !
503  RETURN
504  END
subroutine front2(NFRLIQ, LIHBOR, LIUBOR, X, Y, NBOR, KP1BOR, DEJAVU, NPOIN, NPTFR, KLOG, LISTIN, NUMLIQ, MAXFRO)
Definition: front2.f:8