Télécharger la plaquette Cast3M


fr en

Télécharger assem0.eso

Retour à la liste

Numérotation des lignes : 

assem0

  1. C ASSEM0    SOURCE    CHAT      09/10/09    21:15:50     6519
  2.       SUBROUTINE ASSEM0(MRIGI1,ICLE,INWUIT)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8 (A-H,O-Z)
  5.  
  6.  
  7. -INC PPARAM
  8. -INC CCOPTIO
  9. -INC SMRIGID
  10. -INC SMLMOTS
  11. -INC SMLREEL
  12. -INC SMLENTI
  13.  
  14. ***      SAVE INWAIT
  15.       SEGMENT  INWAIT
  16.          INTEGER IIM(IRM)
  17.       ENDSEGMENT
  18.  
  19. *      write(6,fmt='(''  entree dans assem0 icle'',i6)') icle
  20.  
  21.       IF(ICLE.EQ.1) THEN
  22. C
  23. C   on est en premier passage il faut triturer les matrices
  24. C
  25.          MRIGID=MRIGI1
  26.          SEGACT MRIGID*MOD
  27.          INWAIT=0
  28.  
  29. *      write(6,fmt='(''  entree dans assem0 mrigi1 ichole  ''
  30. *     * ,2i6)')
  31. *     *  mrigid,ichole
  32.  
  33.          IF(ICHOLE.NE.0) THEN
  34.             SEGDES,MRIGID
  35.             RETURN
  36.          ENDIF
  37.          IRM=IRIGEL(/2)
  38.          SEGINI INWAIT
  39.          inWUIT=inwait
  40.          MLMOTS=NORINC
  41.  
  42. *      write( 6,fmt='('' assem0 norinc '', i6)') norinc
  43.  
  44. C    ... Si la normalisation est AUTOmatique ...
  45.          IF( MLMOTS.EQ. -1) THEN
  46. C       ... Taille initiale de la liste de coeff. ...
  47.             JG=20
  48. C       ... Taille initiale de la liste de nom de DDL ...
  49.             JGN=4
  50.             JGM=20
  51. C       ... JFIN = nombre de différents noms de DDL ...
  52.             JFIN=0
  53. C       ... Taille actuelle des MLMOTS et MLREEL ...
  54.             JMAX=20
  55.             SEGINI MLMOTS,MLREEL
  56.             CALL SAVSEG(MLMOTS)
  57.             CALL SAVSEG(MLREEL)
  58. C       ... Boucle sur les zones élémentaires ...
  59.             DO 50 I=1,IRM
  60.                DESCR=IRIGEL(3,I)
  61.                SEGACT DESCR
  62.                NLIGRP=NOELEP(/1)
  63.                NLIGRD=NLIGRP
  64. C          ... On met le premier DDL dans MLMOTS ...
  65.                IF(I.EQ.1) THEN
  66.                   JFIN=1
  67.                   MOTS(JFIN)=LISINC(1)
  68.                ENDIF
  69. C          ... Si les autres n'y sont pas encore, on les rajoute ...
  70.                DO 51 J=1,NLIGRP
  71.                   DO 52 K=1,JFIN
  72.                      IF( LISINC(J).EQ.MOTS(K)) GO TO 51
  73.    52             CONTINUE
  74.                   JFIN=JFIN+1
  75. C             ... En agrandissant (s'il le faut) les MLMOTS et MLREEL ...
  76.                   IF(JFIN.GT.JMAX)THEN
  77.                      JG = JMAX+20
  78.                      JGM=JG
  79.                      SEGADJ MLMOTS,MLREEL
  80.                   ENDIF
  81.                   MOTS(JFIN)=LISINC(J)
  82.    51          CONTINUE
  83. C          ... Initialisation de MLREE1 de taille NLIGRP ...
  84.                JG=NLIGRP
  85.                SEGINI MLREE1
  86.                XMATRI=IRIGEL(4,I)
  87.                SEGACT XMATRI
  88. C          ... Dans lequel on mettra les maxi des NLIGRP termes diagonaux
  89. C              de toutes le matrices élémentaires de la zone ...
  90.                DO 53 J=1,RE(/3)
  91. *                  XMATRI=IMATTT(J)
  92. *                  SEGACT XMATRI
  93.                   DO 54 K=1,NLIGRP
  94.                      IF(ABS(RE(K,K,J)).GT.MLREE1.PROG(K))
  95.      &                    MLREE1.PROG(K)=ABS(RE(K,K,J))
  96.    54             CONTINUE
  97.    53          CONTINUE
  98. C          ... Puis, pour chaque DDL différent, on met le maxi des
  99. C              termes diagonaux concernés dans MLREEL ...
  100.                DO 55 J=1,NLIGRP
  101.                   DO 56 K=1,JFIN
  102.                      IF( LISINC(J).EQ.MOTS(K)) GO TO 57
  103.    56             CONTINUE
  104.    57             CONTINUE
  105.                   IF(MLREE1.PROG(J).GT.PROG(K))
  106.      &                     PROG(K)=MLREE1.PROG(J)
  107.    55          CONTINUE
  108.                SEGSUP MLREE1
  109.    50       CONTINUE
  110. C       ... Toutes les zones étant parcourues, la taille des MLMOTS et MLREEL
  111. C           est définitive, on peut donc l'ajuster ...
  112.             JG=JFIN
  113.             JGM=JG
  114.             SEGADJ MLMOTS,MLREEL
  115. C       ... Les coefficients valent 0.8/sqrt(terme maxi) pour les DDL
  116. C           "physiques" et 1 pour les multiplicateurs de Lagrange ...
  117.             DO 59 K=1,JFIN
  118.                IF(PROG(K).EQ.0.D0) PROG(K)=1.D0
  119.                PROG(K)=0.8D0 / SQRT(PROG(K))
  120.                IF(MOTS(K).EQ.'LX') PROG(K)=1.D0
  121.    59       CONTINUE
  122. C       ... Il n'y a pas de normalisation des variables duales ...
  123.             MLMOT1=0
  124.             MLREE1=0
  125.             LIN=JFIN
  126.             NORINC=MLMOTS
  127.             NORVAL=MLREEL
  128.  
  129. *           write(6,fmt='('' norinc '',4( A4,2x))')(mots(kk),kk=1,
  130. *     *  mots(/2))
  131. *           write(6,fmt='('' norval '',4(e12.5,2x))')(prog(kk),kk=1,
  132. *     *  prog(/1))
  133.  
  134. C    ... Sinon (cad. la normalisation n'est pas automatique) ...
  135.          ELSE
  136.             SEGACT MLMOTS
  137.             LIN=MOTS(/2)
  138.             MLREEL=NORVAL
  139.             SEGACT MLREEL
  140.             MLMOT1=NORIND
  141.             MLREE1=NORVAD
  142.             IF(MLMOT1.NE.0)  THEN
  143.                SEGACT MLMOT1
  144.                SEGACT MLREE1
  145.                LIND=MLREE1.PROG(/1)
  146.             ENDIF
  147.          ENDIF
  148.  
  149. C ... La normalisation proprement dite commence ici ...
  150. C
  151. C   BOUCLE 1 SUR LES SOUS ZONES ELEMENTAIRES
  152. C
  153.          DO 1 I=1,IRM
  154.             DESCR=IRIGEL(3,I)
  155.             SEGACT DESCR
  156.             NLIGRP=NOELEP(/1)
  157.             NLIGRD=NLIGRP
  158.             JG=NLIGRP
  159.             MLREE2=0
  160.             MLREE3=0
  161. C
  162. C   existe-t-il des inconnues a normer dans la matrice
  163. C   si oui création de MLREE2 ET MLREE3 qui serviront de coef
  164. C   multiplicateurs
  165. C
  166.             DO 2 K=1,NLIGRP
  167. C          ... Pour toute variable primale, on cherche si elle
  168. C              est présente dans la liste des variables à normaliser ...
  169.                DO 3 L=1,LIN
  170.                   IF(LISINC(K).EQ.MOTS(L)) THEN
  171. C                ... Si c'est le cas, on vérifie si les listes des
  172. C                    coefficients sont initialisées ...
  173.                      IF(MLREE2.EQ.0)  THEN
  174.                         SEGINI MLREE2
  175.                         DO 33 M=1,NLIGRP
  176.                            MLREE2.PROG(M)=1.D0
  177.    33                   CONTINUE
  178.                         IF(MLMOT1.EQ.0) THEN
  179.                            SEGINI MLREE3
  180.                            DO 34 M=1,NLIGRP
  181.                               MLREE3.PROG(M)=1.D0
  182.    34                      CONTINUE
  183.                         ENDIF
  184.                      ENDIF
  185. C                ... Puis, on met le bon coefficient au bon endroit ...
  186.                      MLREE2.PROG(K)=PROG(L)
  187. C                ... Si la normalisation des variables duales n'a pas été
  188. C                    demandée, les coeff. seront les mêmes que pour les
  189. C                    variables primales ...
  190.                      IF(MLREE3.NE.0) MLREE3.PROG(K)= PROG(L)
  191.                   ENDIF
  192.     3          CONTINUE
  193.     2       CONTINUE
  194.  
  195. C       ... Si au contraire, la normalisation des variables duales a
  196. C           été demandée, on refait la même chose pour les variables
  197. C           duales ...
  198.             IF(MLMOT1.EQ.0) GO TO 6
  199.             NLIGRD=NOELED(/1)
  200.             JG=NLIGRD
  201.             MLREE3=0
  202.             DO 4 K=1,NLIGRD
  203.                DO 5 L=1,LIND
  204.                   IF(LISDUA(K).EQ.MLMOT1.MOTS(L)) THEN
  205.                      IF(MLREE3.EQ.0)  THEN
  206.                         SEGINI MLREE3
  207.                         DO 35 M=1,NLIGRD
  208.                            MLREE3.PROG(M)=1.D0
  209.    35                   CONTINUE
  210.                      ENDIF
  211. C                ... c.a.d. on met les coefficients au bons endroits dans MLREE3 ...
  212.                      MLREE3.PROG(K)=MLREE1.PROG(L)
  213.                   ENDIF
  214.     5          CONTINUE
  215.     4       CONTINUE
  216.     6       CONTINUE
  217. C
  218. C  si MLREE2*MLREE3 NE 0 il faut multiplier les matrices
  219. C  On va le faire et on cree un nouveau segment IMATRI
  220. C  l'ancien etant stocke dans inwait
  221. C
  222.             IF(MLREE2*MLREE3.EQ.0) GO TO 15
  223.             XMATRI=IRIGEL(4,I)
  224.             SEGACT XMATRI
  225.             NELRIG=RE(/3)
  226. C       ... On met le pointeur IMATRI dans INWAIT ...
  227.             IIM(I)=XMATRI
  228. C       ... On créé un nouveau IMATRI et on le met dans MRIGID ...
  229.             SEGINI,XMATR1=XMATRI
  230.             IRIGEL(4,I)=XMATR1
  231. C       ... Puis on parcourt les matrices élémentaires ...
  232.             DO 7 K=1,NELRIG
  233. *               XMATRI=IMATTT(K)
  234. C          ... Chaque nouvelle matrice est égale au début à la précédente ...
  235. *               SEGINI,XMATR1=XMATRI
  236. *               IMATR1.IMATTT(K)=XMATR1
  237. C          ... Boucle sur les variables duales ...
  238.                DO 8 L=1,NLIGRD
  239.                   COE=MLREE3.PROG(L)
  240.                   IF(COE.EQ.1.D0) GO TO 8
  241. C             ... Si le coefficient est différent de 1 ...
  242.                   DO 9 M=1,NLIGRP
  243. C                ... On multiplie la ligne N° L par ce coeff. ...
  244.                      XMATR1.RE(L,M,k)=XMATR1.RE(L,M,k)*COE
  245.    9              CONTINUE
  246.    8           CONTINUE
  247. C          ... Boucle sur les variables primales ...
  248.                DO 10 L=1,NLIGRP
  249.                   COE=MLREE2.PROG(L)
  250.                   IF(COE.EQ.1.D0) GO TO 10
  251. C             ... Si le coefficient est différent de 1 ...
  252.                   DO 11 M=1,NLIGRD
  253. C                ... On multiplie la colonne N° L par ce coeff. ...
  254.                      XMATR1.RE(M,L,k)=XMATR1.RE(M,L,k)*COE
  255.   11              CONTINUE
  256.   10           CONTINUE
  257. *               SEGDES XMATR1
  258.     7       CONTINUE
  259. C       ... Ménage ...
  260.             SEGDES XMATR1,XMATRI
  261.   15        CONTINUE
  262.             SEGDES DESCR
  263.     1    CONTINUE
  264.  
  265. *      write(6,fmt='('' assem0 norinc norval norind norvad'',4i6)
  266. *     * ') norinc,norval,norind,norvad
  267.  
  268.          SEGDES INWAIT
  269.          SEGDES MRIGID
  270.  
  271.       ELSEIF(ICLE.EQ.2) THEN
  272.  
  273. *        write(6,fmt='('' assem0   inwait '' ,i6 )') inwait
  274.  
  275. C    ... Destruction des matrices normalisées et remise dans MRIGID
  276. C        des matrices d'origine conservées dans INWAIT ...
  277.          inwait=inwuit
  278.          IF(INWAIT.EQ.0) RETURN
  279.          SEGACT INWAIT
  280.          MRIGID=MRIGI1
  281.          SEGACT MRIGID*MOD
  282.          DO 20 I=1,IIM(/1)
  283.             IF(IIM(I).EQ.0) GO TO 20
  284.             XMATRI=IRIGEL(4,I)
  285. **            SEGACT IMATRI
  286. *            DO 21 L=1,IMATTT(/1)
  287. *               XMATRI=IMATTT(L)
  288. *               SEGSUP XMATRI
  289. *   21       CONTINUE
  290.             SEGSUP XMATRI
  291.             IRIGEL(4,I) = IIM(I)
  292.    20    CONTINUE
  293.          SEGSUP INWAIT
  294.          SEGDES MRIGID
  295.          INWUIT=0
  296.       ENDIF
  297.  
  298. *      write(6,fmt='(''  sortie assem0  norinc norval'',2i6)')
  299. *     *  norinc,norval
  300. *      mlmots = norinc
  301. *      mlreel = norval
  302. *      segact mlmots,mlreel
  303. *      write(6,fmt='(''  sortie assem0  definitive'')')
  304.       RETURN
  305.       END
  306.  
  307.  
  308.  
  309.  
  310.  
  311.  
  312.  
  313.  
  314.  
  315.  
  316.  

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