33 use element_mod , only : nixtg
34
35
36
37#include "implicit_f.inc"
38#include "comlock.inc"
39
40
41
42#include "com04_c.inc"
43#include "vect01_c.inc"
44#include "param_c.inc"
45#include "tabsiz_c.inc"
46
47
48
49 INTEGER ,INTENT (IN) :: IXTG(NIXTG,NUMELTG)
50 INTEGER ,INTENT (INOUT) :: LGAUGE(3,NBGAUGE)
51 my_real,
INTENT(IN) :: phi(sphi),x(3,numnod)
52 my_real,
INTENT(INOUT) :: gauge(llgauge,nbgauge)
53 TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECT
54
55
56
57 INTEGER I,II,J,N,IG,IAD2
59 . det_12_1p,det_23_2p,det_31_3p,y1,z1,y2,z2,y3,z3
60
61
62
63
64
65
66 DO i=1,llt
67 ii=i+nft
68 iad2 = ale_connect%ee_connect%iad_connect(ii)
69 DO j=1,3
70 n= ale_connect%ee_connect%connected(iad2 + j - 1)
71 IF(n<=0)cycle
72 ig=nint(phi(n))
73 IF(ig==0)cycle
74
75
76 xg = gauge(2,ig)
77 yg = gauge(3,ig)
78 zg = gauge(4,ig)
79
80
81
82
83 y1 = x(2,ixtg(2,ii))
84 z1 = x(3,ixtg(2,ii))
85 y2 = x(2,ixtg(3,ii))
86 z2 = x(3,ixtg(3,ii))
87 y3 = x(2,ixtg(4,ii))
88 z3 = x(3,ixtg(4,ii))
89
90 det_12_1p = (y2-y1)*(zg-z1)-(z2-z1)*(yg-y1)
91 det_23_2p = (y3-y2)*(zg-z2)-(z3-z2)*(yg-y2)
92 det_31_3p = (y1-y3)*(zg-z3)-(z1-z3)*(yg-y3)
93
95 IF ( (det_12_1p >=0 .AND. det_23_2p >=0 .AND.
96 . det_31_3p >=0) .OR.
97 . (det_12_1p <=0 .AND. det_23_2p <=0 .AND.
98 . det_31_3p <=0) )
alpha = one
99
100#include "lockon.inc"
101 IF(
alpha == one .AND.
alpha >= gauge(5,ig))
THEN
102
103
104 lgauge(1,ig)=-ii
106 ENDIF
107#include "lockoff.inc"
108 ENDDO
109 ENDDO
110
111
112 RETURN