Download nif.eso

Back to the list

   1 : C NIF       SOURCE    GOUNAND   06/08/04    21:17:10     5520           
   2 :       SUBROUTINE NIF(LRFVOL,JXCOPG,PNM1,
   3 :      $     FNPG,DFNPG,
   4 :      $     IMPR,IRET)
   5 :       IMPLICIT REAL*8 (A-H,O-Z)
   6 :       IMPLICIT INTEGER (I-N)
   7 : C***********************************************************************
   8 : C NOM         : NIF
   9 : C PROJET      : Noyau linéaire NLIN
  10 : C DESCRIPTION : Calcul des valeurs des fonctions de forme et des
  11 : C               valeurs de leurs dérivées premières aux points de Gauss.
  12 : C     
  13 : C               On utilise la relation (produit) suivante :
  14 : C     
  15 : C               < Ni (point) > = < Pi (point) > [Pn]^{-1}
  16 : C     
  17 : C               avec (cf. CALPN) :
  18 : C               [PN] = ( P1(ksi1)  .....   Pn(ksi1))
  19 : C                      (   ...     .....     ...   )
  20 : C                      ( P1(ksin)  .....   Pn(ksin))
  21 : C               n = nb. ddl sur l'élément (NDFN)
  22 : C               ksii = coords. du ieme noeud d'approximation
  23 : C                      dans l'espace de référence (de dimension
  24 : C                      NDIML)
  25 : C               Pi   = ieme polynome d'interpolation sur
  26 : C                      l'élément de référence.
  27 : C               Ni   = ieme fonction nodale d'interpolation sur
  28 : C                      l'élément de référence.
  29 : C               point= point quelconque sur l'élément de référence
  30 : C                      (donc en particulier les points de Gauss)
  31 : C LANGAGE     : Esope
  32 : C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  33 : C               mél : gounand@semt2.smts.cea.fr
  34 : C***********************************************************************
  35 : C APPELES          : VALPOL, NOMINC
  36 : C APPELE PAR       : KFNREF
  37 : C***********************************************************************
  38 : C ENTREES            : LRFVOL, MYPG, PNM1
  39 : C ENTREES/SORTIES    : -
  40 : C SORTIES            : FNPG, DFNPG
  41 : C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  42 : C***********************************************************************
  43 : C VERSION    : v1, 16/09/99, version initiale
  44 : C HISTORIQUE : v1, 16/09/99, création
  45 : C HISTORIQUE : v2, 10/05/00, modif. du segment ELREF
  46 : C HISTORIQUE :
  47 : C HISTORIQUE :
  48 : C***********************************************************************
  49 : C Prière de PRENDRE LE TEMPS de compléter les commentaires
  50 : C en cas de modification de ce sous-programme afin de faciliter
  51 : C la maintenance !
  52 : C***********************************************************************
  53 : -INC CCOPTIO
  54 : CBEGININCLUDE SELREF
  55 :       SEGMENT ELREF
  56 :       CHARACTER*(LNNOM)  NOMLRF
  57 :       CHARACTER*(LNFORM) FORME 
  58 :       CHARACTER*(LNTYPL) TYPEL 
  59 :       CHARACTER*(LNESP)  ESPACE
  60 :       INTEGER DEGRE
  61 :       REAL*8 XCONOD(NDIMEL,NBNOD)
  62 :       INTEGER NPQUAF(NBDDL)
  63 :       INTEGER NUMCMP(NBDDL)
  64 :       INTEGER QUENOD(NBDDL)
  65 :       INTEGER ORDDER(NDIMEL,NBDDL)
  66 :       POINTEUR MBPOLY.POLYNS
  67 :       ENDSEGMENT
  68 :       SEGMENT ELREFS
  69 :       POINTEUR LISEL(0).ELREF
  70 :       ENDSEGMENT
  71 : CENDINCLUDE SELREF
  72 :       POINTEUR LRFVOL.ELREF
  73 : CBEGININCLUDE SPOLYNO
  74 :       SEGMENT POLYNO
  75 :       REAL*8  COEMON(NBMON)
  76 :       INTEGER EXPMON(NDIML,NBMON)
  77 :       ENDSEGMENT 
  78 :       SEGMENT POLYNS
  79 :       POINTEUR LIPOLY(NBPOLY).POLYNO
  80 :       ENDSEGMENT
  81 : CENDINCLUDE SPOLYNO
  82 :       POINTEUR MYBPOL.POLYNS
  83 :       POINTEUR MYPOLY.POLYNO
  84 : CBEGININCLUDE SMMREEL
  85 :       SEGMENT MMREEL
  86 :       REAL*8 MAT(JLIG,JCOL)
  87 :       ENDSEGMENT
  88 : CENDINCLUDE SMMREEL
  89 :       POINTEUR PNM1.MMREEL
  90 : -INC SMLENTI
  91 :       INTEGER JG
  92 :       POINTEUR ORDER1.MLENTI
  93 : -INC SMLREEL
  94 :       POINTEUR VECTPI.MLREEL
  95 : CBEGININCLUDE SMCHAEL
  96 :       SEGMENT MCHAEL
  97 :       POINTEUR IMACHE(N1).MELEME
  98 :       POINTEUR ICHEVA(N1).MCHEVA
  99 :       ENDSEGMENT
 100 :       SEGMENT MCHEVA
 101 :       REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM)
 102 :       ENDSEGMENT
 103 :       SEGMENT LCHEVA
 104 :       POINTEUR LISCHE(NBCHE).MCHEVA
 105 :       ENDSEGMENT
 106 : CENDINCLUDE SMCHAEL
 107 :       POINTEUR JXCOPG.MCHEVA
 108 :       INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM
 109 :       POINTEUR FNPG.MCHEVA
 110 :       POINTEUR DFNPG.MCHEVA
 111 : *     
 112 :       INTEGER IMPR,IRET
 113 : *     
 114 : * Fonction Blas (produit scalaire)
 115 : * 
 116 :       REAL*8 DDOT
 117 : *
 118 :       INTEGER NBMONO
 119 :       INTEGER NDIMQR,       NPGFAC,NBFN,NBELFV
 120 :       INTEGER IDIMQR,JDIMQR,IPGFAC,IBFN,IBELFV
 121 : *     
 122 : * Executable statements
 123 : *     
 124 :       
 125 :       IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans nif'
 126 : *     
 127 : * Initialisations
 128 : *     
 129 : *      SEGPRT,JXCOPG
 130 : *      SEGPRT,LRFVOL
 131 : *      SEGPRT,PNM1
 132 :       SEGACT PNM1
 133 :       SEGACT JXCOPG
 134 :       NDIMQR=JXCOPG.VELCHE(/4)
 135 :       NPGFAC=JXCOPG.VELCHE(/5)
 136 :       NBELFV=JXCOPG.VELCHE(/6)
 137 :       SEGACT LRFVOL
 138 :       MYBPOL=LRFVOL.MBPOLY
 139 :       SEGDES LRFVOL
 140 :       SEGACT MYBPOL
 141 :       SEGACT MYBPOL.LIPOLY(*)
 142 :       NBFN=MYBPOL.LIPOLY(/1)
 143 :       SEGACT PNM1
 144 :       JG=NDIMQR
 145 :       SEGINI ORDER1
 146 :       JG=NBFN
 147 :       SEGINI VECTPI
 148 : *     
 149 : * On calcule les valeurs des fonctions de forme aux points de Gauss
 150 : *     
 151 :       NBLIG=1
 152 :       NBCOL=NBFN
 153 :       N2LIG=1
 154 :       N2COL=1
 155 :       NBPOI=NPGFAC
 156 :       NBELM=NBELFV
 157 :       SEGINI FNPG
 158 :       DO 1 IDIMQR=1,NDIMQR
 159 :          ORDER1.LECT(IDIMQR)=0
 160 :  1    CONTINUE
 161 :       DO IBELFV=1,NBELFV
 162 :          DO 3 IPGFAC=1,NPGFAC
 163 : * Calcul de < P (pg) > = < P1(pg) ... Pnbfn(pg) > où pg est le
 164 : * IPGFACieme point de Gauss
 165 :             DO 32 IBFN=1,NBFN
 166 :                MYPOLY=MYBPOL.LIPOLY(IBFN)
 167 :                NBMONO=MYPOLY.EXPMON(/2)
 167 :                CALL VALPOL(NDIMQR,NBMONO,
 168 :      $              JXCOPG.VELCHE(1,1,1,1,IPGFAC,IBELFV),
 169 :      $              MYPOLY.COEMON,MYPOLY.EXPMON,
 170 :      $              ORDER1.LECT,
 171 :      $              VECTPI.PROG(IBFN),
 172 :      $              IMPR,IRET)
 173 :                IF (IRET.NE.0) GOTO 9999
 174 :  32         CONTINUE
 175 : * On calcule : < N (pg) > = < P (pg) > [Pn]^{-1}
 176 :             DO 34 IBFN=1,NBFN
 177 :                FNPG.VELCHE(1,IBFN,1,1,IPGFAC,IBELFV)=
 178 :      $              DDOT(NBFN,VECTPI.PROG,1,PNM1.MAT(1,IBFN),1)
 179 :  34         CONTINUE
 180 :  3       CONTINUE
 181 :          IF (IMPR.GT.3) THEN
 182 :             WRITE(IOIMP,*) 'Face n° ',IBELFV
 183 :             WRITE(IOIMP,*) 'Ordre de dérivation / coordonnée de réf. :'
 184 :             WRITE(IOIMP,4003) (ORDER1.LECT(IDIMQR),IDIMQR=1,NDIMQR)
 185 :             DO 5 IPGFAC=1,NPGFAC
 186 :                WRITE(IOIMP,*) 'Noeud de coordonnées :'
 187 :                WRITE(IOIMP,4004)
 188 :      $              (JXCOPG.VELCHE(1,1,1,IDIMQR,IPGFAC,IBELFV),IDIMQR=1
 189 :      $              ,NDIMQR) 
 190 :                WRITE(IOIMP,*) 'FNPG.VELCHE(nb.fns.forme) :'
 191 :                WRITE(IOIMP,4004)
 192 :      $              (FNPG.VELCHE(1,IBFN,1,1,IPGFAC,IBELFV),IBFN=1,NBFN
 193 :      $              )
 194 :  5          CONTINUE
 195 :          ENDIF
 196 :       ENDDO
 197 :       SEGDES FNPG
 198 : *     
 199 : * On calcule les valeurs des dérivées premières des fonctions
 200 : * de forme aux points de Gauss
 201 : *     
 202 :       NBLIG=1
 203 :       NBCOL=NBFN
 204 :       N2LIG=1
 205 :       N2COL=NDIMQR
 206 :       NBPOI=NPGFAC
 207 :       NBELM=NBELFV
 208 :       SEGINI DFNPG
 209 :       DO 7 IDIMQR=1,NDIMQR
 210 :          DO 72 JDIMQR=1,NDIMQR
 211 :             IF (JDIMQR.EQ.IDIMQR) THEN
 212 :                ORDER1.LECT(JDIMQR)=1
 213 :             ELSE
 214 :                ORDER1.LECT(JDIMQR)=0
 215 :             ENDIF
 216 :  72      CONTINUE
 217 :          DO IBELFV=1,NBELFV
 218 :             DO 74 IPGFAC=1,NPGFAC
 219 :                DO 742 IBFN=1,NBFN
 220 : * Calcul de < dP/dksi_idimqr (pg) > où pg est le
 221 : * IPGFACieme point de Gauss
 222 :                   MYPOLY=MYBPOL.LIPOLY(IBFN)
 223 :                   NBMONO=MYPOLY.EXPMON(/2)
 223 :                   CALL VALPOL(NDIMQR,NBMONO,
 224 :      $                 JXCOPG.VELCHE(1,1,1,1,IPGFAC,IBELFV),
 225 :      $                 MYPOLY.COEMON,MYPOLY.EXPMON,
 226 :      $                 ORDER1.LECT,
 227 :      $                 VECTPI.PROG(IBFN),
 228 :      $                 IMPR,IRET)
 229 :                   IF (IRET.NE.0) GOTO 9999
 230 :  742           CONTINUE
 231 : * On calcule : < dN/dksi_idimqr (pg) > = < dP/dksi_idimqr (pg) > [Pn]^{-1}
 232 :                DO 744 IBFN=1,NBFN
 233 :                   DFNPG.VELCHE(1,IBFN,1,IDIMQR,IPGFAC,IBELFV)=
 234 :      $                 DDOT(NBFN,VECTPI.PROG,1,PNM1.MAT(1,IBFN),1)
 235 :  744           CONTINUE
 236 :  74         CONTINUE
 237 :             IF (IMPR.GT.3) THEN
 238 :                WRITE(IOIMP,*)
 239 :      $              'Ordre de dérivation / coordonnée de réf. :'
 240 :                WRITE(IOIMP,4003) (ORDER1.LECT(JDIMQR),JDIMQR=1,NDIMQR)
 241 :                DO 76 IPGFAC=1,NPGFAC
 242 :                   WRITE(IOIMP,*) 'Face n° ',IBELFV
 243 :                   WRITE(IOIMP,*) 'Noeud de coordonnées :'
 244 :                   WRITE(IOIMP,4004) (JXCOPG.VELCHE(1,1,1,JDIMQR,IPGFAC
 245 :      $                 ,IBELFV),JDIMQR=1,NDIMQR)
 246 :                   WRITE(IOIMP,*) 'DFNPG.VELCHE(nb.fns.forme) :'
 247 :                   WRITE(IOIMP,4004) 
 248 :      $                 (DFNPG.VELCHE(1,IBFN,1,IDIMQR,IPGFAC,IBELFV),
 249 :      $                 IBFN=1,NBFN)
 250 :  76            CONTINUE
 251 :             ENDIF
 252 :          ENDDO
 253 :  7    CONTINUE
 254 :       SEGDES DFNPG
 255 :       SEGSUP VECTPI
 256 :       SEGSUP ORDER1
 257 :       SEGDES PNM1
 258 :       SEGDES MYBPOL.LIPOLY(*)
 259 :       SEGDES MYBPOL
 260 :       SEGDES JXCOPG
 261 : *     
 262 : * Normal termination
 263 : *     
 264 :       IRET=0
 265 :       RETURN
 266 : *     
 267 : * Format handling
 268 : *     
 269 :  4003 FORMAT (2X,6(1X,I3))
 270 :  4004 FORMAT (2X,6(1X,1PE13.5))
 271 : *     
 272 : * Error handling
 273 : *     
 274 :  9999 CONTINUE
 275 :       IRET=1
 276 :       WRITE(IOIMP,*) 'An error was detected in subroutine nif'
 277 :       RETURN
 278 : *     
 279 : * End of subroutine NIF
 280 : *     
 281 :       END
 282 : 
 283 :  
 284 :  

© Cast3M 2003 - All rights reserved.
Disclaimer