Télécharger la plaquette Cast3M


fr en

Télécharger yctscl.eso

Retour à la liste

Numérotation des lignes : 

yctscl

  1. C YCTSCL    SOURCE    CB215821  26/08/24    21:18:59     12622          
  2.       SUBROUTINE YCTSCL
  3. C
  4. C     VERSION VECTORISEE
  5. C
  6. C     Les éléments sont groupés en paquets de LRV éléments, LRV étant
  7. C     la longueur des registres vectoriels de la machine cible, i.e
  8. C     64 sur Cray, 128 ou 256 sur IBM 3090VF. On promène une fenêtre
  9. C     de longueur LRV sur la boucle générale de longueur NEL.
  10. C
  11.      & (HR,RPG,DRR,LE,NEL,K0,IES,NP,IAXI,
  12.      & IPADI,IPADS,IPADF,IKOMP,IKAS,
  13.      & ALFE,IND1,UN,INDU,NPTS,TN,NPTD,QE,IKS,
  14.      & HRN,G,NPTI,
  15.      & ALT,SGT,
  16.      & VOLU,COTE,NELZ,IDCEN,IPG,
  17.      & DTM1,DT,DTT1,DTT2,NUEL,DIAEL,FN)
  18.  
  19.       IMPLICIT INTEGER(I-N)
  20.       IMPLICIT REAL*8 (A-H,O-Z)
  21.  
  22. C***********************************************************************
  23. C
  24. C  CE SP DISCRETISE UNE EQUATION GENERALE DE TRANSPORT-DIFFUSION AVEC
  25. C  SOURCE.
  26. C    EN 2D SUR LES ELEMENTS QUA4 ET TRI3      PLAN OU AXI
  27. C    EN 3D SUR LES ELEMENTS CUB8 ET PRI6
  28. C  LES OPERATEURS SONT "SOUS-INTEGRES"
  29. C
  30. C
  31. C                        APPELE PAR YTSCAL
  32. C
  33. C
  34. C***********************************************************************
  35.  
  36. -INC CCVQUA4
  37. -INC CCREEL
  38. C
  39. C Longueur des registres vectoriels de la machine cible
  40. C On prend 64 pour ne pas augmenter la taille des tableaux
  41. C nécessaires à la vectorisation.
  42. C
  43.       PARAMETER(LRV=64)
  44.  
  45.       DIMENSION UN(NPTS,IES),HRN(NPTI),TN(NPTD)
  46.       DIMENSION COTE(NELZ,IES),VOLU(NELZ),QE(*)
  47.       DIMENSION ALFE(*),ALT(*)
  48.  
  49.       DIMENSION IPADI(*),IPADS(*),IPADF(*),LE(NP,1)
  50.       DIMENSION HR(NEL,NP,IES),RPG(1),DRR(NP,NEL)
  51.  
  52.       DIMENSION BF(9,9)
  53.       DIMENSION QGGT(8,8),Q1(8,8),Q2(8,8),Q3(8,8)
  54.  
  55.       DIMENSION AIRE(LRV)
  56.       DIMENSION AL(LRV),AH(LRV),AP(LRV)
  57.       DIMENSION ALFT(LRV),QT(LRV)
  58.       DIMENSION UIX(LRV,9),UIY(LRV,9),UIZ(LRV,9)
  59.       DIMENSION TETAC(LRV,9),TETAD(LRV,9),TETA(LRV,9)
  60.       DIMENSION UMI(LRV,3)
  61.       DIMENSION SBF(LRV,9)
  62.       DIMENSION WT(LRV,9),CHGLD(LRV),CHGLP(LRV)
  63.  
  64.       REAL*8    G(NPTS),FN(NP,*)
  65.  
  66.       SAVE IPAS,QGGT,Q1,Q2,Q3
  67.  
  68.       DATA CD/1.D0/
  69.  
  70.       DATA IPAS/0/
  71. C************************************************************************
  72. C
  73. C           INITIALISATIONS DIVERSES
  74. C
  75.        ZERMA=XPETIT
  76.  
  77.       NK=K0
  78. C                           ********
  79. C                           *  2D  *
  80. C                           ********
  81.  
  82.  
  83.       IF(IES.EQ.3)GO TO 10
  84.  
  85.       IAX1=0
  86.       IAX2=0
  87.       IF(IAXI.EQ.1)IAX2=1
  88.       IF(IAXI.EQ.2)IAX1=1
  89.  
  90. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  91.  
  92. C                           DIFFERENCES TRIANGLE / QUADRANGLE
  93.       QUA4=0.D0
  94.       IF(NP.EQ.4)QUA4=1.D0
  95. C
  96. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  97.  
  98. C
  99. C Calcul du nombre de paquets de LRV éléments
  100. C
  101.       NNN=MOD(NEL,LRV)
  102.       IF(NNN.EQ.0) NPACK=NEL/LRV
  103.       IF(NNN.NE.0) NPACK=1+(NEL-NNN)/LRV
  104.       KPACKD=1
  105.       KPACKF=NPACK
  106. C
  107. C ******* BOUCLE SUR LES PAQUETS DE LRV ELEMENTS **********
  108. C
  109.       DO 7001 KPACK=KPACKD,KPACKF
  110. C
  111. C ======= A L'INTERIEUR DE CHAQUE PAQUET DE LRV ELEMENTS =======
  112. C
  113. C 1. Calcul des limites du paquet courant.
  114.       KDEB=1+(KPACK-1)*LRV
  115.       KFIN=MIN(NEL,KDEB+LRV-1)
  116. C
  117.       DO 7002 K=KDEB,KFIN
  118.       KP=K-KDEB+1
  119.       NK=K+K0
  120.       NK1=(1-IND1)*(NK-1)+1
  121.       ALFT(KP)=ALFE(NK1)+ZERMA
  122.       AIRE(KP)=VOLU(NK)
  123.       AL(KP)=COTE(NK,1)
  124.       AH(KP)=COTE(NK,2)
  125.  7002 CONTINUE
  126.  
  127.       IF((IKOMP.EQ.0.AND.IKAS.EQ.5).OR.
  128.      &(IKOMP.EQ.1.AND.IKAS.EQ.6))THEN
  129.       DO 7003 K=KDEB,KFIN
  130.       KP=K-KDEB+1
  131.       NK=K+K0
  132.       ALFT(KP)=ALFT(KP)+ALT(NK)/SGT
  133.  7003 CONTINUE
  134.       ENDIF
  135.  
  136. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  137.  
  138.       DO 81065 I=1,NP
  139.       DO 7006 K=KDEB,KFIN
  140.       KP=K-KDEB+1
  141.       NF=IPADF(LE(I,K))
  142.       NU=IPADS(LE(I,K))
  143.       NI=IPADI(LE(I,K))
  144.       NFU=(1-INDU)*(NU-1)+1
  145.       UIX(KP,I)=UN(NFU,1)
  146.       UIY(KP,I)=UN(NFU,2)
  147.       TETAC(KP,I)=HRN(NI)
  148.       TETAD(KP,I)=TN(NF)
  149.  7006 CONTINUE
  150. 81065 CONTINUE
  151.  
  152.       CALL KSUPG1(WT,UMI,CHGLD,CHGLP,KDEB,KFIN,LRV,
  153.      &HRN,IPADI,UN,ALFT,NPTS,NEL,NP,DRR,HR,FN,
  154.      &AIRE,AL,AH,AP,IDCEN,IPADS,LE,QUA4,IKOMP,
  155.      &DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  156. C
  157. C Initialisation de la variable d'accumulation SBF au terme source
  158. C
  159.  
  160.       IF(IKOMP.EQ.0)THEN
  161.  
  162.          DO 70021 K=KDEB,KFIN
  163.          KP=K-KDEB+1
  164.          NK=K+K0
  165.          NKS=(1-IKS)*(NK-1)+1
  166.          QT(KP)=QE(NKS)
  167. 70021    CONTINUE
  168.  
  169.          IF(IPG.EQ.0)THEN
  170.          DO 81066 I=1,NP
  171.          DO 70062 K=KDEB,KFIN
  172.          KP=K-KDEB+1
  173.          SBF(KP,I)=-QT(KP)*DRR(I,K)
  174. 70062    CONTINUE
  175. 81066    CONTINUE
  176.          ELSE
  177.          DO 81067 I=1,NP
  178.          DO 71062 K=KDEB,KFIN
  179.          KP=K-KDEB+1
  180.          SBF(KP,I)=-QT(KP)*WT(KP,I)
  181. 71062    CONTINUE
  182. 81067    CONTINUE
  183.          ENDIF
  184.  
  185.       ELSEIF(IKOMP.EQ.1)THEN
  186.  
  187.          DO 70023 K=KDEB,KFIN
  188.          KP=K-KDEB+1
  189.          NK=K+K0
  190.          NKS=(1-IKS)*(NK-1)+1
  191.          QT(KP)=QE(NKS)
  192. 70023    CONTINUE
  193.  
  194.        IF(IPG.EQ.0)THEN
  195.          DO 81068 I=1,NP
  196.          DO 70064 K=KDEB,KFIN
  197.          KP=K-KDEB+1
  198.          SBF(KP,I)=-QT(KP)*DRR(I,K)
  199. 70064    CONTINUE
  200. 81068    CONTINUE
  201.        ELSE
  202.          DO 81069 I=1,NP
  203.          DO 71064 K=KDEB,KFIN
  204.          KP=K-KDEB+1
  205.          SBF(KP,I)=-QT(KP)*WT(KP,I)
  206. 71064    CONTINUE
  207. 81069    CONTINUE
  208.        ENDIF
  209.  
  210.       ENDIF
  211.  
  212. C Le coeur du calcul ...
  213.  
  214.       IF(IKOMP.EQ.0)THEN
  215.  
  216.       DO 81071 I=1,NP
  217.       DO 81070 J= 1,NP
  218.       DO 7014 K=KDEB,KFIN
  219.       KP=K-KDEB+1
  220.       ZVGG=AIRE(KP)*CHGLP(KP)*VGGT(J,I)
  221.  
  222.       ZVGT=AIRE(KP)*(
  223.      &     HR(K,I,1)*HR(K,J,1)*ALFT(KP)
  224.      &+    HR(K,I,2)*HR(K,J,2)*ALFT(KP)
  225.      &+    CHGLD(KP)*VGGT(J,I) )
  226.  
  227.       V2=(UMI(KP,1)*HR(K,J,1)+UMI(KP,2)*HR(K,J,2))*WT(KP,I)
  228.  
  229.       SBF(KP,I)=SBF(KP,I)+TETAC(KP,J)*(ZVGG+V2)+ TETAD(KP,J)*ZVGT
  230.  
  231.  7014 CONTINUE
  232. 81070 CONTINUE
  233. 81071 CONTINUE
  234.  
  235.       ELSEIF(IKOMP.EQ.1)THEN
  236.  
  237.       DO 81073 I=1,NP
  238.       DO 81072 J= 1,NP
  239.       DO 7015 K=KDEB,KFIN
  240.       KP=K-KDEB+1
  241.       ZVGG=AIRE(KP)*CHGLP(KP)*VGGT(J,I)
  242.  
  243.       ZVGT=AIRE(KP)*(
  244.      &     HR(K,I,1)*HR(K,J,1)*ALFT(KP)
  245.      &+    HR(K,I,2)*HR(K,J,2)*ALFT(KP)
  246.      &+    CHGLD(KP)*VGGT(J,I) )
  247.  
  248.       V2=(UIX(KP,J)*HR(K,J,1)+UIY(KP,J)*HR(K,J,2))*WT(KP,I)
  249.  
  250.       SBF(KP,I)=SBF(KP,I)+TETAC(KP,J)*(ZVGG+V2)+ TETAD(KP,J)*ZVGT
  251.  
  252.  7015 CONTINUE
  253. 81072 CONTINUE
  254. 81073 CONTINUE
  255.  
  256.       ENDIF
  257. C
  258. C     Fin de l'accumulation dans SBF.
  259. C     On ajoute ces incréments G.
  260. C
  261.       DO 81074 I=1,NP
  262.       DO 7017 K=KDEB,KFIN
  263.       KP=K-KDEB+1
  264.       NF=IPADS(LE(I,K))
  265.       G(NF) = G(NF)-SBF(KP,I)
  266.  7017 CONTINUE
  267. 81074 CONTINUE
  268.  
  269.  7001 CONTINUE
  270.  
  271. C     WRITE(6,*)' G DANS YCTSCL '
  272. C     WRITE(6,1984)(M,G(M),M=1,NPTS)
  273.  1984 FORMAT(7(1X,I4,2X,1PE11.4))
  274.  
  275. C     CALL ARRET(0)
  276.       IPAS=1
  277.       RETURN
  278.  
  279. C                           ********
  280. C                           *  3D  *
  281. C                           ********
  282.  
  283.  10   CONTINUE
  284.  
  285.       IF(IPAS.EQ.0)CALL CALHRH(QGGT,Q1,Q2,Q3,IES)
  286.       CUB8=0.D0
  287.       IF(NP.EQ.8)CUB8=1.D0
  288.  
  289. C
  290. C Calcul du nombre de paquets de LRV éléments
  291. C
  292.       NNN=MOD(NEL,LRV)
  293.       IF(NNN.EQ.0) NPACK=NEL/LRV
  294.       IF(NNN.NE.0) NPACK=1+(NEL-NNN)/LRV
  295.       KPACKD=1
  296.       KPACKF=NPACK
  297. C
  298. C ******* BOUCLE SUR LES PAQUETS DE LRV ELEMENTS **********
  299. C
  300.       DO 8001 KPACK=KPACKD,KPACKF
  301. C
  302. C ======= A L'INTERIEUR DE CHAQUE PAQUET DE LRV ELEMENTS =======
  303. C
  304. C 1. Calcul des limites du paquet courant.
  305.       KDEB=1+(KPACK-1)*LRV
  306.       KFIN=MIN(NEL,KDEB+LRV-1)
  307. C
  308.       DO 8002 K=KDEB,KFIN
  309.       KP=K-KDEB+1
  310.       NK=K+K0
  311.       NK1=(1-IND1)*(NK-1)+1
  312.       ALFT(KP)=ALFE(NK1)+ZERMA
  313.       AIRE(KP)=VOLU(NK)
  314.       AL(KP)=COTE(NK,1)
  315.       AH(KP)=COTE(NK,2)
  316.       AP(KP)=COTE(NK,3)
  317.  8002 CONTINUE
  318.       IF((IKOMP.EQ.0.AND.IKAS.EQ.5).OR.
  319.      &(IKOMP.EQ.1.AND.IKAS.EQ.6))THEN
  320.       DO 8003 K=KDEB,KFIN
  321.       KP=K-KDEB+1
  322.       NK=K+K0
  323.       ALFT(KP)=ALFT(KP)+ALT(NK)/SGT
  324.  8003 CONTINUE
  325.       ENDIF
  326.  
  327. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  328.  
  329.       DO 81075 I=1,NP
  330.       DO 8006 K=KDEB,KFIN
  331.       KP=K-KDEB+1
  332.       NF=IPADF(LE(I,K))
  333.       NU=IPADS(LE(I,K))
  334.       NI=IPADI(LE(I,K))
  335.       NFU=(1-INDU)*(NU-1)+1
  336.       UIX(KP,I)=UN(NFU,1)
  337.       UIY(KP,I)=UN(NFU,2)
  338.       UIZ(KP,I)=UN(NFU,3)
  339.       TETAC(KP,I)=HRN(NI)
  340.       TETAD(KP,I)=TN(NF)
  341.  8006 CONTINUE
  342. 81075 CONTINUE
  343.  
  344.       CALL KSUPG1(WT,UMI,CHGLD,CHGLP,KDEB,KFIN,LRV,
  345.      &HRN,IPADI,UN,ALFT,NPTS,NEL,NP,DRR,HR,FN,
  346.      &AIRE,AL,AH,AP,IDCEN,IPADS,LE,CUB8,IKOMP,
  347.      &DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  348.  
  349. C
  350. C Initialisation de la variable d'accumulation SBF au terme source
  351. C                                                                      M
  352.       IF(IKOMP.EQ.0)THEN
  353.  
  354.          DO 80021 K=KDEB,KFIN
  355.          KP=K-KDEB+1
  356.          NK=K+K0
  357.          NKS=(1-IKS)*(NK-1)+1
  358.          QT(KP)=QE(NKS)
  359. 80021    CONTINUE
  360.  
  361.        IF(IPG.EQ.0)THEN
  362.          DO 81076 I=1,NP
  363.          DO 80062 K=KDEB,KFIN
  364.          KP=K-KDEB+1
  365.          SBF(KP,I)=-QT(KP)*DRR(I,K)
  366. 80062    CONTINUE
  367. 81076    CONTINUE
  368.        ELSE
  369.          DO 81077 I=1,NP
  370.          DO 81062 K=KDEB,KFIN
  371.          KP=K-KDEB+1
  372.          SBF(KP,I)=-QT(KP)*WT(KP,I)
  373. 81062    CONTINUE
  374. 81077    CONTINUE
  375.        ENDIF
  376.  
  377.       ELSEIF(IKOMP.EQ.1)THEN
  378.  
  379.          DO 80023 K=KDEB,KFIN
  380.          KP=K-KDEB+1
  381.          NK=K+K0
  382.          NKS=(1-IKS)*(NK-1)+1
  383.          QT(KP)=QE(NKS)
  384. 80023    CONTINUE
  385.  
  386.        IF(IPG.EQ.0)THEN
  387.          DO 81078 I=1,NP
  388.          DO 80064 K=KDEB,KFIN
  389.          KP=K-KDEB+1
  390.          SBF(KP,I)=-QT(KP)*DRR(I,K)
  391. 80064    CONTINUE
  392. 81078    CONTINUE
  393.        ELSE
  394.          DO 81079 I=1,NP
  395.          DO 81064 K=KDEB,KFIN
  396.          KP=K-KDEB+1
  397.          SBF(KP,I)=-QT(KP)*WT(KP,I)
  398. 81064    CONTINUE
  399. 81079    CONTINUE
  400.        ENDIF
  401.  
  402.       ENDIF
  403.  
  404. C Le coeur du calcul ...
  405.  
  406.       IF(IKOMP.EQ.0)THEN
  407.  
  408.       DO 81081 I=1,NP
  409.       DO 81080 J= 1,NP
  410.       DO 8014 K=KDEB,KFIN
  411.       KP=K-KDEB+1
  412.  
  413.       ZVGG=AIRE(KP)*CHGLP(KP)*QGGT(J,I)
  414.  
  415.       ZVGT=AIRE(KP)*(
  416.      &     HR(K,I,1)*HR(K,J,1)*ALFT(KP)
  417.      &+    HR(K,I,2)*HR(K,J,2)*ALFT(KP)
  418.      &+    HR(K,I,3)*HR(K,J,3)*ALFT(KP) )
  419.      &+    CHGLD(KP)*QGGT(J,I)
  420.  
  421.       V2=UMI(KP,1)*HR(K,J,1)+UMI(KP,2)*HR(K,J,2)+UMI(KP,3)*HR(K,J,3)
  422.  
  423.       SBF(KP,I)=SBF(KP,I)
  424.      &         +TETAC(KP,J)*(ZVGG+V2*WT(KP,I))+TETAD(KP,J)*ZVGT
  425.  
  426.  8014 CONTINUE
  427. 81080 CONTINUE
  428. 81081 CONTINUE
  429.  
  430.       ELSEIF(IKOMP.EQ.1)THEN
  431.  
  432.       DO 81083 I=1,NP
  433.       DO 81082 J= 1,NP
  434.       DO 8015 K=KDEB,KFIN
  435.  
  436.       KP=K-KDEB+1
  437.  
  438.       ZVGG=AIRE(KP)*CHGLP(KP)*QGGT(J,I)
  439.  
  440.       ZVGT=AIRE(KP)*(
  441.      &     HR(K,I,1)*HR(K,J,1)*ALFT(KP)
  442.      &+    HR(K,I,2)*HR(K,J,2)*ALFT(KP)
  443.      &+    HR(K,I,3)*HR(K,J,3)*ALFT(KP) )
  444.      &+    CHGLD(KP)*QGGT(J,I)
  445.  
  446.       V2=UIX(KP,J)*HR(K,J,1)+UIY(KP,J)*HR(K,J,2)+UIZ(KP,J)*HR(K,J,3)
  447.  
  448.       SBF(KP,I)=SBF(KP,I)
  449.      &         +TETAC(KP,J)*(ZVGG+V2*WT(KP,I))+TETAD(KP,J)*ZVGT
  450.  
  451.  8015 CONTINUE
  452. 81082 CONTINUE
  453. 81083 CONTINUE
  454.  
  455.       ENDIF
  456.  
  457. C
  458. C     Fin de l'accumulation dans SBF.
  459. C     On ajoute ces incréments G.
  460. C
  461.       DO 81084 I=1,NP
  462.       DO 8017 K=KDEB,KFIN
  463.       KP=K-KDEB+1
  464.       NF=IPADS(LE(I,K))
  465.       G(NF) = G(NF)-SBF(KP,I)
  466.  8017 CONTINUE
  467. 81084 CONTINUE
  468.  
  469.  8001 CONTINUE
  470.  
  471. C     WRITE(6,*)' G DANS YCTSCL '
  472. C     WRITE(6,1984)(M,G(M),M=1,NPTS)
  473.  
  474. C     CALL ARRET(0)
  475.       IPAS=1
  476.       RETURN
  477.  1002 FORMAT(10(1X,1PE11.4))
  478.       END
  479.  
  480.  
  481.  
  482.  
  483.  
  484.  
  485.  
  486.  
  487.  
  488.  

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