Télécharger la plaquette Cast3M


fr en

Télécharger zclpls.eso

Retour à la liste

Numérotation des lignes : 

zclpls

  1. C ZCLPLS    SOURCE    CB215821  26/08/24    21:19:02     12622          
  2.       SUBROUTINE ZCLPLS(HR,RPG,LE,NEL,K0,NPTD,NX,IES,NP,IAXI,IPADL,
  3.      & COEFF,IK1,
  4.      & TN,G,
  5.      & VOLU,COTE,NELZ,DME,DIAM,
  6.      & DT,DTT2,NUEL,DIAEL)
  7. C
  8. C VERSION VECTORISEE (CF XCVTIT POUR PLUS DE DETAILS)
  9. C
  10.       IMPLICIT INTEGER(I-N)
  11.       IMPLICIT REAL*8 (A-H,O-Z)
  12. C NOMBRE MAXI DE POINTS PAR ELEMENT : NPX
  13.       PARAMETER(NPX=9)
  14. C LONGUEUR DES REGISTRES VECTORIELS DE LA MACHINE CIBLE
  15.       PARAMETER(LRV=64)
  16. C***********************************************************************
  17. C
  18. C  CE SP DISCRETISE LE LAPLACIEN D UNE VARIABLE SCALAIRE
  19. C
  20. C    EN 1D SUR L ELEMENT SEG2
  21. C    EN 2D SUR L ELEMENT QUA4       PLAN OU AXI
  22. C    EN 3D SUR L ELEMENT CUB8 PRI6 et TET4
  23. C  LES OPERATEURS SONT SOUS INTEGRES
  24. C
  25. C  COEFF(SCAL DOMA)  LE COEFFICIENT DU LAPLACIEN  ( NU )
  26. C       (SCAL ELEM)
  27. C
  28. C     NUELD(NELZ): NUMERO DE L'ELEMENT NK DE LA ZONE DANS LE DOMAINE
  29. C     (NELZ      : NOMBRE D'ELEMENTS DE LA ZONE)
  30. C     ROC(NELD)  : ROC PAR ELEMENT DONNE SUR TOUT LE DOMAINE
  31. C
  32. C
  33. C***********************************************************************
  34. -INC CCVQUA4
  35. -INC CCREEL
  36. C
  37.       DIMENSION TN(NPTD,NX),GG(8)
  38.       DIMENSION COEFF(*)
  39.       DIMENSION COTE(NELZ,IES),DME(*),VOLU(*),DIAM(*)
  40.  
  41.       DIMENSION IPADL(1),LE(NP,1)
  42.       DIMENSION HR(NEL,NP,IES),RPG(1)
  43.  
  44.       REAL*8    G(NPTD,NX)
  45.       REAL*8    ROC(1)
  46.       DIMENSION NUELD(1)
  47.       DIMENSION QGGT(8,8),Q1(8,8),Q2(8,8),Q3(8,8)
  48. C
  49. C
  50.  
  51. C* -TABLEAUX ADDITIONNELS POUR LA VECTORISATION   **********************
  52.       DIMENSION AIRE(LRV),ALF (LRV),COEF(LRV),CLSR(LRV),GINC(LRV)
  53.       DIMENSION AL  (LRV),AH  (LRV),AP  (LRV),CFM (LRV),XMA (LRV)
  54.       DIMENSION XMB(LRV),XMD(LRV),XMH(LRV),XMI(LRV),AHL(LRV)
  55.       DIMENSION ALH(LRV),DPR(LRV),TETA(LRV,NPX,3),F(LRV,NPX),BF(LRV,NPX)
  56. C
  57. C --- TABLEAU POUR L'OPTION RAPIDE ---
  58.       DIMENSION ALP(LRV)
  59. C***
  60.       SAVE IPAS,QGGT,Q1,Q2,Q3
  61.       DATA IPAS/0/
  62.  
  63. C
  64. C           INITIALISATIONS DIVERSES
  65. C
  66. C      WRITE(6,*)' ROC=',(ROC(MM),MM=1,NEL)
  67. C      WRITE(6,*)' NUELD=',(NUELD(MM),MM=1,NEL)
  68.  
  69.       NK=K0
  70.       CALL INITD(BF,64,0.D0)
  71.  
  72.       IF(IES.EQ.2)GO TO 20
  73.       IF(IES.EQ.3)GO TO 30
  74. C**********
  75. C                                          *   1D*
  76. C**********
  77. C
  78. C
  79. C K EST LE NUMERO DE L'ELEMENT
  80. C KP PERMET DE SE SITUER A L'INTERIEUR DES TABLEAUX DE LRV ELEMENTS
  81.       NPACK=INT(NEL/FLOAT(LRV))+1
  82.       KPACKD=1
  83.       KPACKF=NPACK
  84. C
  85. C*** BOUCLE SUR LES ELEMENTS ***
  86. C
  87.       DO 70001 KPACK=KPACKD,KPACKF
  88. C     DO 40 K=1,NEL
  89. C
  90. C POUR CHAQUE PAQUET DE LRV ELEMENTS
  91.       KDEB=1+(KPACK-1)*LRV
  92.       KFIN=MIN(NEL,KDEB+LRV-1)
  93.       DO 70002 K=KDEB,KFIN
  94.       KP=K-KDEB+1
  95.       NK=K+K0
  96.       K1=1+(1-IK1)*(NK-1)
  97.       COEF(KP)=COEFF(K1)
  98.       CLSR(KP)=COEFF(K1)
  99.       AIRE(KP)=VOLU(NK)
  100.       XMA(KP)=AIRE(KP)
  101. 70002 CONTINUE
  102.       DO 90008 I=1,NP
  103.       DO 90007 N=1,NX
  104.       DO 4 K=KDEB,KFIN
  105.       KP=K-KDEB+1
  106.       NF=IPADL(LE(I,K))
  107.       TETA(KP,I,N)=TN(NF,N)
  108.     4 CONTINUE
  109. 90007 CONTINUE
  110. 90008 CONTINUE
  111. C
  112. C
  113. C
  114.       DO 70003 K=KDEB,KFIN
  115.       NK=K+K0
  116.       KP=K-KDEB+1
  117.       DT0=DT
  118.       DT2=0.5*XMA(KP)*XMA(KP)/CLSR(KP)
  119.       IF(DT2.LT.DT)DT=DT2
  120.       IF(DT.EQ.DT0)GO TO 152
  121.       DTT2=DT2
  122.       DIAEL=XMA(KP)
  123.       NUEL=NK
  124.  152  CONTINUE
  125. 70003 CONTINUE
  126. C
  127. C
  128.       N=1
  129. C --- CALCUL DE L'INCREMENT GINC ---------------------------
  130.       DO 70104 K=KDEB,KFIN
  131.       KP=K-KDEB+1
  132.       GINC(KP)=COEF(KP)*(TETA(KP,1,N)-TETA(KP,2,N))/AIRE(KP)
  133. 70104 CONTINUE
  134. C --- ACCUMULATION DANS LE TABLEAU G -----------------------
  135.       DO 70004 K=KDEB,KFIN
  136.       KP=K-KDEB+1
  137.       NF=IPADL(LE(1,K))
  138.       G(NF,N)=G(NF,N) +  GINC(KP)
  139.       NF=IPADL(LE(2,K))
  140.       G(NF,N)=G(NF,N) -  GINC(KP)
  141. 70004 CONTINUE
  142. 70001 CONTINUE
  143.   40  CONTINUE
  144.       IPAS=1
  145.       RETURN
  146. C                           ***********
  147. C                           *   2D    *
  148. C                           ***********
  149.   20  CONTINUE
  150. C
  151.  
  152. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  153.  
  154. C                           DIFFERENCES TRIANGLE / QUADRANGLE
  155.       QUA4=0.D0
  156.       IF(NP.EQ.4)QUA4=1.D0
  157.  
  158. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  159.  
  160.  
  161. C
  162. C K EST LE NUMERO DE L'ELEMENT
  163. C KP PERMET DE SE SITUER A L'INTERIEUR DES TABLEAUX DE LRV ELEMENTS
  164.       NPACK=INT(NEL/FLOAT(LRV))+1
  165.       KPACKD=1
  166.       KPACKF=NPACK
  167. C
  168. C*** BOUCLE SUR LES ELEMENTS ***
  169. C
  170.       DO 80001 KPACK=KPACKD,KPACKF
  171. C
  172. C POUR CHAQUE PAQUET DE LRV ELEMENTS
  173.       KDEB=1+(KPACK-1)*LRV
  174.       KFIN=MIN(NEL,KDEB+LRV-1)
  175.       DO 80002 K=KDEB,KFIN
  176.       KP=K-KDEB+1
  177.       NK=K+K0
  178.       K1=1+(1-IK1)*(NK-1)
  179.       XMI(KP)=DIAM(NK)
  180.       XMA(KP)=DME(NK)
  181.       COEF(KP)=COEFF(K1)
  182.       CLSR(KP)=COEFF(K1)
  183.       AIRE(KP)=VOLU(NK)
  184. CC
  185.       AL(KP)=COTE(NK,1)
  186.       AH(KP)=COTE(NK,2)
  187.       XMB(KP)=(AL(KP)+AH(KP))/2.
  188.       XMD(KP)=((AL(KP)*AH(KP))**2)/(AL(KP)**2+AH(KP)**2)
  189.       AHL(KP)=AH(KP)/AL(KP)
  190.       ALH(KP)=AL(KP)/AH(KP)
  191.       DPR(KP)=1.D0
  192.       IF(IAXI.NE.0)DPR(KP)=2.D0*XPI*RPG(K)
  193. 80002 CONTINUE
  194. CC
  195. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  196. C
  197.       DO 90010 I=1,NP
  198.       DO 90009 N=1,NX
  199.       DO 5 K=KDEB,KFIN
  200.       KP=K-KDEB+1
  201.       NF=IPADL(LE(I,K))
  202.       TETA(KP,I,N)=TN(NF,N)
  203.     5 CONTINUE
  204. 90009 CONTINUE
  205. 90010 CONTINUE
  206. C
  207. C
  208. C
  209.       DO 80004 K=KDEB,KFIN
  210.       KP=K-KDEB+1
  211.       NK=K+K0
  212.       DT0=DT
  213.       DT2=0.5*XMD(KP)/CLSR(KP)
  214.       IF(DT2.LT.DT)DT=DT2
  215.       IF(DT.EQ.DT0)GO TO 252
  216.       DTT2=DT2
  217.       DIAEL=XMB(KP)
  218.       NUEL=NK
  219.  252  CONTINUE
  220. 80004 CONTINUE
  221. C
  222. C --- VERSION RAPIDE -----------------------
  223.       CF12=1.D0/12.D0
  224.       DO 80011 K=KDEB,KFIN
  225.       KP=K-KDEB+1
  226.       ALP(KP)=CF12*(AHL(KP)+ALH(KP))*DPR(KP)*QUA4
  227. 80011 CONTINUE
  228. C ------------------------------------------
  229. C
  230.       DO 90011 N=1,NX
  231.       DO 90012 I= 1,NP
  232.       DO 80010 K=KDEB,KFIN
  233.       KP=K-KDEB+1
  234.       BF(KP,I)=0.D0
  235. 80010 CONTINUE
  236. 90012 CONTINUE
  237.  
  238.       DO 1 I=1,NP
  239.       DO 90013 J= 1,NP
  240.       DO 3 K=KDEB,KFIN
  241.       KP=K-KDEB+1
  242.       BF(KP,I)=BF(KP,I)+TETA(KP,J,N)*
  243.      &       ((HR(K,I,1)*HR(K,J,1)+HR(K,I,2)*HR(K,J,2))*AIRE(KP)
  244.      &        +ALP(KP)*VGGT(J,I))*COEF(KP)
  245.    3  CONTINUE
  246. 90013 CONTINUE
  247.  
  248. C --- ACCUMULATION DANS LE TABLEAU G ------------------------------
  249.       DO 80006 K=KDEB,KFIN
  250.       KP=K-KDEB+1
  251.       NF=IPADL(LE(I,K))
  252.       G(NF,N)=G(NF,N)+BF(KP,I)
  253. 80006 CONTINUE
  254.     1 CONTINUE
  255. 90011 CONTINUE
  256. 80001 CONTINUE
  257.   50  CONTINUE
  258.  
  259. C     WRITE(6,*)' SUB XCLPLS G(1,='
  260. C     WRITE(6,1002) (G(MM,1),MM=1,NPTD)
  261. C     WRITE(6,*)' SUB XCLPLS G(2,='
  262. C     WRITE(6,1002) (G(MM,2),MM=1,NPTD)
  263. C     WRITE(6,*)' FIN **** '
  264.  
  265.       IPAS=1
  266.       RETURN
  267.  30   CONTINUE
  268. C                           ***********
  269. C                           *   3D    *
  270. C                           ***********
  271.  
  272.       IF(IPAS.EQ.0)CALL CALHRH(QGGT,Q1,Q2,Q3,IES)
  273.       CUB8=0.D0
  274.       IF(NP.EQ.8)CUB8=1.D0
  275.  
  276.  1003 FORMAT('   XCVTI ',10I10)
  277. C
  278. C K EST LE NUMERO DE L'ELEMENT
  279. C KP PERMET DE SE SITUER A L'INTERIEUR DES TABLEAUX DE LRV ELEMENTS
  280.       NPACK=INT(NEL/FLOAT(LRV))+1
  281.       KPACKD=1
  282.       KPACKF=NPACK
  283. C
  284. C*** BOUCLE SUR LES ELEMENTS ***
  285. C
  286.       DO 90001 KPACK=KPACKD,KPACKF
  287. C     DO 60 K=1,NEL
  288. C
  289. C POUR CHAQUE PAQUET DE LRV ELEMENTS
  290.       KDEB=1+(KPACK-1)*LRV
  291.       KFIN=MIN(NEL,KDEB+LRV-1)
  292.       DO 90002 K=KDEB,KFIN
  293.       KP=K-KDEB+1
  294.       NK=K+K0
  295.       K1=1+(1-IK1)*(NK-1)
  296.       AL(KP)=COTE(NK,1)
  297.       AH(KP)=COTE(NK,2)
  298.       AP(KP)=COTE(NK,3)
  299.       CFM(KP)=AL(KP)*AH(KP)/AP(KP)+AL(KP)*AP(KP)/AH(KP)+
  300.      &        AP(KP)*AH(KP)/AL(KP)
  301.       XMI(KP)=DIAM(NK)
  302.       XMA(KP)=DME(NK)
  303.       COEF(KP)=COEFF(K1)
  304. C      CLSR(KP)=COEFF(K1)/ROC(NUELD(NK))
  305.       CLSR(KP)=COEFF(K1)
  306.       AIRE(KP)=VOLU(NK)
  307.       XMB(KP)=(XMA(KP)+XMI(KP))*0.5
  308.       XMD(KP)=((XMA(KP)*XMI(KP))**2)/(2.*XMA(KP)**2+XMI(KP)**2)
  309.       XMH(KP)=AIRE(KP)**0.33
  310. 90002 CONTINUE
  311. C
  312.       DO 90015 I=1,NP
  313.       DO 90014 N=1,NX
  314.       DO 15 K=KDEB,KFIN
  315.       KP=K-KDEB+1
  316.       NF=IPADL(LE(I,K))
  317.       TETA(KP,I,N)=TN(NF,N)
  318.    15 CONTINUE
  319. 90014 CONTINUE
  320. 90015 CONTINUE
  321. C     WRITE(6,*)' K ** ',K,' TETA ','    VOLU=',AIRE,'  COEF=',COEF
  322. C    &,'  XMH=',XMH
  323. C     WRITE(6,1002)(TETA(MM,1),MM=1,8)
  324. C
  325. C
  326.       DO 90003 K=KDEB,KFIN
  327.       KP=K-KDEB+1
  328.       NK=K+K0
  329.       DT0=DT
  330.       DT2=0.5*XMD(KP)/CLSR(KP)
  331.       IF(DT2.LT.DT)DT=DT2
  332.       IF(DT.EQ.DT0)GO TO 353
  333.       DTT2=DT2
  334.       DIAEL=XMB(KP)
  335.       NUEL=NK
  336.  353  CONTINUE
  337. 90003 CONTINUE
  338. C
  339.       DO 90016 N=1,NX
  340.       DO 90017 I=1,NP
  341.       DO 90004 K=KDEB,KFIN
  342.       KP=K-KDEB+1
  343.       BF(KP,I)=0.
  344. 90004 CONTINUE
  345. 90017 CONTINUE
  346.       DO 11 I=1,NP
  347.       DO 90018 J=1,NP
  348.       DO 13 K=KDEB,KFIN
  349.       KP=K-KDEB+1
  350.       GEO1=CFM(KP)*QGGT(J,I)*CUB8
  351.       BF(KP,I)=BF(KP,I)+TETA(KP,J,N)*
  352.      &(HR(K,I,1)*HR(K,J,1)+HR(K,I,2)*HR(K,J,2)+HR(K,I,3)*HR(K,J,3))
  353.      & *AIRE(KP)+COEF(KP)*XMH(KP)*GEO1
  354.   13  CONTINUE
  355. 90018 CONTINUE
  356.       DO 90006 K=KDEB,KFIN
  357.       KP=K-KDEB+1
  358.       NF=IPADL(LE(I,K))
  359.       G(NF,N)=G(NF,N)+BF(KP,I)
  360. 90006 CONTINUE
  361.    11 CONTINUE
  362. 90016 CONTINUE
  363. 90001 CONTINUE
  364. C     WRITE(6,*)' GG(I)='
  365. C     WRITE(6,1002)GG
  366.    60 CONTINUE
  367. C     CALL ARRET(0)
  368.       IPAS=1
  369.       RETURN
  370.  1001 FORMAT('  XCVTI',I10,6E12.5)
  371.  1002 FORMAT(10(1X,1PE11.4))
  372.       END
  373.  
  374.  
  375.  
  376.  
  377.  
  378.  
  379.  
  380.  
  381.  
  382.  

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