30
31
32
33#include "implicit_f.inc"
34
35
36
37#include "com04_c.inc"
38#include "param_c.inc"
39
40
41
42 INTEGER ICODR(*),ISKEW(*),ISKWN(LISKN,*),IXTG(NIXTG,*),
43 . IXTG1(4,*),(NNPBY,*)
45 . x(3,*),skew(lskew,*)
46
47
48
49 INTEGER I,J,II,,IC2,IC3,N1,N2,N3,J1(3),IS,NELTG3
51 . lx,ly,lz,ll,ll1,e1x,e1y,e1z,e2x,e2y,e2z,e3x,e3y,e3z
52
53 neltg3 = numeltg-numeltg6
54 DO i=1,numeltg6
55
56 ii = i + neltg3
57 n1=ixtg(2,ii)
58 n2=ixtg(3,ii)
59 n3=ixtg(4,ii)
60 ic1=icodr(n1)
61 ic2=icodr(n2)
62 ic3=icodr(n3)
63
64 IF (ic1>0.AND.ic2>0)THEN
65 lx=x(1,n1)-x(1,n2)
66 ly=x(2,n1)-x(2,n2)
67 lz=x(3,n1)-x(3,n2)
68 ll =sqrt(lx*lx+ly*ly+lz*lz)
69 j1(1)=ic1/4
70 j1(2)=(ic1-4*j1(1))/2
71 j1(3)=(ic1-4*j1(1)-2*j1(2))
72 is=iskew(n1)
73 IF (is==1) THEN
74 e1x=skew(1,is)
75 e1y=skew(2,is)
76 e1z=skew(3,is)
77 e2x=skew(4,is)
78 e2y=skew(5,is)
79 e2z=skew(6,is)
80 e3x=skew(7,is)
81 e3y=skew(8,is)
82 e3z=skew(9,is)
83 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
84 ELSE
85 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
86 ENDIF
87 IF (abs(ll1)/ll>em6) THEN
88 j1(1)=ic2/4
89 j1(2)=(ic2-4*j1(1))/2
90 j1(3)=(ic2-4*j1(1)-2*j1(2))
91 is=iskew(n2)
92 IF (is==1) THEN
93 e1x=skew(1,is)
94 e1y=skew(2,is)
95 e1z=skew(3,is)
96 e2x=skew(4,is)
97 e2y=skew(5,is)
98 e2z=skew(6,is)
99 e3x=skew(7,is)
100 e3y=skew(8,is)
101 e3z=skew(9,is)
102 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
103 ELSE
104 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
105 ENDIF
106 IF (abs(ll1)/ll>em6) ixtg1(1,i)=-1
107 ENDIF
108
109 ELSEIF (ic2>0.AND.ic3>0)THEN
110 lx=x(1,n3)-x(1,n2)
111 ly=x(2,n3)-x(2,n2)
112 lz=x(3,n3)-x(3,n2)
113 ll =sqrt(lx*lx+ly*ly+lz*lz)
114 j1(1)=ic3/4
115 j1(2)=(ic3-4*j1(1))/2
116 j1(3)=(ic3-4*j1(1)-2*j1(2))
117 is=iskew(n3)
118 IF (is==1) THEN
119 e1x=skew(1,is)
120 e1y=skew(2,is)
121 e1z=skew(3,is)
122 e2x=skew(4,is)
123 e2y=skew(5,is)
124 e2z=skew(6,is)
125 e3x=skew(7,is)
126 e3y=skew(8,is)
127 e3z=skew(9,is)
128 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
129 ELSE
130 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
131 ENDIF
132 IF (abs(ll1)/ll>em6) THEN
133 j1(1)=ic2/4
134 j1(2)=(ic2-4*j1(1))/2
135 j1(3)=(ic2-4*j1(1)-2*j1(2))
136 is=iskew(n2)
137 IF (is==1) THEN
138 e1x=skew(1,is)
139 e1y=skew(2,is)
140 e1z=skew(3,is)
141 e2x=skew(4,is)
142 e2y=skew(5,is)
143 e2z=skew(6,is)
144 e3x=skew(7,is)
145 e3y=skew(8,is)
146 e3z=skew(9,is)
147 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
148 ELSE
149 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
150 ENDIF
151 IF (abs(ll1)/ll>em6) ixtg1(2,i)=-1
152 ENDIF
153
154 ELSEIF (ic1>0.AND.ic3>0)THEN
155 lx=x(1,n3)-x(1,n1)
156 ly=x(2,n3)-x(2,n1)
157 lz=x(3,n3)-x(3,n1)
158 ll =sqrt(lx*lx+ly*ly+lz*lz)
159 j1(1)=ic3/4
160 j1(2)=(ic3-4*j1(1))/2
161 j1(3)=(ic3-4*j1(1)-2*j1(2))
162 is=iskew(n3)
163 IF (is==1) THEN
164 e1x=skew(1,is)
165 e1y=skew(2,is)
166 e1z=skew(3,is)
167 e2x=skew(4,is)
168 e2y=skew(5,is)
169 e2z=skew(6,is)
170 e3x=skew(7,is)
171 e3y=skew(8,is)
172 e3z=skew(9,is)
173 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
174 ELSE
175 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
176 ENDIF
177 IF (abs(ll1)/ll>em6) THEN
178 j1(1)=ic1/4
179 j1(2)=(ic1-4*j1(1))/2
180 j1(3)=(ic1-4*j1(1)-2*j1(2))
181 is=iskew(n1)
182 IF (is==1) THEN
183 e1x=skew(1,is)
184 e1y=skew(2,is)
185 e1z=skew(3,is)
186 e2x=skew(4,is)
187 e2y=skew(5,is)
188 e2z=skew(6,is)
189 e3x=skew(7,is)
190 e3y=skew(8,is)
191 e3z=skew(9,is)
192 ll1=lx*(e1x+e2x+e3x)+ly*(e1y+e2y+e3y)+lz*(e1z+e2z+e3z)
193 ELSE
194 ll1=lx*j1(1)+ly*j1(2)+lz*j1(3)
195 ENDIF
196 IF (abs(ll1)/ll>em6) ixtg1(3,i)=-1
197 ENDIF
198 ENDIF
199 ENDDO
200
201 RETURN