46
47
48
49#include "implicit_f.inc"
50
51
52
53#include "comlock.inc"
54#include "com01_c.inc"
55#include "com04_c.inc"
56
57#include "com08_c.inc"
58#include "param_c.inc"
59#include "task_c.inc"
60
61
62
63 INTEGER NSTRF(*), WEIGHT(*), IAD_ELEM(2,*), FR_ELEM(*),
64 . RG_CUT(*), IAD_CUT(NSPMD+2,*), FR_CUT(*),WEIGHT_MD(*),
65 . IOLDSECT, STABSEN,DIMFB,TABS(STABSEN)
67 . d(3,*), dr(3,*), v(3,*), vr(3,*), a(3,*), ar(3,*), ms(*),
68 . fsav(nthvki,*), secfcum(7,numnod,*), secbuf(*), in(*),
69 . fani(3,*), x(3,*), xsec(4,3,*)
70 DOUBLE PRECISION FBSAV6(12,6,DIMFB)
71 DOUBLE PRECISION, INTENT(INOUT) :: WFEXT
72
73
74
75 INTEGER J, I, K, II, I1, I2, N, KR1,KR2,KR3,K0,KR0,K1,K2,
76 . IFRL1, IFRL2, , ID_SEC,TYPE, LREC, NNOD,KR11,KR12, LENR,
77 . KR21,KR22,NBINTER, NN, LEN, KC, NSIZE, NNODG, SIZE, NNODT,
78 . ISECT
79 my_real dw, tt1, tt2, tt3, vi, dd, d1, d2,wfextl, aold, tnext, deltat,err(8), ff, fold,
alpha, aa, dv, wa(10)
80 real*4 r4
81
82
83
84 IF(nsect==0)RETURN
85
86
87
88 k0=nstrf(25)
89 DO i=1,nsect
90 IF(nstrf(k0)+nstrf(k0+14)>0)THEN
91
92
93
94
95 k2 = k0 + 30 + nstrf(k0+14)
96 IF(iroddl/=0)THEN
97 SIZE = 6
98 ELSE
99 SIZE = 3
100 END IF
101 IF (nspmd > 1) THEN
102 lenr = iad_elem(1,nspmd+1)-iad_elem(1,1)
104 1 nstrf(k2),secfcum(1,1,i),iad_elem,fr_elem,SIZE,
105 2 lenr ,nstrf(k0+6),weight)
106 END IF
107
108 k2 = k0 + 30 + nstrf(k0+14)
110 1 nstrf(k0+6),nstrf(k0+3),nstrf(k0+4),nstrf(k0+5),nstrf(k2),x,
111 2 v ,vr ,fsav(1,i),fani(1,1+2*(i-1)),secfcum(1,1,i),ms,
112 3 in ,nstrf(k0+26),xsec(1,1,i) )
113 ENDIF
114 k0=nstrf(k0+24)
115 ENDDO
116 IF(nstrf(1)==0.AND.nstrf(2)==0)RETURN
117
118
119
120 tnext = secbuf(5)
121 deltat = secbuf(1)
122 lrec = nstrf(6)
123 tt1 = secbuf(2)
124 tt2 = secbuf(3)
125 tt3 = secbuf(4)
126 IF(nstrf(1)>=1.AND.tnext<=tt)THEN
127 secbuf(5) = tnext + deltat
128
129 k0 = nstrf(25)
130
131 kc = 1
132 IF(ispmd==0 .AND. ioldsect == 1) THEN
134 r4 = tt
138 ENDIF
139 DO n=1,nsect
140 TYPE=nstrf(k0)
141 IF(ispmd==0 .AND. ioldsect /= 1 .AND. TYPE >= 1) then
143 r4 = tt
147 ENDIF
148 nbinter = nstrf(k0+14)
149 k1 = k0+30
150 k2=k1+nbinter
151 nnod = nstrf(k0+6)
152 TYPE=nstrf(k0)
153 IF(type>=1)THEN
154
155 id_sec=nstrf(k0+23)
156 IF(nspmd==1) THEN
160 ELSEIF(ispmd==0) THEN
163 nnodg = iad_cut(nspmd+2,n)
165 ENDIF
166 IF(iroddl/=0)THEN
167
168
169
170 IF(nspmd>1) THEN
171 IF (ispmd==0) THEN
172 nsize = iad_cut(nspmd+1,n)
173 nnodg = iad_cut(nspmd+2,n)
174 ELSE
175 nsize = 0
176 nnodg = 0
177 ENDIF
179 1 nnod ,nstrf(k2),d ,dr ,rg_cut(kc),
180 2 iad_cut(1,n),nsize ,nnodg,weight,2 )
181 ELSE
182 DO i=1,nnod
183 r4 = d(1,nstrf(k2+i-1))
185 r4 = d(2,nstrf(k2+i-1))
187 r4 = d(3,nstrf(k2+i-1))
189 r4 = dr(1,nstrf(k2+i-1))
191 r4 = dr(2,nstrf(k2+i-1))
193 r4 = dr(3,nstrf(k2+i-1))
195 ENDDO
196 ENDIF
197 ELSE
198
199
200
201 IF(nspmd>1) THEN
202 IF (ispmd==0) THEN
203 nsize = iad_cut(nspmd+1,n)
204 nnodg = iad_cut(nspmd+2,n)
205 ELSE
206 nsize = 0
207 nnodg = 0
208 ENDIF
210 1 nnod ,nstrf(k2),d ,dr ,rg_cut(kc),
211 2 iad_cut(1,n),nsize ,nnodg,weight,1 )
212 ELSE
213 DO i=1,nnod
214 r4 = d(1,nstrf(k2+i-1))
216 r4 = d(2,nstrf(k2+i-1))
218 r4 = d(3,nstrf(k2+i-1))
220 r4 = zero
224 ENDDO
225 ENDIF
226 ENDIF
227 ENDIF
228 IF(type>=2)THEN
229
230 IF(iroddl/=0)THEN
231
232
233
234 IF(nspmd>1) THEN
235 IF (ispmd==0) THEN
236 nsize = iad_cut(nspmd+1,n)
237 nnodg = iad_cut(nspmd+2,n)
238 ELSE
239 nsize = 0
240 nnodg = 0
241 ENDIF
243 1 nnod ,nstrf(k2),secfcum(1,1,n),rg_cut(kc),iad_cut(1,n),
244 2 nsize,nnodg ,weight ,2 )
245 ELSE
246 DO i=1,nnod
247 r4 = secfcum(1,nstrf(k2+i-1),n)
249 r4 = secfcum(2,nstrf(k2+i-1),n)
251 r4 = secfcum(3,nstrf(k2+i-1),n)
253 r4 = secfcum(5,nstrf(k2+i-1),n)
255 r4 = secfcum(6,nstrf(k2+i-1),n)
257 r4 = secfcum(7,nstrf(k2+i-1),n)
259 ENDDO
260 ENDIF
261 ELSE
262
263
264
265 IF(nspmd>1) THEN
266 IF (ispmd==0) THEN
267 nsize = iad_cut(nspmd+1,n)
268 nnodg = iad_cut(nspmd+2,n)
269 ELSE
270 nsize = 0
271 nnodg = 0
272 ENDIF
274 1 nnod ,nstrf(k2),secfcum(1,1,n),rg_cut(kc),iad_cut(1,n),
275 2 nsize,nnodg ,weight ,1 )
276 ELSE
277 DO i=1,nnod
278 r4 = secfcum(1,nstrf(k2+i-1),n)
280 r4 = secfcum(2,nstrf(k2+i-1),n)
282 r4 = secfcum(3,nstrf(k2+i-1),n)
284 r4 = zero
288 ENDDO
289 ENDIF
290 ENDIF
291 ENDIF
292
293 k0 = nstrf(k0+24)
294 IF(type>=1) kc = kc + nnod
295 ENDDO
297 ENDIF
298
299
300
301
302 IF(nstrf(2)>=1)THEN
303
304
305
306 ifrl1=nstrf(7)
307 ifrl2=mod(ifrl1+1,2)
308 k0 = nstrf(25)
309 kr0 = nstrf(26)
310 DO n=1,nsect
311 nnod = nstrf(k0+6)
312 TYPE=nstrf(k0)
313 nbinter = nstrf(k0+14)
314 IF(type>=101)THEN
315 k2 = k0 + 30 + nbinter
316 kr1 = kr0 + 10
317 kr2 = kr1 + 12*nnod
318 kr3 = kr2 + 12*nnod
319 kr21 = kr2 + ifrl2*6*nnod
320 kr22 = kr2 + ifrl1*6*nnod
321 err(4) = zero
322 err(8) = zero
323 DO k=1,3
324 err(k) = zero
325 err(k+4) = zero
326 DO i=1,nnod
327 ii = nstrf(k2+i-1)
328 IF(weight_md(ii)==1)THEN
329 fold = secfcum(k,ii,n)
330 d2 = secbuf(kr22+6*i-7+k)
331 d1 = secbuf(kr21+6*i-7+k)
332 ff = (tt*(d2-d1)+tt2*d1-tt1*d2) / (tt2-tt1)
333 err(k) = err(k) + (ff - fold)
334 err(4) = err(4) + (ff - fold)**2
335 END IF
336 ENDDO
337 IF(iroddl/=0)THEN
338 DO i=1,nnod
339 ii = nstrf(k2+i-1)
340 IF(weight_md(ii)==1)THEN
341 fold = secfcum(k+4,ii,n)
342 d2 = secbuf(kr22+6*i-4+k)
343 d1 = secbuf(kr21+6*i-4+k)
344 ff = (tt*(d2-d1)+tt2*d1-tt1*d2) / (tt2-tt1)
345 err(k+4) = err(k+4) + (ff - fold)
346 err(8) = err(8) + (ff - fold)**2
347 END IF
348 ENDDO
349 ENDIF
350 ENDDO
351 fsav(11,n) = fsav(11,n) + err(1)*dt12
352 fsav(12,n) = fsav(12,n) + err(2)*dt12
353 fsav(13,n) = fsav(13,n) + err(3)*dt12
354 fsav(14,n) = err(4)
355 fsav(16,n) = fsav(16,n) + err(5)*dt12
356 fsav(17,n) = fsav(17,n) + err(6)*dt12
357 fsav(18,n) = fsav(18,n) + err(7)*dt12
358 fsav(19,n) = err(8)
359 ENDIF
360 kr0 = nstrf(k0+25)
361 k0 = nstrf(k0+24)
362 ENDDO
363 ENDIF
364
365
366
367
368 IF(nspmd==1) THEN
370 ELSE
371 nnodt = 0
372 IF(ispmd==0) THEN
373 k0 = nstrf(25)
374 DO i = 1, nsect
375 IF(nstrf(k0)>=100) nnodt = nnodt + iad_cut(nspmd+2,i)
376 k0 = nstrf(k0+24)
377 END DO
378 END IF
379
380
381
383 END IF
384
385
386
387
388 tt1 = secbuf(2)
389 tt2 = secbuf(3)
390 tt3 = secbuf(4)
391 IF(nstrf(2)>=1)THEN
392 ifrl1=nstrf(7)
393 ifrl2=mod(ifrl1+1,2)
394 k0 = nstrf(25)
395 kr0 = nstrf(26)
396 DO n=1,nsect
397 nnod = nstrf(k0+6)
398 TYPE=nstrf(k0)
399 nbinter = nstrf(k0+14)
400 alpha = 1.-secbuf(kr0+2)
401 IF(type>=100.AND.
alpha/=0.0)
THEN
402 k2 = k0 + 30 + nbinter
403 kr1 = kr0 + 10
404 kr2 = kr1 + 12*nnod
405 kr3 = kr2 + 12*nnod
406 kr11 = kr1 + ifrl2*6*nnod
407 kr12 = kr1 + ifrl1*6*nnod
408 dw = secbuf(kr0+1)
409 IF(ispmd==0) THEN
410 wfextl=dw*dt1
411 ELSE
412 wfextl=zero
413 ENDIF
414 wfext=wfext + wfextl
415 dw=zero
416 err(4) = zero
417 err(8) = zero
418 DO k=1,3
419 err(k) = zero
420 err(k+4) = zero
421 DO i=1,nnod
422 ii = nstrf(k2+i-1)
423 d2 = secbuf(kr12+6*i-7+k)
424 d1 = secbuf(kr11+6*i-7+k)
425 dd = ((tt+dt2)*(d2-d1)+tt2*d1-tt1*d2) / (tt2-tt1)
426 vi = (dd-d(k,ii))/dt2
427 aa =
alpha*((vi-v(k,ii))/dt12 - a(k,ii))
428 a(k,ii) = a(k,ii) + aa
429 IF(weight(ii)==1) THEN
430 dv = dt12*a(k,ii)
431 dw = dw + half*(v(k,ii)+half*dv)*ms(ii)*aa
432 err(k)=err(k)+weight_md(ii)*ms(ii)*(vi-v(k,ii)-dv)
433 err(4)=err(4)
434 . + weight_md(ii)*ms(ii)*(vi**2-(v(k,ii)+dv)**2)
435 ENDIF
436 ENDDO
437 IF(iroddl/=0)THEN
438 DO i=1,nnod
439 ii = nstrf(k2+i-1)
440 d2 = secbuf(kr12+6*i-4+k)
441 d1 = secbuf(kr11+6*i-4+k)
442 dd = ((tt+dt2)*(d2-d1)+tt2*d1-tt1*d2) / (tt2-tt1)
443 vi = (dd-dr(k,ii))/dt2
444 aa =
alpha*((vi-vr(k,ii))/dt12 - ar(k,ii))
445 ar(k,ii) = ar(k,ii) + aa
446 IF(weight(ii)==1) THEN
447 dv = dt12*ar(k,ii)
448 dw = dw + half*(vr(k,ii)+half*dv)*in(ii)*aa
449 err(k+4)=err(k+4)
450 . + weight_md(ii)*in(ii)*(vi-vr(k,ii) - dv)
451 err(8)=err(8)
452 . + weight_md(ii)*in(ii)*(vi**2-(vr(k,ii)+dv)**2)
453 ENDIF
454 ENDDO
455 ENDIF
456 ENDDO
457 wfextl=wfextl + dt1*dw
458 wfext=wfext + dt1*dw
459 secbuf(kr0+1) = dw
460
461
462
463 IF(nspmd>1) THEN
464 IF (ispmd==0) THEN
465 wa(1) = secbuf(kr0+1)
466 wa(2) = secbuf(kr0+3)
467 wa(3) = secbuf(kr0+4)
468 len = 3
470 secbuf(kr0+1) = wa(1)
471 secbuf(kr0+3) = wa(2)
472 secbuf(kr0+4) = wa(3)
473
474 ELSE
475 len = 3
476 wa(1) = secbuf(kr0+1)
477 wa(2) = secbuf(kr0+3)
478 wa(3) = secbuf(kr0+4)
480 secbuf(kr0+1) = zero
481 secbuf(kr0+3) = zero
482 secbuf(kr0+4) = zero
483 ENDIF
484 ENDIF
485
486 fsav(22,n) = err(1)
487 fsav(23,n) = err(2)
488 fsav(24,n) = err(3)
489 fsav(25,n) = half*err(4)
490 fsav(26,n) = err(5)
491 fsav(27,n) = err(6)
492 fsav(28,n) = err(7)
493 fsav(29,n) = half*err(8)
494 fsav(30,n) = fsav(30,n) + wfextl + secbuf(kr0+4)
495 isect=0
496 IF(stabsen/=0) isect=tabs(n+1)-tabs(n)
497 IF(isect/=0) THEN
498 fbsav6(7,2:6,isect) = zero
499 fbsav6(7,1,isect)=fsav(30,n)
500 ENDIF
501 ENDIF
502 kr0 = nstrf(k0+25)
503 k0 = nstrf(k0+24)
504 ENDDO
505 ENDIF
506
507 k0=nstrf(25)
508 DO i=1,nsect
509 nnod = nstrf(k0+6)
510 k2 = k0 + 30 + nstrf(k0+14)
511 DO k = 1, nnod
512 n = nstrf(k2+k-1)
513 secfcum(1,n,i)=zero
514 secfcum(2,n,i)=zero
515 secfcum(3,n,i)=zero
516 secfcum(5,n,i)=zero
517 secfcum(6,n,i)=zero
518 secfcum(7,n,i)=zero
519 ENDDO
520 k0=nstrf(k0+24)
521 ENDDO
522
523 RETURN
subroutine section(nnod, n1, n2, n3, nstrf, x, v, vr, fsav, fopta, secfcum, ms, in, ifram, xsec)
subroutine section_read(ttt, nstrf, secbuf)
subroutine section_readp(ttt, nstrf, secbuf, nnodt, iad_cut, fr_cut)
subroutine spmd_wrt_cutf(nnod, nstrf, secfcum, rg_cut, iad_cut, nsize, nnodg, weight, iflg)
subroutine spmd_wrt_cutd(nnod, nstrf, d, dr, rg_cut, iad_cut, nsize, nnodg, weight, iflg)
subroutine spmd_exch_cut(nstrf, secfcum, iad_elem, fr_elem, size, lenr, nnod, weight)
subroutine spmd_glob_dsum9(v, len)
void write_i_c(int *w, int *len)
void write_r_c(float *w, int *len)