36
37
38
42
43
44
45#include "implicit_f.inc"
46
47
48
49#include "mvsiz_p.inc"
50
51
52
53#include "units_c.inc"
54#include "param_c.inc"
55#include "vect01_c.inc"
56#include "scr03_c.inc"
57#include "scr11_c.inc"
58#include "com04_c.inc"
59
60
61
62 INTEGER MAT(*),NIX,IX(NIX,*)
63 INTEGER,INTENT(IN) :: M151_ID
64 my_real pm(npropm,nummat), tb(*),x(3,*),vdet2
65 TYPE(DETONATORS_STRUCT_)::DETONATORS
66
67
68
69 INTEGER I, N, MTL, MT,IOPT,NPE
70 INTEGER NDETPS, NDETSG, NECRAN, , NDETCORD,NDETPS_NO_SHADOW,NDETPS_SHADOW
71 INTEGER :: I_SHADOW_FLAG
72 LOGICAL HAS_DETONATOR(MVSIZ)
73 my_real y1(mvsiz), y2(mvsiz), y3(mvsiz),
74 . z1(mvsiz), z2(mvsiz), z3(mvsiz),
75 . xc(mvsiz), yc(mvsiz), zc(mvsiz), bt(mvsiz),
76 . dl(mvsiz), alt, xlp, ylp, zlp, xlp1,
77 . ylp1, zlp1, xlp2, ylp2, zlp2, xl0, yl0, zl0, xl1, yl1, zl1,
78 . xl2, yl2, zl2, ps1, ps2, dl1, dl2, s1, s2, s3,
79 . nx, ny, nz , nn, vdet
80 INTEGER :: NODE1, NODE2, NODE3, II, GRNOD_ID, INOD, NNOD, NODE_ID
81
82
83
84 ndetps = detonators%N_DET_POINT
85 ndetsg = detonators%N_DET_LINE
86 necran = detonators%N_DET_WAVE_SHAPER
87 ndetpl = detonators%N_DET_PLANE
88 ndetcord = detonators%N_DET_CORD
89
90
91 ndetps_no_shadow = 0
92 ndetps_shadow = 0
93 DO i = 1, detonators%N_DET_POINT
94 i_shadow_flag = detonators%POINT(i)%SHADOW
95 IF(i_shadow_flag == 0)THEN
96 ndetps_no_shadow = ndetps_no_shadow + 1
97 ELSE
98 ndetps_shadow = ndetps_shadow + 1
99 ENDIF
100 ENDDO
101
102
103
104 IF(detonators%N_DET - ndetps_shadow > 0) THEN
105
106 DO i = lft, llt
107 ii = i + nft
108 node1 = ix(2, ii)
109 node2 = ix(3, ii)
110 node3 = ix(4, ii)
111 y1(i) = x(2, node1)
112 y2(i) = x(2, node2)
113 y3(i) = x(2, node3)
114 z1(i) = x(3, node1)
115 z2(i) = x(3, node2)
116 z3(i) = x(3, node3)
117 ENDDO
118
119
120
121
122 DO i=lft,llt
123 tb(i)=ep21
124 has_detonator(i)=.false.
125 xc(i)=zero
126 yc(i)=third*(y1(i)+y2(i)+y3(i))
127 zc(i)=third*(z1(i)+z2(i)+z3(i))
128 END DO
129
130
131
132 !---------------------------------
133 IF(ndetps > 0) THEN
134 DO i=lft,llt
135 DO n=1,ndetps
136 i_shadow_flag = detonators%POINT(n)%SHADOW
137 IF(i_shadow_flag /= 0)cycle
138 mtl=detonators%POINT(n)%MAT
139 grnod_id=detonators%POINT(n)%GRNOD_ID
140 IF(mtl == 0 .OR. mtl == mat(i) .OR. mtl == m151_id) THEN
141
142 IF(grnod_id == 0)THEN
143 alt=detonators%POINT(n)%TDET
144 xlp=detonators%POINT(n)%XDET
145 ylp=detonators%POINT(n)%YDET
146 zlp=detonators%POINT(n)%ZDET
147 dl(i) =(xc(i)-xlp)**2+(yc(i)-ylp)**2+(zc(i)-zlp)**2
148 dl(i)=sqrt(dl(i))
149 bt(i) =alt+dl(i)/pm(38,mat(i))
150 IF(bt(i) < abs(tb(i))) tb(i)=-bt(i)
151 has_detonator(i)= .true.
152
153 ELSE
154 nnod = detonators%POINT(n)%NNOD
155 alt=detonators%POINT(n)%TDET
156 has_detonator(i)=.true.
157 DO inod=1,nnod
158 node_id=detonators%POINT(n)%NODLIST(inod)
159 xlp=x(1,node_id)
160 ylp=x(2,node_id)
161 zlp=x(3,node_id)
162 dl(i) =(xc(i)-xlp)**2+(yc(i)-ylp)**2+(zc(i)-zlp)**2
163 dl(i)=sqrt(dl(i))
164 bt(i) =alt+dl(i)/pm(38,mat(i))
165 IF(bt(i) < abs(tb(i))) tb(i)=-bt(i)
166 ENDDO
167 ENDIF
168 END IF
169 END DO
170 END DO
171 ENDIF
172
173
174
175
176 IF(ndetsg > 0) THEN
177 DO n=1,ndetsg
178 alt=detonators%LINE(n)%TDET
179 mtl=detonators%LINE(n)%MAT
180 xlp1=detonators%LINE(n)%XDET_1
181 ylp1=detonators%LINE(n)%YDET_1
182 zlp1=detonators%LINE(n)%ZDET_1
183 xlp2=detonators%LINE(n)%XDET_2
184 ylp2=detonators%LINE(n)%YDET_2
185 zlp2=detonators%LINE(n)%ZDET_2
186 DO i=lft,llt
187 IF(mtl == 0 .OR. mtl == mat(i) .OR. mtl == m151_id) THEN
188 xl0 =(xlp1-xlp2)
189 yl0 =(ylp1-ylp2)
190 zl0 =(zlp1-zlp2)
191 xl1 =(xc(i)-xlp1)
192 yl1 =(yc(i)-ylp1)
193 zl1 =(zc(i)-zlp1)
194 xl2 =(xc(i)-xlp2)
195 yl2 =(yc(i)-ylp2)
196 zl2 =(zc(i)-zlp2)
197 ps1 =xl1*xl0+yl1*yl0+zl1*zl0
198 ps2 =xl2*xl0+yl2*yl0+zl2*zl0
199 IF(ps1*ps2 > zero)THEN
200 dl1 =sqrt(xl1**2+yl1**2+zl1**2)
201 dl2 =sqrt(xl2**2+yl2**2+zl2**2)
203 ELSE
204 s1 =yl1*zl0 - zl1*yl0
205 s2 =zl1*xl0 - xl1*zl0
206 s3 =xl1*yl0 - yl1*xl0
207 dl(i)=sqrt((s1**2+s2**2+s3**2)/(xl0**2+yl0**2+zl0**2))
208 ENDIF
209 bt(i) =alt+dl(i)/pm(38,mat(i))
210 IF(bt(i) < abs(tb(i))) tb(i)=-bt(i)
211 has_detonator(i)=.true.
212 END IF
213 END DO
214 END DO
215 ENDIF
216
217
218
219
220 IF(necran > 0) THEN
221 DO n=1,necran
222 alt=detonators%WAVE_SHAPER(n)%TDET
223 mtl=detonators%WAVE_SHAPER(n)%MAT
224 vdet =detonators%WAVE_SHAPER(n)%VDET
225 yd =detonators%WAVE_SHAPER(n)%YDET
226 zd =detonators%WAVE_SHAPER(n)%ZDET
227 npe=detonators%WAVE_SHAPER(n)%NUMNOD
228 dto0=alt
229 vdto=pm(38,mat(1))
230 IF(vdet == zero)vdet=pm(38,mat(1))
231
232 CALL ecran1(detonators%WAVE_SHAPER(n),x,vdet)
233
234 DO i=lft,llt
235 IF(mtl /= mat(i) .AND. mtl /= 0 .AND. mtl /= m151_id) cycle
236 yd =detonators%WAVE_SHAPER(n)%YDET
237 zd =detonators%WAVE_SHAPER(n)%ZDET
238 dto0=alt
239 yl=yc(i)
240 zl=zc(i)
241 CALL ecran2(detonators%WAVE_SHAPER(n),x,vdet)
242 bt(i) =dto
243 IF(bt(i) < abs(tb(i))) tb(i)=-bt(i)
244 has_detonator(i)= .true.
245 END DO
246 END DO
247 ENDIF
248
249
250
251
252 IF(ndetpl > 0) THEN
253 DO n=1,ndetpl
254 alt=detonators%PLANE(n)%TDET
255 mtl=detonators%PLANE(n)%MAT
256 xlp=detonators%PLANE(n)%XDET
257 ylp=detonators%PLANE(n)%YDET
258 zlp=detonators%PLANE(n)%ZDET
259 nx=detonators%PLANE(n)%NX
260 ny=detonators%PLANE(n)%NY
261 nz=detonators%PLANE
262 nn=sqrt(nx**2+ny**2+nz**2)
264 dl1=xlp*nx + ylp*ny + zlp*nz
265 dl1 = dl1/nn
266 DO i=lft,llt
267 IF(mtl == 0 .OR. mtl == mat(i) .OR. mtl == m151_id) THEN
268
269
270
271
272 dl(i) = (xc(i)-xlp)*nx + (yc(i)-ylp)*ny + (zc(i)-zlp)*nz
273 dl(i) = abs(dl(i))
274 dl(i) = dl(i)/nn
275 bt(i) =alt+dl(i)/pm(38,mat(i))
276 IF(bt(i) < abs(tb(i))) tb(i)=-bt(i)
277 has_detonator(i)= .true.
278 END IF
279 END DO
280 END DO
281 ENDIF
282
283
284
285
286 IF(ndetcord > 0) THEN
287 DO n=1,ndetcord
288 alt = detonators%CORD(n)%TDET
289 mtl = detonators%CORD(n)%MAT
290 vdet2 = detonators%CORD(n)%VDET
291 iopt = detonators%CORD(n)%IOPT
292 dto0 = alt
293 mt = mtl
294 IF(mt == 0)mt=mat(1)
295 vdto = pm(38,mt)
296 IF(mtl /= mat(1) .AND. mtl /= 0 .AND. mtl /= m151_id) cycle
297 dto0 = alt
298 CALL detcord(detonators%CORD(n),x,xc,yc,zc,vdto,vdet2,alt,bt,tb,has_detonator
299 END do
300 ENDIF
301
302
303 END IF
304
305
306 RETURN
subroutine detcord(detonator_cord, x, xc, yc, zc, vdet, vdet2, alt, bt, tb, has_detonator, iopt)
subroutine ecran1(detonator, x, vdet)
subroutine ecran2(detonator_wave_shaper, x, vdet)