Télécharger la plaquette Cast3M


fr en

Télécharger d2vini.eso

Retour à la liste

Numérotation des lignes : 

d2vini

  1. C D2VINI    SOURCE    BP208322  20/09/18    21:15:13     10718          
  2.        SUBROUTINE D2VINI(ITINIT,KTKAM,KTQ,KTFEX,KTPAS,KTNUM,KTLIAA,
  3.      &                  KTLIAB,KTPHI,KCPR,KOCLFA,KOCLB1,REPRIS,
  4.      &                  RIGIDE,lmodyn)
  5. *--------------------------------------------------------------------*
  6. *                                                                    *
  7. *     Operateur DYNE : algorithme differences centrees               *
  8. *     ________________________________________________               *
  9. *                                                                    *
  10. *     Initialisation de l'algorithme ou reprise de calcul.           *
  11. *                                                                    *
  12. *     Parametres:                                                    *
  13. *                                                                    *
  14. * e   ITINIT  Table contenant les conditions initiales ou les        *
  15. *             champs necessaires a la reprise du calcul              *
  16. * es  KTKAM   Segment contenant les vecteurs XK, XASM et XM          *
  17. * es  KTQ     Segment des variables de mouvement generalisees        *
  18. *             (et des travaux)
  19. * es  KTFEX   Segment contenant les chargements libres               *
  20. * es  KTPAS   Segment des variables au cours d'un pas de temps       *
  21. * es  KTNUM   Segment contenant les parametres numeriques            *
  22. * e   KTLIAA  Segment descriptif des liaisons en base A              *
  23. * e   KTLIAB  Segment descriptif des liaisons en base B              *
  24. * e   KTPHI   Segment contenant les deformees modales                *
  25. * e   KCPR    Segment des points                                     *
  26. * -   KOCLFA  Segment contenant les tableaux locaux de la subroutine *
  27. *             D2VLFA                                                 *
  28. * -   KOCLB1  Segment contenant les tableaux locaux de la subroutine *
  29. *             D2VLB1                                                 *
  30. * e   REPRIS  Vrai si reprise de calcul                              *
  31. * e   RIGIDE  Vrai si corps rigide
  32. *                                                                    *
  33. *     Auteur, date de creation:                                      *
  34. *                                                                    *
  35. *     Denis ROBERT-MOUGIN, le 25 mai 1989.                           *
  36. *                                                                    *
  37. *--------------------------------------------------------------------*
  38.       IMPLICIT INTEGER(I-N)
  39.       IMPLICIT REAL*8(A-H,O-Z)
  40.  
  41. -INC PPARAM
  42. -INC CCOPTIO
  43. -INC CCNOYAU
  44. -INC SMCOORD
  45. -INC SMTABLE
  46. -INC CCASSIS
  47. *
  48.       SEGMENT,MTQ
  49.          REAL*8 Q1(NA1,4),Q2(NA1,4),Q3(NA1,4)
  50.          REAL*8 WEXT(NA1,2),WINT(NA1,2)
  51.       ENDSEGMENT
  52.       SEGMENT,MTFEX
  53.          REAL*8  FEXA(NPFEXA,NPC1,2)
  54.          REAL*8  FEXPSM(NPLB,NPC1,2,IDIMB)
  55.          REAL*8  FTEXB(NPLB,NPC1,2,IDIM)
  56. *         INTEGER IFEXA(NPFEXA),IFEXB(NPFEXB)
  57.       ENDSEGMENT
  58.       SEGMENT,MTPAS
  59.          REAL*8 FTOTA(NA1,4),FTOTB(NPLB,IDIMB),FTOTBA(NA1)
  60.          REAL*8 XPTB(NPLB,2,IDIMB),FINERT(NA1,4)
  61.          REAL*8 XVALA(NLIAA,4,NTVAR),XVALB(NLIAB,4,NTVAR)
  62.          REAL*8 FEXB(NPLB,2,IDIM),XCHPFB(2,NLIAB,4,NPLB)
  63.       ENDSEGMENT
  64.       SEGMENT,MTNUM
  65.          REAL*8 XDT(NPC1),XTEMPS(NPC1)
  66.       ENDSEGMENT
  67.       SEGMENT,MTKAM
  68.          REAL*8 XK(NA1,NB1K),XASM(NA1,NB1C),XM(NA1,NB1M)
  69.          REAL*8 XOPER(NB1,NB1,NOPER)
  70.       ENDSEGMENT
  71.       SEGMENT,MTLIAA
  72.          INTEGER IPALA(NLIAA,NIPALA),IPLIA(NLIAA,NPLAA),JPLIA(NPLA)
  73.          REAL*8  XPALA(NLIAA,NXPALA)
  74.       ENDSEGMENT
  75.       SEGMENT,MTLIAB
  76.          INTEGER IPALB(NLIAB,NIPALB),IPLIB(NLIAB,NPLBB),JPLIB(NPLB)
  77.          REAL*8  XPALB(NLIAB,NXPALB)
  78.          REAL*8  XABSCI(NLIAB,NIP),XORDON(NLIAB,NIP)
  79.       ENDSEGMENT
  80.       SEGMENT,MTPHI
  81.          INTEGER IBASB(NPLB),IPLSB(NPLB),INMSB(NSB),IORSB(NSB)
  82.          INTEGER IAROTA(NSB)
  83.          REAL*8  XPHILB(NSB,NPLSB,NA2,IDIMB)
  84.       ENDSEGMENT
  85.       SEGMENT,ICPR(nbpts)
  86.  
  87. * Segment "local" pour D2VLFA ...
  88.       SEGMENT,LOCLFA
  89.          REAL*8 FTEST(NA1,4)
  90.       ENDSEGMENT
  91. * Segment "local" pour D2VLB1 ...
  92.       SEGMENT,LOCLB1
  93.          REAL*8 FTEST2(NPLB,6)
  94.       ENDSEGMENT
  95.       segment mwinit
  96.         integer jpdep,jpvit,jrepr
  97.       endsegment
  98. c       SEGMENT FAMOR(NA1,4)
  99.       SEGMENT FAMOR(NA1)
  100. *
  101. c *     TODO : KTOTXA et KTOTVA a allouer dans le segment LOCLFA ?
  102. c       REAL*8 KTOTXA(NA1,NA1), KTOTVA(NA1,NA1) 
  103. c       cf. devpas d2vpas d3vpas
  104. c       en attendant, petite bidouille locale
  105.       segment LOCLF2
  106.          REAL*8 KTOTXA(NA1,NA1), KTOTVA(NA1,NA1) 
  107.       endsegment
  108.       segment LOCLF3
  109.       REAL*8 KTOTXB(NPLB,IDIMB,IDIMB),KTOTVB(NPLB,IDIMB,IDIMB) 
  110.       endsegment
  111.  
  112.       LOGICAL L0,L1,REPRIS,RIGIDE,lmodyn
  113.       CHARACTER*8 TYPRET,TYPIND,CHARRE
  114.       CHARACTER*(19) CHAI1
  115.       ILC1 = 19
  116.       DATA CHAI1 /'VARIABLES_LIAISON_A'/
  117. *
  118.       MTFEX  = KTFEX
  119.       MTPAS  = KTPAS
  120.       MTNUM  = KTNUM
  121.       MTKAM  = KTKAM
  122.       MTQ    = KTQ
  123.       MTLIAA = KTLIAA
  124.       MTLIAB = KTLIAB
  125.       MTPHI  = KTPHI
  126.       LOCLFA = KOCLFA
  127.       LOCLB1 = KOCLB1
  128.       IDEPL  = 0
  129.       IVITE  = 0
  130.       ITABL  = 0
  131.       ICH1   = 0
  132.       ICH2   = 0
  133.       ICH3   = 0
  134.       ICH4   = 0
  135.       ICH5   = 0
  136.       ICH6   = 0
  137.       NA1    = Q1(/1)
  138.       NB1K   = XK(/2)
  139.       NB1C   = XASM(/2)
  140.       NB1M   = XM(/2)
  141.       NB1    = XOPER(/1)
  142.       NLIAA  = IPALA(/1)
  143.       NLIAB  = IPALB(/1)
  144.       NPLB   = JPLIB(/1)
  145.       NSB    = XPHILB(/1)
  146.       NPLSB  = XPHILB(/2)
  147.       NA2    = XPHILB(/3)
  148.       IDIMB  = XPHILB(/4)
  149.       NPFEXA = FEXA(/1)
  150.       NPC1   = FEXPSM(/2)
  151.       NIP    = XABSCI(/2)
  152.       IERRD = 0
  153.       segini,LOCLF2,LOCLF3
  154. *
  155. *     S'agit-il d'une initialisation ou d'une reprise de calcul ? ------
  156. *
  157.       IF ( REPRIS ) THEN
  158.          IF (IIMPI.EQ.333) THEN
  159.             WRITE(IOIMP,*)'D2VINI : ---> reprise de calcul'
  160.          ENDIF
  161.          if (lmodyn) then
  162.            mwinit = itinit
  163.            segact mwinit
  164.            ITABL = jrepr
  165.          else
  166.          CALL ACCTAB(ITINIT,'MOT',I0,X0,'REPRISE',L0,IP0,
  167.      &                    'TABLE',I1,X1,' ',L1,ITABL)
  168.          endif
  169. *
  170. *        Reprise du calcul, on remplit: info initiales
  171. *
  172.          CALL ACCTAB(ITABL,'MOT',I0,X0,'DEPLACEMENT',L0,IP0,
  173.      &                 'CHPOINT',I1,X1,' ',L1,ICH1)
  174.          IF (IERR.NE.0) RETURN
  175. *
  176.          CALL ACCTAB(ITABL,'MOT',I0,X0,'VITESSE',L0,IP0,
  177.      &                 'CHPOINT',I1,X1,' ',L1,ICH2)
  178.          IF (IERR.NE.0) RETURN
  179. *
  180.          CALL ACCTAB(ITABL,'MOT',I0,X0,'FORCES_BASE_A',L0,IP0,
  181.      &                 'CHPOINT',I1,X1,' ',L1,ICH5)
  182.          IF (IERR.NE.0) RETURN
  183.  
  184.          CALL ACCTAB(ITABL,'MOT',I0,X0,'ACCELERATION',L0,IP0,
  185.      &                 'CHPOINT',I1,X1,' ',L1,ICH7)
  186.          IF (IERR.NE.0) RETURN
  187.  
  188. *
  189.          CALL ACCTAB(ITABL,'MOT',I0,X0,'TRAVAIL_EXTERIEUR',L0,IP0,
  190.      &                 'CHPOINT',I1,X1,' ',L1,ICH9)
  191.          IF (IERR.NE.0) RETURN
  192. *
  193.          CALL ACCTAB(ITABL,'MOT',I0,X0,'TRAVAIL_INTERIEUR',L0,IP0,
  194.      &                 'CHPOINT',I1,X1,' ',L1,ICH10)
  195.          IF (IERR.NE.0) RETURN
  196. *
  197. c        Q1(2)=q_0   
  198.          IF (ICH1.NE.0) THEN
  199.             CALL DYNE18(ICH1,KTQ,1,2,KCPR)
  200.          ELSE
  201.             CALL ERREUR(487)
  202.             RETURN
  203.          ENDIF
  204. *
  205. c        Q2(2)=\dot{q}_0
  206.          IF (ICH2.NE.0) THEN
  207.             CALL DYNE18(ICH2,KTQ,2,2,KCPR)
  208.          ELSE
  209.             CALL ERREUR(487)
  210.             RETURN
  211.          ENDIF
  212.  
  213. *        FTOTA(2)=F_0
  214.          IF (ICH5.NE.0) THEN
  215.             CALL DYNE23(ICH5,KTPAS,2)
  216.          ELSE
  217.             CALL ERREUR(487)
  218.             RETURN
  219.          ENDIF
  220.  
  221.          IF (ICH7.NE.0) THEN
  222.             CALL DYNE18(ICH7,KTQ,3,2,KCPR)
  223.          ELSE
  224.             CALL ERREUR(487)
  225.             RETURN
  226.          ENDIF
  227.  
  228.          IF (ICH9.NE.0) THEN
  229.             CALL DYNE18(ICH9,KTQ,4,2,KCPR)
  230.          ENDIF
  231.          IF (ICH10.NE.0) THEN
  232.             CALL DYNE18(ICH10,KTQ,5,2,KCPR)
  233.          ENDIF
  234. *
  235. *
  236.          IF (NLIAA.NE.0) THEN
  237. *
  238. *           l'indice VARIABLES_LIAISON_A n'existe que pour
  239. *           les liaisons POLYNOMIALES
  240. *
  241.             IPOLY = 0
  242.             MTABLE = ITABL
  243.             SEGACT,MTABLE
  244.             NIND1 = MLOTAB
  245.             if(nbesc.ne.0)segact ipiloc
  246.             DO 100 I=1,NIND1
  247.                TYPIND = MTABTI(I)
  248.                IF (TYPIND.EQ.'MOT     ') THEN
  249.                   IP = MTABII(I)
  250.                   ID = IPCHAR(IP)
  251.                   IFI = IPCHAR(IP+1)
  252.                   IL1 = IFI - ID
  253.                   IF (IL1.EQ.ILC1) THEN
  254.                      IF (CHAI1.EQ.ICHARA(ID:IFI-1)) THEN
  255.                         IPOLY = 1
  256.                      ENDIF
  257.                   ENDIF
  258.                ENDIF
  259.  100        CONTINUE
  260.             if(nbesc.ne.0)SEGDES,IPILOC
  261.             SEGDES,MTABLE
  262.             IF (IIMPI.EQ.333) THEN
  263.                WRITE(IOIMP,*) 'D2VINI : IPOLY = ',IPOLY
  264.             ENDIF
  265.             IF (IPOLY.NE.0) THEN
  266.                CALL ACCTAB(ITABL,'MOT',I0,X0,'VARIABLES_LIAISON_A',
  267.      &                     L0,IP0,'TABLE',I1,X1,' ',L1,ITREFR)
  268.                CALL DYNA14(ITREFR,KTLIAA)
  269.                IF (IERR.NE.0) RETURN
  270.             ENDIF
  271.          ENDIF
  272.          IF (NLIAB.NE.0) THEN
  273. *old            CALL DEVRCO(Q1,NA1,XPTB,NPLB,XPHILB,NSB,NPLSB,NA2,IDIMB,
  274.             CALL DEVRCO(Q1,Q2,NA1,XPTB,NPLB,XPHILB,NSB,NPLSB,NA2,IDIMB,
  275. *           --> XPTB(:,2,:)=x_0   XPTB(:,3,:)=\dot{q}_0
  276. *               ok car XPTB(:,3,:) utilisé uniquement pour calcul FnlB
  277.      &                  IBASB,IPLSB,INMSB,IORSB,2,IAROTA)
  278.             CALL ACCTAB(ITABL,'MOT',I0,X0,'VARIABLES_LIAISON_B',
  279.      &                  L0,IP0,'TABLE',I1,X1,' ',L1,ITREFR)
  280.             IF (IERR.NE.0) RETURN
  281.             CALL DYNE14(ITREFR,KTLIAB)
  282.             IF (IERR.NE.0) RETURN
  283.          ENDIF
  284. *
  285.       ELSE
  286. *
  287. *        Chpoints initiaux de deplacement et de vitesse:
  288. *
  289.          IF (ITINIT.NE.0) THEN
  290.           if (lmodyn) then
  291.            mwinit = itinit
  292.            segact mwinit
  293.            idepl = jpdep
  294.            ivite = jpvit
  295.           else
  296.             TYPRET=' '
  297.             CALL ACCTAB(ITINIT,'MOT',I0,X0,'DEPLACEMENT',L0,IP0,
  298.      &                        TYPRET,I1,X1,CHARRE,L1,IDEPL)
  299. *
  300.             TYPRET=' '
  301.             CALL ACCTAB(ITINIT,'MOT',I0,X0,'VITESSE',L0,IP0,
  302.      &                        TYPRET,I1,X1,CHARRE,L1,IVITE)
  303.           endif
  304.          ENDIF
  305. *
  306. *        Mise des valeurs initiales au pas m (indice 3)
  307. c        Q1(2)=q_0   
  308. c        Q2(2)=\dot{q}_0
  309.          IF (IDEPL.NE.0) THEN
  310.             CALL DYNE18(IDEPL,KTQ,1,2,KCPR)
  311.             MTQ = KTQ
  312.          ENDIF
  313.          IVINIT = 0
  314.          IF (IVITE.NE.0) THEN
  315.             CALL DYNE18(IVITE,KTQ,2,2,KCPR)
  316.             IVINIT = 1
  317.          ENDIF
  318. *
  319. *        Ajout des forces de liaison
  320. *
  321.          PDT=XDT(1)
  322.          T=XTEMPS(1)
  323.          IF (NLIAA.NE.0) THEN
  324.             CALL D2VLFA(Q1,Q2,FTOTA,NA1,IPALA,IPLIA,XPALA,XVALA,
  325.      &                  NLIAA,PDT,T,1,2,FINERT,IVINIT,FTEST,
  326.      &                  KTOTXA,KTOTVA,.FALSE.)     
  327.          ENDIF
  328.          IF (NLIAB.NE.0) THEN
  329.              IF (RIGIDE) THEN
  330.                 DO 13 IP=1,NPLB
  331.                   DO 15 ID=1,IDIM
  332.                        FEXB(IP,2,ID) = FTEXB(IP,1,2,ID)
  333.  15               CONTINUE
  334.  13             CONTINUE
  335.             ENDIF
  336.             CALL D2VLF2(Q1,Q2,FTOTA,NA1,IPALB,IPLIB,XPALB,XVALB,
  337.      &           NLIAB,XPHILB,JPLIB,NPLB,IDIMB,FTOTB,FTOTBA,XPTB,
  338.      &           PDT,T,1,IBASB,IPLSB,INMSB,IORSB,NSB,NPLSB,NA2,2,
  339.      &           FEXPSM,NPC1,IERRD,FTEST2,XABSCI,XORDON,
  340.      &           NIP,IAROTA,RIGIDE,FEXB,XCHPFB,
  341.      &           KTOTXB,KTOTVB)
  342.  
  343.             IF (IERRD.NE.0) THEN
  344.                IF (IERRD.EQ.1) CALL ERREUR(528)
  345.                RETURN
  346.             ENDIF
  347.          ENDIF
  348. *
  349. *        Calcul des deplacements et des vitesses au pas debut du initial
  350. *
  351.       IF (IIMPI.EQ.333) THEN
  352.          DO 91 I=1,NA1
  353.             WRITE(IOIMP,*)'D2VINI :'
  354.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',3) =',FTOTA(I,3)
  355.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',4) =',FTOTA(I,4)
  356.             WRITE(IOIMP,*)'D2VINI : Fexa (',I,',1) =',Fexa(I,1)
  357.             WRITE(IOIMP,*)'D2VINI : Fexa (',I,',2) =',Fexa(I,2)
  358.  91         CONTINUE
  359. *        end do
  360.       ENDIF
  361.  
  362. c        Cas selon nature (pleine ou diag) des matrices K,C et M :
  363. c------- Cas K,C ou M pleine -----------------------------------
  364.          IF(NB1K.GT.1.OR.NB1C.GT.1.OR.NB1M.GT.1) THEN
  365.  
  366. c          F_0 =  F^liaison_0 + F^ext_0 - K * Q_0
  367.            DO 61 I = 1,NA1
  368.               FTOTA(I,2) = FTOTA(I,2) + FEXA(I,1,2)
  369.  61        CONTINUE
  370.            CALL DEVLK0(Q1,XK,FTOTA,NA1,NB1K,2)
  371.            SEGINI,FAMOR
  372. c            CALL DEVLC0(Q2,XASM,FAMOR,NA1,NB1C,2)
  373.            IF (NB1C.EQ.1) THEN
  374.              DO 611 I=1,NA1
  375.                FAMOR(I) = XASM(I,1)*Q2(I,2)
  376.  611         CONTINUE
  377.            ELSE
  378.              DO 612 I=1,NA1
  379.                FAMOR(I) = 0.D0
  380.              DO 613 IB=1,NB1
  381.                FAMOR(I) = FAMOR(I) + XASM(I,IB)*Q2(IB,2)
  382.  613         CONTINUE
  383.  612         CONTINUE
  384.            ENDIF
  385. c          \ddot Q_0 = [M+dt/2C]-1 * (F_0 - F^AMOR_0)
  386. c         -cas matrice M ou C pleine
  387.            IF(NB1.NE.1) THEN
  388.              DO 63 I=1,NA1
  389.                Q3(I,2) = 0.D0
  390.              DO 64 J=1,NB1
  391.                Q3(I,2) = Q3(I,2) + XOPER(I,J,1)*(FTOTA(J,2)-FAMOR(J))
  392.  64          CONTINUE
  393.  63          CONTINUE
  394.            ELSE
  395.              DO 62 I=1,NA1
  396.                Q3(I,2) = (FTOTA(I,2)-FAMOR(I))/(XM(I,1)-FINERT(I,2))
  397.  62          CONTINUE
  398.            ENDIF
  399.            SEGSUP,FAMOR
  400.  
  401. c       -cas matrices M et C diagonales
  402.  
  403. c------- Cas K,C et M diagonal -----------------------------------
  404.          ELSE
  405.  
  406.            DO 10 I = 1,NA1
  407.              UNSM = 1.D0 / ( XM(I,1) - FINERT(I,2) )
  408.              FTOTA(I,2) = FTOTA(I,2) + FEXA(I,1,2) - XK(I,1)*Q1(I,2)
  409.      &                  - XASM(I,1)*Q2(I,2)
  410.              Q3(I,2) = FTOTA(I,2)*UNSM
  411.  10        CONTINUE
  412.  
  413.          ENDIF
  414. *
  415.       ENDIF
  416.       IF (IIMPI.EQ.333) THEN
  417.          DO 30 I=1,NA1
  418.             WRITE(IOIMP,*)'D2VINI :'
  419.             WRITE(IOIMP,*)'D2VINI :  Q1(',I,',1) =',Q1(I,1)
  420.             WRITE(IOIMP,*)'D2VINI :  Q2(',I,',1) =',Q2(I,1)
  421.             WRITE(IOIMP,*)'D2VINI :  Q3(',I,',1) =',Q3(I,1)
  422.             WRITE(IOIMP,*)'D2VINI :'
  423.             WRITE(IOIMP,*)'D2VINI :'
  424.             WRITE(IOIMP,*)'D2VINI :  Q1(',I,',2) =',Q1(I,2)
  425.             WRITE(IOIMP,*)'D2VINI :  Q2(',I,',2) =',Q2(I,2)
  426.             WRITE(IOIMP,*)'D2VINI :  Q3(',I,',2) =',Q3(I,2)
  427.             WRITE(IOIMP,*)'D2VINI :'
  428.             WRITE(IOIMP,*)'D2VINI :'
  429.             WRITE(IOIMP,*)'D2VINI :  Q1(',I,',3) =',Q1(I,3)
  430.             WRITE(IOIMP,*)'D2VINI :  Q2(',I,',3) =',Q2(I,3)
  431.             WRITE(IOIMP,*)'D2VINI :  Q3(',I,',3) =',Q3(I,3)
  432.             WRITE(IOIMP,*)'D2VINI :'
  433.             WRITE(IOIMP,*)'D2VINI :  Q1(',I,',4) =',Q1(I,4)
  434.             WRITE(IOIMP,*)'D2VINI :  Q2(',I,',4) =',Q2(I,4)
  435.             WRITE(IOIMP,*)'D2VINI :  Q3(',I,',4) =',Q3(I,4)
  436.  30         CONTINUE
  437. *        end do
  438.          DO 31 I=1,NA1
  439.             WRITE(IOIMP,*)'D2VINI :',FTOTA(/1),FTOTA(/2)
  440.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',1) =',FTOTA(I,1)
  441.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',2) =',FTOTA(I,2)
  442.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',3) =',FTOTA(I,3)
  443.             WRITE(IOIMP,*)'D2VINI : FTOTA(',I,',4) =',FTOTA(I,4)
  444.  31         CONTINUE
  445. *        end do
  446.       ENDIF
  447. *
  448.       ICPR = KCPR
  449.       SEGSUP,ICPR
  450. *
  451.       END
  452.  
  453.  
  454.  
  455.  
  456.  
  457.  
  458.  
  459.  
  460.  
  461.  
  462.  
  463.  
  464.  
  465.  

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