Télécharger la plaquette Cast3M


fr en

Télécharger depimp.eso

Retour à la liste

Numérotation des lignes : 

depimp

  1. C DEPIMP    SOURCE    PV        21/03/04    21:15:05     10913          
  2.       SUBROUTINE DEPIMP
  3.  
  4. ************************************************************************
  5. *     CE SUBROUTINE SERT A IMPOSER DES VALEURS DE DEPLACEMENTS
  6. *     IMPOSES NON NULS.
  7. *
  8. *  SYNTAXE  TOTO = DEPIMPOSE   BRIG  FLOT
  9. *     OU    TOTO = DEPIMPOSE   BRIG  CHPOI ( COMPOSANTES PRIMALES)
  10. *     OU    TOTO = DEPIMPOSE   BRIG  'RELA' CHPSCAL
  11. *
  12. *  ENTREE : BRIG = OBJET RIGIDITE DE TYPE  BLOQUAGE
  13. *           FLOT = VALEUR DU DEPLACEMENT A IMPOSER
  14. *           CHPOI   = chpoint AVEC LES DDLS PRIMALS
  15. *           CHPSCAL = CHPOINT DE SCALAIRE QUI PRECISE LA
  16. *                     VALEUR A IMPOSER EN CHAQUE POINT.
  17. *
  18. *  SORTIE : TOTO = OBJET DE TYPE CHPOINT (FLX)
  19. *
  20. ************************************************************************
  21.  
  22. ************************************************************************
  23. *     DECLARATIONS ET INITIALISATIONS
  24. ************************************************************************
  25.  
  26.       IMPLICIT INTEGER(I-N)
  27.       CHARACTER*4 charm,charre
  28.       REAL*8 XXA,vval,X0,X1
  29.       LOGICAL ISCALA
  30.       PARAMETER(NCLE=1)
  31.       CHARACTER*4 MOCLE(NCLE)
  32.       DATA MOCLE  /'RELA'/
  33.  
  34. -INC SMRIGID
  35. -INC SMCHPOI
  36.  
  37. -INC PPARAM
  38. -INC CCOPTIO
  39. -INC SMELEME
  40. -INC SMCOORD
  41. -INC SMTABLE
  42.  
  43.       character*4 cnoha
  44.       integer*4 inoha
  45.       data cnoha/'NOHA'/
  46.       equivalence(inoha,cnoha)
  47.  
  48.       SEGMENT SCOLOR
  49.        CHARACTER*(LOCOMP) COLOR(NBELEM)
  50.       ENDSEGMENT
  51.       POINTEUR SCOL1.SCOLOR,SCOL2.SCOLOR,SCOL3.SCOLOR
  52.       SEGMENT ICPR(NNN)
  53.       segment irelat
  54.          logical lrelat(nnn)
  55.       end segment
  56.       logical lllblo,lllrel
  57.  
  58. C     INITIALISATIONS
  59.       ISCALA = .FALSE.
  60.  
  61.  
  62. ************************************************************************
  63. *     LECTURES ET TESTS PRELIMINAIRES DES ENTREES
  64. ************************************************************************
  65.  
  66. C  ****  LECTURE TABLE LIAISONS STATIQUES
  67.       CALL LIRTAB('LIAISONS_STATIQUES',ipt,0,iretou)
  68.       IF (IRETOU.NE.0) THEN
  69.         CALL DEPIM2(IPT)
  70.         RETURN
  71.       ENDIF
  72. C
  73. C  ****  LECTURE D'UN OBJET DE TYPE RIGIDITE
  74. C
  75.       CALL LIROBJ('RIGIDITE',IPOIRI,1,IRETOU)
  76.       IF(IERR.NE.0) RETURN
  77. C
  78. C  ****  LECTURE D'UN FLOTTANT OU D'UN CHPOINT
  79. C
  80. C     LECTURE D'UNE VALEUR 
  81.       CALL LIRREE(XXA,0,IREVAL)
  82.       VVAL=XXA
  83. c     SI ECHEC LECTURE D'UN CHPOINT DE SCALAIRES OU DE DDL PRIMAL
  84.       IF(IREVAL.EQ.0) THEN
  85.  
  86. *         mot-cle 'RELA' ? ==> ISCALA
  87.           CALL LIRMOT(MOCLE,NCLE,ICLE,0)
  88.           IF(IERR.NE.0) RETURN
  89.           IF(ICLE.EQ.1) ISCALA=.TRUE.
  90.  
  91.           CALL LIROBJ('CHPOINT ',ISCA,1,IRETOU)
  92.           IF(IERR.NE.0) RETURN
  93.           MCHPO1=ISCA
  94. c           SEGACT MCHPO1
  95.           CALL ACTOBJ('CHPOINT ',MCHPO1,1)
  96. C         Si le CHPOINT n'a aucune sous-zone, il est vide, alors erreur
  97.           NBSZCH=MCHPO1.IPCHP(/1)
  98.           IF(NBSZCH.LT.1) THEN
  99.              MOTERR(1:8)='CHPOINT '
  100.              INTERR(1)=ISCA
  101.              CALL ERREUR(356)
  102.              RETURN
  103.           ENDIF
  104.  
  105. c         RELA => cas SCALAIRE : 1 zone et 1 composante nommee 'SCAL'
  106.           IF(ISCALA) THEN
  107. c           verif : 1 seule zone
  108.             IF(NBSZCH.NE.1) THEN
  109.               MOTERR(1:8)='CHPOINT '
  110.               INTERR(1)=ISCA
  111. c             Le %m1:8 de pointeur %i1 n'est pas elementaire (n<>1)
  112.               CALL ERREUR(110)
  113.               RETURN
  114.             ENDIF          
  115.             MSOUP1 = MCHPO1.IPCHP(1)
  116. c           segact MSOUP1
  117. c           verif : 1 seule composante
  118.             NBCOMP = MSOUP1.NOHARM(/1)
  119.             IF(NBCOMP.NE.1) THEN
  120. c             Il faut specifier un champ par point avec une seule composante
  121.               CALL ERREUR(180)
  122.               RETURN
  123.             ENDIF           
  124.             IF(MSOUP1.NOCOMP(1).NE.'SCAL') THEN
  125.               MOTERR(1:4)='SCAL'
  126. c             La composante %m1:4 ne peut etre extraite du champ par point specifie
  127. c             car elle en est absente
  128.               CALL ERREUR(181)
  129.               RETURN
  130.             ENDIF
  131. c           ici ISCALA=TRUE et tout va bien !
  132.           ENDIF
  133.  
  134.       ENDIF
  135. c
  136. c  ... test si la RIGIDITE n'est pas vide, si OUI on cree un CHPOINT
  137. c      vide puis on s'en va ...
  138. c
  139.       MRIGID=IPOIRI
  140.       SEGACT,MRIGID
  141.       NNN=IRIGEL(/2)
  142.       IF (NNN.EQ.0) THEN
  143.          NSOUPO=0
  144.          NAT=1
  145.          SEGINI MCHPOI
  146.          MTYPOI='FLX'
  147.          JATTRI(1) = 2
  148.       MOCHDE=' CE CHAMP PAR POINTS A ETE CREE PAR LE SOUS-PROGRAMME'//
  149.      #  ' DEPIMP'
  150.          IFOPOI = IFOUR
  151.          GO TO 252
  152.       ENDIF
  153.  
  154.  
  155. ************************************************************************
  156. *    TRAVAIL
  157. ************************************************************************
  158.  
  159.       IPT2=0
  160.       NOHA=0
  161. C
  162. ************************************************************************
  163. C BOUCLE SUR LES SOUS RIGIDITES . ON VERIFIE QUE LAMBDA EXISTE ET ON
  164. C CONSTRUIT LE SEGMENT GEOMETRIE LX1 LX2 NNOE, DANS scol1.COLOR ON MET LE
  165. C NOM DE L'INCONNUE
  166. ************************************************************************
  167. C
  168.       segini,irelat
  169.       DO 1 NN=1,NNN
  170.          DESCR=IRIGEL(3,NN)
  171.          MELEME=IRIGEL(1,NN)
  172.          NOHAR=IRIGEL(5,NN)
  173.          IF(NOHA.NE.0.AND.NOHA.NE.NOHAR) THEN
  174.             CALL ERREUR ( 25 )
  175.             RETURN
  176.          ENDIF
  177.          lrelat(nn) = .false.
  178. c    ... on va chercher les multiplicateurs dans DESCR ...
  179.          SEGACT,DESCR
  180.          IA=LISINC(/2)
  181.          if (ia.ne.noelep(/1)) then
  182.           write(6,*) ' descr longueur ',descr,ia
  183.           call erreur(5)
  184.          endif
  185.          DO 2 I=1,IA
  186.             IF(LISINC(I).EQ.'LX  ') GO TO 3
  187.     2    CONTINUE
  188. c    ... on n'a pas trouve de multiplicateurs, donc bye ...
  189.          SEGDES,DESCR
  190.          CALL ERREUR(245)
  191.          RETURN
  192. c    ... on a trouve les multiplicateurs ...
  193.    3     CONTINUE
  194.          SEGACT,MELEME
  195.          NBNN=2
  196.          NBELEM=NUM(/2)
  197.          NBREF=0
  198.          NBSOUS=0
  199.          SEGINI,IPT1,SCOL1
  200. c    ... boucle sur les elements de blocage ...
  201.          DO J=1,NUM(/2)
  202.             JB=0
  203. c       ... JA sert a compter les multiplicateurs dans chaque
  204. c           element, un seul est permis
  205.             JA=0
  206. c       ... boucle sur les noeuds de ces elements ...
  207.             DO K=1,NOELEP(/1)
  208. c          ... si c'est un support de multiplicateur, on met son n°
  209. c              dans IPT1 (position 1 ) ...
  210.                IF(LISINC(K).EQ.'LX  ') THEN
  211.                   JA=JA+1
  212.                   if (ja.gt.1) then
  213.                    write(6,*) ' plus que 1 LX dans la matrice ',descr
  214.                    call erreur(5)
  215.                   endif
  216.                   IPT1.NUM(JA,J)=NUM(NOELEP(K),J)
  217. c          ... sinon ...
  218.                ELSE
  219. c             ... on teste si c'est le premier DDL <<physique>>, si OUI ...
  220.                   IF(JB.EQ.0) THEN
  221. c                ... on met son n° dans IPT1 (position 2) ...
  222.                      JB=2
  223.                      IPT1.NUM(JB,J)=NUM(NOELEP(K),J)
  224. C                ... et le nom du DDL dans SCOL1.COLOR ...
  225.                      SCOL1.COLOR(J)=LISINC(K)
  226. c             ... sinon (c.a d. ceci est une relation et non un blocage) ...
  227.                   ELSE
  228. c                ... on teste si le support n'est pas le même que
  229. c                    celui du premier DDL <<physique>> ...
  230.                      IF(IPT1.NUM(JB,J).NE.NUM(NOELEP(K),J)) THEN
  231. c                   ... si c'est le cas on sert une ERREUR en cas de lecture d'un CHPOINT ...
  232.                         IF(IREVAL.ne.1) then
  233.                            CALL ERREUR(794)
  234.                            RETURN
  235.                         endif
  236.                      ENDIF
  237. c                ... et de toute façon on efface le nom du DDL de SCOL1.COLOR ...
  238.                      SCOL1.COLOR(J)=' '
  239.                      lrelat(nn) = .true.
  240.                   ENDIF
  241.                ENDIF
  242.             ENDDO
  243.          ENDDO
  244.  
  245. C
  246. C  SI NN= 1 IPT2 = IPT1; SINON IPT3 = IPT2 + IPT1, PUIS IPT2 = IPT3
  247. C
  248.          SEGDES,DESCR
  249.          IF(IPT2.NE.0) GO TO 5
  250.          IPT2=IPT1
  251.          SCOL2=SCOL1
  252.          GO TO 1
  253.     5    CONTINUE
  254.          NA=IPT1.NUM(/2)
  255.          NB=IPT2.NUM(/2)
  256.          NBELEM=NA+NB
  257.          SEGINI,IPT3,SCOL3
  258.          DO 71 I=1,NA
  259.             SCOL3.COLOR(I)=SCOL1.COLOR(I)
  260.             DO 72 J=1,2
  261.                IPT3.NUM(J,I)=IPT1.NUM(J,I)
  262.    72       CONTINUE
  263.    71    CONTINUE
  264.          DO 8 I=1,NB
  265.             SCOL3.COLOR(I+NA)=SCOL2.COLOR(I)
  266.             DO 9 J=1,2
  267.                IPT3.NUM(J,I+NA)=IPT2.NUM(J,I)
  268.     9       CONTINUE
  269.     8    CONTINUE
  270.          SEGSUP IPT1,SCOL1
  271.          SEGSUP,IPT2,SCOL2
  272.          IPT2=IPT3
  273.          SCOL2=SCOL3
  274.     1 CONTINUE
  275.       SEGDES,MRIGID
  276. c
  277. c ... on va verifier si la matrice contient seulement les blocages, ou seulement
  278. c     les relations (on en a seulement besoin en cas de lecture d'un CHPOINT) ...
  279. c
  280.       if(ireval.eq.0) then
  281.          lllblo = .true.
  282.          lllrel = .true.
  283.          do 77 i=1,nnn
  284.             lllblo = lllblo .and. (.not.lrelat(i))
  285.             lllrel = lllrel .and.       lrelat(i)
  286.  77      continue
  287.          segsup,irelat
  288. c    ... si les deux sont faux, on a un melange de blocages et de relations ...
  289.          if(.not.lllblo .and. .not.lllrel) then
  290.             call erreur(795)
  291.             return
  292.          endif
  293. c    ... si on n'a que des relations, le CHPOINT doit être scalaire ...
  294.          if(lllrel .and. .not.iscala) then
  295. c           ERREUR : "Impossible d'imposer les valeurs aux relations 
  296. c                     avec un champ non scalaire"
  297.             call erreur(796)
  298.             return
  299.          endif
  300.       else
  301.          segsup,irelat
  302.       endif
  303. C
  304. C ON VIENT DE CREER IPT2 CONTENANT DES ELEMENTS COMPOSES DE LX1  NOE
  305. C DANS COLOR ON A LE NOM DE L'INCONNUE A METTRE EN FACE DE NNOE
  306. C
  307.       NSOUPO=1
  308.       NAT=1
  309.       SEGINI,MCHPOI
  310.       MTYPOI='FLX'
  311.       JATTRI(1) = 2
  312.       MOCHDE=' CE CHAMP PAR POINTS A ETE CREE PAR LE SOUS-PROGRAMME'//
  313.      #  ' DEPIMP'
  314.       IFOPOI=IFOUR
  315.  
  316.       NC=1
  317.       SEGINI,MSOUPO
  318.       IPCHP(1)=MSOUPO
  319.       NOCOMP(1)='FLX'
  320.       NOHARM(1)=NOHAR
  321.       write (charm,fmt='(A4)') nohar
  322.       if (nohar.eq.inoha) noharm(1)=nifour
  323. C
  324. ************************************************************************
  325. C  CREATION DE L'ELEMENT SUPPORT GEOMETRIQUE ET EN MEME TEMPS DES
  326. C  VALEURS VPOCHA
  327. ************************************************************************
  328. C
  329.       NBNN=1
  330.       NBELEM=IPT2.NUM(/2)
  331.       SEGINI MELEME
  332.       IGEOC=MELEME
  333.       ITYPEL=1
  334.  
  335.       N=IPT2.NUM(/2)
  336.       SEGINI,MPOVAL
  337.       IPOVAL=MPOVAL
  338.  
  339. c ... Si on a lu un reel, il n'y a pas grand chose a faire ...
  340.       IF(IREVAL.NE.0) GO TO 250
  341. C
  342. c   
  343. C   + CAS DU CHPOINT SCALAIRE ------------------------------------------
  344. c     (on teste seulement ISCALA car on a deja verifie que cela va 
  345. c      ensemble avec LLLREL)
  346.       IF(ISCALA) THEN
  347. c           write(*,*) '>>> DEPI d un chpoint SCALAIRE <<<'
  348.           MSOUP1=MCHPO1.IPCHP(1)
  349.           SEGACT MSOUP1
  350.           MPOVA1=MSOUP1.IPOVAL
  351.           SEGACT MPOVA1
  352.           NNN=nbpts
  353.           SEGINI ICPR
  354.           IPT3=MSOUP1.IGEOC
  355.           SEGACT IPT3
  356.           NNU=IPT3.NUM(/2)
  357. c         numerotation locale
  358.           DO 25 IUY=1,NNU
  359.              ICPR(IPT3.NUM(1,IUY))=IUY
  360.   25      CONTINUE
  361.           DO 26 IU=1,IPT2.NUM(/2)
  362.              NUM(1,IU)=IPT2.NUM(1,IU)
  363.              INOD2=IPT2.NUM(2,IU)
  364.              ID=ICPR(INOD2)
  365.              IF(ID.EQ.0) THEN
  366. c              ERREUR : "Un point de l'objet rigidite n'est pas  
  367. c                        inclus dans le champ de scalaire"
  368.                CALL ERREUR(244)
  369.                RETURN               
  370.              ELSEIF(ID.EQ.-1) THEN
  371. c              Le noeud apparait dans plusieurs relations --> ERREUR :
  372. c              "On ne peut avoir 2 relations sur un meme ddl  noeud %i1"
  373.                INTERR(1)=INOD2
  374.                CALL ERREUR(886) 
  375.                RETURN 
  376.              ELSE
  377.                XXA=MPOVA1.VPOCHA(ID,1)
  378.                VPOCHA(IU,1)=XXA
  379.                ICPR(INOD2)=-1
  380.              ENDIF
  381.    26     CONTINUE
  382.           SEGSUP ICPR          
  383. C
  384. C    + CAS DU CHPOINT D'INCONNUES PRIMALES -----------------------------
  385.        ELSE
  386.  
  387.          NNN=nbpts
  388.          SEGINI ICPR
  389.          JB=1
  390.          DO 36 J=1,IPT2.NUM(/2)
  391.             NUM(1,JB)=IPT2.NUM(1,J)
  392.             JB=JB+1
  393.    36    CONTINUE
  394.          DO 31 I=1,MCHPO1.IPCHP(/1)
  395.             DO 40 J=1,NNN
  396.                ICPR(J)=0
  397.    40       CONTINUE
  398.             MSOUP1=MCHPO1.IPCHP(I)
  399.             SEGACT MSOUP1
  400.             MPOVA1=MSOUP1.IPOVAL
  401.             SEGACT MPOVA1
  402.             IPT1=MSOUP1.IGEOC
  403.             SEGACT IPT1
  404.             IA=0
  405.             DO 32 J=1,IPT1.NUM(/2)
  406.                ID=IPT1.NUM(1,J)
  407.                IF(ICPR(ID).EQ.0) THEN
  408.                   IA=IA+1
  409.                   ICPR(ID)=IA
  410.                ELSE
  411.                   CALL ERREUR(245)
  412.                   RETURN
  413.                ENDIF
  414.   32        CONTINUE
  415.             DO 33 J=1,IPT2.NUM(/2)
  416.                ID=IPT2.NUM(2,J)
  417.                JB=JB+2
  418.                IF(ICPR(ID).EQ.0) GO TO 33
  419.                DO 34 K=1,MSOUP1.NOCOMP(/2)
  420.                   IF(MSOUP1.NOCOMP(K).EQ.SCOL2.COLOR(J)) GO TO 35
  421.    34          CONTINUE
  422.                GO TO 33
  423.    35          CONTINUE
  424.                JD=ICPR(ID)
  425.                XXA=MPOVA1.VPOCHA(JD,K)
  426.                JA=J
  427.                VPOCHA(JA,1)=XXA
  428.    33       CONTINUE
  429.    31    CONTINUE
  430.          SEGSUP ICPR
  431.  
  432.       ENDIF
  433. C   + FIN CAS DES CHPOINTS SCALAIRE OU PAS -----------------------------
  434. c     le chpoint d'entree est inutile -> segdes
  435.       CALL ACTOBJ('CHPOINT ',MCHPO1,0)
  436.       GO TO 251
  437.  
  438.  
  439. C     CAS DU FLOTTANT --------------------------------------------------
  440. C ... En cas de lecture d'un reel le remplissage du segment MPOVAL est assez simple ...
  441.   250 CONTINUE
  442.       DO 10 I=1,N
  443.          VPOCHA(I,1)=VVAL
  444.    10 CONTINUE
  445. c ... celui du segment MELEME n'est pas plus complique ...
  446.       DO 11 I=1,IPT2.NUM(/2)
  447.          NUM(1,I)=IPT2.NUM(1,I)
  448.    11 CONTINUE
  449.  
  450.  
  451. c     TOUS LES CAS -----------------------------------------------------
  452.   251 CONTINUE
  453.       SEGSUP IPT2,SCOL2
  454.   252 CONTINUE
  455. c     chpoint de sortie -> segact
  456.       CALL ACTOBJ('CHPOINT ',MCHPOI,1)
  457.       CALL ECROBJ('CHPOINT ',MCHPOI)
  458.  
  459.       END
  460.  
  461.  
  462.  
  463.  
  464.  

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