Télécharger la plaquette Cast3M


fr en

Télécharger xlapl.eso

Retour à la liste

Numérotation des lignes : 

xlapl

  1. C XLAPL     SOURCE    CHAT      05/01/13    04:14:37     5004
  2.       SUBROUTINE XLAPL(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,
  3.      & COEF,COG,COES,KITT,KJTT,IK1,
  4.      & LE,NBEL,K0,XCOOR,AIMPL,IKOMP,
  5.      & AF1,AF2,AF3,
  6.      & AS1,AS2,AS3,
  7.      & NINC,IHV,IARG,S2)
  8. C
  9.       IMPLICIT INTEGER(I-N)
  10.       IMPLICIT REAL*8 (A-H,O-Z)
  11. C************************************************************************
  12. C
  13. C     CALCUL MATRICE ELEMENTAIRE DU LAPLACIEN
  14. C
  15. C
  16. C************************************************************************
  17.  
  18. C     DIMENSION FN(NP,NPG),HR(KES,NP,NPG),PGSQ(NPG),RPG(NPG)
  19.       DIMENSION LE(NP,NBEL)
  20.       DIMENSION XYZ(IDIM,NP),FN(NP,NPG),GR(IDIM,NP,NPG),PG(NPG)
  21.       DIMENSION HR(IDIM,NP,NPG),PGSQ(NPG),RPG(NPG)
  22.       DIMENSION AF1(NBEL,NP,NP),AF2(NBEL,NP,NP),AF3(NBEL,NP,NP)
  23.       DIMENSION AS1(NBEL,NP,NP),AS2(NBEL,NP,NP),AS3(NBEL,NP,NP)
  24.       DIMENSION CL(9),COEF(1),COES(IDIM,IDIM,*)
  25.       DIMENSION AF(4,4,2,2),XCOOR(*)
  26. C
  27.       DIMENSION COE(3,3),S2(IDIM,1),COG(2,1)
  28.  
  29. C     write(6,*)' XLAPL  KITT,KJTT,IK1 =',KITT,KJTT,IK1
  30.       NK=K0
  31.       DO 108 KE=1,NBEL
  32.       NK=NK+1
  33.       DO 109 I=1,NP
  34.       J=LE(I,KE)
  35.       DO 109 N=1,IDIM
  36.       XYZ(N,I)=XCOOR((J-1)*(IDIM+1)+N)
  37.  109  CONTINUE
  38.  
  39.       CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,
  40.      *IDIM,NP,NPG,IAXI,AIRE)
  41. C
  42.       IAX=0
  43. C     WRITE(6,*)
  44. C     WRITE(6,*)' KITT ',KITT,' KJTT ',KJTT
  45. C     if(iax.eq.0)go to 108
  46.       IF(IHV.EQ.1)IAX=IAXI
  47.       IF(KITT.EQ.4)GO TO 70
  48.       IF(KITT.EQ.3)GO TO 80
  49.       IF(KITT.EQ.1)GO TO 80
  50.       IF(KJTT.EQ.4)GO TO 5
  51.       IC=1+(1-IK1)*(NK-1)
  52. C     WRITE(6,*)' IC ',IC ,COEF(IC)
  53.       DO 1 I=1,NP
  54.       DO 1 J=1,NP
  55.       U=0.D0
  56.       DO 2 L=1,NPG
  57.       V=0.D0
  58.       DO 3 N=1,IDIM
  59.       V=V+HR(N,I,L)*HR(N,J,L)
  60.  3    CONTINUE
  61.       U=U+V*PGSQ(L)
  62.  2    CONTINUE
  63.       U=U*COEF(IC)
  64.       AF1(KE,J,I)=AF1(KE,J,I)+U*AIMPL
  65.       AS1(KE,J,I)=AS1(KE,J,I)+U*(AIMPL-1.D0)
  66.  1    CONTINUE
  67.  
  68.       IF(NINC.GE.2)THEN
  69.       DO 1762 I=1,NP
  70.       DO 1762 J=1,NP
  71.       AF2(KE,I,J)=AF1(KE,I,J)
  72.       AS2(KE,I,J)=AS1(KE,I,J)
  73.  1762 CONTINUE
  74.       ENDIF
  75.       IF(NINC.GE.3)THEN
  76.       DO 1763 I=1,NP
  77.       DO 1763 J=1,NP
  78.       AF3(KE,I,J)=AF1(KE,I,J)
  79.       AS3(KE,I,J)=AS1(KE,I,J)
  80.  1763 CONTINUE
  81.       ENDIF
  82.  
  83.       IF(IHV.EQ.1)THEN
  84.       IF(IAXI.EQ.1)THEN
  85.       DO 41 I=1,NP
  86.       DO 41 J=1,NP
  87.       U=0.D0
  88.       DO 42 L=1,NPG
  89.       U=U+FN(I,L)*FN(J,L)/RPG(L)/RPG(L)*PGSQ(L)
  90.  42   CONTINUE
  91.       U=U*COEF(IC)
  92.       AF2(KE,J,I)=AF2(KE,J,I)+U*AIMPL
  93.       AS2(KE,J,I)=AS2(KE,J,I)+U*(AIMPL-1.D0)
  94.  41   CONTINUE
  95.       ELSEIF(IAXI.EQ.2)THEN
  96.       DO 43 I=1,NP
  97.       DO 43 J=1,NP
  98.       U=0.D0
  99.       DO 44 L=1,NPG
  100.       U=U+FN(I,L)*FN(J,L)/RPG(L)/RPG(L)*PGSQ(L)
  101.  44   CONTINUE
  102.       U=U*COEF(IC)
  103.       AF1(KE,J,I)=AF1(KE,J,I)+U*AIMPL
  104.       AS1(KE,J,I)=AS1(KE,J,I)+U*(AIMPL-1.D0)
  105.  43   CONTINUE
  106.       ENDIF
  107.       ENDIF
  108.       GO TO 108
  109.  
  110.  5    CONTINUE
  111.       DO 51 L=1,NPG
  112.       C=0.D0
  113.       DO 52 I=1,NP
  114.       IU=LE(I,KE)
  115.       C=C+COEF(IU)*FN(I,L)
  116.  52   CONTINUE
  117.       CL(L)=C
  118.  51   CONTINUE
  119. C
  120. C     write(6,*)' BCL 15 '
  121.       DO 15 I=1,NP
  122.       DO 15 J=I,NP
  123.       U=0.D0
  124.       DO 12 L=1,NPG
  125.       V=0.D0
  126.       DO 13 N=1,IDIM
  127.       V=V+HR(N,I,L)*HR(N,J,L)
  128.  13   CONTINUE
  129.       U=U+V*PGSQ(L)*CL(L)
  130.  12   CONTINUE
  131.       DO 14 K=1,NINC
  132.       AF(J,I,K,K)=AF(J,I,K,K)+U
  133.       IF(I.NE.J)AF(I,J,K,K)=AF(I,J,K,K)+U
  134.  14   CONTINUE
  135.       IF(IAX.EQ.0)GO TO 15
  136.       U=0.D0
  137.       DO 410 L=1,NPG
  138.       U=U+FN(I,L)*FN(J,L)/RPG(L)/RPG(L)*PGSQ(L)*CL(L)
  139.  410  CONTINUE
  140.       KU=3-IAX
  141.       AF(J,I,KU,KU)=AF(J,I,KU,KU)+U
  142.       IF(I.NE.J)AF(I,J,KU,KU)=AF(I,J,KU,KU)+U
  143.  15   CONTINUE
  144.       GO TO 108
  145. C
  146. C
  147.  70   CONTINUE
  148.       IF(NINC.NE.1)CALL ARRET(0)
  149.       IF(KJTT.EQ.4)CALL ARRET(0)
  150.       IC=1+(1-IK1)*(NK-1)
  151.       DO 71 I=1,NP
  152.       DO 71 J=1,NP
  153.       U=0.D0
  154.       DO 72 L=1,NPG
  155.       UL=0.D0
  156.       DO 73 M=1,IDIM
  157.       UN=0.D0
  158.       DO 74 N=1,IDIM
  159.       UN=UN+COES(N,M,IC)*HR(N,J,L)
  160.  74   CONTINUE
  161.       UL=UL+UN*HR(M,I,L)
  162.  73   CONTINUE
  163.       U=U+UL*PGSQ(L)
  164.  72   CONTINUE
  165.       AF(J,I,1,1)=AF(J,I,1,1)+U
  166.  71   CONTINUE
  167.       GO TO 108
  168. C
  169. C
  170.    80 CONTINUE
  171.       IF(KJTT.EQ.4)CALL ARRET(0)
  172.       IF(IARG.EQ.0)CALL ARRET(0)
  173.       IC=1+(1-IK1)*(NK-1)
  174.       ALT=COG(1,IC)-COG(2,IC)
  175.       AT=COG(2,IC)
  176. C     WRITE(6,*)' ALT ',ALT,' AT ',AT
  177.       IF(IDIM.EQ.2) THEN
  178.       CO=SQRT(S2(1,IC)*S2(1,IC)+S2(2,IC)*S2(2,IC))
  179.       ENDIF
  180.       IF(IDIM.EQ.3) THEN
  181.       CO=SQRT(S2(1,IC)*S2(1,IC)+S2(2,IC)*S2(2,IC)+S2(3,IC)*S2(3,IC))
  182.       ENDIF
  183. C     WRITE(6,*)' CO ',CO
  184. C     WRITE(6,*)' S2 ',S2(1,IC),S2(2,IC)
  185.       IF(IDIM.EQ.2) THEN
  186.       COE(1,1)=ALT*S2(1,IC)*S2(1,IC)/CO+AT*CO
  187.       COE(2,2)=ALT*S2(2,IC)*S2(2,IC)/CO+AT*CO
  188.       COE(1,2)=ALT*S2(1,IC)*S2(2,IC)/CO
  189.       COE(2,1)=COE(1,2)
  190.       ENDIF
  191.       IF(IDIM.EQ.3) THEN
  192.       COE(1,1)=ALT*S2(1,IC)*S2(1,IC)/CO+AT*CO
  193.       COE(2,2)=ALT*S2(2,IC)*S2(2,IC)/CO+AT*CO
  194.       COE(3,3)=ALT*S2(3,IC)*S2(3,IC)/CO+AT*CO
  195.       COE(1,2)=ALT*S2(1,IC)*S2(2,IC)/CO
  196.       COE(1,3)=ALT*S2(1,IC)*S2(3,IC)/CO
  197.       COE(2,3)=ALT*S2(2,IC)*S2(3,IC)/CO
  198.       COE(2,1)=COE(1,2)
  199.       COE(3,1)=COE(1,3)
  200.       COE(3,2)=COE(2,3)
  201.       ENDIF
  202. C     WRITE(6,*)' COE ',COE(1,1),COE(1,2),COE(1,3)
  203. C     WRITE(6,*)'     ',COE(2,1),COE(2,2),COE(2,3)
  204. C     WRITE(6,*)'     ',COE(3,1),COE(3,2),COE(3,3)
  205. C
  206.       DO 81 I=1,NP
  207.       DO 81 J=1,NP
  208.       U=0.D0
  209.       DO 82 L=1,NPG
  210.       UL=0.D0
  211.       DO 83 M=1,IDIM
  212.       UN=0.D0
  213.       DO 84 N=1,IDIM
  214.       UN=UN+COE(N,M)*HR(N,J,L)
  215.  84   CONTINUE
  216.       UL=UL+UN*HR(M,I,L)
  217.  83   CONTINUE
  218.       U=U+UL*PGSQ(L)
  219.  82   CONTINUE
  220.       AF(J,I,1,1)=AF(J,I,1,1)+U
  221.  81   CONTINUE
  222.  
  223.  108  CONTINUE
  224.  
  225. C     write(6,*)' FIN XLAPL'
  226.  
  227.       RETURN
  228.  1002 FORMAT(10(1X,1PE11.4))
  229.  1001 FORMAT(20(1X,I5))
  230.       END
  231.  
  232.  
  233.  
  234.  
  235.  
  236.  
  237.  
  238.  

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