34
35
36
37 use element_mod , only : nixt
38
39
40
41#include "implicit_f.inc"
42
43
44
45#include "mvsiz_p.inc"
46
47
48
49#include "param_c.inc"
50#include "vect01_c.inc"
51#include "scr12_c.inc"
52
53
54
55 INTEGER NC(NIXT,*),IPART(*), MXT(*), NC1(*), NC2(*)
57 . x(3,*),geo(npropg,*),pm(npropm,*),ms(*),stifn(*),
58 . v(3,*),partsav(20,*),mst(*),stifint(*),stt(*),
area(*),
59 . x1(*), x2(*), y1(*), y2(*), z1(*), z2(*)
60
61
62
63 INTEGER I, IP, I1, I2
65 . xx,yy,zz,xy,yz,zx
67 . al(mvsiz),
68 . rho(mvsiz),ems(mvsiz),
69 . sti(mvsiz)
70
71
72 DO i=lft,llt
73 rho(i) =pm(89,mxt(i))
74 ENDDO
75
76
77
78 DO i=lft,llt
79 al(i)=sqrt((x2(i)-x1(i))**2
80 . +(y2(i)-y1(i))**2
81 . +(z2(i)-z1(i))**2)
82 ems(i)=rho(i)*al(i)*
area(i)*half
83 ENDDO
84
85
86
87
88 DO i=lft,llt
89 mst(i)=ems(i)
90 ENDDO
91
92
93
94
95 DO i=lft,llt
96 sti(i) = pm(20,mxt(i)) *
area(i) / al(i)
97 stifn(nc1(i))=stifn(nc1(i))+sti(i)
98 stifn(nc2(i))=stifn(nc2(i))+sti(i)
99 ENDDO
100 IF (i7stifs /= 0) THEN
101 DO i=lft,llt
102 stt(i)=sti(i)
103 ENDDO
104 ENDIF
105
106 DO i=lft,llt
107 i1 = nc1(i)
108 i2 = nc2(i)
109
110 ip=ipart(i)
111 partsav(1,ip)=partsav(1,ip) + two*ems(i)
112 partsav(2,ip)=partsav(2,ip) + ems(i)*(x(1,i1)+x(1,i2))
113 partsav(3,ip)=partsav(3,ip) + ems(i)*(x(2,i1)+x(2,i2))
114 partsav(4,ip)=partsav(4,ip) + ems(i)*(x(3,i1)+x(3,i2))
115 xx = (x(1,i1)*x(1,i1)+x(1,i2)*x(1,i2))
116 xy = (x(1,i1)*x(2,i1)+x(1,i2)*x(2,i2))
117 yy = (x(2,i1)*x(2,i1)+x(2,i2)*x(2,i2))
118 yz = (x(2,i1)*x(3,i1)+x(2,i2)*x(3,i2))
119 zz = (x(3,i1)*x(3,i1)+x(3,i2)*x(3,i2))
120 zx = (x(3,i1)*x(1,i1)+x(3,i2)*x(1,i2))
121 partsav(5,ip) =partsav(5,ip) + ems(i) * (yy+zz)
122 partsav(6,ip) =partsav(6,ip) + ems(i) * (zz+xx)
123 partsav(7,ip) =partsav(7,ip) + ems(i) * (xx+yy)
124 partsav(8,ip) =partsav(8,ip) - ems(i) * xy
125 partsav(9,ip) =partsav(9,ip) - ems(i) * yz
126 partsav(10,ip)=partsav(10,ip) - ems(i) * zx
127
128 partsav(11,ip)=partsav(11,ip) + ems(i)*(v(1,i1)+v(1,i2))
129 partsav(12,ip)=partsav(12,ip) + ems(i)*(v(2,i1)+v(2,i2))
130 partsav(13,ip)=partsav(13,ip) + ems(i)*(v(3,i1)+v(3,i2))
131 partsav(14,ip)=partsav(14,ip) + half * ems(i) *
132 . (v(1,i1)*v(1,i1)+v(2,i1)*v(2,i1)+v(3,i1)*v(3,i1)
133 . +v(1,i2)*v(1,i2)+v(2,i2)*v(2,i2)+v(3,i2)*v(3,i2))
134 ENDDO
135
136 RETURN
subroutine area(d1, x, x2, y, y2, eint, stif0)