Télécharger la plaquette Cast3M


fr en

Télécharger ycvi.eso

Retour à la liste

Numérotation des lignes : 

ycvi

  1. C YCVI      SOURCE    CHAT      05/01/13    04:16:42     5004
  2.       SUBROUTINE YCVI
  3. C
  4. C
  5. C     VERSION VECTORISEE
  6. C
  7. C     Les {l{ments sont group{s en paquets de LRV {l{ments, LRV {tant
  8. C     la longueur des registres vectoriels de la machine cible, i.e
  9. C     64 sur Cray, 128 ou 256 sur IBM 3090VF. On prom}ne une fenetre
  10. C     de longueur LRV sur la boucle g{n{rale de longueur NEL.
  11. C
  12. C
  13.      &     (HR,RPG,DRR,LE,NEL,K0,NPTD,IES,NP,IAXI,IKOMP,IKAS,
  14.      &     COEFF,IK1,RGE,IKG,NELG,TN,IKT,TREF,IKREF,IPADS,
  15.      &     UN,IPADU,NPTU,GN,F,IPADI,VF,IPADF,NPTF,
  16.      &     VOLU,COTE,NELZ,IDCEN,IPG,
  17.      &     DTM1,DT,DTT1,DTT2,NUEL,DIAEL,FN)
  18.       IMPLICIT INTEGER(I-N)
  19.       IMPLICIT REAL*8 (A-H,O-Z)
  20. C***********************************************************************
  21. C
  22. C  CE SP DISCRETISE LES EQUATIONS DE NAVIER STOKES
  23. C    EN 2D SUR LES ELEMENTS QUA4 ET TRI3      PLAN OU AXI
  24. C    EN 3D SUR LES ELEMENTS CUB8 ET PRI6
  25. C  LES OPERATEURS SONT "SOUS-INTEGRES"
  26. C
  27. C  SYNTAXE :
  28. C
  29. C         NS(NU,UN,RGE,DE)  INCO GN :
  30. C
  31. C  COEFFICIENTS :
  32. C  --------------
  33. C
  34. C  UN(NPTD,IES)       CHAMPS DE VITESSE TRANSPORTANT
  35. C  COEFF(SCAL DOMA)  VISCOSITE CINEMATIQE MOLECULAIRE( NU )
  36. C       (SCAL ELEM)
  37. C   RGE(NELG,IES)       TERME SOURCE
  38. C
  39. C
  40. C  INCONNUES :
  41. C  -----------
  42. C
  43. C   UN(NPTD,IES)        CHAMPS DE VITESSE TRANSPORTANT
  44. C   GN(NPTD,IES)        CHAMPS DE VITESSE TRANSPORTE
  45. C   VN(NPTD,IES)        CHAMPS DE VITESSE DU FLUIDE
  46. C
  47. C
  48. C
  49. C
  50. C  TABLEAUX DE TRAVAIL :
  51. C  ---------------------        N     N-1
  52. C                     (D + D1) U - D U              N-1    T  N
  53. C                     -------------------  = F - A U    - C  P
  54. C                               DT
  55. C
  56. C                              N-1
  57. C   F(NPTD,IES)   CONTIENT   A U   - F     (VITESSE)
  58. C
  59. C
  60. C
  61. C***********************************************************************
  62.  
  63. -INC CCVQUA4
  64. -INC CCREEL
  65. *-
  66.  
  67. -INC PPARAM
  68. -INC CCOPTIO
  69. -INC SMCOORD
  70. C
  71. C Longueur des registres vectoriels de la machine cible
  72. C On prend 64 pour ne pas augmenter la taille des tableaux
  73. C n{cessaires @ la vectorisation.
  74. C
  75.       PARAMETER(LRV=64)
  76.       DIMENSION UN(NPTU,IES),GN(NPTD,IES),VF(NPTF,IES)
  77.       DIMENSION TN(*),TREF(*)
  78.       DIMENSION COEFF(*),RGE(NELG,IES)
  79.       DIMENSION COTE(NELZ,IES),VOLU(NELZ),KLIP(100)
  80.       DIMENSION IPADI(*),LE(NP,1),IPADU(*),IPADF(*),IPADS(*)
  81.       DIMENSION HR(NEL,NP,IES),RPG(1),DRR(NP,NEL)
  82.       DIMENSION QGGT(8,8),Q1(8,8),Q2(8,8),Q3(8,8)
  83.       DIMENSION COEF(LRV),AIRE(LRV)
  84.       DIMENSION WX(LRV,9),WY(LRV,9),WZ(LRV,9)
  85.       DIMENSION AL(LRV),AH(LRV),AP(LRV)
  86. C UIX,... vitesse transportante
  87.       DIMENSION UIX(LRV,9),UIY(LRV,9),UIZ(LRV,9)
  88. C GIX,... vitesse massique ou transportée ou inconnue du/dt
  89.       DIMENSION GIX(LRV,9),GIY(LRV,9),GIZ(LRV,9)
  90. C VIX,... vitesse du fluide
  91.       DIMENSION VIX(LRV,9),VIY(LRV,9),VIZ(LRV,9)
  92.       DIMENSION UMI(LRV,3),VMI(LRV,3)
  93.       DIMENSION COEFT(LRV),RGX(LRV),RGY(LRV),RGZ(LRV)
  94.       DIMENSION TO1(LRV),TO2X(LRV),TO2Y(LRV)
  95.       DIMENSION SAF1(LRV,9),SAF2(LRV,9),SAF3(LRV,9)
  96.       DIMENSION CHGLD(LRV),CHGLPX(LRV),CHGLPY(LRV),CHGLPZ(LRV)
  97.       DIMENSION F(NPTD,*),FN(NP,*)
  98.       SAVE IPAS,QGGT,Q1,Q2,Q3
  99.       DATA IPAS/0/
  100. C
  101. C           INITIALISATIONS DIVERSES
  102. C
  103. C      WRITE(IOIMP,*)' DEBUT YCVI ',' IKAS=',ikas,' IKOMP=',ikomp,
  104. C     $     'IDCEN=',idcen
  105. C      WRITE(IOIMP,*)' NPTD=',nptd,'IPAS=',IPAS
  106.       IF(IPAS.EQ.0)CALL CALHRH(QGGT,Q1,Q2,Q3,IES)
  107.  
  108. C                           ********
  109. C                           *  2D  *
  110. C                           ********
  111.  
  112.       IF(IES.EQ.3)GO TO 10
  113.  
  114. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  115. C                           DIFFERENCES TRIANGLE / QUADRANGLE
  116.       IF(NP.EQ.4)THEN
  117.          QUA4=1.D0
  118.       ELSE
  119.          QUA4=0.D0
  120.       ENDIF
  121. C
  122. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  123.  
  124. C
  125. C Calcul du nombre de paquets de LRV {l{ments
  126. C
  127.       NNN=MOD(NEL,LRV)
  128.       IF(NNN.EQ.0) NPACK=NEL/LRV
  129.       IF(NNN.NE.0) NPACK=1+(NEL-NNN)/LRV
  130.       KPACKD=1
  131.       KPACKF=NPACK
  132. C
  133. C ******* BOUCLE SUR LES PAQUETS DE LRV ELEMENTS **********
  134. C
  135. C     WRITE(IOIMP,*)' DEBUT YCVI 7001'
  136.       DO 7001 KPACK=KPACKD,KPACKF
  137. C
  138. C ======= A L'INTERIEUR DE CHAQUE PAQUET DE LRV ELEMENTS =======
  139. C
  140. C 1.D0 Calcul des limites du paquet courant.
  141.          KDEB=1+(KPACK-1)*LRV
  142.          KFIN=MIN(NEL,KDEB+LRV-1)
  143. C
  144. C Ben voil@, on peut y aller ... i.e. traiter le paquet courant.
  145. C
  146.          DO 7002 K=KDEB,KFIN
  147.             KP=K-KDEB+1
  148.             NK=K+K0
  149.             K1=1+(1-IK1)*(NK-1)
  150.             COEF(KP)=COEFF(K1)
  151.             AIRE(KP)=VOLU(NK)
  152.             AL(KP)=COTE(NK,1)+XPETIT
  153.             AH(KP)=COTE(NK,2)+XPETIT
  154.  7002    CONTINUE
  155. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  156.          DO 7006 I=1,NP
  157.             DO 7016 K=KDEB,KFIN
  158.                KP=K-KDEB+1
  159.                NU=IPADU(LE(I,K))
  160.                NG=IPADI(LE(I,K))
  161.                NF=IPADF(LE(I,K))
  162.                UIX(KP,I)=UN(NU,1)
  163.                UIY(KP,I)=UN(NU,2)
  164.                GIX(KP,I)=GN(NG,1)
  165.                GIY(KP,I)=GN(NG,2)
  166.                VIX(KP,I)=VF(NF,1)
  167.                VIY(KP,I)=VF(NF,2)
  168.  7016       CONTINUE
  169.  7006    CONTINUE
  170.          CALL KSUPG1(WX,UMI,CHGLD,CHGLPX,KDEB,KFIN,LRV,
  171.      &        GN(1,1),IPADI,UN,COEF,NPTD,NEL,NP,DRR,HR,FN,
  172.      &        AIRE,AL,AH,AP,IDCEN,IPADU,LE,QUA4,IKOMP,
  173.      &        DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  174.          CALL KSUPG1(WY,UMI,CHGLD,CHGLPY,KDEB,KFIN,LRV,
  175.      &        GN(1,2),IPADI,UN,COEF,NPTD,NEL,NP,DRR,HR,FN,
  176.      &        AIRE,AL,AH,AP,IDCEN,IPADU,LE,QUA4,IKOMP,
  177.      &        DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  178. C
  179. C Initialisation des variables d'accumulation SAF1,SAF2,SBF
  180. C
  181.          IF(IAXI.NE.0)THEN
  182.             DO 70050 N=1,IDIM
  183.                DO 70051 K=KDEB,KFIN
  184.                   KP=K-KDEB+1
  185.                   VMI(KP,N)=XPETIT
  186. 70051          CONTINUE
  187. 70050       CONTINUE
  188.             DO 70052 N=1,IDIM
  189.                DO 70053 I=1,NP
  190.                   DO 70054 K=KDEB,KFIN
  191.                      KP=K-KDEB+1
  192.                      NF=IPADF(LE(I,K))
  193.                      VMI(KP,N)=VMI(KP,N)+VF(NF,N)*DRR(I,K)
  194. 70054             CONTINUE
  195. 70053          CONTINUE
  196. 70052       CONTINUE
  197.             DO 70055 K=KDEB,KFIN
  198.                KP=K-KDEB+1
  199.                VMI(KP,1)=VMI(KP,1)/AIRE(KP)
  200.                VMI(KP,2)=VMI(KP,2)/AIRE(KP)
  201. 70055       CONTINUE
  202.          ENDIF
  203.          IF(IKOMP.EQ.0)THEN
  204.             IF(IKAS.EQ.1)THEN
  205.                DO 70061 I=1,NP
  206.                   DO 71161 K=KDEB,KFIN
  207.                      KP=K-KDEB+1
  208.                      SAF1(KP,I)=0.D0
  209.                      SAF2(KP,I)=0.D0
  210. 71161             CONTINUE
  211. 70061          CONTINUE
  212.             ELSEIF(IKAS.EQ.2)THEN
  213.                DO 70021 K=KDEB,KFIN
  214.                   KP=K-KDEB+1
  215.                   NK=K+K0
  216.                   NKG=1+(1-IKG)*(NK-1)
  217.                   RGX(KP)=RGE(NKG,1)
  218.                   RGY(KP)=RGE(NKG,2)
  219. 70021          CONTINUE
  220.  
  221.                IF(IPG.EQ.0)THEN
  222.                   DO 70062 I=1,NP
  223.                      DO 70162 K=KDEB,KFIN
  224.                         KP=K-KDEB+1
  225.                         SAF1(KP,I)=(-RGX(KP))*DRR(I,K)
  226.                         SAF2(KP,I)=(-RGY(KP))*DRR(I,K)
  227. 70162                CONTINUE
  228. 70062             CONTINUE
  229.                ELSE
  230.                   DO 71062 I=1,NP
  231.                      DO 71162 K=KDEB,KFIN
  232.                         KP=K-KDEB+1
  233.                         SAF1(KP,I)=(-RGX(KP))*WX(KP,I)
  234.                         SAF2(KP,I)=(-RGY(KP))*WY(KP,I)
  235. 71162                CONTINUE
  236. 71062             CONTINUE
  237.                ENDIF
  238.             ELSEIF(IKAS.EQ.4)THEN
  239.                DO 70022 K=KDEB,KFIN
  240.                   KP=K-KDEB+1
  241.                   NK=K+K0
  242.                   NKG=1+(1-IKG)*(NK-1)
  243.                   RGX(KP)=RGE(NKG,1)
  244.                   RGY(KP)=RGE(NKG,2)
  245. 70022          CONTINUE
  246.                IF(IPG.EQ.0)THEN
  247.                   DO 70063 I=1,NP
  248.                      DO 70163 K=KDEB,KFIN
  249.                         KP=K-KDEB+1
  250.                         NF=1+(1-IKT)*(IPADS(LE(I,K))-1)
  251.                         NFR=1+(1-IKREF)*(IPADS(LE(I,K))-1)
  252. C?     WRITE(IOIMP,*)' NF=',NF,' NFR=',nfr
  253.                         SAF1(KP,I)=(-RGX(KP)*(TREF(NFR)-TN(NF)))*DRR(I,K
  254.      $                       )
  255.                         SAF2(KP,I)=(-RGY(KP)*(TREF(NFR)-TN(NF)))*DRR(I,K
  256.      $                       )
  257. 70163                CONTINUE
  258. 70063             CONTINUE
  259.                ELSE
  260.                   DO 71063 I=1,NP
  261.                      DO 71163 K=KDEB,KFIN
  262.                         KP=K-KDEB+1
  263.                         NF=1+(1-IKT)*(IPADS(LE(I,K))-1)
  264.                         NFR=1+(1-IKREF)*(IPADS(LE(I,K))-1)
  265.                         SAF1(KP,I)=(-RGX(KP)*(TREF(NFR)-TN(NF)))*WX(KP,I
  266.      $                       )
  267.                         SAF2(KP,I)=(-RGY(KP)*(TREF(NFR)-TN(NF)))*WY(KP,I
  268.      $                       )
  269. 71163                CONTINUE
  270. 71063             CONTINUE
  271.                ENDIF
  272.             ENDIF
  273.          ELSEIF(IKOMP.EQ.1)THEN
  274.             IF(IKAS.EQ.2)THEN
  275.                DO 70064 I=1,NP
  276.                   DO 70164 K=KDEB,KFIN
  277.                      KP=K-KDEB+1
  278.                      SAF1(KP,I)=0.D0
  279.                      SAF2(KP,I)=0.D0
  280. 70164             CONTINUE
  281. 70064          CONTINUE
  282.  
  283.             ELSEIF(IKAS.EQ.3)THEN
  284.                DO 70024 K=KDEB,KFIN
  285.                   KP=K-KDEB+1
  286.                   NK=K+K0
  287.                   NKG=1+(1-IKG)*(NK-1)
  288.                   RGX(KP)=RGE(NKG,1)
  289.                   RGY(KP)=RGE(NKG,2)
  290. 70024          CONTINUE
  291.  
  292.                IF(IPG.EQ.0)THEN
  293.                   DO 70065 I=1,NP
  294.                      DO 70165 K=KDEB,KFIN
  295.                         KP=K-KDEB+1
  296.                         SAF1(KP,I)=-RGX(KP)*DRR(I,K)
  297.                         SAF2(KP,I)=-RGY(KP)*DRR(I,K)
  298. 70165                CONTINUE
  299. 70065             CONTINUE
  300.                ELSE
  301.                   DO 71065 I=1,NP
  302.                      DO 71165 K=KDEB,KFIN
  303.                         KP=K-KDEB+1
  304.                         SAF1(KP,I)=-RGX(KP)*WX(KP,I)
  305.                         SAF2(KP,I)=-RGY(KP)*WY(KP,I)
  306. 71165                CONTINUE
  307. 71065             CONTINUE
  308.                ENDIF
  309.             ENDIF
  310.          ENDIF
  311. C Le coeur du calcul ...
  312.          IF(IKOMP.EQ.0)THEN
  313.             DO 70140 I=1,NP
  314.                DO 70141 J= 1,NP
  315.                   DO 70142 K=KDEB,KFIN
  316.                      KP=K-KDEB+1
  317.                      ZVGGX=AIRE(KP)*CHGLPX(KP)*VGGT(J,I)
  318.                      ZVGGY=AIRE(KP)*CHGLPY(KP)*VGGT(J,I)
  319.                      ZVGT=AIRE(KP)*(
  320.      &                    HR(K,I,1)*HR(K,J,1)*COEF(KP)
  321.      &                    +    HR(K,I,2)*HR(K,J,2)*COEF(KP)
  322.      &                    +    CHGLD(KP)*VGGT(J,I) )
  323.                      V2=UMI(KP,1)*HR(K,J,1)+UMI(KP,2)*HR(K,J,2)
  324.                      SAF1(KP,I)=SAF1(KP,I)+(V2*WX(KP,I)+ZVGGX)*GIX(KP,J)
  325.      $                    +ZVGT*VIX(KP,J)
  326.                      SAF2(KP,I)=SAF2(KP,I)+(V2*WY(KP,I)+ZVGGY)*GIY(KP,J)
  327.      $                    +ZVGT*VIY(KP,J)
  328. 70142             CONTINUE
  329. 70141          CONTINUE
  330. 70140       CONTINUE
  331.          ELSEIF(IKOMP.EQ.1)THEN
  332.             DO 70150 I=1,NP
  333.                DO 70151 J= 1,NP
  334.                   DO 70152 K=KDEB,KFIN
  335.                      KP=K-KDEB+1
  336.                      ZVGGX=AIRE(KP)*CHGLPX(KP)*VGGT(J,I)
  337.                      ZVGGY=AIRE(KP)*CHGLPY(KP)*VGGT(J,I)
  338.                      ZVGT=AIRE(KP)*(
  339.      &                    HR(K,I,1)*HR(K,J,1)*COEF(KP)
  340.      &                    +    HR(K,I,2)*HR(K,J,2)*COEF(KP)
  341.      &                    +    CHGLD(KP)*VGGT(J,I) )
  342.                      COEFT3=COEF(KP)/3.D0
  343.                      ZVGUU=AIRE(KP)*(1.D0/AL(KP)/AL(KP)/12.D0
  344.      &                    *VGGT(J,I)*QUA4  )*COEFT3*UIX(KP,J)
  345.                      ZVGUV=AIRE(KP)*( HR(K,I,1)*HR(K,J,2))*COEFT3*UIY(KP
  346.      $                    ,J)
  347.                      ZVGVU=AIRE(KP)*( HR(K,I,2)*HR(K,J,1))*COEFT3*UIX(KP
  348.      $                    ,J)
  349.                      ZVGVV=AIRE(KP)*(1.D0/AH(KP)/AH(KP)/12.D0
  350.      &                    *VGGT(J,I)*QUA4  )*COEFT3*UIY(KP,J)
  351.                      V2=UIX(KP,J)*HR(K,J,1)+UIY(KP,J)*HR(K,J,2)
  352.                      SAF1(KP,I)=SAF1(KP,I)+(V2*WX(KP,I)+ZVGGX)*GIX(KP,J)
  353.      $                    +ZVGT*UIX(KP,J)+ ZVGUU + ZVGUV
  354.                      SAF2(KP,I)=SAF2(KP,I)+(V2*WY(KP,I)+ZVGGY)*GIY(KP,J)
  355.      $                    +ZVGT*UIY(KP,J)+ ZVGVU + ZVGVV
  356. 70152             CONTINUE
  357. 70151          CONTINUE
  358. 70150       CONTINUE
  359.          ENDIF
  360.          IF(IAXI.NE.0) THEN
  361.             DO 70160 I=1,NP
  362.                DO 70161 K=KDEB,KFIN
  363.                   KP=K-KDEB+1
  364.                   R2=1.D0/RPG(K)/RPG(K)*WX(KP,I)
  365.                   SAF1(KP,I)=SAF1(KP,I)+R2*COEF(KP)*VMI(KP,1)
  366. 70161          CONTINUE
  367. 70160       CONTINUE
  368.             IF(IKOMP.EQ.1)THEN
  369.                DO 7118 I=1,NP
  370.                   DO 71181 K=KDEB,KFIN
  371.                      KP=K-KDEB+1
  372.                      R1=1.D0/RPG(K)*WY(KP,I)
  373.                      SAF1(KP,I)=SAF1(KP,I)+R1*UMI(KP,1)*GIX(KP,I)
  374.                      SAF2(KP,I)=SAF2(KP,I)+R1*UMI(KP,1)*GIY(KP,I)
  375. 71181             CONTINUE
  376.  7118          CONTINUE
  377.             ENDIF
  378.          ENDIF
  379. C
  380. C     Fin de l'accumulation dans SAF1,SAF2.
  381. C     On ajoute ces incr{ments @ F.
  382. C
  383.          DO 7017 I=1,NP
  384.             DO 70171 K=KDEB,KFIN
  385.                KP=K-KDEB+1
  386.                NF=IPADI(LE(I,K))
  387.                F(NF,1)=F(NF,1)-SAF1(KP,I)
  388.                F(NF,2)=F(NF,2)-SAF2(KP,I)
  389. 70171       CONTINUE
  390.  7017    CONTINUE
  391.  1960    FORMAT(/,' ***** SUB XCVTIT : IPAT=',I5,'  K=',I5,' *****')
  392.  1961    FORMAT(2X,I5,' * ',4(1X,I5))
  393.  1962    FORMAT(2X,8(1X,1PE11.4))
  394.  1964    FORMAT(4(1X,1PE11.4))
  395.  7001 CONTINUE
  396. C     WRITE(IOIMP,*)' ********** FIN YCVI 2D *****************'
  397. C     CALL ARRET(0)
  398.       IPAS=1
  399.       RETURN
  400. C                           ********
  401. C                           *  3D  *
  402. C                           ********
  403.  10   CONTINUE
  404. C::::::BENET:::SUPPRESION CORRECTION HOURGLASS POUR LES PRISME::29:01:91
  405.       CUB8=0.D0
  406.       IF(NP.EQ.8)CUB8=1.D0
  407. C
  408. C Calcul du nombre de paquets de LRV {l{ments
  409. C
  410.       NNN=MOD(NEL,LRV)
  411.       IF(NNN.EQ.0) NPACK=NEL/LRV
  412.       IF(NNN.NE.0) NPACK=1+(NEL-NNN)/LRV
  413.       KPACKD=1
  414.       KPACKF=NPACK
  415. C
  416. C ******* BOUCLE SUR LES PAQUETS DE LRV ELEMENTS **********
  417. C
  418.       DO 8001 KPACK=KPACKD,KPACKF
  419. C
  420. C ======= A L'INTERIEUR DE CHAQUE PAQUET DE LRV ELEMENTS =======
  421. C
  422. C 1.D0 Calcul des limites du paquet courant.
  423.          KDEB=1+(KPACK-1)*LRV
  424.          KFIN=MIN(NEL,KDEB+LRV-1)
  425. C
  426. C Ben voil@, on peut y aller ... i.e. traiter le paquet courant.
  427. C
  428.          DO 8002 K=KDEB,KFIN
  429.             KP=K-KDEB+1
  430.             NK=K+K0
  431.             K1=1+(1-IK1)*(NK-1)
  432.             COEF(KP)=COEFF(K1)
  433.             AIRE(KP)=VOLU(NK)
  434.             AL(KP)=COTE(NK,1)+XPETIT
  435.             AH(KP)=COTE(NK,2)+XPETIT
  436.             AP(KP)=COTE(NK,3)+XPETIT
  437.  8002    CONTINUE
  438.  
  439. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  440. C Initialisation des UMI avant accumulation
  441.          DO 8006 I=1,NP
  442.             DO 81161 K=KDEB,KFIN
  443.                KP=K-KDEB+1
  444.                NU=IPADU(LE(I,K))
  445.                NF=IPADF(LE(I,K))
  446.                NG=IPADI(LE(I,K))
  447.                UIX(KP,I)=UN(NU,1)
  448.                UIY(KP,I)=UN(NU,2)
  449.                UIZ(KP,I)=UN(NU,3)
  450.                VIX(KP,I)=VF(NF,1)
  451.                VIY(KP,I)=VF(NF,2)
  452.                VIZ(KP,I)=VF(NF,3)
  453.                GIX(KP,I)=GN(NG,1)
  454.                GIY(KP,I)=GN(NG,2)
  455.                GIZ(KP,I)=GN(NG,3)
  456. 81161       CONTINUE
  457.  8006    CONTINUE
  458. C     WRITE(IOIMP,*)'****************************'
  459. C     WRITE(IOIMP,*)' DT,DTT1,DTT2=',DT,DTT1,DTT2
  460.          CALL KSUPG1(WX,UMI,CHGLD,CHGLPX,KDEB,KFIN,LRV,
  461.      &        GN(1,1),IPADI,UN,COEF,NPTD,NEL,NP,DRR,HR,FN,
  462.      &        AIRE,AL,AH,AP,IDCEN,IPADU,LE,CUB8,IKOMP,
  463.      &        DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  464. C     WRITE(IOIMP,1002)umi
  465. C     WRITE(IOIMP,*)' DT,DTT1,DTT2=',DT,DTT1,DTT2
  466.          CALL KSUPG1(WY,UMI,CHGLD,CHGLPY,KDEB,KFIN,LRV,
  467.      &        GN(1,2),IPADI,UN,COEF,NPTD,NEL,NP,DRR,HR,FN,
  468.      &        AIRE,AL,AH,AP,IDCEN,IPADU,LE,CUB8,IKOMP,
  469.      &        DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  470. C     WRITE(IOIMP,*)' DT,DTT1,DTT2=',DT,DTT1,DTT2
  471.          CALL KSUPG1(WZ,UMI,CHGLD,CHGLPZ,KDEB,KFIN,LRV,
  472.      &        GN(1,3),IPADI,UN,COEF,NPTD,NEL,NP,DRR,HR,FN,
  473.      &        AIRE,AL,AH,AP,IDCEN,IPADU,LE,CUB8,IKOMP,
  474.      &        DTM1,DT,DTT1,DTT2,DIAEL,NUEL)
  475. C     WRITE(IOIMP,*)' DT,DTT1,DTT2=',DT,DTT1,DTT2
  476. C
  477. C Initialisation des variables d'accumulation SAF1,SAF2,SBF
  478. C
  479.          IF(IKOMP.EQ.0)THEN
  480.             IF(IKAS.EQ.1)THEN
  481.                DO 80061 I=1,NP
  482.                   DO 80161 K=KDEB,KFIN
  483.                      KP=K-KDEB+1
  484.                      SAF1(KP,I)=0.D0
  485.                      SAF2(KP,I)=0.D0
  486.                      SAF3(KP,I)=0.D0
  487. 80161             CONTINUE
  488. 80061          CONTINUE
  489.             ELSEIF(IKAS.EQ.2)THEN
  490.                DO 80021 K=KDEB,KFIN
  491.                   KP=K-KDEB+1
  492.                   NK=K+K0
  493.                   NKG=1+(1-IKG)*(NK-1)
  494.                   RGX(KP)=RGE(NKG,1)
  495.                   RGY(KP)=RGE(NKG,2)
  496.                   RGZ(KP)=RGE(NKG,3)
  497. 80021          CONTINUE
  498.                IF(IPG.EQ.0)THEN
  499.                   DO 80062 I=1,NP
  500.                      DO 80162 K=KDEB,KFIN
  501.                         KP=K-KDEB+1
  502.                         SAF1(KP,I)=(-RGX(KP))*DRR(I,K)
  503.                         SAF2(KP,I)=(-RGY(KP))*DRR(I,K)
  504.                         SAF3(KP,I)=(-RGZ(KP))*DRR(I,K)
  505. 80162                CONTINUE
  506. 80062             CONTINUE
  507.                ELSE
  508.                   DO 81062 I=1,NP
  509.                      DO 81162 K=KDEB,KFIN
  510.                         KP=K-KDEB+1
  511.                         SAF1(KP,I)=(-RGX(KP))*WX(KP,I)
  512.                         SAF2(KP,I)=(-RGY(KP))*WY(KP,I)
  513.                         SAF3(KP,I)=(-RGZ(KP))*WZ(KP,I)
  514. 81162                CONTINUE
  515. 81062             CONTINUE
  516.                ENDIF
  517.             ELSEIF(IKAS.EQ.4)THEN
  518.                DO 80022 K=KDEB,KFIN
  519.                   KP=K-KDEB+1
  520.                   NK=K+K0
  521.                   NKG=1+(1-IKG)*(NK-1)
  522.                   RGX(KP)=RGE(NKG,1)
  523.                   RGY(KP)=RGE(NKG,2)
  524.                   RGZ(KP)=RGE(NKG,3)
  525. 80022          CONTINUE
  526.                IF(IPG.EQ.0)THEN
  527.                   DO 80063 I=1,NP
  528.                      DO 80163 K=KDEB,KFIN
  529.                         KP=K-KDEB+1
  530.                         NF=1+(1-IKT)*(IPADS(LE(I,K))-1)
  531.                         NFR=1+(1-IKREF)*(IPADS(LE(I,K))-1)
  532.                         SAF1(KP,I)=(-RGX(KP)*(TREF(NFR)-TN(NF)))*DRR(I,K
  533.      $                       )
  534.                         SAF2(KP,I)=(-RGY(KP)*(TREF(NFR)-TN(NF)))*DRR(I,K
  535.      $                       )
  536.                         SAF3(KP,I)=(-RGZ(KP)*(TREF(NFR)-TN(NF)))*DRR(I,K
  537.      $                       )
  538. 80163                CONTINUE
  539. 80063             CONTINUE
  540.                ELSE
  541.                   DO 81063 I=1,NP
  542.                      DO 81163 K=KDEB,KFIN
  543.                         KP=K-KDEB+1
  544.                         NF=1+(1-IKT)*(IPADS(LE(I,K))-1)
  545.                         NFR=1+(1-IKREF)*(IPADS(LE(I,K))-1)
  546.                         SAF1(KP,I)=(-RGX(KP)*(TREF(NFR)-TN(NF)))*WX(KP,I
  547.      $                       )
  548.                         SAF2(KP,I)=(-RGY(KP)*(TREF(NFR)-TN(NF)))*WY(KP,I
  549.      $                       )
  550.                         SAF3(KP,I)=(-RGZ(KP)*(TREF(NFR)-TN(NF)))*WZ(KP,I
  551.      $                       )
  552. 81163                CONTINUE
  553. 81063             CONTINUE
  554.                ENDIF
  555.             ENDIF
  556.          ELSEIF(IKOMP.EQ.1)THEN
  557.             IF(IKAS.EQ.2)THEN
  558.                DO 80064 I=1,NP
  559.                   DO 80164 K=KDEB,KFIN
  560.                      KP=K-KDEB+1
  561.                      SAF1(KP,I)=0.D0
  562.                      SAF2(KP,I)=0.D0
  563.                      SAF3(KP,I)=0.D0
  564. 80164             CONTINUE
  565. 80064          CONTINUE
  566.             ELSEIF(IKAS.EQ.3)THEN
  567.                DO 80024 K=KDEB,KFIN
  568.                   KP=K-KDEB+1
  569.                   NK=K+K0
  570.                   NKG=1+(1-IKG)*(NK-1)
  571.                   RGX(KP)=RGE(NKG,1)
  572.                   RGY(KP)=RGE(NKG,2)
  573.                   RGZ(KP)=RGE(NKG,3)
  574. 80024          CONTINUE
  575.                IF(IPG.EQ.0)THEN
  576.                   DO 80065 I=1,NP
  577.                      DO 80165 K=KDEB,KFIN
  578.                         KP=K-KDEB+1
  579.                         SAF1(KP,I)=-RGX(KP)*DRR(I,K)
  580.                         SAF2(KP,I)=-RGY(KP)*DRR(I,K)
  581.                         SAF3(KP,I)=-RGZ(KP)*DRR(I,K)
  582. 80165                CONTINUE
  583. 80065             CONTINUE
  584.                ELSE
  585.                   DO 81065 I=1,NP
  586.                      DO 81165 K=KDEB,KFIN
  587.                         KP=K-KDEB+1
  588.                         SAF1(KP,I)=-RGX(KP)*WX(KP,I)
  589.                         SAF2(KP,I)=-RGY(KP)*WY(KP,I)
  590.                         SAF3(KP,I)=-RGZ(KP)*WZ(KP,I)
  591. 81165                CONTINUE
  592. 81065             CONTINUE
  593.                ENDIF
  594.             ENDIF
  595.          ENDIF
  596. C Le coeur du calcul ...
  597.          IF(IKOMP.EQ.0)THEN
  598.             DO 80140 I=1,NP
  599.                DO 80141 J= 1,NP
  600.                   DO 80142 K=KDEB,KFIN
  601.                      KP=K-KDEB+1
  602.                      ZVGGX=AIRE(KP)*CHGLPX(KP)*QGGT(J,I)
  603.                      ZVGGY=AIRE(KP)*CHGLPY(KP)*QGGT(J,I)
  604.                      ZVGGZ=AIRE(KP)*CHGLPZ(KP)*QGGT(J,I)
  605.                      ZVGT=AIRE(KP)*(
  606.      &                    HR(K,I,1)*HR(K,J,1)*COEF(KP)
  607.      &                    +    HR(K,I,2)*HR(K,J,2)*COEF(KP)
  608.      &                    +    HR(K,I,3)*HR(K,J,3)*COEF(KP) )
  609.      &                    +    CHGLD(KP)*QGGT(J,I)
  610.                      V2=UMI(KP,1)*HR(K,J,1)+UMI(KP,2)*HR(K,J,2)
  611.      &                    +UMI(KP,3)*HR(K,J,3)
  612.                      SAF1(KP,I)=SAF1(KP,I)+(V2*WX(KP,I)+ZVGGX)*GIX(KP,J)
  613.      $                    +ZVGT*VIX(KP,J)
  614.                      SAF2(KP,I)=SAF2(KP,I)+(V2*WY(KP,I)+ZVGGY)*GIY(KP,J)
  615.      $                    +ZVGT*VIY(KP,J)
  616.                      SAF3(KP,I)=SAF3(KP,I)+(V2*WZ(KP,I)+ZVGGZ)*GIZ(KP,J)
  617.      $                    +ZVGT*VIZ(KP,J)
  618. 80142             CONTINUE
  619. 80141          CONTINUE
  620. 80140       CONTINUE
  621.          ELSEIF(IKOMP.EQ.1)THEN
  622.             DO 80150 I=1,NP
  623.                DO 80151 J= 1,NP
  624.                   DO 80152 K=KDEB,KFIN
  625.                      KP=K-KDEB+1
  626.                      ZVGGX=AIRE(KP)*CHGLPX(KP)*QGGT(J,I)
  627.                      ZVGGY=AIRE(KP)*CHGLPY(KP)*QGGT(J,I)
  628.                      ZVGGZ=AIRE(KP)*CHGLPZ(KP)*QGGT(J,I)
  629.                      ZVGT=AIRE(KP)*(
  630.      &                    HR(K,I,1)*HR(K,J,1)*COEF(KP)
  631.      &                    +    HR(K,I,2)*HR(K,J,2)*COEF(KP)
  632.      &                    +    HR(K,I,3)*HR(K,J,3)*COEF(KP) )
  633.      &                    +    CHGLD(KP)*QGGT(J,I)
  634.                      COEFT3=COEF(KP)/3.D0
  635.                      GEO1=0.D0
  636.                      ZVGUU=(AIRE(KP)*(HR(K,I,1)*HR(K,J,1))+GEO1)*COEFT3
  637.      $                    *UIX(KP,J)
  638.                      ZVGUV=AIRE(KP)*( HR(K,I,1)*HR(K,J,2))*COEFT3*UIY(KP
  639.      $                    ,J)
  640.                      ZVGUW=AIRE(KP)*( HR(K,I,1)*HR(K,J,3))*COEFT3*UIZ(KP
  641.      $                    ,J)
  642.                      ZVGVU=AIRE(KP)*( HR(K,I,2)*HR(K,J,1))*COEFT3*UIX(KP
  643.      $                    ,J)
  644.                      ZVGVV=(AIRE(KP)*(HR(K,I,2)*HR(K,J,2))+GEO1)*COEFT3
  645.      $                    *UIY(KP,J)
  646.                      ZVGVW=AIRE(KP)*( HR(K,I,2)*HR(K,J,3))*COEFT3*UIZ(KP
  647.      $                    ,J)
  648.                      ZVGWU=AIRE(KP)*( HR(K,I,3)*HR(K,J,1))*COEFT3*UIX(KP
  649.      $                    ,J)
  650.                      ZVGWV=AIRE(KP)*( HR(K,I,3)*HR(K,J,2))*COEFT3*UIY(KP
  651.      $                    ,J)
  652.                      ZVGWW=(AIRE(KP)*(HR(K,I,3)*HR(K,J,3))+GEO1)*COEFT3
  653.      $                    *UIZ(KP,J)
  654.                      V2=(UIX(KP,J)*HR(K,J,1)+UIY(KP,J)*HR(K,J,2)+UIZ(KP
  655.      $                    ,J)*HR(K,J,3))*WZ(KP,I)
  656.                      SAF1(KP,I)=SAF1(KP,I)+(V2+ZVGGX)*GIX(KP,J)+ZVGT
  657.      $                    *UIX(KP,J)+ ZVGUU + ZVGUV + ZVGUW
  658.                      SAF2(KP,I)=SAF2(KP,I)+(V2+ZVGGY)*GIY(KP,J)+ZVGT
  659.      $                    *UIY(KP,J)+ ZVGVU + ZVGUV + ZVGVW
  660.                      SAF3(KP,I)=SAF3(KP,I)+(V2+ZVGGZ)*GIZ(KP,J)+ZVGT
  661.      $                    *UIZ(KP,J)+ ZVGWU + ZVGWV + ZVGWW
  662. 80152             CONTINUE
  663. 80151          CONTINUE
  664. 80150       CONTINUE
  665.          ENDIF
  666. C
  667. C     Fin de l'accumulation dans SAF1,SAF2.
  668. C     On ajoute ces incr{ments @ F.
  669. C
  670.          DO 8017 I=1,NP
  671.             DO 80171 K=KDEB,KFIN
  672.                KP=K-KDEB+1
  673.                NF=IPADI(LE(I,K))
  674.                F(NF,1)=F(NF,1)-SAF1(KP,I)
  675.                F(NF,2)=F(NF,2)-SAF2(KP,I)
  676.                F(NF,3)=F(NF,3)-SAF3(KP,I)
  677. 80171       CONTINUE
  678.  8017    CONTINUE
  679.  8001 CONTINUE
  680. C     WRITE(IOIMP,*)' ********** FIN YCVI 3D *****************'
  681.       IPAS=1
  682.       RETURN
  683.  1002 FORMAT(10(1X,1PE11.4))
  684.  1001 FORMAT(20(1X,I5))
  685.       END
  686.  
  687.  
  688.  
  689.  
  690.  
  691.  
  692.  
  693.  
  694.  
  695.  
  696.  

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