39
40
41
44
45
46
47#include "implicit_f.inc"
48#include "comlock.inc"
49
50
51
52#include "parit_c.inc"
53#include "com01_c.inc"
54#include "com08_c.inc"
55#include "scr18_c.inc"
56#include "scr16_c.inc"
57#include "com06_c.inc"
58#include "scr07_c.inc"
59#include "scr14_c.inc"
60
61
62
63 INTEGER NSN, NMN, NTY, IBC, IMAST
64 INTEGER IRECT(4,*), LMSR(*), MSR(*), NSV(*), ILOC(*), IRTL(*),
65 . IRTLO(*), ICODT(*), ISKY(*)
66
68 . x(3,*), e(*), stf(*), cst(2,*), fric0(3,*), frigap(*),
69 . stfn(*), fsav(*),fskyi(lskyi,nfskyi),ptmax, areas(*),
70 . fcont(3,*),fncont(3,*), ftcont(3,*)
71 TYPE(H3D_DATABASE) :: H3D_DATA
72
73
74
75 INTEGER IX(2), II, I, J, K, L, M, IMP, I3, I2, JJ, J3, J2, LOLD,
76 . NISKYL, IPAS
77
79 . h(2), n2, n3, fric, gap, ym1, zm1, ym2, zm2, ys, zs, t2, t3,
80 . xl, ans, ss, stif, fni, fyi, fzi, ss0, fti, ds, anst, fmax,
81 . stfri, ax, fs, ft
82
83 fric=frigap(1)
84 gap =frigap(2)
85
86 DO 500 ii=1,nsn
87 ipas = 0
88 i=nsv(ii)
89 j=iloc(ii)
90 k=msr(j)
91 l=irtl(ii)
92 m=msr(irect(1,l))
93 ix(1)=m
94 ym1=x(2,m)
95 zm1=x(3,m)
96 m=msr(irect(2,l))
97 ix(2)=m
98 ym2=x(2,m)
99 zm2=x(3,m)
100 ys =x(2,i)
101 zs =x(3,i)
102 IF(n2d==1)THEN
103 ax=ys
104 ELSE
105 ax=one
106 ENDIF
107 t2=ym2-ym1
108 t3=zm2-zm1
109 xl=sqrt(t2**2+t3**2)
110 t2=t2/xl
111 t3=t3/xl
112 n2= t3
113 n3=-t2
114
115 imp=0
116 i3=3*i
117 i2=i3-1
118
119 ans =n2*(ys-ym1)+n3*(zs-zm1)
120 ans =ans-gap
121 IF(ans>zero)GOTO 120
122 h(2)=t2*(ys-ym1)+t3*(zs-zm1)
123 h(2)=h(2)/xl
124 h(1)=one-h(2)
125 ss=h(2)-h(1)
126 IF(ss> onep05)GO TO 120
127 IF(ss<-onep05)GO TO 120
130
131 IF(nty==5)THEN
132
133 IF (stfn(ii)<zero) THEN
134 stif = zero
135 ELSE
136 stif=stf(l)
137 ENDIF
138 ELSE
139 stif=stf(l)*stfn(ii)/
max(em20,(stf(l)+stfn(ii)))
140 ENDIF
141 fni=ans*stif
142 fyi=n2*fni
143 fzi=n3*fni
144 imp=1
145
146
147
148 fsav(2)=fsav(2)+fyi*imast*dt12*ax
149 fsav(3)=fsav(3)+fzi*imast*dt12*ax
150
151 IF(iparit==0)THEN
152 DO 100 jj=1,2
153 j3=3*ix(jj)
154 j2=j3-1
155 e(j2)=e(j2)+fyi*h(jj)
156 e(j3)=e(j3)+fzi*h(jj)
157 100 CONTINUE
158 e(i2)=e(i2)-fyi
159 e(i3)=e(i3)-fzi
160 ELSE
161#include "lockon.inc"
162 niskyl = nisky
163 nisky = nisky + 3
164#include "lockoff.inc"
165 ipas = 1
166
167 IF(kdtint==0)THEN
168 fskyi(niskyl+1,1)= zero
169 fskyi(niskyl+1,2)= fyi*h(1)
170 fskyi(niskyl+1,3)= fzi*h(1)
171 fskyi(niskyl+1,4)= zero
172 isky(niskyl+1) = ix(1)
173
174 fskyi(niskyl+2,1)= zero
175 fskyi(niskyl+2,2)= fyi*h(2)
176 fskyi(niskyl+2,3)= fzi*h(2)
177 fskyi(niskyl+2,4)= zero
178 isky(niskyl+2) = ix(2)
179
180 fskyi(niskyl+3,1)= zero
181 fskyi(niskyl+3,2)= -fyi
182 fskyi(niskyl+3,3)= -fzi
183 fskyi(niskyl+3,4)= zero
184 isky(niskyl+3) = i
185 ELSE
186 fskyi(niskyl+1,1)= zero
187 fskyi(niskyl+1,2)= fyi*h(1)
188 fskyi(niskyl+1,3)= fzi*h(1)
189 fskyi(niskyl+1,4)= zero
190 fskyi(niskyl+1,5)= zero
191 isky(niskyl+1) = ix(1)
192
193 fskyi(niskyl+2,1)= zero
194 fskyi(niskyl+2,2)= fyi*h(2)
195 fskyi(niskyl+2,3)= fzi*h(2)
196 fskyi(niskyl+1,4)= zero
197 fskyi(niskyl+1,5)= zero
198 isky(niskyl+2) = ix(2)
199
200 fskyi(niskyl+3,1)= zero
201 fskyi(niskyl+3,2)= -fyi
202 fskyi(niskyl+3,3)= -fzi
203 fskyi(niskyl+1,4)= zero
204 fskyi(niskyl+1,5)= zero
205 isky(niskyl+3) = i
206 ENDIF
207 ENDIF
208
209 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0.AND.
210 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
211 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
212 fcont(2,ix(1)) =fcont(2,ix(1)) + fyi*h(1)
213 fcont(3,ix(1)) =fcont(3,ix(1)) + fzi*h(1)
214 fcont(2,ix(2)) =fcont(2,ix(2)) + fyi*h(2)
215 fcont(3,ix(2)) =fcont(3,ix(2)) + fzi*h(2)
216
217 fcont(2,i)=fcont(2,i)- fyi
218 fcont(3,i)=fcont(3,i)- fzi
219 ENDIF
220
221 IF(anim_v(12)+outp_v(12)+h3d_data%N_VECT_PCONT>0.AND.
222 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP) .OR.
223 . (manim>=4.AND.manim<=15).OR.h3d_data%MH3D/=0))THEN
224
225 fncont(2,ix(1)) =fncont(2,ix(1)) + fyi*h(1)
226 fncont(3,ix(1)) =fncont(3,ix(1)) + fzi*h(1)
227 fncont(2,ix(2)) =fncont(2,ix(2)) + fyi*h(2)
228 fncont(3,ix(2)) =fncont(3,ix(2)) + fzi*h(2)
229
230 fncont(2,i)=fncont(2,i)- fyi
231 fncont(3,i)=fncont(3,i)- fzi
232 ENDIF
233
234 IF(ibc/=0)
CALL ibcoff(ibc,icodt(i))
235
236 120 CONTINUE
237 IF(fric==zero)GO TO 500
238 IF(imp==0) THEN
239 irtlo(ii)=0
240 fric0(2,ii)=zero
241 fric0(3,ii)=zero
242 GO TO 500
243 ENDIF
244
245 lold=irtlo(ii)
246 IF(lold==0)THEN
247 irtlo(ii)=l
248 cst(1,ii)=ss
249 GO TO 500
250 ENDIF
251
252 ss0=cst(1,ii)
253 fti=fric0(1,ii)
254 ds=ss-ss0
255 anst=half*ds*xl
256 fmax=-
min(fric*fni,zero)
257 stfri=em01*stif
258 fti=fti + anst*stfri
259
260 IF(fti>fmax)THEN
261 fti=fmax
262 ELSE
263 IF(fti<-fmax)THEN
264 fti=-fmax
265 ELSE
266 fric0(1,ii)=fti
267 irtlo(ii)=l
268 cst(1,ii)=ss
269 ENDIF
270 ENDIF
271
272
273 fs = ptmax*areas(ii)/sqrt(three)
274 ft =fti
275 IF(fs>zero) THEN
276 IF(fti>fs)THEN
277 ft=fs
278 ELSEIF(fti<-fs)THEN
279 ft=-fs
280 ENDIF
281 ENDIF
282
283 fyi=t2*ft
284 fzi=t3*ft
285
286
287
288 fsav(5)=fsav(5)+fyi*imast*dt12*ax
289 fsav(6)=fsav(6)+fzi*imast*dt12*ax
290
291 IF(iparit==0)THEN
292 DO 400 jj=1,2
293 j3=3*ix(jj)
294 j2=j3-1
295 e(j2)=e(j2)+fyi*h(jj)
296 400 e(j3)=e(j3)+fzi*h(jj)
297 e(i2)=e(i2)-fyi
298 e(i3)=e(i3)-fzi
299 ELSE
300
301 IF(ipas==0) THEN
302#include "lockon.inc"
303 niskyl = nisky
304 nisky = nisky + 3
305#include "lockoff.inc"
306 IF(kdtint==0)THEN
307 fskyi(niskyl,1)= zero
308 fskyi(niskyl+1,2)= fyi*h(1)
309 fskyi(niskyl+1,3)= fzi*h(1)
310 fskyi(niskyl+1,4)= zero
311 isky(niskyl+1) = ix(1)
312
313 fskyi(niskyl+2,1)= zero
314 fskyi(niskyl+2,2)= fyi*h(2)
315 fskyi(niskyl+2,3)= fzi*h(2)
316 fskyi(niskyl+2,4)= zero
317 isky(niskyl+2) = ix(2)
318
319 fskyi(niskyl+3,1)= zero
320 fskyi(niskyl+3,2)= -fyi
321 fskyi(niskyl+3,3)= -fzi
322 fskyi(niskyl+3,4)= zero
323 isky(niskyl+3) = i
324 ELSE
325 fskyi(niskyl,1)= zero
326 fskyi(niskyl+1,2)= fyi*h(1)
327 fskyi(niskyl+1,3)= fzi*h(1)
328 fskyi(niskyl+1,4)= zero
329 fskyi(niskyl+1,5)= zero
330 isky(niskyl+1) = ix(1)
331
332 fskyi(niskyl+2,1)= zero
333 fskyi(niskyl+2,2)= fyi*h(2)
334 fskyi(niskyl+2,3)= fzi*h(2)
335 fskyi(niskyl+2,4)= zero
336 fskyi(niskyl+1,5)= zero
337 isky(niskyl+2) = ix(2)
338
339 fskyi(niskyl+3,1)= zero
340 fskyi(niskyl+3,2)= -fyi
341 fskyi(niskyl+3,3)= -fzi
342 fskyi(niskyl+3,4)= zero
343 fskyi(niskyl+1,5)= zero
344 isky(niskyl+3) = i
345 ENDIF
346 ELSE
347
348 fskyi(niskyl+1,2)= fskyi(niskyl+1,2)+fyi*h(1)
349 fskyi(niskyl+1,3)= fskyi(niskyl+1,3)+fzi*h(1)
350
351 fskyi(niskyl+2,2)= fskyi(niskyl+2,2)+fyi*h(2)
352 fskyi(niskyl+2,3)= fskyi(niskyl+2,3)+fzi*h(2)
353
354 fskyi(niskyl+3,2)= fskyi(niskyl+3,2)-fyi
355 fskyi(niskyl+3,3)= fskyi(niskyl+3,3)-fzi
356 ENDIF
357 ENDIF
358
359 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0.AND.
360 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP
361 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
362 fcont(2,ix(1)) =fcont(2,ix(1)) + fyi*h(1)
363 fcont(3,ix(1)) =fcont(3,ix(1)) + fzi*h(1)
364 fcont(2,ix(2)) =fcont(2,ix(2)) + fyi*h(2)
365 fcont(3,ix(2)) =fcont(3,ix(2)) + fzi*h(2)
366
367 fcont(2,i)=fcont(2,i)- fyi
368 fcont(3,i)=fcont(3,i)- fzi
369 ENDIF
370 IF(anim_v(12)+outp_v(12)+h3d_data%N_VECT_PCONT>0.AND.
371 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
372 . (manim>=4.AND.manim<=15).OR.h3d_data%MH3D/=0))THEN
373 ftcont(2,ix(1)) =ftcont(2,ix(1)) + fyi*h(1)
374 ftcont(3,ix(1)) =ftcont(3,ix(1)) + fzi*h(1)
375 ftcont(2,ix(2)) =ftcont(2,ix(2)) + fyi*h(2)
376 ftcont(3,ix(2)) =ftcont(3,ix(2)) + fzi*h(2)
377
378 ftcont(2,i)=ftcont(2,i)- fyi
379 ftcont(3,i)=ftcont(3,i)- fzi
380 ENDIF
381
382 500 CONTINUE
383 RETURN
subroutine ibcoff(ibc, icodt)