35
36
37
38
39
40#include "implicit_f.inc"
41
42
43
44#include "mvsiz_p.inc"
45
46
47
48#include "param_c.inc"
49#include "vect01_c.inc"
50#include "scr12_c.inc"
51#include "com04_c.inc"
52
53
54
55 INTEGER IPART(*),IGEO(NPROPGI,*),
56 . NC1(MVSIZ), NC2(MVSIZ), IMAT, MXG(MVSIZ)
57 INTEGER, INTENT(IN) :: NINTEMP
59 . geo(npropg,*), pm(npropm,*),
60 . stifn(*),stifr(*),v(3,*),partsav(20,*),msp(*),inp(*),
61 . stp(*),
62 . x1(mvsiz), x2(mvsiz),
63 . y1(mvsiz), y2(mvsiz),
64 . z1(mvsiz), z2(mvsiz),
area(mvsiz), al(mvsiz), strp(*),
65 . mcpp(*),temp(*)
66
67
68
69 INTEGER :: I,IP,I1,I2,J,IPY,IPZ,IPA,IGTYP,NIP
70 my_real :: xx,yy,zz,xy,yz,zx,sti, e, g, aa, bb, ari,ini,ryi,rzi,stir, dmp, sl2i, ll,temp0
71 my_real :: rho(mvsiz), ems(mvsiz),coefi(mvsiz),tin(mvsiz), tixx(mvsiz), tiyy(mvsiz), tizz(mvsiz),rhocp(mvsiz)
72 my_real :: facdt(mvsiz), phii(mvsiz),kphi(mvsiz),cst,phmax,phmin,fsh(mvsiz)
73
74 ipy = 200
75 ipz = 300
76 ipa = 400
77 temp0 = pm(79,imat)
78 DO i=lft,llt
79 igtyp = nint(geo(12,mxg(i)))
80 IF (igtyp == 18) THEN
81 rho(i) = pm(89,imat)
83 tiyy(i)= zero
84 tizz(i)= zero
85 nip = igeo(3,mxg(i))
86 DO j=1,nip
87 ari = geo(ipa+j,mxg(i))
88 ini = ari*ari*one_over_12
89 ryi = geo(ipy+j,mxg(i))
90 rzi = geo(ipz+j,mxg(i))
92 tiyy(i) = tiyy(i) + ini + ari * ryi*ryi
93 tizz(i) = tizz(i) + ini + ari * rzi*rzi
94 ENDDO
95 tixx(i) = tiyy(i) + tizz(i)
96 geo( 1,mxg(i)) =
area(i)
97 geo( 4,mxg(i)) = tixx(i)
98 geo( 2,mxg(i)) = tiyy(i)
99 geo(18,mxg(i)) = tizz(i)
100 ELSE
101 area(i)=geo(1,mxg(i))
102 tixx(i)=geo(4,mxg(i))
103 tiyy(i)=geo(2,mxg(i))
104 tizz(i)=geo(18,mxg(i))
105 rho(i) =pm(89,imat)
106 ENDIF
107 ENDDO
108
109
110
111 DO i=lft,llt
112 e = pm(20,imat)
113 g = pm(22,imat)
114 cst = six_over_5*e/g
115 bb =
max(tiyy(i),tizz(i),em30)
116 sl2i=
area(i)*al(i)**2 /bb
117 facdt(i) = one_over_12*sl2i
118 phmax = cst/facdt(i)
119 phmin =
min(tiyy(i),tizz(i))*phmax/bb
120 kphi(i) = (four+phmin)/(one+phmin)
121 phii(i) = kphi(i)/(one+facdt(i))
122 phii(i) =
max(one,phii(i))
123 fsh(i) = al(i)/(facdt(i)+cst)
124 fsh(i) =
max(one,fsh(i))
125 coefi(i) =one_over_12
126 ENDDO
127 IF (igtyp == 18) THEN
128 DO i=lft,llt
129 fsh(i) = one
130 kphi(i) =
max(one,twelve*facdt(i))
131 IF (kphi(i) > twelve ) coefi(i) =one
132 ENDDO
133 END IF
134
135
136
137 DO i=lft,llt
138 ll = onep2*al(i)
139 ems(i)=rho(i)*al(i)*
area(i)* half
140 tin(i)=onep2*ems(i)*ll**2*coefi(i) + rho(i)*(al(i)*half)
141 . *
max(tiyy(i),tizz(i))
142 IF (facdt(i)<one) tin(i)=phii(i)*tin(i)
143 tin(i)=
max(tin(i),rho(i)*al(i)/two*tixx(i))
144 ENDDO
145 IF( jthe > 0 ) THEN
146 DO i=lft,llt
147 rhocp(i) = pm(69,imat)
148 mcpp(i) = rhocp(i)*al(i)*
area(i)* half
149 ENDDO
150 ENDIF
151
152
153
154 DO i=lft,llt
155 msp(i) = ems(i)
156 ENDDO
157
158
159
160 DO i=lft,llt
161 inp(i) = tin(i)
162 ENDDO
163
164
165
166 IF(i7stifs/=0)THEN
167 DO i=lft,llt
168 e = pm(20,imat)
169 sti = e *
area(i) / al(i)
170 stp(i) = sti
171 ENDDO
172 ENDIF
173
174
175
176 DO i=lft,llt
177 e = pm(20,imat)
178 g = pm(22,imat)
179
180 dmp =
max(geo(16,mxg(i)),geo(17,mxg(i)))
181 dmp =dmp*sqrt(two)
182 aa =(sqrt(one +dmp*dmp)-dmp)
183 aa = al(i) * aa * aa
184 bb =
max(tiyy(i),tizz(i))
185 stir =
max(g*tixx(i),kphi(i)*e*bb) / aa
186 sti = fsh(i)*
area(i) * e / aa
187
188 stifn(nc1(i))=stifn(nc1(i))+sti
189 stifn(nc2(i))=stifn(nc2(i))+sti
190 stifr(nc1(i))=stifr(nc1(i))+stir
191 stifr(nc2(i))=stifr(nc2(i))+stir
192 strp(i)=stir
193 ENDDO
194
195 DO i=lft,llt
196 i1 = nc1(i)
197 i2 = nc2(i)
198
199 ip=ipart(i)
200 partsav(1,ip)=partsav(1,ip) + two*ems(i)
201 partsav(2,ip)=partsav(2,ip) + ems(i)*(x1(i)+x2(i))
202 partsav(3,ip)=partsav(3,ip) + ems(i)*(y1(i)+y2(i))
203 partsav(4,ip)=partsav(4,ip) + ems(i)*(z1(i)+z2(i))
204 xx = (x1(i)*x1(i)+x2(i)*x2(i))
205 xy = (x1(i)*y1(i)+x2(i)*y2(i))
206 yy = (y1(i)*y1(i)+y2(i)*y2(i))
207 yz = (y1(i)*z1(i)+y2(i)*z2(i))
208 zz = (z1(i)*z1(i)+z2(i)*z2(i))
209 zx = (z1(i)*x1(i)+z2(i)*x2(i))
210 partsav(5,ip) =partsav(5,ip) + two*tin(i) + ems(i) * (yy+zz)
211 partsav(6,ip) =partsav(6,ip) + two*tin(i) + ems(i) * (zz+xx)
212 partsav(7,ip) =partsav(7,ip) + two*tin(i) + ems(i) * (xx+yy)
213 partsav(8,ip) =partsav(8,ip) - ems(i) * xy
214 partsav(9,ip) =partsav(9,ip) - ems(i) * yz
215 partsav(10,ip)=partsav(10,ip) - ems(i) * zx
216
217 partsav(11,ip)=partsav(11,ip) + ems(i)*(v(1,i1)+v(1,i2))
218 partsav(12,ip)=partsav(12,ip) + ems(i)*(v(2,i1)+v(2,i2))
219 partsav(13,ip)=partsav(13,ip) + ems(i)*(v(3,i1)+v(3,i2))
220 partsav(14,ip)=partsav(14,ip) + half * ems(i) *
221 . (v(1,i1)*v(1,i1)+v(2,i1)*v(2,i1)+v(3,i1)*v(3,i1)
222 . +v(1,i2)*v(1,i2)+v(2,i2)*v(2,i2)+v(3,i2)*v(3,i2))
223 ENDDO
224 IF (jthe > 0) THEN
225 IF (nintemp > 0 ) THEN
226 DO i= lft,llt
227 IF(temp(nc1(i))== zero) temp(nc1(i)) = temp0
228 IF(temp(nc2(i))== zero) temp(nc2(i)) = temp0
229 ENDDO
230 ELSE
231 DO i=lft,llt
232 temp(nc1(i)) = pm(79,imat)
233 temp(nc2(i)) = pm(79,imat)
234 ENDDO
235 ENDIF
236 ENDIF
237
238 RETURN
subroutine area(d1, x, x2, y, y2, eint, stif0)