1 : C BCOQ4 SOURCE CHAT 06/03/29 21:15:38 5360
2 : SUBROUTINE BCOQ4
3 : & (IGAU,XEL,SHPTOT,SHP,BGENE,DJAC,EXCEN,NOPLAN,IARR,IFS)
4 : IMPLICIT INTEGER(I-N)
5 : IMPLICIT REAL*8(A-H,O-Z)
6 : ************************************************************************
7 : *
8 : * B C O Q 4
9 : * ---------
10 : *
11 : * FONCTION:
12 : * ---------
13 : *
14 : * Calcul de la matrice "B" pour un COQ4 excentr{ ou non.
15 : *
16 : * PARAMETRES: (E)=ENTREE (S)=SORTIE
17 : * -----------
18 : *
19 : * IGAU (E) Num{ro du point d'int{gration.
20 : * XEL (E) Coordonn{es locales des noeuds de l'{l{ment.
21 : * SHPTOT (E) Fonctions de forme et d{riv{es dans l'espace de
22 : * r{f{rence.
23 : * SHPTOT(1,... = fonction,
24 : * SHPTOT(2,... = d{riv{e par rapport @ Qsi,
25 : * SHPTOT(3,... = d{riv{e par rapport @ Eta.
26 : * EXCEN (E) Excentrement.
27 : * NOPLAN (E) = 1 si {l{ment non plan,
28 : * = 0 sinon.
29 : * IFS (E) = 1 si demande de Matrice "B" compl}te
30 : * (pour SIGMA et EPSILON),
31 : * = 0 sinon.
32 : * SHP (S) Fonctions de forme et d{riv{es dans l'espace
33 : * g{om{trique.
34 : * SHP(1,... = fonction,
35 : * SHP(2,... = d{riv{e par rapport @ "X" local,
36 : * SHP(3,... = d{riv{e par rapport @ "Y" local.
37 : * - au point d'int{gration n.5 si IFS=1,
38 : * - au point d'int{gration IGAU courant sinon.
39 : * DJAC (S) Jacobien au point d'int{gration IGAU.
40 : * BGENE (S) Matrice "B" au point d'int{gration IGAU.
41 : * IARR (S) Diff{rent de 0 si erreur.
42 : *
43 : PARAMETER (LRE=24,NST=8,NBNO=4)
44 : REAL*8 XEL(3,*),BGENE(NST,*),SHP(6,*),SHPTOT(6,NBNO,*)
45 : *
46 : * CONSTANTES:
47 : * -----------
48 : *
49 : * NBNO = Nombre de noeuds.
50 : * LRE = Nombre de colonnes de la Matrice "B".
51 : * NST = Nombre de composantes de contraintes.
52 : *
53 : * VARIABLES:
54 : * ----------
55 : *
56 : * FULL = .TRUE. si on veut une Matrice "B" compl}te.
57 : * GAUSS = .TRUE. si utilisation d'un point d'int{gration n.1 @ 4.
58 : * SHP7 = Fonctions d'interpolation mixte pour le cisaillement X-Z
59 : * SHP8 = Fonctions d'interpolation mixte pour le cisaillement Y-Z
60 : *
61 : REAL*8 SHP7(4:5,NBNO),SHP8(4:5,NBNO)
62 : LOGICAL FULL,GAUSS
63 : *
64 : * MODE DE FONCTIONNEMENT:
65 : * -----------------------
66 : *
67 : * POUR LA RIGIDITE:
68 : *
69 : * Pour la membrane,
70 : * - int{gration r{duite (1 point) pour le cisaillement,
71 : * - int{gration r{duite de tout si {l{ment non plan.
72 : *
73 : * Pour le cisaillement transverse, voir:
74 : * DONEA / LAMAIN (1987) - ISPRA
75 : * ({l{ment de type "A.N.S.")
76 : *
77 : * Pour le couplage membrane-flexion (excentrement), int{gration
78 : * r{duite analogue @ membrane, non justifi{e au niveau de la
79 : * Rigidit{, mais indispensable pour le calcul pratique des forces
80 : * nodales @ partir des contraintes, par "BSIGMA".
81 : *
82 : * POUR LES CONTRAINTES:
83 : *
84 : * Pour {viter des surcontraintes dues @ un "hourglass", il ne faut
85 : * pas calculer aux 4 points de Gauss des contraintes qui ont {t{
86 : * calcul{es au centre lors de la cr{ation de la Matrice de Rigidit{.
87 : * Par cons{quent, il faut prendre
88 : * - toujours : le cisaillement de membrane au centre,
89 : * - quand NOPLAN=1: toutes les contraintes de membrane au centre.
90 : * (centre = point d'int{gration n.5)
91 : * Calcul au centre dans les memes conditions pour le couplage
92 : * membrane-flexion.
93 : *
94 : * Cependant, pour que le champ de contraintes produit soit facile
95 : * @ lire par l'utilisateur et @ traiter par les op{rateurs tels
96 : * que:
97 : * - CHANGER CHPOINT,
98 : * - TRACER, ...
99 : * on inscrit aux points de Gauss 1 @ 4 les valeurs de contrainte
100 : * au point d'int{gration n.5 pour les parties concern{es.
101 : * De plus, on remplit compl}tement le champ de valeurs au point
102 : * d'int{gration n.5 .
103 : *
104 : * AUTEUR, DATE DE CREATION:
105 : * -------------------------
106 : *
107 : * Pascal MANIGOT 17 Juin 1991
108 : *
109 : ************************************************************************
110 : *
111 : IARR=0
112 : FULL = IFS .EQ. 1
113 : GAUSS = IGAU .LE. 4
114 : *
115 : *
116 : * Matrice jacobienne:
117 : * -------------------
118 : *
119 : DO 110 I=1,NBNO
120 : SHP(1,I)=SHPTOT(1,I,IGAU)
121 : SHP(2,I)=SHPTOT(2,I,IGAU)
122 : SHP(3,I)=SHPTOT(3,I,IGAU)
123 : 110 CONTINUE
124 : * END DO
125 : *
126 : DXDQSI = 0.D0
127 : DXDETA = 0.D0
128 : DYDQSI = 0.D0
129 : DYDETA = 0.D0
130 : DO 120 I=1,NBNO
131 : DXDQSI=DXDQSI+SHP(2,I)*XEL(1,I)
132 : DXDETA=DXDETA+SHP(3,I)*XEL(1,I)
133 : DYDQSI=DYDQSI+SHP(2,I)*XEL(2,I)
134 : DYDETA=DYDETA+SHP(3,I)*XEL(2,I)
135 : 120 CONTINUE
136 : * END DO
137 : DJAC=DXDQSI*DYDETA-DXDETA*DYDQSI
138 : IF(DJAC.LE. 0.D0) THEN
139 : IARR=1
140 : RETURN
141 : END IF
142 : *
143 : * Fonctions de forme dans l'espace g{om{trique:
144 : DQSIDX= DYDETA/DJAC
145 : DQSIDY=-DXDETA/DJAC
146 : DETADX=-DYDQSI/DJAC
147 : DETADY= DXDQSI/DJAC
148 : DO 130 I=1,NBNO
149 : TT=SHP(2,I)*DQSIDX+SHP(3,I)*DETADX
150 : SHP(3,I)=SHP(2,I)*DQSIDY+SHP(3,I)*DETADY
151 : SHP(2,I)=TT
152 : 130 CONTINUE
153 : * END DO
154 : *
155 : *
156 : * Matrice "B":
157 : * ------------
158 : *
158 : CALL ZERO(BGENE,NST,LRE)
159 : *
160 : IF (GAUSS .OR. FULL) THEN
161 : * Pr{paration pour cisaillement transverse:
162 : *
163 : N1 = 1
164 : N2 = 2
165 : DX1 = SHPTOT(2,N1,IGAU)*XEL(1,N1)+SHPTOT(2,N2,IGAU)*XEL(1,N2)
166 : DY1 = SHPTOT(2,N1,IGAU)*XEL(2,N1)+SHPTOT(2,N2,IGAU)*XEL(2,N2)
167 : N1 = 3
168 : N2 = 4
169 : DX2 = SHPTOT(2,N1,IGAU)*XEL(1,N1)+SHPTOT(2,N2,IGAU)*XEL(1,N2)
170 : DY2 = SHPTOT(2,N1,IGAU)*XEL(2,N1)+SHPTOT(2,N2,IGAU)*XEL(2,N2)
171 : N1 = 1
172 : N2 = 4
173 : DX3 = SHPTOT(3,N1,IGAU)*XEL(1,N1)+SHPTOT(3,N2,IGAU)*XEL(1,N2)
174 : DY3 = SHPTOT(3,N1,IGAU)*XEL(2,N1)+SHPTOT(3,N2,IGAU)*XEL(2,N2)
175 : N1 = 2
176 : N2 = 3
177 : DX4 = SHPTOT(3,N1,IGAU)*XEL(1,N1)+SHPTOT(3,N2,IGAU)*XEL(1,N2)
178 : DY4 = SHPTOT(3,N1,IGAU)*XEL(2,N1)+SHPTOT(3,N2,IGAU)*XEL(2,N2)
179 : *
180 : D2JAC = 2.D0 * DJAC
181 : *
182 : SHP7(4,1) = (- DYDETA*DY1 + DYDQSI*DY3) / D2JAC
183 : SHP7(4,2) = (- DYDETA*DY1 + DYDQSI*DY4) / D2JAC
184 : SHP7(4,3) = (- DYDETA*DY2 + DYDQSI*DY4) / D2JAC
185 : SHP7(4,4) = (- DYDETA*DY2 + DYDQSI*DY3) / D2JAC
186 : SHP7(5,1) = ( DYDETA*DX1 - DYDQSI*DX3) / D2JAC
187 : SHP7(5,2) = ( DYDETA*DX1 - DYDQSI*DX4) / D2JAC
188 : SHP7(5,3) = ( DYDETA*DX2 - DYDQSI*DX4) / D2JAC
189 : SHP7(5,4) = ( DYDETA*DX2 - DYDQSI*DX3) / D2JAC
190 : *
191 : SHP8(4,1) = ( DXDETA*DY1 - DXDQSI*DY3) / D2JAC
192 : SHP8(4,2) = ( DXDETA*DY1 - DXDQSI*DY4) / D2JAC
193 : SHP8(4,3) = ( DXDETA*DY2 - DXDQSI*DY4) / D2JAC
194 : SHP8(4,4) = ( DXDETA*DY2 - DXDQSI*DY3) / D2JAC
195 : SHP8(5,1) = (- DXDETA*DX1 + DXDQSI*DX3) / D2JAC
196 : SHP8(5,2) = (- DXDETA*DX1 + DXDQSI*DX4) / D2JAC
197 : SHP8(5,3) = (- DXDETA*DX2 + DXDQSI*DX4) / D2JAC
198 : SHP8(5,4) = (- DXDETA*DX2 + DXDQSI*DX3) / D2JAC
199 : *
200 : END IF
201 : *
202 : DO 300 I=1,NBNO
203 : K = 6*(I-1)
204 : *
205 : IF (GAUSS .OR. FULL) THEN
206 : IF (NOPLAN .EQ. 0) THEN
207 : * Membrane:
208 : BGENE(1,K+1)=SHP(2,I)
209 : BGENE(2,K+2)=SHP(3,I)
210 : * Couplage membrane-flexion:
211 : BGENE(1,K+5) = SHP(2,I) * EXCEN
212 : BGENE(2,K+4) = -SHP(3,I) * EXCEN
213 : ENDIF
214 : * Flexion:
215 : BGENE(4,K+5)=SHP(2,I)
216 : BGENE(5,K+4)=-SHP(3,I)
217 : BGENE(6,K+4)=-SHP(2,I)
218 : BGENE(6,K+5)=SHP(3,I)
219 : * Cisaillement transverse:
220 : BGENE(7,K+3)=SHP(2,I)
221 : BGENE(7,K+4) = SHP7(4,I)
222 : BGENE(7,K+5) = SHP7(5,I)
223 : BGENE(8,K+3)=SHP(3,I)
224 : BGENE(8,K+4) = SHP8(4,I)
225 : BGENE(8,K+5) = SHP8(5,I)
226 : ENDIF
227 : *
228 : IF (.NOT. GAUSS) THEN
229 : IF (NOPLAN.EQ.1 .OR. FULL) THEN
230 : * Membrane:
231 : BGENE(1,K+1)=SHP(2,I)
232 : BGENE(2,K+2)=SHP(3,I)
233 : * Couplage membrane-flexion:
234 : BGENE(1,K+5) = SHP(2,I) * EXCEN
235 : BGENE(2,K+4) = -SHP(3,I) * EXCEN
236 : ENDIF
237 : * Cisaillement de membrane et de couplage membrane-flexion:
238 : BGENE(3,K+1)=SHP(3,I)
239 : BGENE(3,K+2)=SHP(2,I)
240 : BGENE(3,K+4) = -SHP(2,I) * EXCEN
241 : BGENE(3,K+5) = SHP(3,I) * EXCEN
242 : ENDIF
243 : *
244 : 300 CONTINUE
245 : * END DO
246 : *
247 : *
248 : IF (GAUSS .AND. FULL) THEN
249 : * On finit de compl{ter la Matrice "B", pour les contraintes ou
250 : * les d{formations, par des valeurs calcul{es au centre:
251 : *
252 : DO 500 I=1,NBNO
253 : SHP(1,I)=SHPTOT(1,I, 5 )
254 : SHP(2,I)=SHPTOT(2,I, 5 )
255 : SHP(3,I)=SHPTOT(3,I, 5 )
256 : 500 CONTINUE
257 : * END DO
258 : * FONCTIONS DE FORME DANS LA G{OM{TRIE R{ELLE:
258 : CALL JACOBI(XEL,SHP,2,NBNO,DJAC5)
259 : IF(DJAC5 .LE. 0.D0) THEN
260 : IARR=1
261 : RETURN
262 : END IF
263 : *
264 : DO 520 I=1,NBNO
265 : K = 6*(I-1)
266 : IF (NOPLAN .EQ. 1) THEN
267 : * Membrane:
268 : BGENE(1,K+1)=SHP(2,I)
269 : BGENE(2,K+2)=SHP(3,I)
270 : * Couplage membrane-flexion:
271 : BGENE(1,K+5) = SHP(2,I) * EXCEN
272 : BGENE(2,K+4) = -SHP(3,I) * EXCEN
273 : ENDIF
274 : * Cisaillement de membrane et de couplage membrane-flexion:
275 : BGENE(3,K+1)=SHP(3,I)
276 : BGENE(3,K+2)=SHP(2,I)
277 : BGENE(3,K+4) = -SHP(2,I) * EXCEN
278 : BGENE(3,K+5) = SHP(3,I) * EXCEN
279 : 520 CONTINUE
280 : * END DO
281 : *
282 : ENDIF
283 : *
284 : END
285 :
286 :
© Cast3M 2003 - All rights reserved.
Disclaimer