28 SUBROUTINE cbadef(JFT,JLT,NG,VCORE,AREA,CDET,VQN,VQ,VJFI,
29 1 VXYZ,RXYZ,VDEF,VNRM,VASTN,
30 2 HX,HY,VETA,VKSI,BM,BMF,BF,BC,TC,NPLAT,IPLAT,
40#include "implicit_f.inc"
48 INTEGER NPG,NG,JFT,JLT,NPLAT,IPLAT(*),ISROT
51 . rxyz(mvsiz,8),vcore(mvsiz,12),vxyz(mvsiz,12),
52 . vqn(mvsiz,9,4),vksi(4,4),veta(4,4),
53 . bm(mvsiz,36),bmf(mvsiz,36),bf(mvsiz,24),bc(mvsiz,40),hx(mvsiz,4),hy(mvsiz,4)
55 . vnrm(mvsiz,12),vastn(mvsiz,16),vjfi(mvsiz,3,2,4),
56 . vq(mvsiz,3,3,4),vdef(mvsiz,8),dt1,
area(*)
58 . cdet(*),tc(mvsiz,2,2),brz(mvsiz,4,4)
62 INTEGER PT,PT00,PT0,I,J,EP,IUN,NG1,K
64 . VPG(2,NPG),PG1,PG,THK,DETJ,DET,
65 . tfn(3,2),bcx,bcy,bxy(3),byx(3),
66 . v1(2),v2(2),rv1,rv2,
67 . c1,c2,vt1,vt2,bc1,bc2,vb1,
68 . vjf1(2,3),vjf(3,3),tbi(2,2),tbc(2,2),v11(4)
70 . a_1,c11,c12,c21,c22,cc,beta1,ksi1,ksiy1,beta2,ksi2,ksiy2
71 parameter(pg=.577350269189626)
72 parameter(pg1=-.577350269189626)
77 DATA vpg/pg1,pg1,pg,pg1,pg,pg,pg1,pg/
79#include "vectorize.inc"
84 bm(ep,1)=vcore(ep,1)+hx(ep,ng)*vcore(ep,5)
85 bm(ep,2)=vcore(ep,2)+hx(ep,ng)*vcore(ep,6)
86 bm(ep,3)=hx(ep,ng)*fourth
88 bm(ep,5)=vcore(ep,3)+hy(ep,ng)*vcore(ep,5)
89 bm(ep,6)=vcore(ep,4)+hy(ep,ng)*vcore(ep,6)
90 bm(ep,7)=hy(ep,ng)*fourth
93 vdef(ep,1)=bm(ep,1)*vxyz(ep,1)+bm(ep,2)*vxyz(ep,4)
94 1 +bm(ep,3)*vxyz(ep,7)
95 vdef(ep,2)=bm(ep,5)*vxyz(ep,2)+bm(ep,6)*vxyz(ep,5)
96 1 +bm(ep,7)*vxyz(ep,8)
98 vdef(ep,6)=bm(ep,1)*rxyz(ep,2)+bm(ep,2)*rxyz(ep,4)
99 1 +bm(ep,3)*rxyz(ep,6)
100 vdef(ep,7)=-(bm(ep,5)*rxyz(ep,1)+bm(ep,6)*rxyz(ep,3)
101 1 +bm(ep,7)*rxyz(ep,5))
102 vdef(ep,8)=-(bm(ep,1)*rxyz(ep,1)+bm(ep,2)*rxyz(ep,3)
103 1 +bm(ep,3)*rxyz(ep,5))
104 2 +bm(ep,5)*rxyz(ep,2)+bm(ep,6)*rxyz(ep,4)
105 3 +bm(ep,7)*rxyz(ep,6)
108 a_1 = 0.25/
max(cdet(ep),em20)
109 c11=(vcore(ep,10)+vcore(ep,12)*vpg(1,ng))*a_1
110 c12=-(vcore(ep,8)+vcore(ep,12)*vpg(2,ng))*a_1
111 c21=-(vcore(ep,9)+vcore(ep,11)*vpg(1,ng))*a_1
112 c22=( vcore(ep,7)+vcore(ep,11)*vpg(2,ng))*a_1
113 beta1=vcore(ep,12)+vcore(ep,8)*vpg(2,ng)
114 ksi1=vcore(ep,12)+vcore(ep,10)*vpg(1,ng)
115 beta2=vcore(ep,11)+vcore(ep,7)*vpg(2,ng)
116 ksi2=vcore(ep,11)+vcore(ep,9)*vpg(1,ng)
120 bc(ep,3)= beta1*c11+ksi1*c12
121 bc(ep,4)= beta1*c21+ksi1*c22
122 bc(ep,5)= -beta2*c11-ksi2*c12
123 bc(ep,6)= -beta2*c21-ksi2*c22
127 bc(ep,9)= beta1*c11-ksi1*c12
128 bc(ep,10)=beta1*c21-ksi1*c22
129 bc(ep,11)=-beta2*c11+ksi2*c12
130 bc(ep,12)=-beta2*c21+ksi2*c22
132 beta1=vcore(ep,8)+vcore(ep,12)*vpg(2,ng)
133 ksi1=vcore(ep,10)+vcore(ep,12)*vpg(1,ng)
134 beta2=vcore(ep,7)+vcore(ep,11)*vpg(2,ng)
135 ksi2=vcore(ep,9)+vcore(ep,11)*vpg(1,ng)
138 bc(ep,13)=c11*vpg(2,ng)+c12*vpg(1,ng)
139 bc(ep,14)=c21*vpg(2,ng)+c22*vpg(1,ng)
140 bc(ep,15)=-beta1*c11-ksi1*c12
141 bc(ep,16)=-beta1*c21-ksi1*c22
142 bc(ep,17)=c11*beta2+c12*ksi2
143 bc(ep,18)=c21*beta2+c22*ksi2
152 vdef(ep,4)=bc(ep,1)*vxyz(ep,3)+bc(ep,7)*vxyz(ep,6)
153 1 +bc(ep,13)*vxyz(ep,9)+bc(ep,3)*rxyz(ep,1)
154 2 +bc(ep,9)*rxyz(ep,3)+bc(ep,15)*rxyz(ep,7)
155 3 +bc(ep,5)*rxyz(ep,2)+bc(ep,11)*rxyz(ep,4)
156 4 +bc(ep,17)*rxyz(ep,8)
157 vdef(ep,5)=bc(ep,2)*vxyz(ep,3)+bc(ep,8)*vxyz(ep,6)
158 1 +bc(ep,14)*vxyz(ep,9)+bc(ep,4)*rxyz(ep,1)
159 2 +bc(ep,10)*rxyz(ep,3)+bc(ep,16)*rxyz(ep,7)
160 3 +bc(ep,6)*rxyz(ep,2)+bc(ep,12)*rxyz(ep,4)
161 4 +bc(ep,18)*rxyz(ep,8)
164#include "vectorize.inc"
168 vdef(ep,3)=bm(ep,1)*vxyz(ep,2)+bm(ep,2)*vxyz(ep,5)
169 1 +bm(ep,3)*vxyz(ep,8)
170 1 +bm(ep,5)*vxyz(ep,1)+bm(ep,6)*vxyz(ep,4)
171 1 +bm(ep,7)*vxyz(ep,7)
175#include "vectorize.inc"
181 tfn(1,1)=vksi(1,ng)*vqn(ep,7,1)+vksi(2,ng)*vqn(ep,7,2)
182 1 +vksi(3,ng)*vqn(ep,7,3)+vksi(4,ng)*vqn(ep,7,4)
183 tfn(2,1)=vksi(1,ng)*vqn(ep,8,1)+vksi(2,ng)*vqn(ep,8,2)
184 1 +vksi(3,ng)*vqn(ep,8,3)+vksi(4,ng)*vqn(ep,8,4)
185 tfn(3,1)=vksi(1,ng)*vqn(ep,9,1)+vksi(2,ng)*vqn(ep,9,2)
186 1 +vksi(3,ng)*vqn(ep,9,3)+vksi(4,ng)*vqn(ep,9,4)
187 tfn(1,2)=veta(1,ng)*vqn(ep,7,1)+veta(2,ng)*vqn(ep,7,2)
188 1 +veta(3,ng)*vqn(ep,7,3)+veta(4,ng)*vqn(ep,7,4)
189 tfn(2,2)=veta(1,ng)*vqn(ep,8,1)+veta(2,ng)*vqn(ep,8,2)
190 1 +veta(3,ng)*vqn(ep,8,3)+veta(4,ng)*vqn(ep,8,4)
191 tfn(3,2)=veta(1,ng)*vqn(ep,9,1)+veta(2,ng)*vqn(ep,9,2)
192 1 +veta(3,ng)*vqn(ep,9,3)+veta(4,ng)*vqn(ep,9,4)
196 tbi(2,2)=vjfi(ep,1,1,ng)*tfn(1,1)+vjfi(ep,2,1,ng)*tfn(2,1)
197 1 + vjfi(ep,3,1,ng)*tfn(3,1)
198 tbi(2,1)=vjfi(ep,1,2,ng)*tfn(1,1)+vjfi(ep,2,2,ng)*tfn(2,1)
199 1 + vjfi(ep,3,2,ng)*tfn(3,1)
200 tbi(1,2)=vjfi(ep,1,1,ng)*tfn(1,2)+vjfi(ep,2,1,ng)*tfn(2,2)
201 1 + vjfi(ep,3,1,ng)*tfn(3,2)
202 tbi(1,1)=vjfi(ep,1,2,ng)*tfn(1,2)+vjfi(ep,2,2,ng)*tfn(2,2)
203 1 + vjfi(ep,3,2,ng)*tfn(3,2)
205 thk =-(tbi(1,1)+tbi(2,2))
212 tc(ep,1,1)=vjfi(ep,1,1,ng)*vq(ep,1,1,ng)+vjfi(ep,2,1,ng)
213 + *vq(ep,2,1,ng)+ vjfi(ep,3,1,ng)*vq(ep,3,1,ng)
214 tc(ep,2,1)=vjfi(ep,1,2,ng)*vq(ep,1,1,ng)+vjfi(ep,2,
215 + *vq(ep,2,1,ng)+ vjfi(ep,3,2,ng)*vq(ep,3,1,ng)
216 tc(ep,1,2)=vjfi(ep,1,1,ng)*vq(ep,1,2,ng)+vjfi(ep,2,1,ng)
217 + *vq(ep,2,2,ng)+ vjfi(ep,3,1,ng)*vq(ep,3,2,ng)
218 tc(ep,2,2)=vjfi(ep,1,2,ng)*vq(ep,1,2,ng)+vjfi(ep,2,2,ng)
219 + *vq(ep,2,2,ng)+ vjfi(ep,3,2,ng)*vq(ep,3,2,ng)
223 tbc(1,1)=tbi(1,1)*tc(ep,1,1)+tbi(1,2)*tc(ep,2,1)
224 tbc(2,1)=tbi(2,1)*tc(ep,1,1)+tbi(2,2)*tc(ep,2,1)
225 tbc(1,2)=tbi(1,1)*tc(ep,1,2)+tbi(1,2)*tc(ep,2,2)
226 tbc(2,2)=tbi(2,1)*tc(ep,1,2)+tbi(2,2)*tc(ep,2,2)
240 vt1=vq(ep,1,1,ng)*vxyz(ep,1)+vq(ep,2,1,ng)*vxyz(ep,2)
241 1 +vq(ep,3,1,ng)*vxyz(ep,3)
243 vt2=vq(ep,1,2,ng)*vxyz(ep,1)+vq(ep,2,2,ng)*vxyz(ep,2)
244 1 +vq(ep,3,2,ng)*vxyz(ep,3)
246 c1=vksi(1,ng)*tc(ep,1,1)+veta(1,ng)*tc(ep,2,1)
247 c2=vksi(1,ng)*tc(ep,1,2)+veta(1,ng)*tc(ep,2,2)
249 bc1=vksi(1,ng)*tbc(1,1)+veta(1,ng)*tbc(2,1)
250 bc2=vksi(1,ng)*tbc(1,2)+veta(1,ng)*tbc(2,2)
253 bxy(1)=vq(ep,1,1,ng)*c2
254 bxy(2)=vq(ep,2,1,ng)*c2
255 bxy(3)=vq(ep,3,1,ng)*c2
256 byx(1)=vq(ep,1,2,ng)*c1
257 byx(2)=vq(ep,2,2,ng)*c1
258 byx(3)=vq(ep,3,2,ng)*c1
259 brz(i,1,1)=-bxy(1)+byx(1)
260 brz(i,2,1)=-bxy(2)+byx(2)
261 brz(i,3,1)=-bxy(3)+byx(3)
262 bm(ep,3)=bxy(1)+byx(1)
263 bm(ep,6)=bxy(2)+byx(2)
264 bm(ep,9)=bxy(3)+byx(3)
265 vdef(ep,3)= c1*vt2+c2*vt1
270 bm(ep,1)=vq(ep,1,1,ng)*c1
271 bmf(ep,1)=thk*bm(ep,1)+vq(ep,1,1,ng)*bc1
273 bm(ep,2)=vq(ep,1,2,ng)*c2
274 bmf(ep,2)=thk*bm(ep,2)+vq(ep,1,2,ng)*bc2
276 bmf(ep,3)=thk*bm(ep,3)+vq(ep,1,1,ng)*bc2+vq(ep,1,2,ng)*bc1
279 bm(ep,4)=vq(ep,2,1,ng)*c1
280 bmf(ep,4)=thk*bm(ep,4)+vq(ep,2,1,ng)*bc1
281 bm(ep,5)=vq(ep,2,2,ng)*c2
282 bmf(ep,5)=thk*bm(ep,5)+vq(ep,2,2,ng)*bc2
283 bmf(ep,6)=thk*bm(ep,6)+vq(ep,2,1,ng)*bc2+vq(ep,2,2,ng)*bc1
286 bm(ep,7)=vq(ep,3,1,ng)*c1
287 bmf(ep,7)=thk*bm(ep,7)+vq(ep,3,1,ng)*bc1
288 bm(ep,8)=vq(ep,3,2,ng)*c2
289 bmf(ep,8)=thk*bm(ep,8)+vq(ep,3,2,ng)*bc2
290 bmf(ep,9)=thk*bm(ep,9)+vq(ep,3,1,ng)*bc2+vq(ep,3,2,ng)*bc1
300 vdef(ep,8)= bc1*vt2+bc2*vt1
303 v1(2)=vq(ep,1,1,ng)*vqn(ep,1,1)+vq(ep,2,1,ng)*vqn(ep,2,1)
304 1 +vq(ep,3,1,ng)*vqn(ep,3,1)
305 v1(1)=-(vq(ep,1,1,ng)*vqn(ep,4,1)+vq(ep,2,1,ng)*vqn
306 1 +vq(ep,3,1,ng)*vqn(ep,6,1))
307 v2(2)=vq(ep,1,2,ng)*vqn(ep,1,1)+vq(ep,2,2,ng)*vqn(ep,2,1)
308 1 +vq(ep,3,2,ng)*vqn(ep,3,1
309 v2(1)=-(vq(ep,1,2,ng)*vqn(ep,4,1)+vq(ep,2,2,ng)*vqn(ep,5,1)
310 1 +vq(ep,3,2,ng)*vqn(ep,6,1))
312 rv1=v1(1)*rxyz(ep,1)+v1(2)*rxyz(ep,2)
314 rv2=v2(1)*rxyz(ep,1)+v2(2)*rxyz(ep,2)
319 bf(ep,3)=v1(1)*c2+v2(1)*c1
322 bf(ep,6)=v1(2)*c2+v2(2)*c1
324 vdef(ep,6)= vdef(ep,6)+c1*rv1
325 vdef(ep,7)= vdef(ep,7)+c2*rv2
326 vdef(ep,8)= vdef(ep,8)+c1*rv2+c2*rv1
330 vt1=vq(ep,1,1,ng)*vxyz(ep,4)+vq(ep,2,1,ng)*vxyz(ep,5)
331 1 +vq(ep,3,1,ng)*vxyz(ep,6)
333 vt2=vq(ep,1,2,ng)*vxyz(ep,4)+vq(ep,2,2,ng)*vxyz(ep,5)
334 1 +vq(ep,3,2,ng)*vxyz(ep,6)
337 c1=vksi(2,ng)*tc(ep,1,1)+veta(2,ng)*tc(ep,2,1)
338 c2=vksi(2,ng)*tc(ep,1,2)+veta(2,ng)*tc(ep,2,2)
340 bc1=vksi(2,ng)*tbc(1,1)+veta(2,ng)*tbc(2,1)
341 bc2=vksi(2,ng)*tbc(1,2)+veta(2,ng)*tbc(2,2)
345 bxy(1)=vq(ep,1,1,ng)*c2
346 bxy(2)=vq(ep,2,1,ng)*c2
347 bxy(3)=vq(ep,3,1,ng)*c2
348 byx(1)=vq(ep,1,2,ng)*c1
349 byx(2)=vq(ep,2,2,ng)*c1
350 byx(3)=vq(ep,3,2,ng)*c1
351 brz(i,1,2)=-bxy(1)+byx(1)
352 brz(i,2,2)=-bxy(2)+byx(2)
353 brz(i,3,2)=-bxy(3)+byx(3)
354 bm(ep,12)=bxy(1)+byx(1)
355 bm(ep,15)=bxy(2)+byx(2)
356 bm(ep,18)=bxy(3)+byx(3)
357 vdef(ep,3)= vdef(ep,3)+ c1*vt2+c2*vt1
362 bm(ep,10)=vq(ep,1,1,ng)*c1
363 bmf(ep,10)=thk*bm(ep,10)+vq(ep,1,1,ng)*bc1
364 bm(ep,11)=vq(ep,1,2,ng)*c2
365 bmf(ep,11)=thk*bm(ep,11)+vq(ep,1,2,ng)*bc2
366 bmf(ep,12)=thk*bm(ep,12)+vq(ep,1,1,ng)*bc2+vq(ep,1,2,ng)*bc1
370 bm(ep,13)=vq(ep,2,1,ng)*c1
371 bmf(ep,13)=thk*bm(ep,13)+vq(ep,2,1,ng)*bc1
372 bm(ep,14)=vq(ep,2,2,ng)*c2
373 bmf(ep,14)=thk*bm(ep,14)+vq(ep,2,2,ng)*bc2
374 bmf(ep,15)=thk*bm(ep,15)+vq(ep,2,1,ng)*bc2+vq(ep,2,2,ng)*bc1
378 bm(ep,16)=vq(ep,3,1,ng)*c1
379 bmf(ep,16)=thk*bm(ep,16)+vq(ep,3,1,ng)*bc1
380 bm(ep,17)=vq(ep,3,2,ng)*c2
381 bmf(ep,17)=thk*bm(ep,17)+vq(ep,3,2,ng)*bc2
382 bmf(ep,18)=thk*bm(ep,18)+vq(ep,3,1,ng)*bc2+vq(ep,3,2,ng)*bc1
385 vdef(ep,1)= vdef(ep,1)+ c1*vt1
386 vdef(ep,2)= vdef(ep,2)+ c2*vt2
389 vdef(ep,6)= vdef(ep,6)+ bc1*vt1
390 vdef(ep,7)= vdef(ep,7)+ bc2*vt2
391 vdef(ep,8)= vdef(ep,8)+ bc1*vt2+bc2*vt1
394 v1(2)=vq(ep,1,1,ng)*vqn(ep,1,2)+vq(ep,2,1,ng)*vqn(ep,2,2)
395 1 +vq(ep,3,1,ng)*vqn(ep,3,2)
396 v1(1)=-(vq(ep,1,1,ng)*vqn(ep,4,2)+vq(ep,2,1,ng)*vqn(ep,5,2)
397 1 +vq(ep,3,1,ng)*vqn(ep,6,2))
398 v2(2)=vq(ep,1,2,ng)*vqn(ep,1,2)+vq(ep,2,2,ng)*vqn(ep,2,2)
399 1 +vq(ep,3,2,ng)*vqn(ep,3,2)
400 v2(1)=-(vq(ep,1,2,ng)*vqn(ep,4,2)+vq(ep,2,2,ng)*vqn(ep,5,2)
401 1 +vq(ep,3,2,ng)*vqn(ep,6,2))
403 rv1=v1(1)*rxyz(ep,2+1)+v1(2)*rxyz(ep,2+2)
405 rv2=v2(1)*rxyz(ep,2+1)+v2(2)*rxyz(ep,2+2)
410 bf(ep,9)=v1(1)*c2+v2(1)*c1
413 bf(ep,12)=v1(2)*c2+v2(2)*c1
415 vdef(ep,6)= vdef(ep,6)+c1*rv1
416 vdef(ep,7)= vdef(ep,7)+c2*rv2
417 vdef(ep,8)= vdef(ep,8)+c1*rv2+c2*rv1
421 vt1=vq(ep,1,1,ng)*vxyz(ep,7)+vq(ep,2,1,ng)*vxyz(ep,8)
422 1 +vq(ep,3,1,ng)*vxyz(ep,9)
424 vt2=vq(ep,1,2,ng)*vxyz(ep,7)+vq(ep,2,2,ng)*vxyz(ep,8)
425 1 +vq(ep,3,2,ng)*vxyz(ep,9)
428 c1=vksi(3,ng)*tc(ep,1,1)+veta(3,ng)*tc(ep,2,1)
429 c2=vksi(3,ng)*tc(ep,1,2)+veta(3,ng)*tc(ep,2,2)
431 bc1=vksi(3,ng)*tbc(1,1)+veta(3,ng)*tbc(2,1)
432 bc2=vksi(3,ng)*tbc(1,2)+veta(3,ng)*tbc(2,2)
435 bxy(1)=vq(ep,1,1,ng)*c2
436 bxy(2)=vq(ep,2,1,ng)*c2
437 bxy(3)=vq(ep,3,1,ng)*c2
438 byx(1)=vq(ep,1,2,ng)*c1
439 byx(2)=vq(ep,2,2,ng)*c1
440 byx(3)=vq(ep,3,2,ng)*c1
441 brz(i,1,3)=-bxy(1)+byx(1)
442 brz(i,2,3)=-bxy(2)+byx(2)
443 brz(i,3,3)=-bxy(3)+byx(3)
444 bm(ep,21)=bxy(1)+byx(1)
445 bm(ep,24)=bxy(2)+byx(2)
446 bm(ep,27)=bxy(3)+byx(3)
447 vdef(ep,3)= vdef(ep,3)+ c1*vt2+c2*vt1
452 bm(ep,19)=vq(ep,1,1,ng)*c1
453 bmf(ep,19)=thk*bm(ep,19)+vq(ep,1,1,ng)*bc1
454 bm(ep,20)=vq(ep,1,2,ng)*c2
455 bmf(ep,20)=thk*bm(ep,20)+vq(ep,1,2,ng)*bc2
456 bmf(ep,21)=thk*bm(ep,21)+vq(ep,1,1,ng)*bc2+vq(ep,1,2,ng)*bc1
459 bm(ep,22)=vq(ep,2,1,ng)*c1
460 bmf(ep,22)=thk*bm(ep,22)+vq(ep,2,1,ng)*bc1
461 bm(ep,23)=vq(ep,2,2,ng)*c2
462 bmf(ep,23)=thk*bm(ep,23)+vq(ep,2,2,ng)*bc2
463 bmf(ep,24)=thk*bm(ep,24)+vq(ep,2,1,ng)*bc2+vq(ep,2,2,ng)*bc1
466 bm(ep,25)=vq(ep,3,1,ng)*c1
467 bmf(ep,25)=thk*bm(ep,25)+vq(ep,3,1,ng)*bc1
468 bm(ep,26)=vq(ep,3,2,ng)*c2
469 bmf(ep,26)=thk*bm(ep,26)+vq(ep,3,2,ng)*bc2
470 bmf(ep,27)=thk*bm(ep,27)+vq(ep,3,1,ng)*bc2+vq(ep,3,2,ng)*bc1
473 vdef(ep,1)= vdef(ep,1)+ c1*vt1
474 vdef(ep,2)= vdef(ep,2)+ c2*vt2
477 vdef(ep,6)= vdef(ep,6)+ bc1*vt1
478 vdef(ep,7)= vdef(ep,7)+ bc2*vt2
479 vdef(ep,8)= vdef(ep,8)+ bc1*vt2+bc2*vt1
482 v1(2)=vq(ep,1,1,ng)*vqn(ep,1,3)+vq(ep,2,1,ng)*vqn(ep,2,3)
483 1 +vq(ep,3,1,ng)*vqn(ep,3,3)
484 v1(1)=-(vq(ep,1,1,ng)*vqn(ep,4,3)+vq(ep,2,1,ng)*vqn(ep,5,3)
485 1 +vq(ep,3,1,ng)*vqn(ep,6,3))
486 v2(2)=vq(ep,1,2,ng)*vqn(ep,1,3)+vq(ep,2,2,ng)*vqn(ep,2,3)
487 1 +vq(ep,3,2,ng)*vqn(ep,3,3)
488 v2(1)=-(vq(ep,1,2,ng)*vqn(ep,4,3)+vq(ep,2,2,ng)*vqn(ep,5,3)
489 1 +vq(ep,3,2,ng)*vqn(ep,6,3))
491 rv1=v1(1)*rxyz(ep,4+1)+v1(2)*rxyz(ep,4+2)
493 rv2=v2(1)*rxyz(ep,4+1)+v2(2)*rxyz(ep,4+2)
498 bf(ep,15)=v1(1)*c2+v2(1)*c1
501 bf(ep,18)=v1(2)*c2+v2(2)*c1
503 vdef(ep,6)= vdef(ep,6)+c1*rv1
504 vdef(ep,7)= vdef(ep,7)+c2*rv2
505 vdef(ep,8)= vdef(ep,8)+c1*rv2+c2*rv1
508 vt1=vq(ep,1,1,ng)*vxyz(ep,10)+vq(ep,2,1,ng)*vxyz(ep,11)
509 1 +vq(ep,3,1,ng)*vxyz(ep,12)
511 vt2=vq(ep,1,2,ng)*vxyz(ep,10)+vq(ep,2,2,ng)*vxyz(ep,11)
512 1 +vq(ep,3,2,ng)*vxyz(ep,12)
515 c1=vksi(4,ng)*tc(ep,1,1)+veta(4,ng)*tc(ep,2,1)
516 c2=vksi(4,ng)*tc(ep,1,2)+veta(4,ng)*tc(ep,2,2)
518 bc1=vksi(4,ng)*tbc(1,1)+veta(4,ng)*tbc(2,1)
519 bc2=vksi(4,ng)*tbc(1,2)+veta(4,ng)*tbc(2,2)
522 bxy(1)=vq(ep,1,1,ng)*c2
523 bxy(2)=vq(ep,2,1,ng)*c2
524 bxy(3)=vq(ep,3,1,ng)*c2
525 byx(1)=vq(ep,1,2,ng)*c1
526 byx(2)=vq(ep,2,2,ng)*c1
527 byx(3)=vq(ep,3,2,ng)*c1
528 brz(i,1,4)=-bxy(1)+byx(1)
529 brz(i,2,4)=-bxy(2)+byx(2)
530 brz(i,3,4)=-bxy(3)+byx(3)
531 bm(ep,30)=bxy(1)+byx(1)
532 bm(ep,33)=bxy(2)+byx(2)
533 bm(ep,36)=bxy(3)+byx(3)
534 vdef(ep,3)= vdef(ep,3)+ c1*vt2+c2*vt1
538 bm(ep,28)=vq(ep,1,1,ng)*c1
539 bmf(ep,28)=thk*bm(ep,28)+vq(ep,1,1,ng)*bc1
540 bm(ep,29)=vq(ep,1,2,ng)*c2
541 bmf(ep,29)=thk*bm(ep,29)+vq(ep,1,2,ng)*bc2
543 bmf(ep,30)=thk*bm(ep,30)+vq(ep,1,1,ng)*bc2+vq(ep,1,2,ng)*bc1
546 bm(ep,31)=vq(ep,2,1,ng)*c1
547 bmf(ep,31)=thk*bm(ep,31)+vq(ep,2,1,ng)*bc1
548 bm(ep,32)=vq(ep,2,2,ng)*c2
549 bmf(ep,32)=thk*bm(ep,32)+vq(ep,2,2,ng)*bc2
550 bmf(ep,33)=thk*bm(ep,33)+vq(ep,2,1,ng)*bc2+vq(ep,2,2,ng)*bc1
553 bm(ep,34)=vq(ep,3,1,ng)*c1
554 bmf(ep,34)=thk*bm(ep,34)+vq(ep,3,1,ng)*bc1
555 bm(ep,35)=vq(ep,3,2,ng)*c2
556 bmf(ep,35)=thk*bm(ep,35)+vq(ep,3,2,ng)*bc2
557 bmf(ep,36)=thk*bm(ep,36)+vq(ep,3,1,ng)*bc2+vq(ep,3,2,ng)*bc1
560 vdef(ep,1)= vdef(ep,1)+ c1*vt1
561 vdef(ep,2)= vdef(ep,2)+ c2*vt2
564 vdef(ep,6)= vdef(ep,6)+ bc1*vt1
565 vdef(ep,7)= vdef(ep,7)+ bc2*vt2
566 vdef(ep,8)= vdef(ep,8)+ bc1*vt2+bc2*vt1
569 v1(2)=vq(ep,1,1,ng)*vqn(ep,1,4)+vq(ep,2,1,ng)*vqn(ep,2,4)
570 1 +vq(ep,3,1,ng)*vqn(ep,3,4)
571 v1(1)=-(vq(ep,1,1,ng)*vqn(ep,4,4)+vq(ep,2,1,ng)*vqn(ep,5,4)
572 1 +vq(ep,3,1,ng)*vqn(ep,6,4))
573 v2(2)=vq(ep,1,2,ng)*vqn(ep,1,4)+vq(ep,2,2,ng)*vqn(ep,2,4)
574 1 +vq(ep,3,2,ng)*vqn(ep,3,4)
575 v2(1)=-(vq(ep,1,2,ng)*vqn(ep,4,4)+vq(ep,2,2,ng)*vqn(ep,5,4)
576 1 +vq(ep,3,2,ng)*vqn(ep,6,4))
578 rv1=v1(1)*rxyz(ep,6+1)+v1(2)*rxyz(ep,6+2)
580 rv2=v2(1)*rxyz(ep,6+1)+v2(2)*rxyz(ep,6+2)
585 bf(ep,21)=v1(1)*c2+v2(1)*c1
590 vdef(ep,6)= vdef(ep,6)+thk*vdef(ep,1
591 vdef(ep,7)= vdef(ep,7)+thk*vdef(ep,2)+c2*rv2
592 vdef(ep,8)= vdef(ep,8)+thk*vdef(ep,3)+c1*rv2+c2*rv1
623 bc(ep,7)=v11(1)*vastn(ep,1)
624 bc(ep,9)=v11(1)*vastn(ep,2)
625 bc(ep,17)=v11(1)*vastn(ep,3)
626 bc(ep,19)=v11(1)*vastn(ep,4)
646 bc(ep,37)=v11(2)*vastn(ep,5)
647 bc(ep,39)=v11(2)*vastn(ep,6)
648 bc(ep,27)=v11(2)*vastn(ep,7)
649 bc(ep,29)=v11(2)*vastn(ep,8)
669 bc(ep,8)=v11(3)*vastn(ep,9)
670 bc(ep,10)=v11(3)*vastn(ep,10)
671 bc(ep,38)=v11(3)*vastn(ep,11)
672 bc(ep,40)=v11(3)*vastn(ep,12)
680 c1=v11(4)*vnrm(ep,10)
684 c1=v11(4)*vnrm(ep,11)
688 c1=v11(4)*vnrm(ep,12)
692 bc(ep,18)=v11(4)*vastn(ep,13)
693 bc(ep,20)=v11(4)*vastn(ep,14)
694 bc(ep,28)=v11(4)*vastn(ep,15)
695 bc(ep,30)=v11(4)*vastn(ep,16)
700 1 bc(ep,1)*vxyz(ep,1)+bc(ep,3)*vxyz(ep,2)
701 2 +bc(ep,5)*vxyz(ep,3)+bc(ep,7)*rxyz(ep,1)
702 3 +bc(ep,9)*rxyz(ep,2)
704 1 bc(ep,2)*vxyz(ep,1)+bc(ep,4)*vxyz(ep,2)
705 2 +bc(ep,6)*vxyz(ep,3)+bc(ep,8)*rxyz(ep,1)
706 3 +bc(ep,10)*rxyz(ep,2)
709 1 +bc(ep,11)*vxyz(ep,4)+bc(ep,13)*vxyz(ep,5)
710 2 +bc(ep,15)*vxyz(ep,6)+bc(ep,17)*rxyz(ep,3)
711 3 +bc(ep,19)*rxyz(ep,4)
713 1 +bc(ep,12)*vxyz(ep,4)+bc(ep,14)*vxyz(ep,5)
714 2 +bc(ep,16)*vxyz(ep,6)+bc(ep,18)*rxyz(ep,3)
715 3 +bc(ep,20)*rxyz(ep,4)
718 1 +bc(ep,21)*vxyz(ep,7)+bc(ep,23)*vxyz(ep,8)
719 2 +bc(ep,25)*vxyz(ep,9)+bc(ep,27)*rxyz(ep,5)
720 3 +bc(ep,29)*rxyz(ep,6)
722 1 +bc(ep,22)*vxyz(ep,7)+bc(ep,24)*vxyz(ep,8)
723 2 +bc(ep,26)*vxyz(ep,9)+bc(ep,28)*rxyz(ep,5)
724 3 +bc(ep,30)*rxyz(ep,6)
727 1 +bc(ep,31)*vxyz(ep,10)+bc(ep,33)*vxyz(ep,11)
728 2 +bc(ep,35)*vxyz(ep,12)+bc(ep,37)*rxyz(ep,7)
729 3 +bc(ep,39)*rxyz(ep,8)
731 1 +bc(ep,32)*vxyz(ep,10)+bc(ep,34)*vxyz(ep,11)
732 2 +bc(ep,36)*vxyz(ep,12)+bc(ep,38)*rxyz(ep,7)
733 3 +bc(ep,40)*rxyz(ep,8)
734 vdef(ep,4)=tc(ep,1,1)*bcx+tc(ep,2,1)*bcy
735 vdef(ep,5)=tc(ep,1,2)*bcx+tc(ep,2,2)*bcy
1544 2 HX,HY,VETA,VKSI,BM,NPLAT,IPLAT,
1554#include "implicit_f.inc"
1555#include "mvsiz_p.inc"
1559 INTEGER NPG,NG,JFT,JLT,NPLAT,IPLAT(*),ISROT
1562 . vcore(mvsiz,12),vksi(4,4),veta(4,4),
1563 . bm(mvsiz,36),hx(mvsiz,4),hy(mvsiz,4)
1565 . vjfi(mvsiz,3,2,4),vq(mvsiz,3,3,4)
1569 INTEGER PT,PT00,PT0,I,J,EP,IUN,NG1,K
1571 . VPG(2,NPG),PG1,PG,THK,DETJ,DET,
1572 . (3,2),BCX,BCY,BXY(3),BYX(3),
1573 . V1(2),V2(2),RV1,RV2,
1574 . C1,C2,VT1,VT2,BC1,BC2,VB1,TC(MVSIZ,2,2),
1575 . VJF1(2,3),VJF(3,3),TBI(2,2),TBC(2,2),V11(4)
1577 . A_1,C11,C12,C21,C22,CC,BETA1,KSI1,KSIY1,BETA2,KSI2,KSIY2
1578 parameter(pg=.577350269189626)
1579 parameter(pg1=-.577350269189626)
1584 DATA vpg/pg1,pg1,pg,pg1,pg,pg,pg1,pg/
1586#include "vectorize.inc"
1590 bm(ep,1)=vcore(ep,1)+hx(ep,ng)*vcore(ep,5)
1591 bm(ep,2)=vcore(ep,2)+hx(ep,ng)*vcore(ep,6)
1592 bm(ep,3)=hx(ep,ng)*fourth
1594 bm(ep,5)=vcore(ep,3)+hy(ep,ng)*vcore(ep,5)
1595 bm(ep,6)=vcore(ep,4)+hy(ep,ng)*vcore(ep,6)
1596 bm(ep,7)=hy(ep,ng)*fourth
1600#include "vectorize.inc"
1601 DO 150 i=nplat+1,jlt
1606 tc(ep,1,1)=vjfi(ep,1,1,ng)*vq(ep,1,1,ng)+vjfi(ep,2,1,ng)
1607 + *vq(ep,2,1,ng)+ vjfi(ep,3,1,ng)*vq(ep,3,1,ng)
1608 tc(ep,2,1)=vjfi(ep,1,2,ng)*vq(ep,1,1,ng)+vjfi(ep,2,2,ng)
1609 + *vq(ep,2,1,ng)+ vjfi(ep,3,2,ng)*vq(ep,3,1,ng)
1610 tc(ep,1,2)=vjfi(ep,1,1,ng)*vq(ep,1,2,ng)+vjfi(ep,2,1,ng)
1611 + *vq(ep,2,2,ng)+ vjfi(ep,3,1,ng)*vq(ep,3,2,ng)
1612 tc(ep,2,2)=vjfi(ep,1,2,ng)*vq(ep,1,2,ng)+vjfi(ep,2,2,ng)
1613 + *vq(ep,2,2,ng)+ vjfi(ep,3,2,ng)*vq(ep,3,2,ng)
1623 c1=vksi(1,ng)*tc(ep,1,1)+veta(1,ng)*tc(ep,2,1)
1624 c2=vksi(1,ng)*tc(ep,1,2)+veta(1,ng)*tc(ep,2,2)
1627 bxy(1)=vq(ep,1,1,ng)*c2
1628 bxy(2)=vq(ep,2,1,ng)*c2
1629 bxy(3)=vq(ep,3,1,ng)*c2
1630 byx(1)=vq(ep,1,2,ng)*c1
1631 byx(2)=vq(ep,2,2,ng)*c1
1632 byx(3)=vq(ep,3,2,ng)*c1
1633 bm(ep,3)=bxy(1)+byx(1)
1634 bm(ep,6)=bxy(2)+byx(2)
1635 bm(ep,9)=bxy(3)+byx(3)
1637 bm(ep,1)=vq(ep,1,1,ng)*c1
1639 bm(ep,2)=vq(ep,1,2,ng)*c2
1641 bm(ep,4)=vq(ep,2,1,ng)*c1
1642 bm(ep,5)=vq(ep,2,2,ng)*c2
1644 bm(ep,7)=vq(ep,3,1,ng)*c1
1645 bm(ep,8)=vq(ep,3,2,ng)*c2
1650 c1=vksi(2,ng)*tc(ep,1,1)+veta(2,ng)*tc(ep,2,1)
1651 c2=vksi(2,ng)*tc(ep,1,2)+veta(2,ng)*tc(ep,2,2)
1654 bxy(1)=vq(ep,1,1,ng)*c2
1655 bxy(2)=vq(ep,2,1,ng)*c2
1656 bxy(3)=vq(ep,3,1,ng)*c2
1657 byx(1)=vq(ep,1,2,ng)*c1
1658 byx(2)=vq(ep,2,2,ng)*c1
1659 byx(3)=vq(ep,3,2,ng)*c1
1660 bm(ep,12)=bxy(1)+byx(1)
1661 bm(ep,15)=bxy(2)+byx(2)
1662 bm(ep,18)=bxy(3)+byx(3)
1664 bm(ep,10)=vq(ep,1,1,ng)*c1
1665 bm(ep,11)=vq(ep,1,2,ng)*c2
1668 bm(ep,13)=vq(ep,2,1,ng)*c1
1669 bm(ep,14)=vq(ep,2,2,ng)*c2
1672 bm(ep,16)=vq(ep,3,1,ng)*c1
1673 bm(ep,17)=vq(ep,3,2,ng)*c2
1678 c1=vksi(3,ng)*tc(ep,1,1)+veta(3,ng)*tc(ep,2,1)
1679 c2=vksi(3,ng)*tc(ep,1,2)+veta(3,ng)*tc(ep,2,2)
1681 bxy(1)=vq(ep,1,1,ng)*c2
1682 bxy(2)=vq(ep,2,1,ng)*c2
1683 bxy(3)=vq(ep,3,1,ng)*c2
1684 byx(1)=vq(ep,1,2,ng)*c1
1685 byx(2)=vq(ep,2,2,ng)*c1
1686 byx(3)=vq(ep,3,2,ng)*c1
1687 bm(ep,21)=bxy(1)+byx(1)
1688 bm(ep,24)=bxy(2)+byx(2)
1689 bm(ep,27)=bxy(3)+byx(3)
1692 bm(ep,20)=vq(ep,1,2,ng)*c2
1694 bm(ep,22)=vq(ep,2,1,ng)*c1
1695 bm(ep,23)=vq(ep,2,2,ng)*c2
1697 bm(ep,25)=vq(ep,3,1,ng)*c1
1698 bm(ep,26)=vq(ep,3,2,ng)*c2
1701 c1=vksi(4,ng)*tc(ep,1,1)+veta(4,ng)*tc(ep,2,1)
1702 c2=vksi(4,ng)*tc(ep,1,2)+veta(4,ng)*tc
1704 bxy(1)=vq(ep,1,1,ng)*c2
1705 bxy(2)=vq(ep,2,1,ng)*c2
1706 bxy(3)=vq(ep,3,1,ng)*c2
1707 byx(1)=vq(ep,1,2,ng)*c1
1708 byx(2)=vq(ep,2,2,ng)*c1
1709 byx(3)=vq(ep,3,2,ng)*c1
1710 bm(ep,30)=bxy(1)+byx(1)
1711 bm(ep,33)=bxy(2)+byx(2)
1712 bm(ep,36)=bxy(3)+byx(3)
1713 bm(ep,28)=vq(ep,1,1,ng)*c1
1714 bm(ep,29)=vq(ep,1,2,ng)*c2
1716 bm(ep,31)=vq(ep,2,1,ng)*c1
1717 bm(ep,32)=vq(ep,2,2,ng)*c2
1719 bm(ep,34)=vq(ep,3,1,ng)*c1
1720 bm(ep,35)=vq(ep,3,2,ng)*c2