36
37
38
40
41
42
43#include "implicit_f.inc"
44
45
46
47#include "mvsiz_p.inc"
48#include "param_c.inc"
49
50
51
52 INTEGER, INTENT(IN) :: NEL
53 INTEGER, INTENT(IN) :: NPT
54 INTEGER, INTENT(IN) :: ISMSTR
55 INTEGER, INTENT(IN) :: ISROT
56 INTEGER NGL(*)
57
58 my_real,
DIMENSION(NPROPG) ,
INTENT(IN) :: geo
59 my_real,
DIMENSION(MVSIZ) ,
INTENT(OUT) :: lc
61 . offg(*),volg(*),deltax(*), volg0(*),
62 . rx(*), ry(*), rz(*), sx(*), sy(*), sz(*),tx(*), ty(*), tz(*)
63 TYPE(DT_), INTENT(INOUT) :: DT
64
65
66
67 INTEGER I,ITET,IDELM,ICST
69 . a1,a2,a3,a4,
70 . a1x,a2x,a3x,a4x,a1y,a2y,a3y,a4y,a1z,a2z,a3z,a4z,
71 . vdefmin,vdefmax,aspmin,asptet,min_aspect,min_defv,delm
72
73 min_aspect = dt%BRICK_CST_COL_MIN
74 min_defv = dt%BRICK_CST_DEFV_MIN
75 vdefmin = geo(190)
76 vdefmax = geo(191)
77 aspmin = geo(192)
78 asptet = geo(193)
79 delm = vdefmin +vdefmax+aspmin+asptet
80 icst = 1
81 IF((min_aspect+min_defv)==zero.AND.
82 . (ismstr==1.OR.ismstr==3.OR.ismstr==11)) icst=0
83 idelm = 1
84 idelm = dt%IDEL_BRICK
85 IF(idelm ==0.AND.delm==zero) idelm = 0
86 IF((idelm+icst) >0) THEN
87 IF(isrot /= 1)THEN
88 DO i=1,nel
89 a1x = ry(i)*sz(i)-rz(i)*sy(i)
90 a1y = rz(i)*sx(i)-rx(i)*sz(i)
91 a1z = rx(i)*sy(i)-ry(i)*sx(i)
92 a1 = a1x*a1x+a1y*a1y+a1z*a1z
93
94 a2x = sy(i)*tz(i)-sz(i)*ty(i)
95 a2y = sz(i)*tx(i)-sx(i)*tz(i)
96 a2z = sx(i)*ty(i)-sy(i)*tx(i)
97 a2 = a2x*a2x+a2y*a2y+a2z*a2z
98
99 a3x = ty(i)*rz(i)-tz(i)*ry(i)
100 a3y = tz(i)*rx(i)-tx(i)*rz(i)
101 a3z = tx(i)*ry(i)-ty(i)*rx(i)
102 a3 = a3x*a3x+a3y*a3y+a3z*a3z
103
104 a4x = a1x+a2x+a3x
105 a4y = a1y+a2y+a3y
106 a4z = a1z+a2z+a3z
107 a4 = a4x*a4x+a4y*a4y+a4z*a4z
108
109 lc(i) = six*volg(i)/sqrt(
max(a1,a2,a3,a4))
110 END DO
111 ELSE
112 DO i=1,nel
113 lc(i) = deltax(i)
114 END DO
115 END IF
116 END IF
117
118 RETURN