FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
precheck_LIB_3d.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7contains
8 !***********************************************************************
9 ! 3D SOLID Element: PreCheck
10 ! PRE_341(XX,YY,ZZ,vol,almax,almin)
11 ! PRE_351(XX,YY,ZZ,vol,almax,almin)
12 ! PRE_361(XX,YY,ZZ,vol,almax,almin)
13 ! PRE_342(XX,YY,ZZ,vol,almax,almin)
14 ! PRE_352(XX,YY,ZZ,vol,almax,almin)
15 ! PRE_362(XX,YY,ZZ,vol,almax,almin)
16 !----------------------------------------------------------------------*
17 subroutine pre_341( XX,YY,ZZ,vol,almax,almin )
18 !----------------------------------------------------------------------*
19 !
20 ! CALCULATION 3D 4-NODE SOLID ELEMENT
21 !
22 use hecmw
24 implicit none
25 ! I/F VARIABLES
26 real(kind=kreal) xx(*),yy(*),zz(*),vol,almax,almin
27 ! LOCAL VARIABLES
28 integer(kind=kint) NN
29 integer(kind=kint) NG
30 parameter(nn=4,ng=2)
31 real(kind=kreal) h(nn),hl1(nn),hl2(nn),hl3(nn),hl4(nn)
32 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
33 integer(kind=kint) I,L,L1,L2,L3
34 real(kind=kreal) xl1,xl2,xl3
35 real(kind=kreal) x1,x2,x3,x4
36 real(kind=kreal) a1,a2,a3,a4,a5,a6
37 !C
38 vol = 0.0
39 ! LOOP FOR INTEGRATION POINTS
40 do l3=1,ng
41 xl3=xg(ng,l3)
42 x3 =(xl3+1.0)*0.5
43 do l2=1,ng
44 xl2=xg(ng,l2)
45 x2 =(1.0-x3)*(xl2+1.0)*0.5
46 do l1=1,ng
47 xl1=xg(ng,l1)
48 x1=(1.0-x2-x3)*(xl1+1.0)*0.5
49 ! INTERPOLATION FUNCTION
50 x4=1.0-x1-x2-x3
51 h(1)=x1
52 h(2)=x2
53 h(3)=x3
54 h(4)=x4
55 ! DERIVATIVE OF INTERPOLATION FUNCTION
56 ! FOR L1-COORDINATE
57 hl1(1)= 1.0
58 hl1(2)= 0.0
59 hl1(3)= 0.0
60 hl1(4)= 0.0
61 ! FOR L2-COORDINATE
62 hl2(1)= 0.0
63 hl2(2)= 1.0
64 hl2(3)= 0.0
65 hl2(4)= 0.0
66 ! FOR L3-COORDINATE
67 hl3(1)= 0.0
68 hl3(2)= 0.0
69 hl3(3)= 1.0
70 hl3(4)= 0.0
71 ! FOR L4-COORDINATE
72 hl4(1)= 0.0
73 hl4(2)= 0.0
74 hl4(3)= 0.0
75 hl4(4)= 1.0
76 ! JACOBI MATRIX
77 xj11=0.0
78 xj21=0.0
79 xj31=0.0
80 xj12=0.0
81 xj22=0.0
82 xj32=0.0
83 xj13=0.0
84 xj23=0.0
85 xj33=0.0
86 do i=1,nn
87 xj11=xj11+(hl4(i)-hl1(i))*xx(i)
88 xj21=xj21+(hl4(i)-hl2(i))*xx(i)
89 xj31=xj31+(hl4(i)-hl3(i))*xx(i)
90 xj12=xj12+(hl4(i)-hl1(i))*yy(i)
91 xj22=xj22+(hl4(i)-hl2(i))*yy(i)
92 xj32=xj32+(hl4(i)-hl3(i))*yy(i)
93 xj13=xj13+(hl4(i)-hl1(i))*zz(i)
94 xj23=xj23+(hl4(i)-hl2(i))*zz(i)
95 xj33=xj33+(hl4(i)-hl3(i))*zz(i)
96 enddo
97 ! DETERMINANT OF JACOBIAN
98 det=xj11*xj22*xj33 &
99 +xj12*xj23*xj31 &
100 +xj13*xj21*xj32 &
101 -xj13*xj22*xj31 &
102 -xj12*xj21*xj33 &
103 -xj11*xj23*xj32
104 ! WEIGT VALUE AT GAUSSIAN POINT
105 wg=wgt(ng,l1)*wgt(ng,l2)*wgt(ng,l3)*det*(1.0-x3)*(1.0-x2-x3)*0.125
106 do i = 1, nn
107 vol = vol + h(i)*wg
108 enddo
109 enddo
110 enddo
111 enddo
112
113 a1 = sqrt( (xx(2)-xx(1))**2+(yy(2)-yy(1))**2+(zz(2)-zz(1))**2 )
114 a2 = sqrt( (xx(3)-xx(2))**2+(yy(3)-yy(2))**2+(zz(3)-zz(2))**2 )
115 a3 = sqrt( (xx(1)-xx(3))**2+(yy(1)-yy(3))**2+(zz(1)-zz(3))**2 )
116 a4 = sqrt( (xx(4)-xx(1))**2+(yy(4)-yy(1))**2+(zz(4)-zz(1))**2 )
117 a5 = sqrt( (xx(4)-xx(2))**2+(yy(4)-yy(2))**2+(zz(4)-zz(2))**2 )
118 a6 = sqrt( (xx(4)-xx(3))**2+(yy(4)-yy(3))**2+(zz(4)-zz(3))**2 )
119 almax = dmax1( a1,a2,a3,a4,a5,a6 )
120 almin = dmin1( a1,a2,a3,a4,a5,a6 )
121
122 end subroutine pre_341
123 !----------------------------------------------------------------------*
124 subroutine pre_351( XX,YY,ZZ,vol,almax,almin )
125 !----------------------------------------------------------------------*
126 !
127 ! CALCULATION 3D 6-NODE SOLID ELEMENT
128 !
129 use hecmw
131 implicit none
132 ! I/F VARIABLES
133 integer(kind=kint) nline
134 real(kind=kreal) xx(*),yy(*),zz(*),vol,tline,almax,almin
135 ! LOCAL VARIABLES
136 integer(kind=kint) NN
137 integer(kind=kint) NG
138 parameter(nn=6,ng=2)
139 real(kind=kreal) h(nn),hl1(nn),hl2(nn),hl3(nn),hz(nn)
140 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
141 real(kind=kreal) xji11,xji21,xji31,xji12,xji22,xji32,xji13,xji23,xji33
142 integer(kind=kint) I,L1,L2,LZ
143 real(kind=kreal) x1,x2,x3,xl1,xl2,zi
144 real(kind=kreal) a1,a2,a3,a4,a5,a6,a7,a8,a9
145 !C
146 vol = 0.0
147 ! LOOP FOR INTEGRATION POINTS
148 do lz=1,ng
149 zi=xg(ng,lz)
150 do l2=1,ng
151 xl2=xg(ng,l2)
152 x2 =(xl2+1.0)*0.5
153 do l1=1,ng
154 xl1=xg(ng,l1)
155 x1=0.5*(1.0-x2)*(xl1+1.0)
156 ! INTERPOLATION FUNCTION
157 x3=1.0-x1-x2
158 h(1)=x1*(1.0-zi)*0.5
159 h(2)=x2*(1.0-zi)*0.5
160 h(3)=x3*(1.0-zi)*0.5
161 h(4)=x1*(1.0+zi)*0.5
162 h(5)=x2*(1.0+zi)*0.5
163 h(6)=x3*(1.0+zi)*0.5
164 ! DERIVATIVE OF INTERPOLATION FUNCTION
165 ! FOR L1-COORDINATE
166 hl1(1)=(1.0-zi)*0.5
167 hl1(2)=0.0
168 hl1(3)=0.0
169 hl1(4)=(1.0+zi)*0.5
170 hl1(5)=0.0
171 hl1(6)=0.0
172 ! FOR L2-COORDINATE
173 hl2(1)=0.0
174 hl2(2)=(1.0-zi)*0.5
175 hl2(3)=0.0
176 hl2(4)=0.0
177 hl2(5)=(1.0+zi)*0.5
178 hl2(6)=0.0
179 ! FOR L3-COORDINATE
180 hl3(1)=0.0
181 hl3(2)=0.0
182 hl3(3)=(1.0-zi)*0.5
183 hl3(4)=0.0
184 hl3(5)=0.0
185 hl3(6)=(1.0+zi)*0.5
186 ! FOR Z-COORDINATE
187 hz(1)=-x1*0.5
188 hz(2)=-x2*0.5
189 hz(3)=-x3*0.5
190 hz(4)= x1*0.5
191 hz(5)= x2*0.5
192 hz(6)= x3*0.5
193 ! JACOBI MATRIX
194 xj11=0.0
195 xj21=0.0
196 xj31=0.0
197 xj12=0.0
198 xj22=0.0
199 xj32=0.0
200 xj13=0.0
201 xj23=0.0
202 xj33=0.0
203 do i=1,nn
204 xj11=xj11+(hl1(i)-hl3(i))*xx(i)
205 xj21=xj21+(hl2(i)-hl3(i))*xx(i)
206 xj31=xj31+hz(i)*xx(i)
207 xj12=xj12+(hl1(i)-hl3(i))*yy(i)
208 xj22=xj22+(hl2(i)-hl3(i))*yy(i)
209 xj32=xj32+hz(i)*yy(i)
210 xj13=xj13+(hl1(i)-hl3(i))*zz(i)
211 xj23=xj23+(hl2(i)-hl3(i))*zz(i)
212 xj33=xj33+hz(i)*zz(i)
213 enddo
214 ! DETERMINANT OF JACOBIAN
215 det=xj11*xj22*xj33 &
216 +xj12*xj23*xj31 &
217 +xj13*xj21*xj32 &
218 -xj13*xj22*xj31 &
219 -xj12*xj21*xj33 &
220 -xj11*xj23*xj32
221 ! WEIGHT VALUE AT GAUSSIAN POINT
222 wg=wgt(ng,l1)*wgt(ng,l2)*wgt(ng,lz)*det*(1.0-x2)*0.25
223 do i = 1, nn
224 vol = vol + h(i)*wg
225 enddo
226 enddo
227 enddo
228 enddo
229
230 a1 = sqrt( (xx(2)-xx(1))**2+(yy(2)-yy(1))**2+(zz(2)-zz(1))**2 )
231 a2 = sqrt( (xx(3)-xx(2))**2+(yy(3)-yy(2))**2+(zz(3)-zz(2))**2 )
232 a3 = sqrt( (xx(1)-xx(3))**2+(yy(1)-yy(3))**2+(zz(1)-zz(3))**2 )
233 a4 = sqrt( (xx(5)-xx(4))**2+(yy(5)-yy(4))**2+(zz(5)-zz(4))**2 )
234 a5 = sqrt( (xx(6)-xx(5))**2+(yy(6)-yy(5))**2+(zz(6)-zz(5))**2 )
235 a6 = sqrt( (xx(4)-xx(6))**2+(yy(4)-yy(6))**2+(zz(4)-zz(6))**2 )
236 a7 = sqrt( (xx(4)-xx(1))**2+(yy(4)-yy(1))**2+(zz(4)-zz(1))**2 )
237 a8 = sqrt( (xx(5)-xx(2))**2+(yy(5)-yy(2))**2+(zz(5)-zz(2))**2 )
238 a9 = sqrt( (xx(6)-xx(3))**2+(yy(6)-yy(3))**2+(zz(6)-zz(3))**2 )
239 almax = dmax1( a1,a2,a3,a4,a5,a6,a7,a8,a9 )
240 almin = dmin1( a1,a2,a3,a4,a5,a6,a7,a8,a9 )
241
242 end subroutine pre_351
243 !----------------------------------------------------------------------*
244 subroutine pre_361( XX,YY,ZZ,vol,almax,almin )
245 !----------------------------------------------------------------------*
246 !
247 ! CALCULATION 3D 8-NODE SOLID ELEMENT
248 !
249 use hecmw
251 implicit none
252 ! I/F VARIABLES
253 real(kind=kreal) xx(*),yy(*),zz(*),vol,almax,almin
254 ! LOCAL VARIABLES
255 integer(kind=kint) NN
256 integer(kind=kint) NG
257 parameter(nn=8,ng=2)
258 real(kind=kreal) h(nn),hr(nn),hs(nn),ht(nn)
259 real(kind=kreal) ri,si,ti,rp,sp,tp,rm,sm,tm
260 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
261 integer(kind=kint) I,LX,LY,LZ
262 real(kind=kreal) a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12
263 !C
264 vol = 0.0
265 ! LOOP FOR INTEGRATION POINTS
266 do lx=1,ng
267 ri=xg(ng,lx)
268 do ly=1,ng
269 si=xg(ng,ly)
270 do lz=1,ng
271 ti=xg(ng,lz)
272 rp=1.0+ri
273 sp=1.0+si
274 tp=1.0+ti
275 rm=1.0-ri
276 sm=1.0-si
277 tm=1.0-ti
278 ! INTERPOLATION FUNCTION
279 h(1)=0.125*rm*sm*tm
280 h(2)=0.125*rp*sm*tm
281 h(3)=0.125*rp*sp*tm
282 h(4)=0.125*rm*sp*tm
283 h(5)=0.125*rm*sm*tp
284 h(6)=0.125*rp*sm*tp
285 h(7)=0.125*rp*sp*tp
286 h(8)=0.125*rm*sp*tp
287 ! DERIVATIVE OF INTERPOLATION FUNCTION
288 ! FOR R-COORDINATE
289 hr(1)=-.125*sm*tm
290 hr(2)= .125*sm*tm
291 hr(3)= .125*sp*tm
292 hr(4)=-.125*sp*tm
293 hr(5)=-.125*sm*tp
294 hr(6)= .125*sm*tp
295 hr(7)= .125*sp*tp
296 hr(8)=-.125*sp*tp
297 ! FOR S-COORDINATE
298 hs(1)=-.125*rm*tm
299 hs(2)=-.125*rp*tm
300 hs(3)= .125*rp*tm
301 hs(4)= .125*rm*tm
302 hs(5)=-.125*rm*tp
303 hs(6)=-.125*rp*tp
304 hs(7)= .125*rp*tp
305 hs(8)= .125*rm*tp
306 ! FOR T-COORDINATE
307 ht(1)=-.125*rm*sm
308 ht(2)=-.125*rp*sm
309 ht(3)=-.125*rp*sp
310 ht(4)=-.125*rm*sp
311 ht(5)= .125*rm*sm
312 ht(6)= .125*rp*sm
313 ht(7)= .125*rp*sp
314 ht(8)= .125*rm*sp
315 ! JACOBI MATRIX
316 xj11=0.0
317 xj21=0.0
318 xj31=0.0
319 xj12=0.0
320 xj22=0.0
321 xj32=0.0
322 xj13=0.0
323 xj23=0.0
324 xj33=0.0
325 do i=1,nn
326 xj11=xj11+hr(i)*xx(i)
327 xj21=xj21+hs(i)*xx(i)
328 xj31=xj31+ht(i)*xx(i)
329 xj12=xj12+hr(i)*yy(i)
330 xj22=xj22+hs(i)*yy(i)
331 xj32=xj32+ht(i)*yy(i)
332 xj13=xj13+hr(i)*zz(i)
333 xj23=xj23+hs(i)*zz(i)
334 xj33=xj33+ht(i)*zz(i)
335 enddo
336 ! DETERMINANT OF JACOBIAN
337 det=xj11*xj22*xj33 &
338 +xj12*xj23*xj31 &
339 +xj13*xj21*xj32 &
340 -xj13*xj22*xj31 &
341 -xj12*xj21*xj33 &
342 -xj11*xj23*xj32
343 ! WEIGHT VALUE AT GAUSSIAN POINT
344 wg=wgt(ng,lx)*wgt(ng,ly)*wgt(ng,lz)*det
345 do i=1,nn
346 vol = vol + h(i)*wg
347 enddo
348 enddo
349 enddo
350 enddo
351
352 a1 = sqrt( (xx(2)-xx(1))**2+(yy(2)-yy(1))**2+(zz(2)-zz(1))**2 )
353 a2 = sqrt( (xx(3)-xx(2))**2+(yy(3)-yy(2))**2+(zz(3)-zz(2))**2 )
354 a3 = sqrt( (xx(4)-xx(3))**2+(yy(4)-yy(3))**2+(zz(4)-zz(3))**2 )
355 a4 = sqrt( (xx(1)-xx(4))**2+(yy(1)-yy(4))**2+(zz(1)-zz(4))**2 )
356
357 a5 = sqrt( (xx(6)-xx(5))**2+(yy(6)-yy(5))**2+(zz(6)-zz(5))**2 )
358 a6 = sqrt( (xx(7)-xx(6))**2+(yy(7)-yy(6))**2+(zz(7)-zz(6))**2 )
359 a7 = sqrt( (xx(8)-xx(7))**2+(yy(8)-yy(7))**2+(zz(8)-zz(7))**2 )
360 a8 = sqrt( (xx(5)-xx(8))**2+(yy(5)-yy(8))**2+(zz(5)-zz(8))**2 )
361
362 a9 = sqrt( (xx(5)-xx(1))**2+(yy(5)-yy(1))**2+(zz(5)-zz(1))**2 )
363 a10 = sqrt( (xx(6)-xx(2))**2+(yy(6)-yy(2))**2+(zz(6)-zz(2))**2 )
364 a11 = sqrt( (xx(7)-xx(3))**2+(yy(7)-yy(3))**2+(zz(7)-zz(3))**2 )
365 a12 = sqrt( (xx(8)-xx(4))**2+(yy(8)-yy(4))**2+(zz(8)-zz(4))**2 )
366 almax = dmax1( a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12 )
367 almin = dmin1( a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12 )
368
369 end subroutine pre_361
370 !----------------------------------------------------------------------*
371 subroutine pre_342( XX,YY,ZZ,vol,almax,almin )
372 !----------------------------------------------------------------------*
373 !
374 ! CALCULATION 3D 10-NODE SOLID ELEMENT
375 !
376 use hecmw
378 implicit none
379 ! I/F VARIABLES
380 real(kind=kreal) xx(*),yy(*),zz(*),vol,almax,almin
381 ! LOCAL VARIABLES
382 integer(kind=kint) NN
383 integer(kind=kint) NG
384 parameter(nn=10,ng=3)
385 real(kind=kreal) h(nn),hl1(nn),hl2(nn),hl3(nn),hl4(nn)
386 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
387 integer(kind=kint) I,L1,L2,L3
388 real(kind=kreal) xl1,xl2,xl3
389 real(kind=kreal) x1,x2,x3,x4
390 real(kind=kreal) a1,a2,al1,al2,al3,al4,al5,al6
391 !
392 vol = 0.0
393 ! LOOP FOR INTEGRATION POINTS
394 do l3=1,ng
395 xl3=xg(ng,l3)
396 x3 =(xl3+1.0)*0.5
397 do l2=1,ng
398 xl2=xg(ng,l2)
399 x2 =(1.0-x3)*(xl2+1.0)*0.5
400 do l1=1,ng
401 xl1=xg(ng,l1)
402 x1=(1.0-x2-x3)*(xl1+1.0)*0.5
403 ! INTERPOLATION FUNCTION
404 x4=1.0-x1-x2-x3
405 h(1)=x1*(2.0*x1-1.0)
406 h(2)=x2*(2.0*x2-1.0)
407 h(3)=x3*(2.0*x3-1.0)
408 h(4)=x4*(2.0*x4-1.0)
409 h(5)=4.0*x1*x2
410 h(6)=4.0*x2*x3
411 h(7)=4.0*x1*x3
412 h(8)=4.0*x1*x4
413 h(9)=4.0*x2*x4
414 h(10)=4.0*x3*x4
415 ! DERIVATIVE OF INTERPOLATION FUNCTION
416 ! FOR L1-COORDINATE
417 hl1(1)= 4.0*x1-1.0
418 hl1(2)= 0.0
419 hl1(3)= 0.0
420 hl1(4)= 0.0
421 hl1(5)= 4.0*x2
422 hl1(6)= 0.0
423 hl1(7)= 4.0*x3
424 hl1(8)= 4.0*x4
425 hl1(9)= 0.0
426 hl1(10)= 0.0
427 ! FOR L2-COORDINATE
428 hl2(1)= 0.0
429 hl2(2)= 4.0*x2-1.0
430 hl2(3)= 0.0
431 hl2(4)= 0.0
432 hl2(5)= 4.0*x1
433 hl2(6)= 4.0*x3
434 hl2(7)= 0.0
435 hl2(8)= 0.0
436 hl2(9)= 4.0*x4
437 hl2(10)= 0.0
438 ! FOR L3-COORDINATE
439 hl3(1)= 0.0
440 hl3(2)= 0.0
441 hl3(3)= 4.0*x3-1.0
442 hl3(4)= 0.0
443 hl3(5)= 0.0
444 hl3(6)= 4.0*x2
445 hl3(7)= 4.0*x1
446 hl3(8)= 0.0
447 hl3(9)= 0.0
448 hl3(10)= 4.0*x4
449 ! FOR L4-COORDINATE
450 hl4(1)= 0.0
451 hl4(2)= 0.0
452 hl4(3)= 0.0
453 hl4(4)= 4.0*x4-1.0
454 hl4(5)= 0.0
455 hl4(6)= 0.0
456 hl4(7)= 0.0
457 hl4(8)= 4.0*x1
458 hl4(9)= 4.0*x2
459 hl4(10)= 4.0*x3
460 ! JACOBI MATRIX
461 xj11=0.0
462 xj21=0.0
463 xj31=0.0
464 xj12=0.0
465 xj22=0.0
466 xj32=0.0
467 xj13=0.0
468 xj23=0.0
469 xj33=0.0
470 do i=1,nn
471 xj11=xj11+(hl4(i)-hl1(i))*xx(i)
472 xj21=xj21+(hl4(i)-hl2(i))*xx(i)
473 xj31=xj31+(hl4(i)-hl3(i))*xx(i)
474 xj12=xj12+(hl4(i)-hl1(i))*yy(i)
475 xj22=xj22+(hl4(i)-hl2(i))*yy(i)
476 xj32=xj32+(hl4(i)-hl3(i))*yy(i)
477 xj13=xj13+(hl4(i)-hl1(i))*zz(i)
478 xj23=xj23+(hl4(i)-hl2(i))*zz(i)
479 xj33=xj33+(hl4(i)-hl3(i))*zz(i)
480 enddo
481 ! DETERMINANT OF JACOBIAN
482 det=xj11*xj22*xj33 &
483 +xj12*xj23*xj31 &
484 +xj13*xj21*xj32 &
485 -xj13*xj22*xj31 &
486 -xj12*xj21*xj33 &
487 -xj11*xj23*xj32
488 ! WEIGHT VALUE AT GAUSSIAN POINT
489 wg=wgt(ng,l1)*wgt(ng,l2)*wgt(ng,l3)*det*(1.0-x3)*(1.0-x2-x3)*0.125
490 do i = 1, nn
491 vol = vol + h(i)*wg
492 enddo
493 enddo
494 enddo
495 enddo
496
497 a1 = sqrt( (xx(5)-xx(1))**2+(yy(5)-yy(1))**2+(zz(5)-zz(1))**2 )
498 a2 = sqrt( (xx(2)-xx(5))**2+(yy(2)-yy(5))**2+(zz(2)-zz(5))**2 )
499 al1 = a1 + a2
500 a1 = sqrt( (xx(6)-xx(2))**2+(yy(6)-yy(2))**2+(zz(6)-zz(2))**2 )
501 a2 = sqrt( (xx(3)-xx(6))**2+(yy(3)-yy(6))**2+(zz(3)-zz(6))**2 )
502 al2 = a1 + a2
503 a1 = sqrt( (xx(7)-xx(3))**2+(yy(7)-yy(3))**2+(zz(7)-zz(3))**2 )
504 a2 = sqrt( (xx(1)-xx(7))**2+(yy(1)-yy(7))**2+(zz(1)-zz(7))**2 )
505 al3 = a1 + a2
506 a1 = sqrt( (xx(8)-xx(1))**2+(yy(8)-yy(1))**2+(zz(8)-zz(1))**2 )
507 a2 = sqrt( (xx(4)-xx(8))**2+(yy(4)-yy(8))**2+(zz(4)-zz(8))**2 )
508 al4 = a1 + a2
509 a1 = sqrt( (xx(9)-xx(2))**2+(yy(9)-yy(2))**2+(zz(9)-zz(2))**2 )
510 a2 = sqrt( (xx(4)-xx(9))**2+(yy(4)-yy(9))**2+(zz(4)-zz(9))**2 )
511 al5 = a1 + a2
512 a1 = sqrt( (xx(10)-xx(3))**2+(yy(10)-yy(3))**2+(zz(10)-zz(3))**2 )
513 a2 = sqrt( (xx(4)-xx(10))**2+(yy(4)-yy(10))**2+(zz(4)-zz(10))**2 )
514 al6 = a1 + a2
515 almax = dmax1( al1,al2,al3,al4,al5,al6 )
516 almin = dmin1( al1,al2,al3,al4,al5,al6 )
517
518 end subroutine pre_342
519 !----------------------------------------------------------------------*
520 subroutine pre_352( XX,YY,ZZ,vol,almax,almin )
521 !----------------------------------------------------------------------*
522 !
523 ! CALCULATION 3D 15-NODE SOLID ELEMENT
524 !
525 use hecmw
527 implicit none
528 ! I/F VARIABLES
529 real(kind=kreal) xx(*),yy(*),zz(*),vol,tline,almax,almin
530 ! LOCAL VARIABLES
531 integer(kind=kint) NN
532 integer(kind=kint) NG
533 parameter(nn=15,ng=3)
534 real(kind=kreal) h(nn),hl1(nn),hl2(nn),hl3(nn),hz(nn)
535 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
536 integer(kind=kint) I,L1,L2,LZ
537 real(kind=kreal) x1,x2,x3,xl1,xl2,zi
538 real(kind=kreal) a1,a2,al1,al2,al3,al4,al5,al6,al7,al8,al9
539 !C
540 vol = 0.0
541 ! LOOP FOR INTEGRATION POINTS
542 do lz=1,ng
543 zi=xg(ng,lz)
544 do l2=1,ng
545 xl2=xg(ng,l2)
546 x2 =(xl2+1.0)*0.5
547 do l1=1,ng
548 xl1=xg(ng,l1)
549 x1=0.5*(1.0-x2)*(xl1+1.0)
550 ! INTERPOLATION FUNCTION
551 x3=1.0-x1-x2
552 h(1)= 0.5*x1*(2.0*x1-2.-zi)*(1.0-zi)
553 h(2)= 0.5*x2*(2.0*x2-2.-zi)*(1.0-zi)
554 h(3)= 0.5*x3*(2.0*x3-2.-zi)*(1.0-zi)
555 h(4)= 0.5*x1*(2.0*x1-2.+zi)*(1.0+zi)
556 h(5)= 0.5*x2*(2.0*x2-2.+zi)*(1.0+zi)
557 h(6)= 0.5*x3*(2.0*x3-2.+zi)*(1.0+zi)
558 h(7)= 2.0*x1*x2*(1.0-zi)
559 h(8)= 2.0*x2*x3*(1.0-zi)
560 h(9)= 2.0*x1*x3*(1.0-zi)
561 h(10)=2.0*x1*x2*(1.0+zi)
562 h(11)=2.0*x2*x3*(1.0+zi)
563 h(12)=2.0*x1*x3*(1.0+zi)
564 h(13)=x1*(1.0-zi**2)
565 h(14)=x2*(1.0-zi**2)
566 h(15)=x3*(1.0-zi**2)
567 ! DERIVATIVE OF INTERPOLATION FUNCTION
568 ! FOR L1-COORDINATE
569 hl1(1)= 0.5*(4.0*x1-2.-zi)*(1.0-zi)
570 hl1(2)= 0.0
571 hl1(3)= 0.0
572 hl1(4)= 0.5*(4.0*x1-2.+zi)*(1.0+zi)
573 hl1(5)= 0.0
574 hl1(6)= 0.0
575 hl1(7)= 2.0*x2*(1.0-zi)
576 hl1(8)= 0.0
577 hl1(9)= 2.0*x3*(1.0-zi)
578 hl1(10)=2.0*x2*(1.0+zi)
579 hl1(11)=0.0
580 hl1(12)=2.0*x3*(1.0+zi)
581 hl1(13)=1.0-zi**2
582 hl1(14)=0.0
583 hl1(15)=0.0
584 ! FOR L2-COORDINATE
585 hl2(1)= 0.0
586 hl2(2)= 0.5*(4.0*x2-2.-zi)*(1.0-zi)
587 hl2(3)= 0.0
588 hl2(4)= 0.0
589 hl2(5)= 0.5*(4.0*x2-2.+zi)*(1.0+zi)
590 hl2(6)= 0.0
591 hl2(7)= 2.0*x1*(1.0-zi)
592 hl2(8)= 2.0*x3*(1.0-zi)
593 hl2(9)= 0.0
594 hl2(10)=2.0*x1*(1.0+zi)
595 hl2(11)=2.0*x3*(1.0+zi)
596 hl2(12)=0.0
597 hl2(13)=0.0
598 hl2(14)=1.0-zi**2
599 hl2(15)=0.0
600 ! FOR L3-COORDINATE
601 hl3(1)= 0.0
602 hl3(2)= 0.0
603 hl3(3)= 0.5*(4.0*x3-2.-zi)*(1.0-zi)
604 hl3(4)= 0.0
605 hl3(5)= 0.0
606 hl3(6)= 0.5*(4.0*x3-2.+zi)*(1.0+zi)
607 hl3(7)= 0.0
608 hl3(8)= 2.0*x2*(1.0-zi)
609 hl3(9)= 2.0*x1*(1.0-zi)
610 hl3(10)=0.0
611 hl3(11)=2.0*x2*(1.0+zi)
612 hl3(12)=2.0*x1*(1.0+zi)
613 hl3(13)=0.0
614 hl3(14)=0.0
615 hl3(15)=1.0-zi**2
616 ! FOR Z-COORDINATE
617 hz(1)= 0.5*x1*(-2.0*x1+1.0+2.0*zi)
618 hz(2)= 0.5*x2*(-2.0*x2+1.0+2.0*zi)
619 hz(3)= 0.5*x3*(-2.0*x3+1.0+2.0*zi)
620 hz(4)= 0.5*x1*( 2.0*x1-1.0+2.0*zi)
621 hz(5)= 0.5*x2*( 2.0*x2-1.0+2.0*zi)
622 hz(6)= 0.5*x3*( 2.0*x3-1.0+2.0*zi)
623 hz(7)= -2.0*x1*x2
624 hz(8)= -2.0*x2*x3
625 hz(9)= -2.0*x1*x3
626 hz(10)= 2.0*x1*x2
627 hz(11)= 2.0*x2*x3
628 hz(12)= 2.0*x1*x3
629 hz(13)=-2.0*x1*zi
630 hz(14)=-2.0*x2*zi
631 hz(15)=-2.0*x3*zi
632 ! JACOBI MATRIX
633 xj11=0.0
634 xj21=0.0
635 xj31=0.0
636 xj12=0.0
637 xj22=0.0
638 xj32=0.0
639 xj13=0.0
640 xj23=0.0
641 xj33=0.0
642 do i=1,nn
643 xj11=xj11+(hl1(i)-hl3(i))*xx(i)
644 xj21=xj21+(hl2(i)-hl3(i))*xx(i)
645 xj31=xj31+hz(i)*xx(i)
646 xj12=xj12+(hl1(i)-hl3(i))*yy(i)
647 xj22=xj22+(hl2(i)-hl3(i))*yy(i)
648 xj32=xj32+hz(i)*yy(i)
649 xj13=xj13+(hl1(i)-hl3(i))*zz(i)
650 xj23=xj23+(hl2(i)-hl3(i))*zz(i)
651 xj33=xj33+hz(i)*zz(i)
652 enddo
653 ! DETERMINANT OF JACOBIAN
654 det=xj11*xj22*xj33 &
655 +xj12*xj23*xj31 &
656 +xj13*xj21*xj32 &
657 -xj13*xj22*xj31 &
658 -xj12*xj21*xj33 &
659 -xj11*xj23*xj32
660 ! WEIGHT VALUE AT GAUSSIAN POINT
661 wg=wgt(ng,l1)*wgt(ng,l2)*wgt(ng,lz)*det*(1.0-x2)*0.25
662 do i = 1, nn
663 vol = vol + h(i)*wg
664 enddo
665 enddo
666 enddo
667 enddo
668
669 a1 = sqrt( (xx(7)-xx(1))**2+(yy(7)-yy(1))**2+(zz(7)-zz(1))**2 )
670 a2 = sqrt( (xx(2)-xx(7))**2+(yy(2)-yy(7))**2+(zz(2)-zz(7))**2 )
671 al1 = a1 + a2
672 a1 = sqrt( (xx(8)-xx(2))**2+(yy(8)-yy(2))**2+(zz(8)-zz(2))**2 )
673 a2 = sqrt( (xx(3)-xx(8))**2+(yy(3)-yy(8))**2+(zz(3)-zz(8))**2 )
674 al2 = a1 + a2
675 a1 = sqrt( (xx(9)-xx(3))**2+(yy(9)-yy(3))**2+(zz(9)-zz(3))**2 )
676 a2 = sqrt( (xx(1)-xx(9))**2+(yy(1)-yy(9))**2+(zz(1)-zz(9))**2 )
677 al3 = a1 + a2
678
679 a1 = sqrt( (xx(10)-xx(4))**2+(yy(10)-yy(4))**2+(zz(10)-zz(4))**2 )
680 a2 = sqrt( (xx(5)-xx(10))**2+(yy(5)-yy(10))**2+(zz(5)-zz(10))**2 )
681 al4 = a1 + a2
682 a1 = sqrt( (xx(11)-xx(5))**2+(yy(11)-yy(5))**2+(zz(11)-zz(5))**2 )
683 a2 = sqrt( (xx(6)-xx(11))**2+(yy(6)-yy(11))**2+(zz(6)-zz(11))**2 )
684 al5 = a1 + a2
685 a1 = sqrt( (xx(12)-xx(6))**2+(yy(12)-yy(6))**2+(zz(12)-zz(6))**2 )
686 a2 = sqrt( (xx(4)-xx(12))**2+(yy(4)-yy(12))**2+(zz(4)-zz(12))**2 )
687 al6 = a1 + a2
688
689 a1 = sqrt( (xx(13)-xx(1))**2+(yy(13)-yy(1))**2+(zz(13)-zz(1))**2 )
690 a2 = sqrt( (xx(4)-xx(13))**2+(yy(4)-yy(13))**2+(zz(4)-zz(13))**2 )
691 al7 = a1 + a2
692 a1 = sqrt( (xx(14)-xx(2))**2+(yy(14)-yy(2))**2+(zz(14)-zz(2))**2 )
693 a2 = sqrt( (xx(5)-xx(14))**2+(yy(5)-yy(14))**2+(zz(5)-zz(14))**2 )
694 al8 = a1 + a2
695 a1 = sqrt( (xx(15)-xx(3))**2+(yy(15)-yy(3))**2+(zz(15)-zz(3))**2 )
696 a2 = sqrt( (xx(6)-xx(15))**2+(yy(6)-yy(15))**2+(zz(6)-zz(15))**2 )
697 al9 = a1 + a2
698
699 almax = dmax1( al1,al2,al3,al4,al5,al6,al7,al8,al9 )
700 almin = dmin1( al1,al2,al3,al4,al5,al6,al7,al8,al9 )
701
702 end subroutine pre_352
703 !----------------------------------------------------------------------*
704 subroutine pre_362( XX,YY,ZZ,vol,almax,almin )
705 !----------------------------------------------------------------------*
706 !
707 ! CALCULATION 3D 20-NODE SOLID ELEMENT
708 !
709 use hecmw
711 implicit none
712 ! I/F VARIABLES
713 real(kind=kreal) xx(*),yy(*),zz(*),vol,almax,almin
714 ! LOCAL VARIABLES
715 integer(kind=kint) NN
716 integer(kind=kint) NG
717 parameter(nn=20,ng=3)
718 real(kind=kreal) h(nn),hr(nn),hs(nn),ht(nn)
719 real(kind=kreal) ri,si,ti,rp,sp,tp,rm,sm,tm
720 real(kind=kreal) xj11,xj21,xj31,xj12,xj22,xj32,xj13,xj23,xj33,det,wg
721 integer(kind=kint) I,LX,LY,LZ
722 real(kind=kreal) a1,a2,al1,al2,al3,al4,al5,al6,al7,al8,al9,al10,al11,al12
723 !C
724 vol = 0.0
725 ! LOOP FOR INTEGRATION POINTS
726 do lx=1,ng
727 ri=xg(ng,lx)
728 do ly=1,ng
729 si=xg(ng,ly)
730 do lz=1,ng
731 ti=xg(ng,lz)
732 rp=1.0+ri
733 sp=1.0+si
734 tp=1.0+ti
735 rm=1.0-ri
736 sm=1.0-si
737 tm=1.0-ti
738 ! INTERPOLATION FUNCTION
739 h(1)=-0.125*rm*sm*tm*(2.0+ri+si+ti)
740 h(2)=-0.125*rp*sm*tm*(2.0-ri+si+ti)
741 h(3)=-0.125*rp*sp*tm*(2.0-ri-si+ti)
742 h(4)=-0.125*rm*sp*tm*(2.0+ri-si+ti)
743 h(5)=-0.125*rm*sm*tp*(2.0+ri+si-ti)
744 h(6)=-0.125*rp*sm*tp*(2.0-ri+si-ti)
745 h(7)=-0.125*rp*sp*tp*(2.0-ri-si-ti)
746 h(8)=-0.125*rm*sp*tp*(2.0+ri-si-ti)
747 h(9)=0.25*(1.0-ri**2)*sm*tm
748 h(10)=0.25*rp*(1.0-si**2)*tm
749 h(11)=0.25*(1.0-ri**2)*sp*tm
750 h(12)=0.25*rm*(1.0-si**2)*tm
751 h(13)=0.25*(1.0-ri**2)*sm*tp
752 h(14)=0.25*rp*(1.0-si**2)*tp
753 h(15)=0.25*(1.0-ri**2)*sp*tp
754 h(16)=0.25*rm*(1.0-si**2)*tp
755 h(17)=0.25*rm*sm*(1.0-ti**2)
756 h(18)=0.25*rp*sm*(1.0-ti**2)
757 h(19)=0.25*rp*sp*(1.0-ti**2)
758 h(20)=0.25*rm*sp*(1.0-ti**2)
759 ! DERIVATIVE OF INTERPOLATION FUNCTION
760 ! FOR R-COORDINATE
761 hr(1)=-0.125*rm*sm*tm+0.125*sm*tm*(2.0+ri+si+ti)
762 hr(2)=+0.125*rp*sm*tm-0.125*sm*tm*(2.0-ri+si+ti)
763 hr(3)=+0.125*rp*sp*tm-0.125*sp*tm*(2.0-ri-si+ti)
764 hr(4)=-0.125*rm*sp*tm+0.125*sp*tm*(2.0+ri-si+ti)
765 hr(5)=-0.125*rm*sm*tp+0.125*sm*tp*(2.0+ri+si-ti)
766 hr(6)=+0.125*rp*sm*tp-0.125*sm*tp*(2.0-ri+si-ti)
767 hr(7)=+0.125*rp*sp*tp-0.125*sp*tp*(2.0-ri-si-ti)
768 hr(8)=-0.125*rm*sp*tp+0.125*sp*tp*(2.0+ri-si-ti)
769 hr(9 )=-0.50*ri*sm*tm
770 hr(10)=+0.25*(1.0-si**2)*tm
771 hr(11)=-0.50*ri*sp*tm
772 hr(12)=-0.25*(1.0-si**2)*tm
773 hr(13)=-0.50*ri*sm*tp
774 hr(14)=+0.25*(1.0-si**2)*tp
775 hr(15)=-0.50*ri*sp*tp
776 hr(16)=-0.25*(1.0-si**2)*tp
777 hr(17)=-0.25*sm*(1.0-ti**2)
778 hr(18)=+0.25*sm*(1.0-ti**2)
779 hr(19)=+0.25*sp*(1.0-ti**2)
780 hr(20)=-0.25*sp*(1.0-ti**2)
781 ! FOR S-COORDINATE
782 hs(1)=-0.125*rm*sm*tm+0.125*rm*tm*(2.0+ri+si+ti)
783 hs(2)=-0.125*rp*sm*tm+0.125*rp*tm*(2.0-ri+si+ti)
784 hs(3)=+0.125*rp*sp*tm-0.125*rp*tm*(2.0-ri-si+ti)
785 hs(4)=+0.125*rm*sp*tm-0.125*rm*tm*(2.0+ri-si+ti)
786 hs(5)=-0.125*rm*sm*tp+0.125*rm*tp*(2.0+ri+si-ti)
787 hs(6)=-0.125*rp*sm*tp+0.125*rp*tp*(2.0-ri+si-ti)
788 hs(7)=+0.125*rp*sp*tp-0.125*rp*tp*(2.0-ri-si-ti)
789 hs(8)=+0.125*rm*sp*tp-0.125*rm*tp*(2.0+ri-si-ti)
790 hs(9)=-0.25*(1.0-ri**2)*tm
791 hs(10)=-0.50*rp*si*tm
792 hs(11)=+0.25*(1.0-ri**2)*tm
793 hs(12)=-0.50*rm*si*tm
794 hs(13)=-0.25*(1.0-ri**2)*tp
795 hs(14)=-0.50*rp*si*tp
796 hs(15)=+0.25*(1.0-ri**2)*tp
797 hs(16)=-0.50*rm*si*tp
798 hs(17)=-0.25*rm*(1.0-ti**2)
799 hs(18)=-0.25*rp*(1.0-ti**2)
800 hs(19)=+0.25*rp*(1.0-ti**2)
801 hs(20)=+0.25*rm*(1.0-ti**2)
802 ! FOR T-COORDINATE
803 ht(1)=-0.125*rm*sm*tm+0.125*rm*sm*(2.0+ri+si+ti)
804 ht(2)=-0.125*rp*sm*tm+0.125*rp*sm*(2.0-ri+si+ti)
805 ht(3)=-0.125*rp*sp*tm+0.125*rp*sp*(2.0-ri-si+ti)
806 ht(4)=-0.125*rm*sp*tm+0.125*rm*sp*(2.0+ri-si+ti)
807 ht(5)=+0.125*rm*sm*tp-0.125*rm*sm*(2.0+ri+si-ti)
808 ht(6)=+0.125*rp*sm*tp-0.125*rp*sm*(2.0-ri+si-ti)
809 ht(7)=+0.125*rp*sp*tp-0.125*rp*sp*(2.0-ri-si-ti)
810 ht(8)=+0.125*rm*sp*tp-0.125*rm*sp*(2.0+ri-si-ti)
811 ht(9)=-0.25*(1.0-ri**2)*sm
812 ht(10)=-0.25*rp*(1.0-si**2)
813 ht(11)=-0.25*(1.0-ri**2)*sp
814 ht(12)=-0.25*rm*(1.0-si**2)
815 ht(13)=0.25*(1.0-ri**2)*sm
816 ht(14)=0.25*rp*(1.0-si**2)
817 ht(15)=0.25*(1.0-ri**2)*sp
818 ht(16)=0.25*rm*(1.0-si**2)
819 ht(17)=-0.5*rm*sm*ti
820 ht(18)=-0.5*rp*sm*ti
821 ht(19)=-0.5*rp*sp*ti
822 ht(20)=-0.5*rm*sp*ti
823 ! JACOBI MATRIX
824 xj11=0.0
825 xj21=0.0
826 xj31=0.0
827 xj12=0.0
828 xj22=0.0
829 xj32=0.0
830 xj13=0.0
831 xj23=0.0
832 xj33=0.0
833 do i=1,nn
834 xj11=xj11+hr(i)*xx(i)
835 xj21=xj21+hs(i)*xx(i)
836 xj31=xj31+ht(i)*xx(i)
837 xj12=xj12+hr(i)*yy(i)
838 xj22=xj22+hs(i)*yy(i)
839 xj32=xj32+ht(i)*yy(i)
840 xj13=xj13+hr(i)*zz(i)
841 xj23=xj23+hs(i)*zz(i)
842 xj33=xj33+ht(i)*zz(i)
843 enddo
844 ! DETERMINANT OF JACOBIAN
845 det=xj11*xj22*xj33 &
846 +xj12*xj23*xj31 &
847 +xj13*xj21*xj32 &
848 -xj13*xj22*xj31 &
849 -xj12*xj21*xj33 &
850 -xj11*xj23*xj32
851 ! WEIGHT VALUE AT GAUSSIAN POINT
852 wg=wgt(ng,lx)*wgt(ng,ly)*wgt(ng,lz)*det
853 do i=1,nn
854 vol = vol + h(i)*wg
855 enddo
856 enddo
857 enddo
858 enddo
859
860 a1 = sqrt( (xx(9)-xx(1))**2+(yy(9)-yy(1))**2+(zz(9)-zz(1))**2 )
861 a2 = sqrt( (xx(2)-xx(9))**2+(yy(2)-yy(9))**2+(zz(2)-zz(9))**2 )
862 al1 = a1 + a2
863 a1 = sqrt( (xx(10)-xx(2))**2+(yy(10)-yy(2))**2+(zz(10)-zz(2))**2 )
864 a2 = sqrt( (xx(3)-xx(10))**2+(yy(3)-yy(10))**2+(zz(3)-zz(10))**2 )
865 al2 = a1 + a2
866 a1 = sqrt( (xx(11)-xx(3))**2+(yy(11)-yy(3))**2+(zz(11)-zz(3))**2 )
867 a2 = sqrt( (xx(4)-xx(11))**2+(yy(4)-yy(11))**2+(zz(4)-zz(11))**2 )
868 al3 = a1 + a2
869 a1 = sqrt( (xx(12)-xx(4))**2+(yy(12)-yy(4))**2+(zz(12)-zz(4))**2 )
870 a2 = sqrt( (xx(1)-xx(12))**2+(yy(1)-yy(12))**2+(zz(1)-zz(12))**2 )
871 al4 = a1 + a2
872
873 a1 = sqrt( (xx(13)-xx(5))**2+(yy(13)-yy(5))**2+(zz(13)-zz(5))**2 )
874 a2 = sqrt( (xx(6)-xx(13))**2+(yy(6)-yy(13))**2+(zz(6)-zz(13))**2 )
875 al5 = a1 + a2
876 a1 = sqrt( (xx(14)-xx(6))**2+(yy(14)-yy(6))**2+(zz(14)-zz(6))**2 )
877 a2 = sqrt( (xx(7)-xx(14))**2+(yy(7)-yy(14))**2+(zz(7)-zz(14))**2 )
878 al6 = a1 + a2
879 a1 = sqrt( (xx(15)-xx(7))**2+(yy(15)-yy(7))**2+(zz(15)-zz(7))**2 )
880 a2 = sqrt( (xx(8)-xx(15))**2+(yy(8)-yy(15))**2+(zz(8)-zz(15))**2 )
881 al7 = a1 + a2
882 a1 = sqrt( (xx(16)-xx(8))**2+(yy(16)-yy(8))**2+(zz(16)-zz(8))**2 )
883 a2 = sqrt( (xx(5)-xx(16))**2+(yy(5)-yy(16))**2+(zz(5)-zz(16))**2 )
884 al8 = a1 + a2
885
886 a1 = sqrt( (xx(17)-xx(1))**2+(yy(17)-yy(1))**2+(zz(17)-zz(1))**2 )
887 a2 = sqrt( (xx(5)-xx(17))**2+(yy(5)-yy(17))**2+(zz(5)-zz(17))**2 )
888 al9 = a1 + a2
889 a1 = sqrt( (xx(18)-xx(2))**2+(yy(18)-yy(2))**2+(zz(18)-zz(2))**2 )
890 a2 = sqrt( (xx(6)-xx(18))**2+(yy(6)-yy(18))**2+(zz(6)-zz(18))**2 )
891 al10 = a1 + a2
892 a1 = sqrt( (xx(19)-xx(3))**2+(yy(19)-yy(3))**2+(zz(19)-zz(3))**2 )
893 a2 = sqrt( (xx(7)-xx(19))**2+(yy(7)-yy(19))**2+(zz(7)-zz(19))**2 )
894 al11 = a1 + a2
895 a1 = sqrt( (xx(20)-xx(4))**2+(yy(20)-yy(4))**2+(zz(20)-zz(4))**2 )
896 a2 = sqrt( (xx(8)-xx(20))**2+(yy(8)-yy(20))**2+(zz(8)-zz(20))**2 )
897 al12 = a1 + a2
898
899 almax = dmax1( al1,al2,al3,al4,al5,al6,al7,al8,al9,al10,al11,al12 )
900 almin = dmin1( al1,al2,al3,al4,al5,al6,al7,al8,al9 ,al10,al11,al12)
901
902 end subroutine pre_362
903end module m_precheck_lib_3d
This module provides data for gauss quadrature.
Definition: GaussM.f90:6
real(kind=kreal), dimension(3, 3) wgt
wieght of gauss points
Definition: GaussM.f90:10
real(kind=kreal), dimension(3, 3) xg
abscissa of gauss points
Definition: GaussM.f90:9
Definition: hecmw.f90:6
This module provides function to check input data of 3d static analysis.
subroutine pre_352(xx, yy, zz, vol, almax, almin)
subroutine pre_342(xx, yy, zz, vol, almax, almin)
subroutine pre_362(xx, yy, zz, vol, almax, almin)
subroutine pre_361(xx, yy, zz, vol, almax, almin)
subroutine pre_341(xx, yy, zz, vol, almax, almin)
subroutine pre_351(xx, yy, zz, vol, almax, almin)