Télécharger la plaquette Cast3M


fr en

Télécharger ntape6.eso

Retour à la liste

Numérotation des lignes : 

ntape6

  1. C NTAPE6    SOURCE    CHAT      05/01/13    02:02:22     5004
  2.       SUBROUTINE NTAPE6(MCP,MCQ,IVMINU,IVMINL,IVMAXU,IVMAXL,IVLAMB,
  3.      * M,N,NVD,IVFP,IVFQ,MVDU,MVDL,IVB,IVD,IVN,II,KK,IVDR,IDVD,
  4.      * NDR,TERMIN,IVLL,IVUL,IPBASE)
  5.       IMPLICIT INTEGER(I-N)
  6.       IMPLICIT REAL*8(A-H,O-Z)
  7.       LOGICAL TERMIN
  8. -INC TMXMAT
  9. -INC SMLENTI
  10.  
  11. -INC PPARAM
  12. -INC CCOPTIO
  13. -INC SMLREEL
  14.       POINTEUR MLREE4.MLREEL,MLREE5.MLREEL,MLREE6.MLREEL,MLREE7.MLREEL
  15.       IVXL1=0
  16.       IVXU1=0
  17.       IVU1=0
  18.       IVN1=0
  19.       IVD1=0
  20.       IVGM1=0
  21.       IVGE1=0
  22.       IVLAM1=0
  23.       N11 = N + 1
  24.       II=-1
  25.       KK=-1
  26.       MLREEL=IVDR
  27.       JG=N11
  28.       SEGINI MLREE1,MLREE2
  29.       IVPZ=MLREE2
  30.       IVQZ=MLREE1
  31.       MXMAT=MCQ
  32.       CALL MATVE1(XMAT,PROG,M,N11,MLREE1.PROG,1)
  33.       MXMAT=MCP
  34.       CALL MATVE1(XMAT,PROG,M,N11,MLREE2.PROG,1)
  35.       ALPMAX=1.E25
  36.       MLREE1=IVLAMB
  37.       IPB=0
  38.       DO 1 I=1,M
  39.       IF(PROG(I).LT.-1.D-20)   THEN
  40.           IF(ALPMAX+(MLREE1.PROG(I)/PROG(I)).GT.0.D0) THEN
  41.             ALPMAX= -( MLREE1.PROG(I)/PROG(I))
  42.             IPB=I
  43.           ENDIF
  44.       ENDIF
  45.     1 CONTINUE
  46.       IF(IIMPI.EQ.1799)
  47.      *WRITE(IOIMP,FMT='('' ALPHAMAX  IPB = '',E12.5,2X,I3)')ALPMAX,IPB
  48.       NDIS=0
  49.       MLENTI=IDVD
  50.       DO 7 I=1,NVD
  51.          NDIS=NDIS+LECT(I)-1
  52.     7 CONTINUE
  53.       JG=(NDIS-NVD)+2*(N11-NVD)+2
  54.       SEGINI MLREE6,MLENT1,MLENT2
  55.       MLREE6.PROG(1)=0.D0
  56.       MLENT1.LECT(1)=-3
  57.       MLENT2.LECT(1)=0
  58. *
  59. *  CALCUL DES ALPHA CRITIQUES
  60. *  ON COMMENCE PAR LES VARIABLES CONTINUES
  61. *
  62.       IP=1
  63.       MLREEL=IVD
  64.       MLREE1=IVN
  65.       MLREE2=IVMINU
  66.       MLREE3=IVMINL
  67.       MLREE4=IVQZ
  68.       MLREE5=IVPZ
  69.       DO 2 I=NVD+1,N11
  70.          CN=PROG(I)*(MLREE2.PROG(I)**2)
  71.          CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2))
  72.          CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I)
  73.          CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I))
  74.          IP=IP+1
  75.          IF(CD.EQ.0.D0) CD=1.D-35
  76.          MLREE6.PROG(IP)=CN/CD
  77.          MLENT1.LECT(IP)=I-NVD
  78.          MLENT2.LECT(IP)=-1
  79.     2 CONTINUE
  80.       MLREE2=IVMAXU
  81.       MLREE3=IVMAXL
  82.       DO 3 I=NVD+1,N11
  83.          CN=PROG(I)*(MLREE2.PROG(I)**2)
  84.          CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2))
  85.          CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I)
  86.          CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I))
  87.          IP=IP+1
  88.          IF(CD.EQ.0.D0) CD=1.D-35
  89.          MLREE6.PROG(IP)=CN/CD
  90.          MLENT1.LECT(IP)=I-NVD
  91.          MLENT2.LECT(IP)=-1
  92.     3 CONTINUE
  93. *
  94. * POUR LES VARIABLES DISCRETES
  95. *
  96.       IF(NVD.NE.0) THEN
  97.        IF(IIMPI.EQ.1799)
  98.      * WRITE(IOIMP,FMT='('' IL Y A'',I4,'' VARIABLES DISCRETES'')')NVD
  99.         MXMAT=MVDU
  100.         MXMA1=MVDL
  101.         DO 4 I=1,NVD
  102.           DO 4 J = 3,LECT(I)
  103.             CN=PROG(I)*XMAT(I,J)*XMAT(I,J-1)
  104.             CN=CN-(MLREE1.PROG(I)*MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1))
  105.             CD=MLREE4.PROG(I)* MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1)
  106.             CD=CD-(MLREE5.PROG(I)*XMAT(I,J)*XMAT(I,J-1))
  107.             IP=IP+1
  108.             IF(CD.EQ.0.D0) CD=1.D-35
  109.             MLREE6.PROG(IP)= CN/CD
  110.             IF(ABS(MLREE6.PROG(IP)).LT.1D-10)MLREE6.PROG(IP)=0.D0
  111.             MLENT1.LECT(IP)=I
  112.             MLENT2.LECT(IP)=J
  113.     4   CONTINUE
  114.       ENDIF
  115.       MAC2=MLENT1
  116.       MAC3=MLENT2
  117.       SEGSUP MLREE4,MLREE5
  118. *
  119. * ON ENTRE ALPHA MAX DANS LE TABLEAU DES ALPHA
  120. *
  121.       IP=IP+1
  122.       MLREE6.PROG(IP)=ALPMAX
  123. *
  124.       MLREEL=MLREE6
  125. *
  126. *  TRI DES ALPHA
  127. *
  128.       JG=PROG(/1)
  129.       SEGINI MLREE1
  130.       SEGINI MLENT1,MLENT2
  131.       MORDRE=MLENT1
  132.       DO 5 I=1,JG
  133.          MLENT1.LECT(I)=I
  134.    5  CONTINUE
  135.       CALL TRIFLO(PROG,MLREE1.PROG,MLENT1.LECT,MLENT2.LECT,JG)
  136.       IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI
  137.      * '',/,(1X,5E12.5))')(PROG(I),I=1,IP)
  138.       SEGSUP MLREE1,MLENT2
  139. *
  140. *   ELIMINATION DES ALPHA NON ADMISSIBLES
  141. *
  142.       IDL=0
  143.       IFL=IP+2
  144.       DO 6 I=1,IP
  145.          IF(PROG(I).LT.1.D-20) IDL=I
  146.          IF(PROG(I).GT.ALPMAX)THEN
  147.              IF(IFL.EQ.IP+2) IFL=I-1
  148.              MLENT1.LECT(I)=-2
  149.          ENDIF
  150.    6  CONTINUE
  151.       IFL=MIN(IFL,IP)
  152.       IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI
  153.      * ET ELIMINATION '',/,(1X,5E12.5))')(PROG(I),I=IDL,IFL)
  154. *
  155. *   RECHERCHE DU PAS ALPHA OPTIMUM
  156. *
  157.       MLREE1=IVLAMB
  158.       MLREE2=IVDR
  159.       SEGINI,MLREE3=MLREE1
  160.       IVLAM1=MLREE3
  161. *
  162. *  CALCUL DU PREMIER POINT
  163. *
  164.       I=IDL
  165.       INPINF=IDL
  166.       QQZ=PROG(I)
  167.       DO 8 J=1,MLREE1.PROG(/1)
  168.          MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
  169.     8 CONTINUE
  170.       IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' VALEUR DU LAMBDA POUR ALPHA
  171.      * =0 '',/,(1X,5E12.5))')(MLREE3.PROG(I),I=1,M)
  172.       CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
  173.      * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
  174.       CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
  175.      *IVGM1,IVLAM1,IPBASE)
  176.       XLP1=0.D0
  177.       MLREE4=IVGE1
  178.       DO 9 J=1,MLREE1.PROG(/1)
  179.       XLP1=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP1
  180.     9 CONTINUE
  181.       IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='(''  1ER ALPHA , PROD SCAL :''
  182.      *  ,1X,E12.5,'' , '',E12.5)')QQZ,XLP1
  183.          MLREE4=IVXU1
  184.          MLREE5=IVU1
  185.          MLREE6=IVN1
  186.          SEGSUP MLREE4,MLREE5,MLREE6
  187.          MLREE4=IVD1
  188.          MLREE5=IVGM1
  189.          MLREE6=IVGE1
  190.          MLREE7=IVXL1
  191.          SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
  192.       IVXL1=0
  193.       IVXU1=0
  194.       IVU1=0
  195.       IVN1=0
  196.       IVD1=0
  197.       IVGM1=0
  198.       IVGE1=0
  199. *
  200. * CALCUL SUR LE DERNIER POINT
  201. *
  202.       INPSUP=IFL
  203.       I=IFL
  204.       QQZ=PROG(I)
  205.       DO 10 J=1,MLREE1.PROG(/1)
  206.          MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
  207.    10 CONTINUE
  208.       CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
  209.      * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
  210.       CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
  211.      *IVGM1,IVLAM1,IPBASE)
  212.       XLP2=0.D0
  213.       MLREE4=IVGE1
  214.       DO 11 J=1,MLREE1.PROG(/1)
  215.          XLP2=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP2
  216.    11 CONTINUE
  217.       IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' DERNIER ALPHA , PROD SCAL
  218.      * :'',1X,E12.5,'' , '',E12.5)')QQZ,XLP2
  219. *
  220. * RECHERCHE  PAR DICHOTOMIE SUR LES INDICES pour trouver l'interval
  221. * de alpha interessant
  222. *
  223.       IF(XLP1 * XLP2.GE.0.D0) THEN
  224.           ALPHA=ALPMAX
  225.           IF (ALPMAX.LT.1E25) THEN
  226.              KK = -3
  227.              II=IPB
  228.           ENDIF
  229.           GO TO 20
  230.       ENDIF
  231.    12 CONTINUE
  232.       IK=(INPINF+INPSUP)/2
  233.       IF(IK.EQ.INPINF) GO TO 15
  234.          MLREE4=IVXU1
  235.          MLREE5=IVU1
  236.          MLREE6=IVN1
  237.          SEGSUP MLREE4,MLREE5,MLREE6
  238.          MLREE4=IVD1
  239.          MLREE5=IVGM1
  240.          MLREE6=IVGE1
  241.          MLREE7=IVXL1
  242.          SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
  243.       IVXL1=0
  244.       IVXU1=0
  245.       IVU1=0
  246.       IVN1=0
  247.       IVD1=0
  248.       IVGM1=0
  249.       IVGE1=0
  250.       I=IK
  251.       QQZ=PROG(I)
  252.       DO 13 J=1,MLREE1.PROG(/1)
  253.          MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
  254.    13 CONTINUE
  255.       CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
  256.      * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
  257.       CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
  258.      *IVGM1,IVLAM1,IPBASE)
  259.       XLP3=0.D0
  260.       MLREE4=IVGE1
  261.       DO 14 J=1,MLREE1.PROG(/1)
  262.          XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3
  263.    14 CONTINUE
  264.       IF(XLP1*XLP3.LE.0.D0) THEN
  265.          XLP2=XLP3
  266.          INPSUP=I
  267.       ELSE
  268.          XLP1=XLP3
  269.          INPINF=I
  270.       ENDIF
  271.       GO TO 12
  272.    15 CONTINUE
  273. *
  274. * RECHERCHE PAR DICHOTOMIE pour determiner la valeur de alpha donnant le
  275. * min de l(lambda + alpha * grad)
  276. *
  277.       IPP=0
  278.       ALPE=PROG(INPINF)
  279.       ALGR=PROG(INPSUP)
  280.       ALPHA=(ALGR + ALPE)/2.D0
  281.    60 CONTINUE
  282.       IPP=IPP + 1
  283.       IF(IPP.GT.100) THEN
  284.         CALL ERREUR ( 603)
  285.         RETURN
  286.        ENDIF
  287.       ALPHA0=ALPHA
  288.          MLREE4=IVXU1
  289.          MLREE5=IVU1
  290.          MLREE6=IVN1
  291.          SEGSUP MLREE4,MLREE5,MLREE6
  292.          MLREE4=IVD1
  293.          MLREE5=IVGM1
  294.          MLREE6=IVGE1
  295.          MLREE7=IVXL1
  296.          SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
  297.       IVXL1=0
  298.       IVXU1=0
  299.       IVU1=0
  300.       IVN1=0
  301.       IVD1=0
  302.       IVGM1=0
  303.       IVGE1=0
  304. *      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA
  305.       QQZ=ALPHA
  306.       DO 17 J=1,MLREE1.PROG(/1)
  307.          MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
  308.   17  CONTINUE
  309.       CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
  310.      *    IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
  311.       CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
  312.      *IVGM1,IVLAM1,IPBASE)
  313.       XLP3=0.D0
  314.       MLREE4=IVGE1
  315.       DO 50 J=1,MLREE1.PROG(/1)
  316.          XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3
  317.    50 CONTINUE
  318.       IF(XLP1*XLP3.LE.0.D0) THEN
  319.          XLP2=XLP3
  320.          ALGR=ALPHA
  321.       ELSE
  322.          XLP1=XLP3
  323.          ALPE=ALPHA
  324.       ENDIF
  325.       ALPHA =(ALGR + ALPE)/2
  326.       EPSI=1.D-10
  327.       IF(ALPHA.LT.1.D-12.AND.NDR.GE.1)THEN
  328.          TERMIN=.TRUE.
  329.          GO TO 20
  330.       ENDIF
  331.       IF((ABS(ALPHA-ALPHA0)/ALPHA).GE.EPSI)   GO TO 60
  332.    20 CONTINUE
  333.          MLENT1=MAC2
  334.          MLENT2=MAC3
  335.          MLENT3=MORDRE
  336.        IF(IIMPI.EQ.1799)WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA
  337.       IF (KK.NE.-3)THEN
  338. *
  339. ****************** TEST SI ARRET SUR PLAN DE DISCONTINUITE ************
  340. *
  341.          ALPHA1=PROG(INPSUP)
  342.          IF((ABS(ALPHA1-ALPHA)/ALPHA).LE.EPSI)THEN
  343.               KK=MLENT2.LECT(MLENT3.LECT(INPSUP))
  344.               II=MLENT1.LECT(MLENT3.LECT(INPSUP))
  345.               ALPHA=ALPHA1
  346.          ENDIF
  347.          ALPHA2=PROG(INPINF)
  348.          IF(((ABS(ALPHA2-ALPHA)/ALPHA).LE.EPSI.AND.ALPHA2.NE.0.D0)
  349.      *     .OR.ALPHA.LT.1.D-12)THEN
  350.               KK=MLENT2.LECT(MLENT3.LECT(INPINF))
  351.               II=MLENT1.LECT(MLENT3.LECT(INPINF))
  352.               ALPHA=ALPHA2
  353.          ENDIF
  354.       ENDIF
  355.       DO 21 J=1,MLREE1.PROG(/1)
  356.          MLREE3.PROG(J)=ALPHA * MLREE2.PROG(J)+MLREE1.PROG(J)
  357.    21 CONTINUE
  358.       IF (KK.EQ.-3)THEN
  359.          MLREE3.PROG(II)=0.
  360.       ENDIF
  361.       ALPHA0=ALPHA
  362.          MLREE4=IVXU1
  363.          MLREE5=IVU1
  364.          MLREE6=IVN1
  365.          SEGSUP MLREE4,MLREE5,MLREE6
  366.          MLREE4=IVD1
  367.          MLREE5=IVGM1
  368.          MLREE6=IVGE1
  369.          MLREE7=IVXL1
  370.          SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
  371.       IVLAMB=IVLAM1
  372.       SEGSUP MLENT1,MLENT2,MLREEL,MLENT3
  373.       RETURN
  374.       END
  375.  
  376.  

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