1 : C INGJ10 SOURCE GOUNAND 05/12/21 21:32:32 5281 2 : SUBROUTINE INGJ10(MYPGS,IMPR,IRET) 3 : IMPLICIT REAL*8 (A-H,O-Z) 4 : IMPLICIT INTEGER (I-N) 5 : C*********************************************************************** 6 : C NOM : INGJ10 7 : C PROJET : Noyau linéaire NLIN 8 : C DESCRIPTION : Remplit le segment des méthodes d'intégration 9 : C avec des méthodes d'intégration numérique de 10 : C Gauss-Jacobi (\alpha=1, \beta=0) à une dimension sur 11 : C l'intervalle [0,1] (ordre 1 à 11). 12 : C 13 : C On intègre donc \int_0^1 (1-x) f(x) dx de manière 14 : C approchée. 15 : C 16 : C Ces méthodes sont utilisés pour générer des formules 17 : C produits pour les éléments de type cônes : triangles. 18 : C 19 : C REFERENCES : Numerical recipes (sous-programme gaujac modifié) 20 : C on a recalculé les poids et points de Gauss en REAL*16 21 : C donc avec environ 32 (plutôt 31) chiffres significatifs 22 : C 23 : C LANGAGE : ESOPE 24 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) 25 : C mél : gounand@semt2.smts.cea.fr 26 : C*********************************************************************** 27 : C APPELES : INIPG, GCSINO 28 : C APPELE PAR : INPGS 29 : C*********************************************************************** 30 : C ENTREES : - 31 : C ENTREES/SORTIES : MYPGS (actif en *MOD) 32 : C SORTIES : - 33 : C CODE RETOUR (IRET) : = 0 si tout s'est bien passé 34 : C*********************************************************************** 35 : C VERSION : v1, 31/05/00, version initiale 36 : C HISTORIQUE : v1, 31/05/00, création 37 : C HISTORIQUE : 38 : C HISTORIQUE : 39 : C*********************************************************************** 40 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 41 : C en cas de modification de ce sous-programme afin de faciliter 42 : C la maintenance ! 43 : C*********************************************************************** 44 : -INC CCOPTIO 45 : CBEGININCLUDE SPOGAU 46 : SEGMENT POGAU 47 : CHARACTER*(LNNPG) NOMPG 48 : CHARACTER*(LNTPG) TYPMPG 49 : CHARACTER*(LNFPG) FORLPG 50 : INTEGER NORDPG 51 : REAL*8 XCOPG(NDLPG,NBPG) 52 : REAL*8 XPOPG(NBPG) 53 : ENDSEGMENT 54 : SEGMENT POGAUS 55 : POINTEUR LISPG(0).POGAU 56 : ENDSEGMENT 57 : CENDINCLUDE SPOGAU 58 : POINTEUR MYPGS.POGAUS 59 : POINTEUR PGCOUR.POGAU 60 : * 61 : INTEGER IMPR,IRET 62 : * 63 : INTEGER DIMSRF 64 : PARAMETER(DIMSRF=1) 65 : REAL*8 XCOR(DIMSRF) 66 : * 67 : * Générateurs pour la méthode de Gauss-Jacobi de degré 1 : GJ10-1-1 : 68 : * 69 : REAL*8 X1D1,P1D1 70 : PARAMETER (X1D1=1.D0/3.D0) 71 : PARAMETER (P1D1=0.5D0) 72 : * 73 : * Générateurs pour la méthode de Gauss-Jacobi de degré 3 : GJ10-3-2 : 74 : * 75 : REAL*8 X1D3,P1D3,X2D3,P2D3 76 : PARAMETER (X1D3=0.15505102572168219018027159252941D0) 77 : PARAMETER (P1D3=0.31804138174397716939436900207515D0) 78 : PARAMETER (X2D3=0.64494897427831780981972840747060D0) 79 : PARAMETER (P2D3=0.18195861825602283060563099792484D0) 80 : * 81 : * Générateurs pour la méthode de Gauss-Jacobi de degré 5 : GJ10-5-3 : 82 : * 83 : REAL*8 X1D5,P1D5,X2D5,P2D5,X3D5,P3D5 84 : PARAMETER (X1D5=0.88587959512703947395546143769455D-1) 85 : PARAMETER (P1D5=0.20093191373895963077219813326462D0) 86 : PARAMETER (X2D5=0.40946686444073471086492625206883D0) 87 : PARAMETER (P2D5=0.22924110635958624669392059455632D0) 88 : PARAMETER (X3D5=0.78765946176084705602524188987600D0) 89 : PARAMETER (P3D5=0.69826979901454122533881272179078D-1) 90 : * 91 : * Générateurs pour la méthode de Gauss-Jacobi de degré 7 : GJ10-7-4 : 92 : * 93 : REAL*8 X1D7,P1D7,X2D7,P2D7,X3D7,P3D7,X4D7,P4D7 94 : PARAMETER (X1D7=0.57104196114517682193121192554116D-1) 95 : PARAMETER (P1D7=0.13550691343148811620826417407793D0) 96 : PARAMETER (X2D7=0.27684301363812382768004599768562D0) 97 : PARAMETER (P2D7=0.20346456801027136079140447593585D0) 98 : PARAMETER (X3D7=0.58359043236891682005669766866292D0) 99 : PARAMETER (P3D7=0.12984754760823244082645620288963D0) 100 : PARAMETER (X4D7=0.86024013565621944784791291887512D0) 101 : PARAMETER (P4D7=0.31180970950008082173875147096569D-1) 102 : * 103 : * Générateurs pour la méthode de Gauss-Jacobi de degré 9 : GJ10-9-5 : 104 : * 105 : REAL*8 X1D9,P1D9,X2D9,P2D9,X3D9,P3D9,X4D9,P4D9,X5D9,P5D9 106 : PARAMETER (X1D9=0.39809857051468742340806690093331D-1) 107 : PARAMETER (P1D9=0.96781590226651679274360971636169D-1) 108 : PARAMETER (X2D9=0.19801341787360817253579213679530D0) 109 : PARAMETER (P2D9=0.16717463809436956549167562309770D0) 110 : PARAMETER (X3D9=0.43797481024738614400501252000523D0) 111 : PARAMETER (P3D9=0.14638698708466980869803786935596D0) 112 : PARAMETER (X4D9=0.69546427335363609451461482372117D0) 113 : PARAMETER (P4D9=0.73908870072616670350633219341705D-1) 114 : PARAMETER (X5D9=0.90146491420117357387650110211225D0) 115 : PARAMETER (P5D9=0.15747914521692276185292316568490D-1) 116 : * 117 : * Générateurs pour la méthode de Gauss-Jacobi de degré 11 : GJ10-11-6 : 118 : * 119 : REAL*8 X1D11,P1D11,X2D11,P2D11,X3D11,P3D11 120 : REAL*8 X4D11,P4D11,X5D11,P5D11,X6D11,P6D11 121 : PARAMETER (X1D11=0.29316427159784891972050276913165D-1) 122 : PARAMETER (P1D11=0.72310330725508683655454326124839D-1) 123 : PARAMETER (X2D11=0.14807859966848429184997685249598D0) 124 : PARAMETER (P2D11=0.13554249723151861684069039663804D0) 125 : PARAMETER (X3D11=0.33698469028115429909705297208078D0) 126 : PARAMETER (P3D11=0.14079255378819892811907683907092D0) 127 : PARAMETER (X4D11=0.55867151877155013208139334180552D0) 128 : PARAMETER (P4D11=0.98661150890655264120584510548346D-1) 129 : PARAMETER (X5D11=0.76923386203005450091688336011564D0) 130 : PARAMETER (P5D11=0.43955165550508975508176624305422D-1) 131 : PARAMETER (X6D11=0.92694567131974111485187396581969D0) 132 : PARAMETER (P6D11=0.87383018136095317560173033123964D-2) 133 : * 134 : INTEGER NOPG 135 : * 136 : * Executable statements 137 : * 138 : IF (IMPR.GT.6) WRITE(IOIMP,*) 'Entrée dans ingj10' 139 : * 140 : * Méthode de nom : GJ10-1-1 141 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 142 : * d'ordre 1 à 1 point 143 : * espace de référence de dimension 1 144 : * 145 : * In INIPG : SEGINI PGCOUR 145 : CALL INIPG('GJ10-1-1','GAUSS-JACOBI10','SEGMENT', 146 : $ 1,1,1, 147 : $ PGCOUR, 148 : $ IMPR,IRET) 149 : IF (IRET.NE.0) GOTO 9999 150 : NOPG=0 151 : XCOR(1)=X1D1 151 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D1,IMPR,IRET) 152 : IF (IRET.NE.0) GOTO 9999 153 : SEGDES PGCOUR 154 : MYPGS.LISPG(**)=PGCOUR 155 : * 156 : * Méthode de nom : GJ10-3-2 157 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 158 : * d'ordre 3 à 2 points 159 : * espace de référence de dimension 1 160 : * 161 : * In INIPG : SEGINI PGCOUR 161 : CALL INIPG('GJ10-3-2','GAUSS-JACOBI10','SEGMENT', 162 : $ 3,2,1, 163 : $ PGCOUR, 164 : $ IMPR,IRET) 165 : IF (IRET.NE.0) GOTO 9999 166 : NOPG=0 167 : XCOR(1)=X1D3 167 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D3,IMPR,IRET) 168 : IF (IRET.NE.0) GOTO 9999 169 : XCOR(1)=X2D3 169 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P2D3,IMPR,IRET) 170 : IF (IRET.NE.0) GOTO 9999 171 : SEGDES PGCOUR 172 : MYPGS.LISPG(**)=PGCOUR 173 : * 174 : * Méthode de nom : GJ10-5-3 175 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 176 : * d'ordre 5 à 3 points 177 : * espace de référence de dimension 1 178 : * 179 : * In INIPG : SEGINI PGCOUR 179 : CALL INIPG('GJ10-5-3','GAUSS-JACOBI10','SEGMENT', 180 : $ 5,3,1, 181 : $ PGCOUR, 182 : $ IMPR,IRET) 183 : IF (IRET.NE.0) GOTO 9999 184 : NOPG=0 185 : XCOR(1)=X1D5 185 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D5,IMPR,IRET) 186 : IF (IRET.NE.0) GOTO 9999 187 : XCOR(1)=X2D5 187 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P2D5,IMPR,IRET) 188 : IF (IRET.NE.0) GOTO 9999 189 : XCOR(1)=X3D5 189 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P3D5,IMPR,IRET) 190 : IF (IRET.NE.0) GOTO 9999 191 : SEGDES PGCOUR 192 : MYPGS.LISPG(**)=PGCOUR 193 : * 194 : * Méthode de nom : GJ10-7-4 195 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 196 : * d'ordre 7 à 4 points 197 : * espace de référence de dimension 1 198 : * 199 : * In INIPG : SEGINI PGCOUR 199 : CALL INIPG('GJ10-7-4','GAUSS-JACOBI10','SEGMENT', 200 : $ 7,4,1, 201 : $ PGCOUR, 202 : $ IMPR,IRET) 203 : IF (IRET.NE.0) GOTO 9999 204 : NOPG=0 205 : XCOR(1)=X1D7 205 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D7,IMPR,IRET) 206 : IF (IRET.NE.0) GOTO 9999 207 : XCOR(1)=X2D7 207 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P2D7,IMPR,IRET) 208 : IF (IRET.NE.0) GOTO 9999 209 : XCOR(1)=X3D7 209 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P3D7,IMPR,IRET) 210 : IF (IRET.NE.0) GOTO 9999 211 : XCOR(1)=X4D7 211 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P4D7,IMPR,IRET) 212 : IF (IRET.NE.0) GOTO 9999 213 : SEGDES PGCOUR 214 : MYPGS.LISPG(**)=PGCOUR 215 : * 216 : * Méthode de nom : GJ10-9-5 217 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 218 : * d'ordre 9 à 5 points 219 : * espace de référence de dimension 1 220 : * 221 : * In INIPG : SEGINI PGCOUR 221 : CALL INIPG('GJ10-9-5','GAUSS-JACOBI10','SEGMENT', 222 : $ 9,5,1, 223 : $ PGCOUR, 224 : $ IMPR,IRET) 225 : IF (IRET.NE.0) GOTO 9999 226 : NOPG=0 227 : XCOR(1)=X1D9 227 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D9,IMPR,IRET) 228 : IF (IRET.NE.0) GOTO 9999 229 : XCOR(1)=X2D9 229 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P2D9,IMPR,IRET) 230 : IF (IRET.NE.0) GOTO 9999 231 : XCOR(1)=X3D9 231 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P3D9,IMPR,IRET) 232 : IF (IRET.NE.0) GOTO 9999 233 : XCOR(1)=X4D9 233 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P4D9,IMPR,IRET) 234 : IF (IRET.NE.0) GOTO 9999 235 : XCOR(1)=X5D9 235 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P5D9,IMPR,IRET) 236 : IF (IRET.NE.0) GOTO 9999 237 : SEGDES PGCOUR 238 : MYPGS.LISPG(**)=PGCOUR 239 : * 240 : * Méthode de nom : GJ10-11-6 241 : * Sur un segment : méthode de Gauss-Jacobi \alpha=1 \beta=0 242 : * d'ordre 11 à 6 points 243 : * espace de référence de dimension 1 244 : * 245 : * In INIPG : SEGINI PGCOUR 245 : CALL INIPG('GJ10-11-6','GAUSS-JACOBI10','SEGMENT', 246 : $ 11,6,1, 247 : $ PGCOUR, 248 : $ IMPR,IRET) 249 : IF (IRET.NE.0) GOTO 9999 250 : NOPG=0 251 : XCOR(1)=X1D11 251 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P1D11,IMPR,IRET) 252 : IF (IRET.NE.0) GOTO 9999 253 : XCOR(1)=X2D11 253 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P2D11,IMPR,IRET) 254 : IF (IRET.NE.0) GOTO 9999 255 : XCOR(1)=X3D11 255 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P3D11,IMPR,IRET) 256 : IF (IRET.NE.0) GOTO 9999 257 : XCOR(1)=X4D11 257 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P4D11,IMPR,IRET) 258 : IF (IRET.NE.0) GOTO 9999 259 : XCOR(1)=X5D11 259 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P5D11,IMPR,IRET) 260 : IF (IRET.NE.0) GOTO 9999 261 : XCOR(1)=X6D11 261 : CALL GCSINO(PGCOUR,NOPG,DIMSRF,XCOR,P6D11,IMPR,IRET) 262 : IF (IRET.NE.0) GOTO 9999 263 : SEGDES PGCOUR 264 : MYPGS.LISPG(**)=PGCOUR 265 : * 266 : * Normal termination 267 : * 268 : IRET=0 269 : RETURN 270 : * 271 : * Format handling 272 : * 273 : * 274 : * Error handling 275 : * 276 : 9999 CONTINUE 277 : IRET=1 278 : WRITE(IOIMP,*) 'An error was detected in subroutine ingj10' 279 : RETURN 280 : * 281 : * End of subroutine INGJ10 282 : * 283 : END 284 : 285 : 286 :
© Cast3M 2003 - All rights reserved.
Disclaimer