Télécharger la plaquette Cast3M


fr en

Télécharger trjpar.eso

Retour à la liste

Numérotation des lignes : 

trjpar

  1. C TRJPAR    SOURCE    CB215821  24/04/12    21:17:23     11897          
  2.       SUBROUTINE TRJPAR(MELEME,IZVIT,IZCENT,IELTFA,IFACEL,MTAB2,
  3.      *                  IICAL,IZPOR,IZDIFF,IZDISP,MCHELM,MMODEL)
  4. C
  5. C
  6. C
  7. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  8. C
  9. C    calcul de trajectoires appelé par TRAJEC decode la table de
  10. C       lacher et pilote le calcul
  11. C
  12. C    sous programmes appelés:
  13. C               TRJCN2 controle des données
  14. C               TRJINI premiere localisation et initialisation des
  15. C                      segments de travail
  16. C               TRJAVA fait avancer une particule
  17. C               TRJRAL en fin de calcul passe des coordonnées de
  18. C                      références aux coordonnées réelles
  19. C    Les tables résultats sont chargées par ce sous programme
  20. C
  21. C     ARGUMENTS
  22. C     E   MELEME pointeur du maillage de calcul ( DOMAINE.MAILLAGE)
  23. C     E   IZVIT  pointeur du segment décrivant les vitesses .
  24. C               Ce segment est créé dans TRJVIT et TRJFLU
  25. C     E   IZCENT pointeur de DOMAINE.CENTRE
  26. C     E   IELTFA pointeur de DOMAINE.ELTFA
  27. C     E   IFACEL pointeur de DOMAINE.FACEL
  28. C     E   MTAB2  pointeur de la table des lachers
  29. C
  30. C     S   MCHELM CHAMELEM resultat
  31. C     S   MMODEL modele sur lequel s'appuie MCHELM
  32. C
  33. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  34. C
  35. C
  36.       IMPLICIT INTEGER(I-N)
  37.       IMPLICIT REAL*8 (A-H,O-Z)
  38.       PARAMETER(MPNO9=9)
  39. C              MPNO9 NOMBRE MAX DE NOEUDS POUR LES ELEMENTS CONSIDERES
  40. C
  41. -INC SMCOORD
  42. -INC SMELEME
  43. -INC SMTABLE
  44. -INC SMMODEL
  45. -INC SMCHAML
  46. -INC SMCHPOI
  47.  
  48. -INC PPARAM
  49. -INC CCOPTIO
  50. -INC CCREEL
  51. -INC SMLREEL
  52.       POINTEUR MLREE4.MLREEL,MLREE5.MLREEL,MLREE6.MLREEL
  53. ***   pteurs sur les MAILLAGES
  54.       POINTEUR IZCENT.MELEME,IELTFA.MELEME,IZLAC.MELEME,IFACEL.MELEME
  55. ***   variable de travail pour sp 'ACCTAB' % table de lacher 'MTAB2'
  56.       CHARACTER*20 MTYPI,CHARI,MTYPR,CHARR
  57.       LOGICAL LOGRE
  58. ***   'IZPART' segment contenant position ref des particules au départ du lacher :
  59. ***   'NPART' nbre de particules lachées
  60. ***   'NLEPA' = n° global de l'element contenant ieme particule lachée,
  61. ***   'NUMPA' = n° ieme particule lachée, 'COORPA' = coords ieme particule lachée
  62.       SEGMENT IZPART
  63.          INTEGER NLEPA(NPART),NUMPA(NPART)
  64.          REAL*8  COORPA(NDIM,NPART)
  65.       ENDSEGMENT
  66.       POINTEUR IZN.IZPART,IZL.IZPART
  67. ***   'IZCOU' segment content pas de tps des mailles traversées (selon critere convectif)
  68.       SEGMENT  IZCOU
  69.          REAL*8 DTCO(NEL),COU
  70.       ENDSEGMENT
  71. ***   'IZSH' segment de travail (pour passage réel -> réf)
  72. ***   'SHP' fonctions de forme et leurs dérivées
  73. ***   'SHY' fonctions de base et leurs dérivées
  74. ***   'XYZL(i, j)' = ieme coord reelle du jeme noeud d'1 elem donne
  75.       SEGMENT IZSH
  76.          REAL*8 SHP(6,MNO9),SHY(12,MNO9),XYZL(3,MNO9)
  77.       ENDSEGMENT
  78. ***   'IZVIT' segment décrivant les vitesses ou flux
  79. ***   'TEMTRA' = tableau des divers tps de palier (cas transitoire)
  80. ***   'NVIPT' = nbre de pas de tps (chgt) (cas transitoire)
  81. ***   'NBS' = nbre de sous-maillages
  82. ***   'IPUN' = pteurs sur les 'IZUN'
  83. ***   'IDUN(I)' = nbre d'éléments avant le sous-maillage I
  84. ***   'IPVPT(J)' = pteur sur les 'IZVPT' pour le pas de tps n° J (cas transitoire)
  85.       SEGMENT IZVIT
  86.          REAL*8 TEMTRA(NVIPT)
  87.          INTEGER IPUN(NBS),IDUN(NBS),IPVPT(NVIPT),IFORML
  88.       ENDSEGMENT
  89. ***   'IPARPO' segment permettant stockage resultats pdt 'avancee particule
  90. ***   'NPOS' = nbre d'elements traversés, 'NAPAR' = n° des élémts traversés
  91. ***   'NUMP' = inutile, 'CREF' = coord sauvegardées de chaque elemt traversee
  92. ***   'TPAR' = tps sauvegardé associé aux coords sauvegardées
  93.       SEGMENT  IPARPO
  94.          INTEGER NAPAR(NPOS),NUMP(NPOS)
  95.          REAL*8 CREF(NDIM,NPOS),TPAR(NPOS)
  96.       ENDSEGMENT
  97.       POINTEUR IZN3.IPARPO
  98. ***   tableau contenant tps et coordonnées de sortie (et non sauvegarde)
  99.       DIMENSION SORTIE(4)
  100.  
  101. ****************************************************************
  102. C
  103. C
  104.       SEGACT MTAB2
  105.       SEGACT IZCENT
  106.       NEL=IZCENT.NUM(/2)
  107.       SEGINI IZCOU
  108.       CALL INITD(DTCO,NEL,0.D0)
  109.       MNO9=MPNO9
  110.       SEGINI IZSH
  111. C
  112. C   RECUPERATION D'ELEMENTS PROVENANT DE LA TABLE DE LACHER
  113.       IVALI=1
  114.       XVALI=0.D0
  115.       IRETI=0
  116.       IVALR=0
  117.       XVALR=0.D0
  118.       IRETR=0
  119.       MTYPI='MOT     '
  120. C     CHARI='CFL'
  121. C     MTYPR='FLOTTANT'
  122.       MTYPR='        '
  123.       CHARR='        '
  124.       COU=0.05D0
  125.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'CFL',.TRUE.,IRETI,
  126.      *            MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  127.       IF(IERR.NE.0)RETURN
  128.       IF(MTYPR.NE.'        ')THEN
  129.           IF(MTYPR.NE.'FLOTTANT')THEN
  130.           MOTERR(1:11)='CFL        '
  131.           MOTERR(12:20)='FLOTTANT'
  132.           CALL ERREUR(627)
  133.           RETURN
  134.           ENDIF
  135.           COU=XVALR
  136.           IF(COU.LT.1.D-8)THEN
  137.               MOTERR(1:8)='CFL     '
  138.               REAERR(1)=REAL(COU)
  139.               REAERR(2)=1.E-8
  140.               CALL ERREUR(41)
  141.               RETURN
  142.           ENDIF
  143.       ENDIF
  144.       IVALI=1
  145.       XVALI=0.D0
  146.       IRETI=0
  147.       IVALR=0
  148.       XVALR=0.D0
  149.       IRETR=0
  150.       MTYPI='MOT     '
  151. C     MTYPR='FLOTTANT'
  152.       MTYPR='        '
  153.       CHARR='        '
  154.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'DELTAT_SAUVE',.TRUE.,
  155.      *            IRETI,MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  156.       IF(MTYPR.NE.'        ')THEN
  157.              IF(MTYPR.NE.'FLOTTANT')THEN
  158.              MOTERR(1:11)='DELTAT_SAUV'
  159.              MOTERR(12:20)='FLOTTANT'
  160.              CALL ERREUR(627)
  161.              RETURN
  162.              ENDIF
  163.       ENDIF
  164.       IF(IERR.NE.0)RETURN
  165.       DELTAT=XVALR
  166.       IF ((IICAL.EQ.3).AND.(DELTAT.GT.0.D0)) THEN
  167.          JG=1
  168.          SEGINI MLREE6
  169.          MLREE6.PROG(1)=DELTAT
  170.          NSAUV=0
  171.       ENDIF
  172.       IVALI=1
  173.       XVALI=0.D0
  174.       IRETI=0
  175.       IVALR=0
  176.       XVALR=0.D0
  177.       IRETR=0
  178.       MTYPI='MOT     '
  179.       MTYPR='FLOTTANT'
  180.       CHARR='        '
  181.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'TEMPS_LIMITE',.TRUE.,
  182.      *            IRETI,MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  183.       IF(IERR.NE.0)RETURN
  184.       TMAX=XVALR
  185.       CALL TRJCN7(TMAX,IZVIT)
  186.       IF(IERR.NE.0)RETURN
  187.  
  188. C     LECTURE DU TABLEAU DES TEMPS DE LACHER
  189.       IVALI=1
  190.       XVALI=0.D0
  191.       IRETI=0
  192.       IVALR=0
  193.       XVALR=0.D0
  194.       IRETR=0
  195.       MTYPI='MOT     '
  196. C     CHARI='BORNES'
  197.       MTYPR='LISTREEL'
  198.       CHARR='        '
  199.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'TEMPS_LACHER',.TRUE.,IRETI,
  200.      *            MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  201.       IF(IERR.NE.0)RETURN
  202.       MLREE3=IRETR
  203.       SEGACT MLREE3
  204.       NLAC=MLREE3.PROG(/1)
  205.  
  206. C     DONNEES SUPLEMENTAIRES POUR CONVECTION_DIFFUSION
  207.       IVALI=1
  208.       XVALI=0.D0
  209.       IRETI=0
  210.       IVALR=0
  211.       XVALR=0.D0
  212.       IRETR=0
  213.       MTYPI='MOT     '
  214.       MTYPR='        '
  215.       CHARR='        '
  216.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'IMPERMEABLE',.TRUE.,IRETI,
  217.      *            MTYPR,IVALR,XVALR,CHARR,LOGRE,IRETR)
  218.       IF(MTYPR.NE.'        ')THEN
  219.            IF(MTYPR.NE.'MAILLAGE')THEN
  220.                 MOTERR(1:11)='IMPERMEABLE'
  221.                 MOTERR(12:20)='MAILLAGE'
  222.                 CALL ERREUR(627)
  223.                 RETURN
  224.            ENDIF
  225.       ENDIF
  226.       IF(IERR.NE.0)RETURN
  227.       IZFACE=IRETR
  228.       IVALI=1
  229.       XVALI=0.D0
  230.       IRETI=0
  231.       IVALR=0
  232.       XVALR=0.D0
  233.       IRETR=0
  234.       MTYPI='MOT         '
  235.       MTYPR='            '
  236.       CHARR='            '
  237. ***   recuperation de la liste des tps de sauvegarde des resultats
  238.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,'TEMPS_SAUVES',.TRUE.,
  239.      $           IRETI,MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  240.       IF (MTYPR.NE.'        ') THEN
  241.          IF (MTYPR.NE.'LISTREEL') THEN
  242.             MOTERR(1:12)='TEMPS_SAUVES'
  243.             MOTERR(13:20)='LISTREEL'
  244.             CALL ERREUR(627)
  245.             RETURN
  246.          ENDIF
  247.       ENDIF
  248.       IF (IERR.NE.0) RETURN
  249.       IZSAUV=IRETR
  250.       IF (IZSAUV.GT.0) THEN
  251.          IF (DELTAT.GT.0.D0) THEN
  252.             MOTERR(1:8)= 'DELTAT_S'
  253.             MOTERR(9:16)='TEMPS_SA'
  254.             CALL ERREUR(135)
  255.             RETURN
  256.          ELSE
  257.             MLREE6=IZSAUV
  258.             SEGACT MLREE6
  259.             NSAUV=MLREE6.PROG(/1)
  260.          ENDIF
  261.       ELSE
  262.          IF (DELTAT.LE.0.D0) NSAUV=-1
  263.       ENDIF
  264.  
  265. C     ON ACTIVE LES TABLEAUX
  266. ***   activation du chamelem de dispersion si nécessaire
  267.       IF (IZDISP.GT.0) THEN
  268.          MCHEL1=IZDISP
  269.          SEGACT MCHEL1
  270.          MCHAML=MCHEL1.ICHAML(1)
  271.          SEGACT MCHAML
  272.          MELVA1=IELVAL(1)
  273.          SEGACT MELVA1
  274.          MELVA2=IELVAL(2)
  275.          SEGACT MELVA2
  276.       ENDIF
  277. ***   activation du chamelem de diffusivité si nécessaire
  278.       IF (IZDIFF.GT.0) THEN
  279.          MCHEL1=IZDIFF
  280.          SEGACT MCHEL1
  281.          MCHAML=MCHEL1.ICHAML(1)
  282.          SEGACT MCHAML
  283.          MELVA3=IELVAL(1)
  284.          SEGACT MELVA3
  285.       ENDIF
  286. ***   activation du maillage des faces impermeables si nécessaire
  287.       IF (IZFACE.GT.0) THEN
  288.          IPT9=IZFACE
  289.          SEGACT IPT9
  290.       ENDIF
  291.  
  292. C     PRISE EN COMPTE DE LA POROSITE
  293. C
  294.       IF(IZPOR.NE.0)THEN
  295.          CALL TRJPOR(IZPOR,IZVIT,MELEME)
  296.          IF(IERR.NE.0)RETURN
  297.       ENDIF
  298. C
  299. C     STRUCTURES DE RESULTATS
  300. C
  301.       L1=5
  302.       N1=1
  303.       N3=6
  304.       SEGINI MCHELM
  305.       TITCHE(1:5)='NOEUD'
  306.       SEGINI MMODEL
  307.  
  308. C     ------------------------------
  309. C     BOUCLE SUR LES TEMPS DE LACHER
  310. C     ------------------------------
  311.  
  312.       EPSILO=10*SQRT(XZPREC)
  313.       NPART0=0
  314.       DO 80 KLAC=1,NLAC
  315.  
  316.       TMIN=MLREE3.PROG(KLAC)
  317.       CALL TRJCN7(TMIN,IZVIT)
  318.       IF(ABS(DELTAT).GT.1.D-15)CALL TRJCN4(TMIN,TMAX,DELTAT)
  319.  
  320. C     LECTURE D UN MAILLAGE DE LACHER
  321.       IVALI=KLAC
  322.       XVALI=0.D0
  323.       IRETI=0
  324.       IVALR=0
  325.       XVALR=0.D0
  326.       IRETR=0
  327.       MTYPI='ENTIER  '
  328.       CHARI='        '
  329.       MTYPR='MAILLAGE'
  330.       CHARR='MAILLAGE'
  331.       CALL ACCTAB(MTAB2,MTYPI,IVALI,XVALI,CHARI,.TRUE.,IRETI,
  332.      *            MTYPR,IVALR,XVALR,CHARR,.TRUE.,IRETR)
  333.       IF(IERR.NE.0)RETURN
  334.       IZLAC=IRETR
  335.       CALL TRJCN2(IZLAC)
  336.       IF(IERR.NE.0)RETURN
  337. C
  338.       CALL TRJINI(IZPART,IZCOU,IZLAC,MELEME,IZVIT,IZCENT,IELTFA,IZSH,
  339.      *            TMIN)
  340.       NPARTI=NUMPA(/1)
  341.       N1=NPART0+NPARTI
  342.       SEGADJ MCHELM
  343.       SEGADJ MMODEL
  344.  
  345. C     ........................................
  346. C
  347. C     BOUCLE SUR LES PARTICULES DU LACHER KLAC
  348. C     ........................................
  349.  
  350.       SEGACT,MCOORD
  351.       DO 50 IPART=1,NPARTI
  352.       IIPOS=0
  353.       IREP=0
  354.       NDIM=COORPA(/1)
  355.       TMINI=TMIN
  356.       DT1=0.D0
  357.       IDD1=1
  358.    53 CONTINUE
  359.  
  360.       IF(IICAL.EQ.1) THEN
  361.           CALL TRJAVA(IZVIT,IZPART,IZN3,IPART,TMINI,TMAX,IZCOU,MELEME,
  362.      *                IELTFA,IZCENT,IFACEL,DELTAT,IZSH,IEROR,DT1)
  363.           IF(IERR.NE.0)RETURN
  364.           CALL TRJRAL(IZN3,IPARPO,MELEME,IZSH)
  365.           SEGSUP IZN3
  366.       ELSEIF(IICAL.EQ.2) THEN
  367. C         convection seule avec méthode analytique
  368.           CALL TRJAV2(IZVIT,IZPART,IZN3,IPART,TMINI,TMAX,MELEME,
  369.      *                IELTFA,IZCENT,IFACEL,DELTAT,IZSH)
  370.           IF(IERR.NE.0)RETURN
  371.           CALL TRJRAL(IZN3,IPARPO,MELEME,IZSH)
  372.           SEGSUP IZN3
  373.       ELSE
  374. C**       convection + dispersion + diffusion
  375.           CALL TRJAV3(EPSILO,NSAUV,MLREE6,MELVA1,MELVA2,
  376.      $                MELVA3,IPT9,IZVIT,IZPART,IPART,TMINI,TMAX,MELEME,
  377.      $                IELTFA,IZCENT,IFACEL,IPARPO,IZSH,SORTIE)
  378.           IF (IERR.NE.0) RETURN
  379. C**       cas ou particule i est eliminee du probleme
  380.           IF (SORTIE(1).LT.0.D0) THEN
  381.               NPART0=NPART0-1
  382.               N1=NPART0+NPARTI
  383.               SEGADJ MCHELM
  384.               SEGADJ MMODEL
  385.               GOTO 50
  386.           ENDIF
  387.       ENDIF
  388.  
  389.       NPOS=NAPAR(/1)
  390.       IF(IIMPI.GT.0)THEN
  391.             WRITE(6,*)' PARTICULE ',IPART,NPOS
  392.             WRITE(6,*)NAPAR(1),NUMP(1),(CREF(J,1),J=1,NDIM),TPAR(1)
  393.             I2=NPOS
  394.             IF(IIMPI.GT.1)I2=2
  395.             DO 10 I=I2,NPOS
  396.                WRITE(6,*)NAPAR(I),NUMP(I),(CREF(J,I),J=1,NDIM),TPAR(I)
  397.    10       CONTINUE
  398.       ENDIF
  399.  
  400. C     CAS OU TOUS LES ELEMENTS N'ONT PAS LA MEME ORIENTATION
  401. C     ON MET LES RESULTATS DANS LE MEME SEGMENT
  402.       NPOS1=NPOS
  403.       IF(IREP.NE.0)THEN
  404.             IZN3=IREP
  405.             NPOS=IIPOS+NPOS+1-IDD1
  406.             SEGADJ IZN3
  407.             DO 54 II=IDD1,NPOS1
  408.                  II1=II+1-IDD1
  409.                  IZN3.CREF(1,IIPOS+II1)=CREF(1,II)
  410.                  IZN3.CREF(2,IIPOS+II1)=CREF(2,II)
  411.                  IZN3.TPAR(IIPOS+II1)=TPAR(II)
  412.    54       CONTINUE
  413.             IF(NDIM.EQ.3)THEN
  414.                  DO 55 II=1,NPOS1
  415.                       IZN3.CREF(3,IIPOS+II+1-IDD1)=CREF(3,II)
  416.    55            CONTINUE
  417.             ENDIF
  418.             SEGSUP IPARPO
  419.             IPARPO=IZN3
  420.       ENDIF
  421. C
  422.       IF(IEROR.EQ.1)THEN
  423.          TMINI=TPAR(NPOS)
  424.          TTEMP=TPAR(NPOS)
  425.          NPART=1
  426.          SEGINI IZN,IZL
  427.          DO 56 I=1,NDIM
  428.             IZL.COORPA(I,1)=CREF(I,NPOS)
  429.    56    CONTINUE
  430.  
  431.          IZN.NLEPA(1)=0
  432.          IZN.NUMPA(1)=1
  433.          IZL.NLEPA(1)=0
  434.          IZL.NUMPA(1)=1
  435.          IIPOS=IIPOS+NPOS1-IDD1
  436.          CALL TRJPEL(IZL,IZN,MELEME,IZVIT,IZCOU,IZCENT,IELTFA,IZSH,
  437.      &               TTEMP)
  438.          IF(IZN.NLEPA(1).NE.0)THEN
  439.             IREP=IPARPO
  440.             DO 57 I=1,NDIM
  441.                COORPA(I,IPART)=IZN.COORPA(I,1)
  442.    57       CONTINUE
  443.             NLEPA(IPART)=IZN.NLEPA(1)
  444.             DT1=0.D0
  445.             IDD1=1
  446. C           CAS D'UN PAS CONSTANT, ON RECALE LE CALCUL AVEC DT1
  447.             IF((ABS(DELTAT).GT.1.D-15).AND.
  448.      *        ((TPAR(NPOS)-TPAR(NPOS-1)).NE.DELTAT))THEN
  449.                 DT1=TPAR(NPOS)-TPAR(NPOS-1)
  450.                 IDD1=2
  451.             ENDIF
  452.             SEGSUP IZL,IZN
  453.             GO TO 53
  454.          ENDIF
  455.          SEGSUP IZL,IZN
  456.       ENDIF
  457. C
  458.       IVALI=IPART+NPART0
  459.       CALL TRJCHA(IPARPO,MCHELM,MMODEL,IVALI)
  460.       SEGSUP IPARPO
  461.    50 CONTINUE
  462.  
  463. C     ............................................
  464. C
  465. C     FIN BOUCLE SUR LES PARTICULES DU LACHER KLAC
  466. C     ............................................
  467.  
  468.       NPART0=NPART0+NPARTI
  469.       SEGSUP IZPART
  470.    80 CONTINUE
  471.  
  472. C     ----------------------------------
  473. C     FIN BOUCLE SUR LES TEMPS DE LACHER
  474. C     ----------------------------------
  475.  
  476.       SEGDES MTAB2
  477.       SEGSUP IZCOU
  478.       SEGSUP IZSH
  479. C
  480. C
  481.       RETURN
  482.       END
  483.  
  484.  
  485.  
  486.  
  487.  
  488.  
  489.  
  490.  

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