Télécharger la plaquette Cast3M


fr en

Télécharger ccgnor.eso

Retour à la liste

Numérotation des lignes : 

ccgnor

  1. C CCGNOR    SOURCE    PV        22/04/22    21:15:03     11344          
  2.       SUBROUTINE CCGNOR(LCOF,NOMLOI,
  3.      $     FC,
  4.      $     IMPR,IRET)
  5.       IMPLICIT REAL*8 (A-H,O-Z)
  6.       IMPLICIT INTEGER (I-N)
  7. C***********************************************************************
  8. C NOM         : CCGNOR
  9. C DESCRIPTION : Calcul des composantes d'un vecteur normal
  10. C
  11. C
  12. C LANGAGE     : ESOPE
  13. C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  14. C               mél : gounand@semt2.smts.cea.fr
  15. C***********************************************************************
  16. C APPELES          :
  17. C APPELE PAR       :
  18. C***********************************************************************
  19. C ENTREES            :
  20. C ENTREES/SORTIES    :
  21. C SORTIES            : -
  22. C TRAVAIL            :
  23. C***********************************************************************
  24. C VERSION    : v1, 09/03/07, version initiale
  25. C HISTORIQUE : v1, 09/03/07, création
  26. C HISTORIQUE :
  27. C HISTORIQUE :
  28. C***********************************************************************
  29. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  30. C en cas de modification de ce sous-programme afin de faciliter
  31. C la maintenance !
  32. C***********************************************************************
  33.  
  34. -INC PPARAM
  35. -INC CCOPTIO
  36. -INC CCREEL
  37. -INC TNLIN      
  38. *-INC SMCHAEL
  39.       INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM,N1
  40.       POINTEUR FC.MCHEVA
  41.       POINTEUR LCOF.LCHEVA
  42.       POINTEUR JMAJAC.MCHEVA
  43.       POINTEUR JMIJAC.MCHEVA
  44.       POINTEUR JDTJAC.MCHEVA
  45.       CHARACTER*8 NOMLOI
  46.       INTEGER ICOF
  47. *
  48. -INC TMXMAT
  49. * Objets temporaires
  50.       POINTEUR JAC.MXMAT,JT.MXMAT,JP.MXMAT
  51.       POINTEUR G.MXMAT,IG.MXMAT
  52. *
  53.       SEGMENT MIMAT2
  54.       INTEGER IMAT2(2,2)
  55.       ENDSEGMENT
  56.       POINTEUR EIJ.MIMAT2
  57.       SEGMENT MIMAT3
  58.       INTEGER IMAT3(3,3,3)
  59.       ENDSEGMENT
  60.       POINTEUR EIJK.MIMAT3
  61. *
  62.       LOGICAL LBID
  63.       INTEGER LAXSP
  64. *
  65.       INTEGER IMPR,IRET
  66.       REAL*8 DETG(1)
  67. *
  68. * Executable statements
  69. *
  70. *      IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans ccgnor'
  71. *      WRITE(IOIMP,*) 'Entrée dans ccgnor'
  72. *      WRITE(IOIMP,*) 'NOMLOI=',NOMLOI
  73. C      IF (.NOT.(IDIM.EQ.1)) THEN
  74. C         WRITE(IOIMP,*) 'IDIM=',IDIM,' ?'
  75. C         GOTO 9999
  76. C      ENDIF
  77.       NLFC=FC.WELCHE(/6)
  78.       NPFC=FC.WELCHE(/5)
  79.       ICOF=0
  80. *
  81.       ICOF=ICOF+1
  82.       JMAJAC=LCOF.LISCHE(ICOF)
  83.       NLJA=JMAJAC.WELCHE(/6)
  84.       NPJA=JMAJAC.WELCHE(/5)
  85.       IREF=JMAJAC.WELCHE(/4)
  86.       IREL=JMAJAC.WELCHE(/3)
  87. *      SEGPRT,JMAJAC
  88. *
  89. *      WRITE(IOIMP,*) 'IREL=',IREL
  90. *      WRITE(IOIMP,*) 'IREF=',IREF
  91.       IF (IREL.NE.IDIM) THEN
  92.          WRITE(IOIMP,*) 'Erreur dims JMAJAC'
  93.          GOTO 9999
  94.       ENDIF
  95.       IF (IREL.NE.IREF+1) THEN
  96.          WRITE(IOIMP,*) 'Le maillage donne nest pas une surface'
  97.          GOTO 9999
  98.       ENDIF
  99.       IF ((IREL.NE.2).AND.(IREL.NE.3)) THEN
  100.          WRITE(IOIMP,*) 'Ne marche quen dimension despace 2 ou 3'
  101.          GOTO 9999
  102.       ENDIF
  103. *
  104. * Objets temporaires et à préconditionner
  105. *
  106.       LDIM1=IREL
  107.       LDIM2=IREF
  108.       SEGINI,JAC
  109.       LDIM1=IREF
  110.       LDIM2=IREL
  111.       SEGINI,JT
  112. *      SEGINI,JP
  113.       LDIM1=IREF
  114.       LDIM2=IREF
  115.       SEGINI,G
  116.       SEGINI,IG
  117. *
  118. * Initialisation des tenseurs de permutation
  119. *
  120.       SEGINI,EIJ
  121.       IMULT=1
  122.       ICPT=0
  123.       DO I=1,2
  124.          DO J=1,2
  125.             IF (I.NE.J) THEN
  126.                ICPT=ICPT+1
  127.                IF (ICPT.EQ.2) THEN
  128.                   ICPT=0
  129.                   IMULT=IMULT*(-1)
  130.                ENDIF
  131.                EIJ.IMAT2(I,J)=IMULT
  132.             ENDIF
  133.          ENDDO
  134.       ENDDO
  135.       SEGINI,EIJK
  136.       IMULT=1
  137.       ICPT=0
  138.       DO I=1,3
  139.          DO J=1,3
  140.             IF (I.NE.J) THEN
  141.                DO K=1,3
  142.                   IF ((K.NE.I).AND.(K.NE.J)) THEN
  143.                      ICPT=ICPT+1
  144.                      IF (ICPT.EQ.2) THEN
  145.                         ICPT=0
  146.                         IMULT=IMULT*(-1)
  147.                      ENDIF
  148.                      EIJK.IMAT3(I,J,K)=IMULT
  149.                   ENDIF
  150.                ENDDO
  151.             ENDIF
  152.          ENDDO
  153.       ENDDO
  154. *      SEGPRT,EIJ
  155. *      SEGPRT,EIJK
  156. *
  157.       DO ILFC=1,NLFC
  158.          IF (NLJA.EQ.1) THEN
  159.             ILJA=1
  160.          ELSE
  161.             ILJA=ILFC
  162.          ENDIF
  163.          DO IPFC=1,NPFC
  164.             IF (NPJA.EQ.1) THEN
  165.                IPJA=1
  166.             ELSE
  167.                IPJA=IPFC
  168.             ENDIF
  169. *
  170. * Copie du jacobien
  171. *
  172.             CALL MAMA(JMAJAC.WELCHE(1,1,1,1,IPJA,ILJA),IREL,IREF,
  173.      $           'COPIE   ',
  174.      $           JAC.XMAT,IREL,IREF,
  175.      $           IMPR,IRET)
  176.             IF (IRET.NE.0) GOTO 9999
  177. *     Calcul de Jt
  178.             CALL MAMA(JAC.XMAT,IREL,IREF,
  179.      $           'TRANSPOS',JT.XMAT,IREF,IREL,
  180.      $           IMPR,IRET)
  181.             IF (IRET.NE.0) GOTO 9999
  182. *     Calcul de G=JtJ
  183.             CALL MAMAMA(JT.XMAT,IREF,IREL,JAC.XMAT,IREL,IREF,
  184.      $           'FOIS    ',G.XMAT,IREF,IREF,
  185.      $           IMPR,IRET)
  186.             IF (IRET.NE.0) GOTO 9999
  187. * Calcul de l'inverse, du déterminant et trace de l'inverse de g
  188.             LBID=.FALSE.
  189.             CALL GEOLI2(IREF,1,1,G.XMAT,IG.XMAT,DETG,LBID,IMPR,IRET)
  190.             IF (IRET.NE.0) GOTO 9999
  191.             XM=SQRT(DETG(1))
  192. **
  193. ** Calcul de la pseudo-inverse J+ = g-1 Jt
  194. **
  195. *               CALL MAMAMA(IG.XMAT,IREF,IREF,JT.XMAT,IREF,IREL,
  196. *     $              'FOIS    ',JP.XMAT,IREF,IREL,
  197. *     $              IMPR,IRET)
  198. *               IF (IRET.NE.0) GOTO 9999
  199.  
  200.             IF (NOMLOI(1:4).EQ.'VNOR') THEN
  201.                CALL CH2INT(NOMLOI(5:5),I,IMPR,IRET)
  202.                IF (IRET.NE.0) GOTO 9999
  203. **
  204. ** Calcul de la pseudo-inverse J+ = g-1 Jt
  205. **
  206. *               CALL MAMAMA(IG.XMAT,IREF,IREF,JT.XMAT,IREF,IREL,
  207. *     $              'FOIS    ',JP.XMAT,IREF,IREL,
  208. *     $              IMPR,IRET)
  209. *               IF (IRET.NE.0) GOTO 9999
  210.  
  211. *
  212.                IF (IREL.EQ.2) THEN
  213.                   CONTRI=0.D0
  214.                   DO J=1,2
  215.                      CONTRI=CONTRI-(EIJ.IMAT2(I,J)*
  216.      $                    JMAJAC.WELCHE(1,1,J,1,IPJA,ILJA))
  217.                   ENDDO
  218.                   CONTRI=CONTRI/XM
  219.                ELSEIF (IREL.EQ.3) THEN
  220. C                  XNUM=XM
  221. C                  XDENO=0.D0
  222. C* C'est louche parce que II ne varie pas de 1 à IREL
  223. C                  DO II=1,IREF
  224. C                     DO IJ=1,IREF
  225. C*                        DO IN=1,IREF
  226. C*                           DO IO=1,IREF
  227. C                              XDENO=XDENO+
  228. C     $                             (EIJ.IMAT2(II,IJ)*
  229. C     $                       JAC.XMAT(II,1)*JAC.XMAT(IJ,2))**2
  230. CC     $                             (EIJ.IMAT2(II,IJ)*EIJ.IMAT2(IN,IO)*
  231. CC     $                             JAC.XMAT(II,1)*JAC.XMAT(IJ,2)*
  232. CC     $                             JAC.XMAT(IN,1)*JAC.XMAT(IO,2))
  233. CC                           ENDDO
  234. CC                        ENDDO
  235. C                     ENDDO
  236. C                  ENDDO
  237. C                  XAL=XNUM/XDENO
  238. C                  CONTRI=0.D0
  239. C                  DO J=1,IREL
  240. C                     DO K=1,IREL
  241. C                        CONTRI=CONTRI+
  242. C     $                       (EIJK.IMAT3(I,J,K)*
  243. C     $                       JAC.XMAT(J,1)*JAC.XMAT(K,2))
  244. C                     ENDDO
  245. C                  ENDDO
  246. C                  CONTRI=CONTRI*XAL
  247.                   CONTRI=0.D0
  248.                   DO IA=1,IREF
  249.                      DO IB=1,IREF
  250.                         DO J=1,IREL
  251.                            DO K=1,IREL
  252.                               CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)*
  253.      $                             EIJK.IMAT3(I,J,K)*
  254.      $                             JMAJAC.WELCHE(1,1,J,IA,IPJA,ILJA)*
  255.      $                             JMAJAC.WELCHE(1,1,K,IB,IPJA,ILJA))
  256. C     $                             JP.XMAT(IA,J)*JP.XMAT(IB,K))
  257.                            ENDDO
  258.                         ENDDO
  259.                      ENDDO
  260.                   ENDDO
  261. *                  SEGPRT,JAC
  262. *                  SEGPRT,JP
  263. *                  WRITE(IOIMP,*) 'XM=',XM
  264.                   CONTRI=CONTRI/(2.D0*XM)
  265. *                  WRITE(IOIMP,*) 'CONTRI=',CONTRI
  266.                ELSE
  267.                   WRITE(IOIMP,*) 'Erreur grave IREL=',IREL
  268.                   GOTO 9999
  269.                ENDIF
  270.             ELSEIF (NOMLOI(1:4).EQ.'VNOJ') THEN
  271.                CALL CH2INT(NOMLOI(5:5),I,IMPR,IRET)
  272.                IF (IRET.NE.0) GOTO 9999
  273.                CALL CH2INT(NOMLOI(6:6),J,IMPR,IRET)
  274.                IF (IRET.NE.0) GOTO 9999
  275.                CALL CH2INT(NOMLOI(7:7),L,IMPR,IRET)
  276.                IF (IRET.NE.0) GOTO 9999
  277.                IF (IREL.EQ.2) THEN
  278.                   CONTRI=0.D0
  279.                   CONTRI=CONTRI-(EIJ.IMAT2(I,J)*
  280.      $                 JMAJAC.WELCHE(1,1,L,1,IPJA,ILJA))
  281.                   CONTRI=CONTRI/XM
  282.                ELSEIF (IREL.EQ.3) THEN
  283.                   CONTRI=0.D0
  284.                   DO IA=1,IREF
  285.                      DO IB=1,IREF
  286. *                        DO J=1,IREL
  287.                            DO K=1,IREL
  288.                               CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)*
  289.      $                             EIJK.IMAT3(I,J,K)*
  290.      $                             JMAJAC.WELCHE(1,1,L,IA,IPJA,ILJA)*
  291.      $                             JMAJAC.WELCHE(1,1,K,IB,IPJA,ILJA))
  292.                               CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)*
  293.      $                             EIJK.IMAT3(I,K,J)*
  294.      $                             JMAJAC.WELCHE(1,1,K,IA,IPJA,ILJA)*
  295.      $                             JMAJAC.WELCHE(1,1,L,IB,IPJA,ILJA))
  296.                            ENDDO
  297. *                        ENDDO
  298.                      ENDDO
  299.                   ENDDO
  300. *                  SEGPRT,JAC
  301. *                  SEGPRT,JP
  302. *                  WRITE(IOIMP,*) 'XM=',XM
  303.                   CONTRI=CONTRI/(2.D0*XM)
  304. *                  WRITE(IOIMP,*) 'CONTRI=',CONTRI
  305.                ELSE
  306.                   WRITE(IOIMP,*) 'Erreur grave IREL=',IREL
  307.                   GOTO 9999
  308.                ENDIF
  309.             ELSE
  310.                WRITE(IOIMP,*) 'Erreur grave NOMLOI=',NOMLOI
  311.                GOTO 9999
  312.             ENDIF
  313.             FC.WELCHE(1,1,1,1,IPFC,ILFC)=
  314.      $           FC.WELCHE(1,1,1,1,IPFC,ILFC)+
  315.      $           CONTRI
  316.          ENDDO
  317.       ENDDO
  318.       SEGSUP,EIJK
  319.       SEGSUP,EIJ
  320. *      SEGPRT,FC
  321.       SEGSUP,JAC
  322. *      SEGSUP,JP
  323.       SEGSUP,JT
  324.       SEGSUP,G
  325.       SEGSUP,IG
  326. *
  327. * Normal termination
  328. *
  329.       IRET=0
  330.       RETURN
  331. *
  332. * Format handling
  333. *
  334. *
  335. * Error handling
  336. *
  337.  9999 CONTINUE
  338.       IRET=1
  339.       WRITE(IOIMP,*) 'An error was detected in subroutine ccgnor'
  340.       RETURN
  341. *
  342. * End of subroutine CCGNOR
  343. *
  344.       END
  345.  
  346.  
  347.  
  348.  
  349.  
  350.  
  351.  
  352.  
  353.  
  354.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales