Télécharger la plaquette Cast3M


fr en

Télécharger intgev.eso

Retour à la liste

Numérotation des lignes : 

intgev

  1. C INTGEV    SOURCE    SP204843  23/01/18    21:15:03     11560          
  2. C
  3. C=======================================================================
  4. C
  5. C   INTEGRATION DE L'ORDONNEE SUR LES ABSCISSES D'UN OBJET DE TYPE
  6. C   EVOLUTION
  7. C
  8. C   APPELEE PAR INTGRA
  9. C
  10. C   L'INTEGRATION EST EFFECTUEE PAR LA METHODE DES TRAPEZES,
  11. C   LE PAS D'INTEGRATION EST CALCULE A CHAQUE INSTANT
  12. C
  13. C   SI IABSO=1, ON INTEGRE LA VALEUR ABSOLUE DES ORDONNEES
  14. C
  15. C=======================================================================
  16. C
  17.       SUBROUTINE INTGEV(IPEVO,IABSO,ILENTA,ILENTB,ILREEA,ILREEB,
  18.      &                  XINT,IPINT,IK)
  19. C
  20.       IMPLICIT INTEGER(I-N)
  21.       IMPLICIT REAL*8(A-H,O-Z)
  22. C
  23.       CHARACTER*6 NOMC
  24. C
  25.       SEGMENT LEVINT(NINT),LOPP(NINT),LRESU(NEVO)
  26. C
  27.       POINTEUR MLREE4.MLREEL,MLREE5.MLREEL,MLREE6.MLREEL
  28.  
  29. -INC PPARAM
  30. -INC CCOPTIO
  31. -INC SMLENTI
  32. -INC SMLREEL
  33. -INC SMEVOLL
  34. -INC SMNUAGE
  35. -INC CCREEL
  36. C
  37.       IPEVOL = IPEVO
  38.       IPINT  = 0
  39.       XINT   = 0.D0
  40.       IK     = 0 
  41.  
  42. C=======================================================================
  43. C     QUELQUES VERIFICATIONS SUR EVOLUTION EN ENTREE
  44. C=======================================================================
  45. C
  46.       CALL ACTOBJ('EVOLUTIO',IPEVOL,1)
  47.       MEVOLL= IPEVOL
  48.       NEVOLL= IEVOLL(/1)
  49.  
  50. c     on n'accepte pas les evolutions VIDES (c'est un choix discutable)
  51.       IF(NEVOLL.EQ.0) THEN
  52.         MOTERR(1:8)='EVOLUTIO'
  53.         CALL ERREUR(555)
  54.         RETURN
  55.       ENDIF
  56.  
  57. C     RESTRICTION AUX EVOLUTIONS AVEC LISTREEL EN ABSC ET ORDO
  58.       SEGINI,MEVOL1=MEVOLL
  59.       N1 = 0 
  60.       DO 10 IE=1,NEVOLL
  61.         KEVOLL=IEVOLL(IE)
  62.         IF (TYPX(1:8).EQ.'LISTREEL'.AND. TYPY(1:8).EQ.'LISTREEL') THEN
  63.           MLREE1 = IPROGX
  64.           IF (MLREE1.PROG(/1).LE.1) THEN
  65.             CALL ERREUR(897)
  66.             RETURN
  67.           ELSE
  68.             N1 = N1 + 1
  69.             MEVOL1.IEVOLL(N1) = KEVOLL
  70.           ENDIF
  71.         ENDIF
  72.  10   CONTINUE
  73.       IF (N1.EQ.0) THEN
  74.         CALL ERREUR(1116)
  75.         RETURN
  76.       ENDIF 
  77.       IF (N1.NE.NEVOLL) THEN
  78.         N = N1
  79.         SEGADJ,MEVOL1
  80.         IPEVOL = MEVOL1
  81.         NEVOLL = N1
  82.       ENDIF
  83.  
  84. C=======================================================================
  85. C     RESTRICTION DE L'EVOLUTION AUX INTERVALLES D'INTEGRATION
  86. C=======================================================================
  87. C
  88. C---- CAS DES LISTENTI -------------------------------------------------
  89. C
  90.       IF (ILENTA.NE.0) THEN
  91.         CALL ACTOBJ('LISTENTI',ILENTA,1)
  92.         CALL ACTOBJ('LISTENTI',ILENTB,1)
  93.         MLENT1 = ILENTA
  94.         MLENT2 = ILENTB
  95.  
  96.         IF (MLENT1.LECT(/1).NE.MLENT2.LECT(/1)) THEN
  97.           CALL ERREUR(1115)
  98.           RETURN
  99.         ENDIF
  100.  
  101. C       LEVINT(i) : POINTEUR SUR EVOLUTION RESTREINTE AU ie INTERVALLE
  102. C       LOPP(i)   : VAUT 1 SI BORNES INTERVALLE INTEGRATION INVERSEE
  103. C                   => PRENDRE VALEYR OPPOSEE INTG DANS CE CAS
  104.         NINT = MLENT1.LECT(/1)
  105.         SEGINI,LEVINT,LOPP
  106.         DO 20 IL=1,NINT
  107.           IA = MLENT1.LECT(IL)
  108.           IB = MLENT2.LECT(IL)
  109.           IF (IA.GT.IB) THEN
  110.             LOPP(IL) = 1
  111.             IX = IA
  112.             IA = IB
  113.             IB = IX
  114.           ELSEIF (IA.EQ.IB) THEN
  115.             LOPP(IL) = 2
  116.             IB = IA + 1
  117.           ENDIF
  118.           CALL ECRENT(IB)
  119.           CALL ECRENT(IA)
  120.           CALL ECRCHA('INDI')
  121.           CALL TEVOLU(IPEVOL,'COMP')
  122. C         CALL PRLIST
  123.           CALL LIROBJ('EVOLUTION',IPEVI,1,IRET)
  124.           IF (IERR.NE.0) RETURN
  125.           LEVINT(IL) = IPEVI
  126.  20     CONTINUE
  127.  
  128. C---- CAS DES LISTREEL -------------------------------------------------
  129. C
  130.       ELSEIF (ILREEA.NE.0) THEN
  131.         CALL ACTOBJ('LISTREEL',ILREEA,1)
  132.         CALL ACTOBJ('LISTREEL',ILREEB,1)
  133.         MLREE1 = ILREEA
  134.         MLREE2 = ILREEB
  135.  
  136.         IF (MLREE1.PROG(/1).NE.MLREE2.PROG(/1)) THEN
  137.           CALL ERREUR(577)
  138.           RETURN
  139.         ENDIF
  140.  
  141.         MEVOLL = IPEVOL
  142.         NINT = MLREE1.PROG(/1)
  143.         SEGINI,LEVINT,LOPP
  144.         DO 25 IL=1,NINT
  145.           XA = MLREE1.PROG(IL)
  146.           XB = MLREE2.PROG(IL)
  147.           IF (XB.LT.XA) THEN
  148.             LOPP(IL) = 1
  149.             XX = XA
  150.             XA = XB
  151.             XB = XX
  152.           ELSEIF (XA.EQ.XB) THEN
  153.             LOPP(IL) = 2
  154.             XB = XA + 1.D0
  155.           ENDIF
  156. C
  157.           N = NEVOLL
  158.           SEGINI, MEVOL1
  159.           DO 250 IE=1,NEVOLL
  160.             KEVOLL = MEVOLL.IEVOLL(IE)
  161.             MLREE3 = KEVOLL.IPROGX
  162.             MLREE4 = KEVOLL.IPROGY
  163.  
  164. C         RESTRICTION KEVOLL A [XA;XB]
  165. C         Je parcours la liste :
  166. C           - quand je suis entre les bornes, je copie
  167. C           - quand je passe les bornes, j'interpolle
  168. C           - quand je suis en dehors, je saute
  169.             JG0    = MLREE3.PROG(/1)
  170.             JG     = 3 * JG0
  171.             SEGINI, MLREE5, MLREE6
  172.             NVAL = 0 
  173.             DO 251 IX=1,JG0
  174.               XI = MLREE3.PROG(IX)
  175.  
  176.               IF (XI.LT.XA) THEN
  177.                 IF (IX.LT.JG0) THEN
  178.                   XII = MLREE3.PROG(IX+1)
  179.                   IF (XII.GT.XA) THEN
  180. C     write(6,*) 'cas 1 : IX , XI', IX , XI
  181.                     NVAL = NVAL+1 
  182.                     MLREE5.PROG(NVAL) = XA
  183.                     FXI = (XII - XA) / (XII - XI)
  184.                     YI1 = MLREE4.PROG(IX)
  185.                     YI2 = MLREE4.PROG(IX+1)
  186.                     MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  187.                     IF (XII.GT.XB) THEN
  188. C     write(6,*) 'cas 1b : IX , XI', IX , XI
  189.                       NVAL = NVAL+1 
  190.                       MLREE5.PROG(NVAL) = XB
  191.                       FXI = (XII - XB) / (XII - XI)
  192.                       YI1 = MLREE4.PROG(IX)
  193.                       YI2 = MLREE4.PROG(IX+1)
  194.                       MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  195.                     ENDIF
  196.                   ENDIF
  197.                 ELSE
  198. C                 Rien a faire !
  199.                 ENDIF
  200.  
  201.               ELSEIF (A_EGALE_B(XI,XA)) THEN
  202. C     write(6,*) 'cas 2 : IX , XI', IX , XI
  203.                   NVAL = NVAL+1 
  204.                   MLREE5.PROG(NVAL) = XI
  205.                   MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  206.  
  207.               ELSE
  208.                 IF (IX.EQ.1) THEN
  209.                   NVAL = NVAL+1
  210.                   MLREE5.PROG(NVAL) = XA
  211.                   MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  212.                 ENDIF
  213.  
  214.                 IF (XI.LT.XB) THEN
  215.                   IF (IX.LT.JG0) THEN
  216.                     XII = MLREE3.PROG(IX+1)
  217.                     IF (XII.GT.XB) THEN
  218. C     write(6,*) 'cas 3 : IX , XI', IX , XI
  219.                       NVAL = NVAL+1 
  220.                       MLREE5.PROG(NVAL) = XI
  221.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  222.                       NVAL = NVAL+1 
  223.                       MLREE5.PROG(NVAL) = XB
  224.                       FXI = (XII - XB) / (XII - XI)
  225.                       YI1 = MLREE4.PROG(IX)
  226.                       YI2 = MLREE4.PROG(IX+1)
  227.                       MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  228.                     ELSEIF (XII.LT.XA) THEN
  229. C                   Xi+1 < XA
  230. C     write(6,*) 'cas 3b : IX , XI', IX , XI
  231.                       NVAL = NVAL+1 
  232.                       MLREE5.PROG(NVAL) = XI
  233.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  234.                       NVAL = NVAL+1 
  235.                       MLREE5.PROG(NVAL) = XA
  236.                       FXI = (XII - XA) / (XII - XI)
  237.                       YI1 = MLREE4.PROG(IX)
  238.                       YI2 = MLREE4.PROG(IX+1)
  239.                       MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  240.                     ELSE
  241. C     write(6,*) 'cas 3c : IX , XI', IX , XI
  242.                       NVAL = NVAL+1 
  243.                       MLREE5.PROG(NVAL) = XI
  244.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  245.                     ENDIF
  246.                   ELSE
  247.                     IF (NVAL.EQ.0) THEN
  248.                       CALL ERREUR(21)
  249.                       RETURN
  250.                     ENDIF
  251.                     IF (XI.GE.MLREE5.PROG(NVAL)) THEN
  252. C     write(6,*) 'cas 5 : IX , XI', IX , XI
  253.                       NVAL = NVAL+1 
  254.                       MLREE5.PROG(NVAL) = XI
  255.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  256.                       NVAL = NVAL+1 
  257.                       MLREE5.PROG(NVAL) = XB
  258.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  259.                     ELSE
  260. C     write(6,*) 'cas 5b : IX , XI', IX , XI
  261.                       NVAL = NVAL+1 
  262.                       MLREE5.PROG(NVAL) = XI
  263.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  264.                       NVAL = NVAL+1 
  265.                       MLREE5.PROG(NVAL) = XA
  266.                       MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  267.                     ENDIF
  268.                   ENDIF
  269.  
  270.                 ELSEIF (A_EGALE_B(XI,XB)) THEN
  271. C     write(6,*) 'cas 6 : IX , XI', IX , XI
  272.                     NVAL = NVAL+1 
  273.                     MLREE5.PROG(NVAL) = XI
  274.                     MLREE6.PROG(NVAL) = MLREE4.PROG(IX)
  275.  
  276.                 ELSEIF (XI.GT.XB) THEN
  277. C     write(6,*) 'XI GT XB',XI,XB
  278.                   IF (IX.GT.1) THEN
  279.                     XII = MLREE3.PROG(IX-1)
  280.                     IF (XII.LT.XB) THEN
  281. C     write(6,*) 'cas 7 : IX , XI', IX , XI
  282.                       NVAL = NVAL+1 
  283.                       MLREE5.PROG(NVAL) = XB
  284.                       FXI = (XII - XB) / (XII - XI)
  285.                       YI1 = MLREE4.PROG(IX)
  286.                       YI2 = MLREE4.PROG(IX-1)
  287.                       MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  288.                       IF (XII.LT.XA) THEN
  289. C     write(6,*) 'cas 7b : IX , XI', IX , XI
  290.                         NVAL = NVAL+1 
  291.                         MLREE5.PROG(NVAL) = XA
  292.                         FXI = (XII - XA) / (XII - XI)
  293.                         YI1 = MLREE4.PROG(IX)
  294.                         YI2 = MLREE4.PROG(IX-1)
  295.                         MLREE6.PROG(NVAL) = FXI*YI1+(1.d0-FXI)*YI2
  296.                       ENDIF
  297.                     ENDIF
  298.                   ELSE
  299. C     write(6,*) 'cas 7c : IX , XI', IX , XI
  300.                     NVAL = NVAL+1
  301.                     MLREE5.PROG(NVAL) = XA
  302.                     MLREE6.PROG(NVAL) = MLREE4.PROG(1)
  303.                     NVAL = NVAL+1
  304.                     MLREE5.PROG(NVAL) = XB
  305.                     MLREE6.PROG(NVAL) = MLREE4.PROG(1)
  306.  
  307.                   ENDIF
  308.                 ENDIF
  309.               ENDIF
  310.  251        CONTINUE
  311.             JG = NVAL
  312.             SEGADJ, MLREE5, MLREE6
  313. C     write(6,*) 'MLREE5 =',(MLREE5.PROG(II),II=1,NVAL)
  314. C     write(6,*) 'MLREE6 =',(MLREE6.PROG(II),II=1,NVAL)
  315.  
  316. C           MAJ IEeme EVOLUTION 
  317.             SEGINI, KEVOL1 
  318.             KEVOL1.IPROGX = MLREE5
  319.             KEVOL1.IPROGY = MLREE6
  320.             KEVOL1.TYPX = KEVOLL.TYPX
  321.             KEVOL1.TYPY = KEVOLL.TYPY
  322.             MEVOL1.IEVOLL(IE) = KEVOL1
  323.  250      CONTINUE
  324. C         Fin de boucle sur les evolutions 
  325.  
  326.           LEVINT(IL) = MEVOL1
  327.  25     CONTINUE
  328. C       Fin de boucle sur les intervalles
  329.  
  330.       ELSE
  331.         NINT = 1
  332.         SEGINI,LEVINT,LOPP
  333.         LEVINT(1) = IPEVOL
  334.       ENDIF
  335.  
  336. C=======================================================================
  337. C     CALCUL DES INTEGRALES 
  338. C=======================================================================
  339.  
  340. C     LRESU(i) : POINTEUR SUR LISTREEL RESULTAT INTEGRATION ie COURBE
  341.       NEVO = NEVOLL
  342.       SEGINI,LRESU
  343.  
  344. C --- BOUCLE SUR LES COURBES KEVOLL ---
  345.       DO 31 I1=1,NEVOLL
  346.  
  347. C       Initialisation LISTREEL SOLUTION POUR I1e COURBE
  348.         JG = NINT
  349.         SEGINI,MLREEL
  350.  
  351. C       BOUCLE SUR LES INTERVALLES D'INTEGRATION
  352.         DO 32 J2=1,NINT
  353.  
  354.           MEVOLL = LEVINT(J2)
  355. C         SEGACT,MEVOLL
  356.           KEVOLL = IEVOLL(I1)
  357.  
  358. C         INTEGRATION I1e COURBE SUR J2e INTERVALLE
  359.           MLREE1=IPROGX
  360. C         SEGACT MLREE1
  361.           MLREE2=IPROGY
  362. C         SEGACT MLREE2
  363.           L1=MLREE1.PROG(/1)
  364.           XINT = 0.D0
  365.           DO 33 K3=1,(L1-1)
  366.             XPAS = MLREE1.PROG(K3+1)-MLREE1.PROG(K3)
  367.             IF (IABSO.EQ.1) THEN
  368.               XVAL = ABS(MLREE2.PROG(K3))+ABS(MLREE2.PROG(K3+1))
  369.             ELSE
  370.               XVAL = MLREE2.PROG(K3)+MLREE2.PROG(K3+1)
  371.             ENDIF
  372.             XINT = 0.5D0*XVAL*XPAS + XINT
  373.  33       CONTINUE
  374.  
  375. C         ON PREND L'OPPOSEE SI BORNES INVERSEE
  376.           IF (LOPP(J2).EQ.1) XINT = -1.D0*XINT
  377.           IF (LOPP(J2).EQ.2) XINT = -0.D0*XINT
  378.           MLREEL.PROG(J2) = XINT
  379.  
  380.  32     CONTINUE
  381.  
  382. C       ON STOCKE LE RESULTAT POUR LA I1e COURBE
  383.         LRESU(I1) = MLREEL
  384.  
  385.  31    CONTINUE
  386.  
  387. C=======================================================================
  388. C     MISE EN FORME DES RESULTATS
  389. C=======================================================================
  390.       IF (NEVOLL.GT.1) THEN
  391.         IF (NINT.GT.1) THEN
  392.           IK = 3
  393.           NVAR = NEVOLL
  394.           NBCOUP = 1
  395.           SEGINI, MNUAGE
  396.           DO 40 I4=1,NVAR
  397.             IF (                  I4.LT.10) WRITE(NOMC,51) I4
  398.             IF (     I4.GE.10.AND.I4.LT.100) WRITE(NOMC,52) I4
  399.             IF (    I4.GE.100.AND.I4.LT.1000) WRITE(NOMC,53) I4
  400.             IF (   I4.GE.1000.AND.I4.LT.10000) WRITE(NOMC,54) I4
  401.             IF (  I4.GE.10000.AND.I4.LT.100000) WRITE(NOMC,55) I4
  402.             IF ( I4.GE.100000.AND.I4.LT.1000000) WRITE(NOMC,56) I4
  403.             IF (I4.GE.1000000) THEN
  404.               CALL ERREUR(21)
  405.               RETURN
  406.             ENDIF
  407.             NUANOM(I4) = 'IE'//NOMC
  408.             NUATYP(I4) = 'LISTREEL'
  409.             SEGINI, NUAVIN
  410.             NUAPOI(I4) = NUAVIN
  411.             NUAINT(1) = LRESU(I4)
  412.  40       CONTINUE
  413.           IPINT = MNUAGE
  414.         ELSE
  415.           IK = 2
  416.           JG = NEVOLL
  417.           SEGINI, MLREE1
  418.           DO 45 I4=1,NEVOLL
  419.             MLREEL = LRESU(I4)
  420.             MLREE1.PROG(I4) = MLREEL.PROG(1)
  421.  45       CONTINUE
  422.           IPINT = MLREE1
  423.         ENDIF
  424.       ELSE
  425.         IF (NINT.GT.1) THEN
  426.           IK = 2
  427.           IPINT = MLREEL
  428.         ELSE
  429.           IK = 1
  430.         ENDIF
  431.       ENDIF
  432.  
  433.  51   FORMAT (I1)
  434.  52   FORMAT (I2)
  435.  53   FORMAT (I3)
  436.  54   FORMAT (I4)
  437.  55   FORMAT (I5)
  438.  56   FORMAT (I6)
  439.  
  440.       RETURN
  441.       END
  442.  
  443.  
  444.  
  445.  
  446.  

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