40
41
42
43
44
45
46
47
48
49#include "implicit_f.inc"
50
51
52
53#include "param_c.inc"
54#include "com01_c.inc"
55#include "com04_c.inc"
56#include "com06_c.inc"
57#include "com08_c.inc"
58#include "spmd_c.inc"
59
60
61
62 INTEGER NODPOR(*), WEIGHT(*), NPORGEO(*)
63 my_real geo(npropg,numgeo),ms(*),x(3,numnod),v(3,numnod),w(3,numnod),
64 . af(3,*),am(3,*),skew(lskew,*), afm(20)
65
66
67
68 INTEGER NAD,IG,I,K,JP,N,,NG, FR_POR(NSPMD+2)
69 my_real g1,g2,g3,r11,r12,r13,r21,r22,r23,r31,r32,r33,r(3,3),
70 . vx,vy,vz,
71 . rfm(3,numporl),f1(numporl),f2(numporl),
72 . f3(numporl),f4(numporl),f5(numporl),f6(numporl)
73 DOUBLE PRECISION FPOR6(6,6)
74
75
76
77 nad=0
78
79 DO ig=1,numgeo
80 IF(int(geo(12,ig)) /= 15)cycle
81 ng= int(geo(31,ig))
82 n = nporgeo(ig)
83 IF(ng == 0)cycle
84 IF(int(geo(20,ig)) == 1)GOTO 146
85
86 g1=geo(24,ig)
87 IF(int(geo(30,ig)) == 0)THEN
88 g2=geo(25,ig)
89 g3=geo(26,ig)
90 ELSEIF(dt1 > 0)THEN
91
92 g2=one/dt1
93 g3=one/dt1
94 ELSE
95 g2=zero
96 g3=zero
97 ENDIF
98 k=int(geo(27,ig))
99 r11=skew(1,k)
100 r12=skew(2,k)
101 r13=skew(3,k)
102 r21=skew(4,k)
103 r22=skew(5,k)
104 r23=skew(6,k)
105 r31=skew(7,k)
106 r32=skew(8,k)
107 r33=skew(9,k)
108 r(1,1)=r11*r11*g1+r21*r21*g2+r31*r31*g3
109 r(1,2)=r11*r12*g1+r21*r22*g2+r31*r32*g3
110 r(1,3)=r11*r13*g1+r21*r23*g2+r31*r33*g3
111 r(2,1)=r(1,2)
112 r(2,2)=r12*r12*g1+r22*r22*g2+r32*r32*g3
113 r(2,3)=r12*r13*g1+r22*r23*g2+r32*r33*g3
114 r(3,1)=r(1,3)
115 r(3,2)=r(2,3)
116 r(3,3)=r13*r13*g1+r23*r23*g2+r33*r33*g3
117 jrb=int(geo(29,ig))
118
119 IF(iale == 1)THEN
120 DO i=nad+1,nad+n
121 jp=nodpor(i)
122 vx=ms(jp)*(v(1,jp)-w(1,jp))+af(1,jp)*dt1
123 vy=ms(jp)*(v(2,jp)-w(2,jp))+af(2,jp)*dt1
124 vz=ms(jp)*(v(3,jp)-w(3,jp))+af(3,jp)*dt1
125 rfm(1,i)=r(1,1)*vx+r(1,2)*vy+r(1,3)*vz
126 rfm(2,i)=r(2,1)*vx+r(2,2)*vy+r(2,3)*vz
127 rfm(3,i)=r(3,1)*vx+r(3,2)*vy+r(3,3)*vz
128 ENDDO
129 ELSE
130 DO i=nad+1,nad+n
131 jp=nodpor(i)
132 vx=ms(jp)*v(1,jp)+af(1,jp)*dt1
133 vy=ms(jp)*v(2,jp)+af(2,jp)*dt1
134 vz=ms(jp)*v(3,jp)+af(3,jp)*dt1
135 rfm(1,i)=r(1,1)*vx+r(1,2)*vy+r(1,3)*vz
136 rfm(2,i)=r(2,1)*vx+r(2,2)*vy+r(2,3)*vz
137 rfm(3,i)=r(3,1)*vx+r(3,2)*vy+r(3,3)*vz
138 ENDDO
139 ENDIF
140
141 DO i=nad+1,nad+n
142 jp=nodpor(i)
143 af(1,jp)=af(1,jp)-rfm(1,i)
144 af(2,jp)=af(2,jp)-rfm(2,i)
145 af(3,jp)=af(3,jp)-rfm(3,i)
146 ENDDO
147
148
149
150 IF(jrb /= 0)THEN
151 DO i=nad+1,nad+n
152 jp=nodpor(i)
153 IF(weight(jp) == 1) THEN
154 f1(i-nad)=rfm(1,i)
155 f2(i-nad)=rfm(2,i)
156 f3(i-nad)=rfm(3,i)
157 f4(i-nad)=rfm(3,i)*(x(2,jp)-x(2,jrb))
158 + -rfm(2,i)*(x(3,jp)-x(3,jrb))
159 f5(i-nad)=rfm(1,i)*(x(3,jp)-x(3,jrb))
160 + -rfm(3,i)*(x(1,jp)-x(1,jrb))
161 f6(i-nad)=rfm(2,i)*(x(1,jp)-x(1,jrb))
162 + -rfm(1,i)*(x(2,jp)-x(2,jrb))
163 ELSE
164 f1(i-nad)=zero
165 f2(i-nad)=zero
166 f3(i-nad)=zero
167 f4(i-nad)=zero
168 f5(i-nad)=zero
169 f6(i-nad)=zero
170 ENDIF
171 ENDDO
172
173
174
175 DO k = 1, 6
176 fpor6(1,k) = zero
177 fpor6(2,k) = zero
178 fpor6(3,k) = zero
179 fpor6(4,k) = zero
180 fpor6(5,k) = zero
181 fpor6(6,k) = zero
182 END DO
189
190 IF(nspmd > 1) THEN
191
192 DO i = 1, nspmd
193 fr_por(i) = 1
194 END DO
195 fr_por(nspmd+1)=ng
196 fr_por(nspmd+2)=1
198 END IF
199
200 afm(1) = fpor6(1,1)+fpor6(1,2)+fpor6(1,3)+
201 + fpor6(1,4)+fpor6(1,5)+fpor6(1,6)
202 afm(2) = fpor6(2,1)+fpor6(2,2)+fpor6(2,3)+
203 + fpor6(2,4)+fpor6(2,5)+fpor6(2,6)
204 afm(3) = fpor6(3,1)+fpor6(3,2)+fpor6(3,3)+
205 + fpor6(3,4)+fpor6(3,5)+fpor6(3,6)
206 afm(4) = fpor6(4,1)+fpor6(4,2)+fpor6(4,3)+
207 + fpor6(4,4)+fpor6(4,5)+fpor6(4,6)
208 afm(5) = fpor6(5,1)+fpor6(5,2)+fpor6(5,3)+
209 + fpor6(5,4)+fpor6(5,5)+fpor6(5,6)
210 afm(6) = fpor6(6,1)+fpor6(6,2)+fpor6(6,3)+
211 + fpor6(6,4)+fpor6(6,5)+fpor6(6,6)
212
213 af(1,jrb)=af(1,jrb)+afm(1)
214 af(2,jrb)=af(2,jrb)+afm(2)
215 af(3,jrb)=af(3,jrb)+afm(3)
216 am(1,jrb)=am(1,jrb)+afm(4)
217 am(2,jrb)=am(2,jrb)+afm(5)
218 am(3,jrb)=am(3,jrb)+afm(6)
219
220
221
222 DO i=nad+1,nad+n
223 jp=nodpor(i)
224 IF(weight(jp) == 1) THEN
225 epor=epor+dt1*
226 + (rfm(1,i)*(v(1,jp)-w(1,jp))
227 + +rfm(2,i)*(v(2,jp)-w(2,jp))
228 + +rfm(3,i)*(v(3,jp)-w(3,jp)))
229 ENDIF
230 ENDDO
231
232 ELSE
233 DO i=nad+1,nad+n
234 jp=nodpor(i)
235 IF(weight(jp) == 1) THEN
236 epor=epor+dt1*
237 + (rfm(1,i)*v(1,jp)+rfm(2,i)*v(2,jp)+rfm(3,i)*v(3,jp))
238 ENDIF
239 ENDDO
240 ENDIF
241 146 nad=nad+n
242 enddo
243
244 RETURN
subroutine sum_6_float(jft, jlt, f, f6, n)
subroutine spmd_exch_fr6(fr, fs6, len)