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