Télécharger la plaquette Cast3M


fr en

Télécharger trjava.eso

Retour à la liste

Numérotation des lignes : 

trjava

  1. C TRJAVA    SOURCE    CB215821  20/11/25    13:41:27     10792          
  2.       SUBROUTINE TRJAVA(IZVIT,IZPART,IZN3,IPART,TMIN,TMAX,IZCOU,
  3.      *      MELEME,IELTFA,IZCENT,IFACEL,DELTAT,IZSH,IEROR,DT1)
  4. C
  5. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  6. C        Issu du sp AVATAR de TRIO-EF ( adapté par F Auriol)
  7. C
  8. C    FAIT AVANCER UNE PARTICULE ( COORDONNEES DE REFERENCES )
  9. C                 JUSQU'A L'INSTANT TMAX OU A L'INSTANT DE SORTIE
  10. C    IZVIT SEGMENT DECRIVANT LES VITESSES ( <--- TRJVIT TRJFLU)
  11. C    IZPART POSITIONS INITIALES DES PARTICULES
  12. C    IZN3  SEGMENT RESULTANT (AJUSTE ICI)
  13. C    IPART NUMERO DE LA PARTICULE
  14. C    TMIN  INSTANT DE DEPART
  15. C    TMAX  INSTANT A NE PAS DEPASSER
  16. C    IZCOU SEGMENT DES DT DONNANT UN NOMBRE DE COURANT DE 1
  17. C                  COU VALEUR DU NOMBRE DE COURANT VOULU
  18. C                       2 PERM
  19. C                       3 TRANS
  20. C    MELEME POINTEUR DU MAILLAGE
  21. C    IELTFA POINTEUR DE DOMAINE.ELTFA
  22. C    IZCENT POINTEUR DE DOMAINE.CENTRE
  23. C    IFACEL POINTEUR DE DOMAINE.FACEL
  24. C    DELTAT PAS DE TEMPS AVEC LEQUEL ON CONSERVE LES RESULTATS
  25. C
  26. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  27. C
  28.       IMPLICIT INTEGER(I-N)
  29.       IMPLICIT REAL*8 (A-H,O-Z)
  30. C
  31.  
  32. -INC PPARAM
  33. -INC CCOPTIO
  34. -INC SMCOORD
  35. -INC SMELEME
  36. -INC SMCHPOI
  37. C
  38.       POINTEUR IZCENT.MELEME,IELTFA.MELEME,IZFAC1.MELEME
  39.       POINTEUR IFACEL.MELEME
  40. C
  41.       SEGMENT IZPART
  42.            INTEGER NLEPA(NPART),NUMPA(NPART)
  43.            REAL*8  COORPA(NDIM,NPART)
  44.       ENDSEGMENT
  45. C
  46.       SEGMENT  IZN3
  47.            INTEGER NAPAR3(NPOS),NUM3(NPOS)
  48.            REAL*8 CREF3(NDIM,NPOS),TPAR(NPOS)
  49.       ENDSEGMENT
  50. C
  51.       SEGMENT  IZCOU
  52.            REAL*8 DTCO(NEL),COU
  53.       ENDSEGMENT
  54. C
  55.       SEGMENT IZSH
  56.             REAL*8 SHP(6,MNO9),SHY(12,MNO9),XYZL(3,MNO9)
  57.       ENDSEGMENT
  58. C
  59.       SEGMENT IZVIT
  60.            REAL*8 TEMTRA(NVIPT)
  61.            INTEGER IPUN(NBS),IDUN(NBS),IPVPT(NVIPT),IFORML
  62.       ENDSEGMENT
  63. C             IDUN(I) nombre d elements avant le sous maillage I
  64. C             IPVPT  pointeurs de izvpt pour chaque pas de temps
  65.       SEGMENT IZVPT
  66.            INTEGER IPUN1(NBS),IPUMAX
  67.       ENDSEGMENT
  68. C
  69.       SEGMENT IZUN
  70.            REAL*8 UN(I1,I2,I3)
  71.       ENDSEGMENT
  72.       POINTEUR IZUN1.IZUN,IZUN2.IZUN
  73.       SEGMENT IZUMAX
  74.            REAL*8 UMAX(NBREL)
  75.       ENDSEGMENT
  76. C
  77.       DIMENSION UELEM(3),XARI(3),XDEP(3),XINT(3)
  78. C
  79. C***
  80. C
  81.       ISENS=0
  82.       IEROR=0
  83.       KIO=0
  84.       NDIM=IDIM
  85.       NPOS=50
  86.       SEGINI IZN3
  87.       SEGACT IZVIT
  88.       IFORMU=IFORML
  89.       IEL1=NLEPA(IPART)
  90.       DO 2 ID1=1,NDIM
  91.         XDEP(ID1)=COORPA(ID1,IPART)
  92.    2  CONTINUE
  93.         ITER=2
  94.         NAPAR3(1)=IEL1
  95.         DO 3 ID1=1,NDIM
  96.           CREF3(ID1,1)=XDEP(ID1)
  97.    3    CONTINUE
  98.         TCOUR=TMIN
  99.         TPAR(1)=TMIN
  100.         NUAPAR=IEL1
  101. C       DT1=0.D0
  102.         SENS=SIGN(1.D0,(TMAX-TMIN))
  103. C
  104.        SEGACT IELTFA,IZCENT,IFACEL
  105.        SEGACT IZSH
  106.       NVIPT=TEMTRA(/1)
  107.    1  CONTINUE
  108. C       IVPT        1 EN PERMANENT
  109. C              EN TRANSITOIRE NUMERO DU PAS DE TEMPS DIRECTEMENT
  110. C              SUPERIEUR AU TEMPS CONSIDERE
  111.  
  112.       IF(NVIPT.EQ.1)THEN
  113.             IVPT=1
  114.       ELSE
  115.             IVPT=2
  116.             CALL TRJTPT(IZVIT,TCOUR,IVPT)
  117.       ENDIF
  118.       CALL MELNEL(IEL1,MELEME,IPT1,NEL0,0)
  119.       NOEL1=IPT1.NUM(/1)
  120.       IELL=IEL1-NEL0
  121.       CALL DOXE(XCOOR,IDIM,NOEL1,IPT1.NUM,IELL,XYZL)
  122.       ITY1=IPT1.ITYPEL
  123. C
  124. C***  NOMBRE DE COURANT
  125. C
  126. C     IF(DTCO(IEL1).EQ.0.D0) THEN
  127. C         WRITE(6,*) 'APPEL A TRJCN5 DANS TRJAVA '
  128.           CALL TRJCN5(ITY1,IZSH,IJK)
  129. C         IF(IERR.GT.0) RETURN
  130.           IF(ISENS.NE.IJK)THEN
  131.                 IF(ISENS.EQ.0)THEN
  132.                       ISENS=IJK
  133.                 ELSE
  134.                       IEROR=1
  135.                       GO TO 10
  136.                 ENDIF
  137.           ENDIF
  138.           IZVPT=IPVPT(IVPT)
  139.           SEGACT IZVPT
  140.           IZUMAX=IPUMAX
  141.           SEGACT IZUMAX
  142.           UEM=UMAX(IEL1)
  143.           SEGDES IZUMAX
  144.           SEGDES IZVPT
  145.           IF(IVPT.NE.1)THEN
  146.                 IZVPT=IPVPT(IVPT-1)
  147.                 SEGACT IZVPT
  148.                 IZUMAX=IPUMAX
  149.                 SEGACT IZUMAX
  150.                 UEM=MAX(UEM,UMAX(IEL1))
  151.                 SEGDES IZUMAX
  152.                 SEGDES IZVPT
  153.           ENDIF
  154.           NUCENT=IZCENT.NUM(1,IEL1)
  155.           CALL MELNEL(IEL1,IELTFA,IZFAC1,NEL0,1)
  156.           CALL TRJCOU(UEM,IZCOU,IZFAC1,IELL,IEL1,NUCENT)
  157. C     ENDIF
  158. C
  159. C***  PAS DE TEMPS
  160. C
  161.         TPREC=TCOUR
  162.       IF(ABS(DELTAT).GT.1.D-15)THEN
  163.           DTT=DELTAT
  164.       ELSE
  165.           DTT=10.D0*COU*DTCO(IEL1)*SENS
  166.           IF (ABS(TMAX-TPREC).LE.ABS(DTT))DTT=TMAX-TPREC
  167.           DT1=0.D0
  168.       ENDIF
  169.       DTI=DTCO(IEL1)*COU*SENS
  170. C     write(6,*) ' dtco ',DTCO(IEL1),IEL1,COU,DTT
  171.       IF(ABS(DTT-DT1).LT.ABS(DTI*1.1D0))THEN
  172.         DTI=DTT-DT1
  173.         NI=1
  174.       ELSE
  175.         NI=INT((DTT-DT1)/DTI)+1
  176.         DTI=(DTT-DT1)/NI
  177.       ENDIF
  178. C     write(6,*)iel1,ni,dti,dt1,dtt,li
  179.       LI=1
  180. C     DT1=0.D0
  181. C
  182.     6 CONTINUE
  183.       DTC=DTI*0.5D0
  184.       DTC1=DTC
  185. C      write(6,*)' trjava li ni ',li,ni,iel1,dti,dtc,dtt
  186.       TTEMP=TPREC+(LI-1)*DTI
  187.       IF(IVPT.NE.1)CALL TRJTPT(IZVIT,TTEMP,IVPT)
  188.       CALL TRJVEL(IZVIT,IZUN,IEL1,IVPT,TTEMP)
  189.       NOUN=UN(/2)
  190.       CALL TRJULO(XYZL,UN(1,1,IELL),UELEM,XDEP,ITY1,NDIM,NOEL1,NOUN,
  191.      * IFORMU,SHP,SHY)
  192.       DO 9 ID1=1,NDIM
  193.         XARI(ID1)=XDEP(ID1)+UELEM(ID1)*DTC
  194.    9  CONTINUE
  195.       ITYG=NUMGEO(ITY1)
  196.       CALL TRJINT(XARI,XDEP,XINT,UELEM,DTT,DTINT,NDIM,ICONT,ITYG,IO)
  197.       IF(IO.NE.0)THEN
  198. C
  199. C***  RATTRAPAGE
  200. C
  201.         KIO=KIO+1
  202. C       IF(MOD(KIO,2).EQ.0)WRITE(6,999)IPART,IEL1,IEL2
  203.         IF(MOD(KIO,2).EQ.0)THEN
  204.                   INTERR(1)=IPART
  205.                   INTERR(2)=IEL1
  206.                   CALL ERREUR(-298)
  207.         ENDIF
  208. C       WRITE(6,*)'TRJRAT 1 ',IPART,IEL1,IEL2,
  209. C    *  XARI,XDEP,XINT,UELEM,DTT,DTINT,NDIM,ICONT,ITYG,IO
  210. C       IF(MOD(KIO,10).EQ.0) STOP
  211.  
  212.         INTERR(1)=IPART
  213.         CALL TRJRAT(IEL1,IEL2,XDEP,XARI,IZCOU,MELEME,IZVIT,DTC1,
  214.      *          IZCENT,IELTFA,XINT,DTINT,ICONT,IZSH,TTEMP)
  215. C        write(6,*) ' apres trjrat ', icont ,iel1,iel2
  216.          IF(ICONT.NE.0)GO TO 16
  217.         IF(IEL2.EQ.0)THEN
  218.           DO 99 ID1=1,NDIM
  219.             XINT(ID1)=XDEP(ID1)
  220.   99      CONTINUE
  221.           DTINT=0.D0
  222.           GO TO 10
  223.         ENDIF
  224.         DO 12 ID1=1,NDIM
  225.           XDEP(ID1)=XARI(ID1)
  226.   12    CONTINUE
  227.         TCOUR=TPREC+(LI-1)*DTI+DTC1
  228.         NUAPAR=IEL2
  229.         DT1=0.D0
  230.         GO TO 13
  231.       ENDIF
  232.    16 CONTINUE
  233.       IF(ICONT.NE.0) THEN
  234. C
  235. C***  LA PARTICULE I EST SORTIE DE L'ELEMENT COURANT
  236. C
  237.         CALL TRJIEL(IEL1,IEL2,ICONT,NF,IFACEL,IZCENT,IELTFA)
  238. C     write(6,*)' apres TRJIEL 1 ',IEL1,IEL2,' nuapar' ,NUAPAR ,ITER
  239.         IF(IEL2 .EQ.0)GO TO 10
  240.         GO TO 14
  241.       ENDIF
  242.       DTC=DTI
  243.       DTC1=DTI
  244.       TTEMP=TPREC+(LI-0.5D0)*DTI
  245.       IF(IVPT.NE.1)CALL TRJTPT(IZVIT,TTEMP,IVPT)
  246.       CALL TRJVEL(IZVIT,IZUN,IEL1,IVPT,TTEMP)
  247.       NOUN=UN(/2)
  248.       CALL TRJULO(XYZL,UN(1,1,IELL),UELEM,XARI,ITY1,NDIM,NOEL1,NOUN,
  249.      * IFORMU,SHP,SHY)
  250.       DO 17 ID1=1,NDIM
  251.         XARI(ID1)=XDEP(ID1)+DTC*UELEM(ID1)
  252.   17  CONTINUE
  253.       ITYG=NUMGEO(ITY1)
  254.       CALL TRJINT(XARI,XDEP,XINT,UELEM,DTC,DTINT,NDIM,ICONT,ITYG,IO)
  255. C     write(6,*)' apres trjint ',io
  256.       IF(IO.NE.0)THEN
  257. C
  258. C***  RATTRAPAGE
  259. C
  260.         KIO=KIO+1
  261. C       IF(MOD(KIO,2).EQ.0)WRITE(6,999)IPART,IEL1,IEL2
  262.         IF(MOD(KIO,2).EQ.0)THEN
  263.                   INTERR(1)=IPART
  264.                   INTERR(2)=IEL1
  265.                   CALL ERREUR(-298)
  266.         ENDIF
  267. C       WRITE(6,*)'TRJRAT 2 ',IPART,IEL1,IEL2,
  268. C    *  XARI,XDEP,XINT,UELEM,DTT,DTINT,NDIM,ICONT,ITYG,IO
  269. C       IF(MOD(KIO,10).EQ.0) STOP
  270.         INTERR(1)=IPART
  271.         CALL TRJRAT(IEL1,IEL2,XDEP,XARI,IZCOU,MELEME,IZVIT,DTC1,
  272.      *                IZCENT,IELTFA,XINT,DTINT,ICONT,IZSH,TTEMP)
  273.         IF(ICONT.NE.0)GO TO 28
  274.         IF(IEL2.EQ.0)THEN
  275.           DO 98 ID1=1,NDIM
  276.             XINT(ID1)=XDEP(ID1)
  277.   98      CONTINUE
  278.           DTINT=0.D0
  279.           GO TO 10
  280.         ENDIF
  281.         IF(IEL1.EQ.IEL2)GO TO 11
  282.         DO 18 ID1=1,NDIM
  283.           XDEP(ID1)=XARI(ID1)
  284.   18    CONTINUE
  285.         TCOUR=TPREC+(LI-1)*DTI+DTC1
  286.         NUAPAR=IEL2
  287.         DT1=0.D0
  288.         GO TO 13
  289.       ENDIF
  290.    28 CONTINUE
  291.       IF(ICONT.NE.0) THEN
  292. C
  293. C***  LA PARTICULE I EST SORTIE DE L'ELEMENT COURANT
  294. C
  295.         CALL TRJIEL(IEL1,IEL2,ICONT,NF,IFACEL,IZCENT,IELTFA)
  296.         IF(IEL2 .EQ.0)GO TO 10
  297.         GO TO 14
  298.       ENDIF
  299. C
  300.   11  CONTINUE
  301.       LI=LI+1
  302.       DO 24 ID1=1,NDIM
  303.           IF(ABS(XARI(ID1)-XDEP(ID1)).GT.1.D-15)GO TO 25
  304.    24 CONTINUE
  305.       IF(IIMPI.GT.0)WRITE(6,*)' LA PARTICULE ',IPART ,
  306.      * ' N AVANCE PLUS ON ARETE LE CALCUL '
  307.       NUAPAR=IEL1
  308.       TCOUR=TPREC+NI*DTI
  309.       GO TO 21
  310.    25 CONTINUE
  311.       DO 19 ID1=1,NDIM
  312.         XDEP(ID1)=XARI(ID1)
  313.   19  CONTINUE
  314.       IF(LI.LE.NI)GO TO 6
  315. C
  316. C***  LA PARTICULE I EST RESTE DANS L'ELEMENT COURANT
  317. C
  318.       NUAPAR=IEL1
  319.       TCOUR=TPREC+NI*DTI
  320.       IF (ABS((TMAX-TCOUR)/TMAX).LE.1.D-04)GO TO 21
  321.         DT1=0.D0
  322.       GO TO 13
  323. C
  324. C***  LA PARTICULE I RESTE DANS LE DOMAINE DE CALCUL
  325. C
  326.   14  CONTINUE
  327.       CALL MELNEL(IEL2,MELEME,IPT2,NEL2,1)
  328.       IELL2=IEL2-NEL2
  329.       JTY=IPT2.ITYPEL
  330. CCCCCCCCCCC
  331. C         write(6,*) 'iel2 ',iel2,'iel1' ,iel1,'ITER',ITER
  332.  
  333.       NOEL2=IPT2.NUM(/1)
  334.       CALL DOXE(XCOOR,IDIM,NOEL2,IPT2.NUM,IELL2,XYZL)
  335.           CALL TRJCN5(JTY,IZSH,IJK)
  336. C         IF(IERR.GT.0) RETURN
  337.           IF(ISENS.NE.IJK) THEN
  338.           IEROR= 1
  339.           IEL2= 0
  340. C         WRITE (6,*) ' JE VAIS EN 10 ', iel1,ITER,ieror
  341.           GO TO 10
  342.           ENDIF
  343. CCCCCCCCCCC
  344.       CALL TRJFAC(IEL2,NF,JFA,IELTFA)
  345.       IOR=1
  346.       IF(NDIM.EQ.3)CALL TRJIOR(IELL,IELL2,IPT1,IPT2,ICONT,JFA,IOR)
  347. C     WRITE(6,*)'  IOR  ',IOR,TCOUR
  348.       ITYG=NUMGEO(ITY1)
  349.       JTYG=NUMGEO(JTY)
  350.       CALL TRJTRJ(XDEP,XINT,ITYG,JTYG,JFA,ICONT,IOR)
  351.       TCOUR=TPREC+(LI-1)*DTI+DTINT
  352.       NUAPAR=IEL2
  353. C     DT1=DTINT+(LI-1)*DTI
  354.       IF(ABS(DELTAT).GT.1.D-15)THEN
  355.            DT1=DTINT+(LI-1)*DTI+DT1
  356.            IEL1=NUAPAR
  357.            GO TO 1
  358.       ENDIF
  359.   13  CONTINUE
  360.       DO 22 ID1=1,NDIM
  361.         CREF3(ID1,ITER)=XDEP(ID1)
  362.   22  CONTINUE
  363.       TPAR(ITER)=TCOUR
  364.       NAPAR3(ITER)=NUAPAR
  365.       IEL1=NUAPAR
  366.       IF ((NPOS-ITER).LE.1) THEN
  367.         NPOS=NPOS+50
  368.         SEGADJ IZN3
  369.       ENDIF
  370.       ITER=ITER+1
  371.       GO TO 1
  372. C
  373. C***  LA PARTICULE I EST SORTIE DU DOMAINE DE CALCUL
  374. C
  375.   10  CONTINUE
  376.       DO 23 ID1=1,NDIM
  377.         XDEP(ID1)=XINT(ID1)
  378.   23  CONTINUE
  379.       TCOUR=TPREC+(LI-1)*DTI+DTINT
  380.       NUAPAR=IEL1
  381.   21  CONTINUE
  382.       DO 20 ID1=1,NDIM
  383.         CREF3(ID1,ITER)=XDEP(ID1)
  384.   20  CONTINUE
  385. C     write(6,*)' IEL1 NUAPAR ITER ', IEL1,NUAPAR,ITER
  386.       TPAR(ITER)=TCOUR
  387.       NAPAR3(ITER)=NUAPAR
  388.       NPOS=ITER
  389.       SEGADJ IZN3
  390. C
  391.  999  FORMAT('PARTICULE',I4,': PROBLEMES DANS LE COIN D''UN ELEMENT !!',
  392.      *        2I7 )
  393. C
  394.       END
  395.  
  396.  
  397.  
  398.  
  399.  
  400.  
  401.  

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