Télécharger la plaquette Cast3M


fr en

Télécharger zctscl.eso

Retour à la liste

Numérotation des lignes : 

zctscl

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

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