1 : C INELQU SOURCE GOUNAND 06/12/19 21:15:22 5612 2 : SUBROUTINE INELQU(MYLRFS,IMPR,IRET) 3 : IMPLICIT REAL*8 (A-H,O-Z) 4 : IMPLICIT INTEGER (I-N) 5 : C*********************************************************************** 6 : C NOM : INELQU 7 : C PROJET : Noyau linéaire NLIN 8 : C DESCRIPTION : Remplit le segment des éléments de référence 9 : C avec les éléments de référence de dimension 2, 10 : C de forme géométrique carrée. 11 : C 12 : C REFERENCES : 13 : C LANGAGE : ESOPE 14 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) 15 : C mél : gounand@semt2.smts.cea.fr 16 : C*********************************************************************** 17 : C APPELES : INILRF, INILAG, GBAPCO, PROBAP 18 : C APPELE PAR : INLRFS 19 : C*********************************************************************** 20 : C ENTREES : - 21 : C ENTREES/SORTIES : MYLRFS 22 : C SORTIES : - 23 : C CODE RETOUR (IRET) : = 0 si tout s'est bien passé 24 : C*********************************************************************** 25 : C VERSION : v1, 23/03/00, version initiale 26 : C HISTORIQUE : v1, 23/03/00, création 27 : C HISTORIQUE : v2, 10/05/00, modif. du segment ELREF 28 : C HISTORIQUE : 29 : C HISTORIQUE : 30 : C*********************************************************************** 31 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 32 : C en cas de modification de ce sous-programme afin de faciliter 33 : C la maintenance ! 34 : C*********************************************************************** 35 : -INC CCOPTIO 36 : CBEGININCLUDE SELREF 37 : SEGMENT ELREF 38 : CHARACTER*(LNNOM) NOMLRF 39 : CHARACTER*(LNFORM) FORME 40 : CHARACTER*(LNTYPL) TYPEL 41 : CHARACTER*(LNESP) ESPACE 42 : INTEGER DEGRE 43 : REAL*8 XCONOD(NDIMEL,NBNOD) 44 : INTEGER NPQUAF(NBDDL) 45 : INTEGER NUMCMP(NBDDL) 46 : INTEGER QUENOD(NBDDL) 47 : INTEGER ORDDER(NDIMEL,NBDDL) 48 : POINTEUR MBPOLY.POLYNS 49 : ENDSEGMENT 50 : SEGMENT ELREFS 51 : POINTEUR LISEL(0).ELREF 52 : ENDSEGMENT 53 : CENDINCLUDE SELREF 54 : POINTEUR MYLRFS.ELREFS 55 : POINTEUR ELCOUR.ELREF 56 : CBEGININCLUDE SPOLYNO 57 : SEGMENT POLYNO 58 : REAL*8 COEMON(NBMON) 59 : INTEGER EXPMON(NDIML,NBMON) 60 : ENDSEGMENT 61 : SEGMENT POLYNS 62 : POINTEUR LIPOLY(NBPOLY).POLYNO 63 : ENDSEGMENT 64 : CENDINCLUDE SPOLYNO 65 : POINTEUR MYBPOL.POLYNS 66 : POINTEUR MBPTMP.POLYNS 67 : INTEGER NBMON,NDIML 68 : POINTEUR MYPOL.POLYNO 69 : * 70 : INTEGER IMPR,IRET 71 : * Elément de nom : L2D0QU1 72 : REAL*8 ZERO 73 : PARAMETER (ZERO=0.D0) 74 : * Elément de nom : L2D1QU3 75 : REAL*8 UNS2 76 : PARAMETER (UNS2=0.5D0) 77 : * Elément de nom : H1D1QU4 78 : REAL*8 UN 79 : PARAMETER (UN=1.D0) 80 : * 81 : INTEGER INDDL 82 : * 83 : * Executable statements 84 : * 85 : IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans inelqu' 86 : * 87 : * Elément de nom : L2D0QU1 88 : * Sur un carré : élément de Lagrange, fonction L2, approximation 89 : * nodale, espace de référence de dimension 2, 1 noeud, 1 degré de 90 : * liberté, degré de l'approximation : 0 91 : * 92 : * In INILRF : SEGINI ELCOUR 92 : CALL INILRF('L2D0QU1','CARRE','LAGRANGE','L2', 93 : $ 2,1,1,0, 94 : $ ELCOUR, 95 : $ IMPR,IRET) 96 : IF (IRET.NE.0) GOTO 9999 97 : ELCOUR.XCONOD(1,1)=ZERO 98 : ELCOUR.XCONOD(2,1)=ZERO 99 : ELCOUR.NPQUAF(1)=9 100 : ELCOUR.NUMCMP(1)=1 101 : * Initialise la correspondance ddl-noeud+ord.der 101 : CALL INILAG(ELCOUR,IMPR,IRET) 102 : IF (IRET.NE.0) GOTO 9999 103 : * Génère une base polynômiale complète (dimension 2, degré 0) 103 : CALL GBAPCO(2,0,MYBPOL,IMPR,IRET) 104 : IF (IRET.NE.0) GOTO 9999 105 : ELCOUR.MBPOLY=MYBPOL 106 : SEGDES ELCOUR 107 : MYLRFS.LISEL(**)=ELCOUR 108 : * 109 : * Elément de nom : L2D1QU3 110 : * Sur un carré : élément de Lagrange, fonction L2, approximation 111 : * nodale, espace de référence de dimension 2, 3 noeuds, 3 degrés de 112 : * liberté, degré de l'approximation : 1 113 : * 114 : * In INILRF : SEGINI ELCOUR 114 : CALL INILRF('L2D1QU3','CARRE','LAGRANGE','L2', 115 : $ 2,3,3,1, 116 : $ ELCOUR, 117 : $ IMPR,IRET) 118 : IF (IRET.NE.0) GOTO 9999 119 : ELCOUR.XCONOD(1,1)=UNS2 120 : ELCOUR.XCONOD(2,1)=UNS2 121 : ELCOUR.XCONOD(1,2)=-UNS2 122 : ELCOUR.XCONOD(2,2)=ZERO 123 : ELCOUR.XCONOD(1,3)=ZERO 124 : ELCOUR.XCONOD(2,3)=-UNS2 125 : ELCOUR.NPQUAF(1)=9 126 : ELCOUR.NUMCMP(1)=1 127 : ELCOUR.NPQUAF(2)=9 128 : ELCOUR.NUMCMP(2)=2 129 : ELCOUR.NPQUAF(3)=9 130 : ELCOUR.NUMCMP(3)=3 131 : * Initialise la correspondance ddl-noeud+ord.der 131 : CALL INILAG(ELCOUR,IMPR,IRET) 132 : IF (IRET.NE.0) GOTO 9999 133 : * Génère une base polynômiale complète (dimension 2, degré 1) 133 : CALL GBAPCO(2,1,MYBPOL,IMPR,IRET) 134 : IF (IRET.NE.0) GOTO 9999 135 : ELCOUR.MBPOLY=MYBPOL 136 : SEGDES ELCOUR 137 : MYLRFS.LISEL(**)=ELCOUR 138 : * 139 : * Elément de nom : H1D1QU4 140 : * Sur un carré : élément de Lagrange, fonction H1, approximation 141 : * nodale, espace de référence de dimension 2, 4 noeuds, 4 degrés de 142 : * liberté, degré de l'approximation : 1 143 : * 144 : * In INILRF : SEGINI ELCOUR 144 : CALL INILRF('H1D1QU4','CARRE','LAGRANGE','H1', 145 : $ 2,4,4,1, 146 : $ ELCOUR, 147 : $ IMPR,IRET) 148 : IF (IRET.NE.0) GOTO 9999 149 : ELCOUR.XCONOD(1,1)=-UN 150 : ELCOUR.XCONOD(2,1)=-UN 151 : ELCOUR.XCONOD(1,2)=UN 152 : ELCOUR.XCONOD(2,2)=-UN 153 : ELCOUR.XCONOD(1,3)=UN 154 : ELCOUR.XCONOD(2,3)=UN 155 : ELCOUR.XCONOD(1,4)=-UN 156 : ELCOUR.XCONOD(2,4)=UN 157 : * Les d.d.l. sont aux noeuds 1,3,5,7 158 : DO 203 INDDL=1,4 159 : ELCOUR.NPQUAF(INDDL)=(2*INDDL)-1 160 : ELCOUR.NUMCMP(INDDL)=1 161 : 203 CONTINUE 162 : * Initialise la correspondance ddl-noeud+ord.der 162 : CALL INILAG(ELCOUR,IMPR,IRET) 163 : IF (IRET.NE.0) GOTO 9999 164 : * Génère une base polynômiale complète (dimension 1, degré 1) 165 : * In GBAPCO : SEGINI MBPTMP.LIPOLY(*) 165 : CALL GBAPCO(1,1,MBPTMP,IMPR,IRET) 166 : IF (IRET.NE.0) GOTO 9999 167 : * Puis la base produit 167 : CALL PROBAP(MBPTMP,MBPTMP,MYBPOL,IMPR,IRET) 168 : IF (IRET.NE.0) GOTO 9999 169 : SEGACT MBPTMP*MOD 170 : SEGSUP MBPTMP.LIPOLY(*) 171 : SEGSUP MBPTMP 172 : ELCOUR.MBPOLY=MYBPOL 173 : SEGDES ELCOUR 174 : MYLRFS.LISEL(**)=ELCOUR 175 : * 176 : * Elément de nom : CRD1QU4 177 : * Sur un carré : élément de Lagrange, fonction continue au centre 178 : * des faces, approximation 179 : * nodale, espace de référence de dimension 2, 4 noeuds, 4 degrés de 180 : * liberté, degré de l'approximation : 1 181 : * 182 : * In INILRF : SEGINI ELCOUR 182 : CALL INILRF('CRD1QU4','CARRE','LAGRANGE','HFAC', 183 : $ 2,4,4,1, 184 : $ ELCOUR, 185 : $ IMPR,IRET) 186 : IF (IRET.NE.0) GOTO 9999 187 : ELCOUR.XCONOD(1,1)=ZERO 188 : ELCOUR.XCONOD(2,1)=-UN 189 : ELCOUR.XCONOD(1,2)=UN 190 : ELCOUR.XCONOD(2,2)=ZERO 191 : ELCOUR.XCONOD(1,3)=ZERO 192 : ELCOUR.XCONOD(2,3)=UN 193 : ELCOUR.XCONOD(1,4)=-UN 194 : ELCOUR.XCONOD(2,4)=ZERO 195 : * Les d.d.l. sont aux noeuds 2,4,6,8 196 : DO INDDL=1,4 197 : ELCOUR.NPQUAF(INDDL)=2*INDDL 198 : ELCOUR.NUMCMP(INDDL)=1 199 : ENDDO 200 : * Initialise la correspondance ddl-noeud+ord.der 200 : CALL INILAG(ELCOUR,IMPR,IRET) 201 : IF (IRET.NE.0) GOTO 9999 202 : * Génère une base polynômiale complète (dimension 2, degré 1) 202 : CALL GBAPCO(2,1,MYBPOL,IMPR,IRET) 203 : IF (IRET.NE.0) GOTO 9999 204 : * On rajoute les polynômes spécifiques à crouzeix-raviart quadrilatère 204 : CALL GPOCRQ(2,MYBPOL,IMPR,IRET) 205 : IF (IRET.NE.0) GOTO 9999 206 : ELCOUR.MBPOLY=MYBPOL 207 : SEGDES ELCOUR 208 : MYLRFS.LISEL(**)=ELCOUR 209 : * 210 : * Elément de nom : H1D2QU9 211 : * Sur un carré : élément de Lagrange, fonction H1, approximation 212 : * nodale, espace de référence de dimension 2, 9 noeuds, 9 degrés de 213 : * liberté, degré de l'approximation : 2 214 : * 215 : * In INILRF : SEGINI ELCOUR 215 : CALL INILRF('H1D2QU9','CARRE','LAGRANGE','H1', 216 : $ 2,9,9,2, 217 : $ ELCOUR, 218 : $ IMPR,IRET) 219 : IF (IRET.NE.0) GOTO 9999 220 : ELCOUR.XCONOD(1,1)=-UN 221 : ELCOUR.XCONOD(2,1)=-UN 222 : ELCOUR.XCONOD(1,2)=ZERO 223 : ELCOUR.XCONOD(2,2)=-UN 224 : ELCOUR.XCONOD(1,3)=UN 225 : ELCOUR.XCONOD(2,3)=-UN 226 : ELCOUR.XCONOD(1,4)=UN 227 : ELCOUR.XCONOD(2,4)=ZERO 228 : ELCOUR.XCONOD(1,5)=UN 229 : ELCOUR.XCONOD(2,5)=UN 230 : ELCOUR.XCONOD(1,6)=ZERO 231 : ELCOUR.XCONOD(2,6)=UN 232 : ELCOUR.XCONOD(1,7)=-UN 233 : ELCOUR.XCONOD(2,7)=UN 234 : ELCOUR.XCONOD(1,8)=-UN 235 : ELCOUR.XCONOD(2,8)=ZERO 236 : ELCOUR.XCONOD(1,9)=ZERO 237 : ELCOUR.XCONOD(2,9)=ZERO 238 : * Les d.d.l. sont aux noeuds 1,2,3,4,5,6,7,8,9 239 : DO 209 INDDL=1,9 240 : ELCOUR.NPQUAF(INDDL)=INDDL 241 : ELCOUR.NUMCMP(INDDL)=1 242 : 209 CONTINUE 243 : * Initialise la correspondance ddl-noeud+ord.der 243 : CALL INILAG(ELCOUR,IMPR,IRET) 244 : IF (IRET.NE.0) GOTO 9999 245 : * Génère une base polynômiale complète (dimension 1, degré 2) 246 : * In GBAPCO : SEGINI MBPTMP.LIPOLY(*) 246 : CALL GBAPCO(1,2,MBPTMP,IMPR,IRET) 247 : IF (IRET.NE.0) GOTO 9999 248 : * Puis la base produit 248 : CALL PROBAP(MBPTMP,MBPTMP,MYBPOL,IMPR,IRET) 249 : IF (IRET.NE.0) GOTO 9999 250 : SEGACT MBPTMP*MOD 251 : SEGSUP MBPTMP.LIPOLY(*) 252 : SEGSUP MBPTMP 253 : ELCOUR.MBPOLY=MYBPOL 254 : SEGDES ELCOUR 255 : MYLRFS.LISEL(**)=ELCOUR 256 : * 257 : * Elément de nom : H1D2QU8 258 : * Sur un carré : élément de Lagrange incomplet (Serendip), 259 : * fonction H1, approximation nodale, espace de référence de dimension 2, 260 : * 8 noeuds, 8 degrés de liberté, degré de l'approximation : 2 261 : * 262 : * In INILRF : SEGINI ELCOUR 262 : CALL INILRF('H1D2QU8','CARRE','LAGRANGE','H1', 263 : $ 2,8,8,2, 264 : $ ELCOUR, 265 : $ IMPR,IRET) 266 : IF (IRET.NE.0) GOTO 9999 267 : ELCOUR.XCONOD(1,1)=-UN 268 : ELCOUR.XCONOD(2,1)=-UN 269 : ELCOUR.XCONOD(1,2)=ZERO 270 : ELCOUR.XCONOD(2,2)=-UN 271 : ELCOUR.XCONOD(1,3)=UN 272 : ELCOUR.XCONOD(2,3)=-UN 273 : ELCOUR.XCONOD(1,4)=UN 274 : ELCOUR.XCONOD(2,4)=ZERO 275 : ELCOUR.XCONOD(1,5)=UN 276 : ELCOUR.XCONOD(2,5)=UN 277 : ELCOUR.XCONOD(1,6)=ZERO 278 : ELCOUR.XCONOD(2,6)=UN 279 : ELCOUR.XCONOD(1,7)=-UN 280 : ELCOUR.XCONOD(2,7)=UN 281 : ELCOUR.XCONOD(1,8)=-UN 282 : ELCOUR.XCONOD(2,8)=ZERO 283 : * Les d.d.l. sont aux noeuds 1,2,3,4,5,6,7,8 284 : DO 211 INDDL=1,8 285 : ELCOUR.NPQUAF(INDDL)=INDDL 286 : ELCOUR.NUMCMP(INDDL)=1 287 : 211 CONTINUE 288 : * Initialise la correspondance ddl-noeud+ord.der 288 : CALL INILAG(ELCOUR,IMPR,IRET) 289 : IF (IRET.NE.0) GOTO 9999 290 : * Génère une base polynômiale complète (dimension 2, degré 2) 290 : CALL GBAPCO(2,2,MYBPOL,IMPR,IRET) 291 : IF (IRET.NE.0) GOTO 9999 292 : * On rajoute (\ksi^2 \eta et \ksi \eta^2) 293 : SEGACT MYBPOL*MOD 294 : NDIML=2 295 : NBMON=1 296 : SEGINI MYPOL 297 : MYPOL.COEMON(1)=UN 298 : MYPOL.EXPMON(1,1)=2 299 : MYPOL.EXPMON(2,1)=1 300 : SEGDES MYPOL 301 : MYBPOL.LIPOLY(**)=MYPOL 302 : NDIML=2 303 : NBMON=1 304 : SEGINI MYPOL 305 : MYPOL.COEMON(1)=UN 306 : MYPOL.EXPMON(1,1)=1 307 : MYPOL.EXPMON(2,1)=2 308 : SEGDES MYPOL 309 : MYBPOL.LIPOLY(**)=MYPOL 310 : * Voilà c'est fait 311 : SEGDES MYBPOL 312 : ELCOUR.MBPOLY=MYBPOL 313 : SEGDES ELCOUR 314 : MYLRFS.LISEL(**)=ELCOUR 315 : * 316 : * Normal termination 317 : * 318 : IRET=0 319 : RETURN 320 : * 321 : * Format handling 322 : * 323 : * 324 : * Error handling 325 : * 326 : 9999 CONTINUE 327 : IRET=1 328 : WRITE(IOIMP,*) 'An error was detected in subroutine inelqu' 329 : RETURN 330 : * 331 : * End of subroutine INELQU 332 : * 333 : END 334 : 335 : 336 : 337 :
© Cast3M 2003 - All rights reserved.
Disclaimer