Télécharger la plaquette Cast3M


fr en

Télécharger kaxk.eso

Retour à la liste

Numérotation des lignes : 

kaxk

  1. C KAXK      SOURCE    MB234859  26/09/04    21:15:02     12638          
  2.       SUBROUTINE KAXK(A1,A2,OBS,NOBS,NT,NP0,NG0,FF,KIMP,EXTINC,RAD)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8 (A-H,O-Z)
  5. C*********************************************************************
  6. C     CALCUL DE  S1.F12   EN TENANT COMPTE DES OBSTRUCTEURS
  7. C    entree
  8. C     A1    :  COORDONNEES FACE 1
  9. C     A2    :  COORDONNEES FACE 2
  10. C     OBS   :  COORDONNES DES OBSTRUCTEURS POTENTIELS
  11. C     NOBS  :  NOMBRE D'OBSTRUCTEURS POTENTIELS
  12. C     NG0   :  NOMBRE DE POINTS DE GAUSS (cas standard)
  13. C     NP0   :  NOMBRE DE POINTS D'INTEGRATION (elements proches)
  14. C     KIMP  :  parametre d'impression
  15. C     EXTINC:  coefficient d'extinction de la cavite si absorbante
  16. C     RAD   :  dimension du pb (le calcul est fait en coor. reduites)
  17. C    resultat
  18. C     FF    : S1.F12
  19. C*********************************************************************
  20.       DIMENSION A1(2,2),A2(2,2),OBS(2,NT)
  21.       DIMENSION AL(2,20),BL(2,20)
  22.       DIMENSION AG(11,10),YA(10),HA(10),YB(10),HB(10)
  23.  
  24. C    LES INTERVALLES D INTEGRATION SONT AL(I,I+1),I=1,NAL
  25. C
  26.       NM=20
  27. C
  28. C     Initialisation des tableaux
  29.       AL=0.D0
  30.       BL=0.D0
  31.       YA=0.D0
  32.       HA=0.D0
  33.       YB=0.D0
  34.       HB=0.D0
  35.  
  36. C     estimation du mode d'integration
  37.       NS=2
  38.       CALL KAXDIS(A1,A2,NS,KIMP,NG0,NP0,NG,NP)
  39.  
  40.       IF(KIMP.GE.3) write(6,*) ' kaxk NG NP ',NG,NP
  41.  
  42.       RI1=A1(1,1)
  43.       ZI1=A1(2,1)
  44.       RI2=A1(1,2)
  45.       ZI2=A1(2,2)
  46.  
  47.       RJ1=A2(1,1)
  48.       ZJ1=A2(2,1)
  49.       RJ2=A2(1,2)
  50.       ZJ2=A2(2,2)
  51.  
  52.       DRI=RI2-RI1
  53.       DRJ=RJ2-RJ1
  54.       DZI=ZI2-ZI1
  55.       DZJ=ZJ2-ZJ1
  56.  
  57.  
  58. C>> MODE D INTEGRATION
  59.       IF(NG.EQ.0) THEN
  60.  
  61.       NA = NP
  62.       NB = NP
  63.  
  64. C>>   INTEGRATION SUR I : A
  65.  
  66.       FF=0.D0
  67.       DA=1./NA
  68.  
  69.       DO 3  IA=1,NA
  70.  
  71.       A = DA/2. + DA*(IA-1)
  72.       RI=(1.-A)*RI1+A*RI2
  73.       ZI=(1.-A)*ZI1+A*ZI2
  74.       DA=1./NA
  75.  
  76. C>>   INTEGRATION SUR J : B
  77.  
  78.       F=0.D0
  79.       DB=1./NB
  80.  
  81.       DO 30 IB=1,NB
  82.  
  83.       G=0.D0
  84.       B = DB/2. + DB*(IB-1)
  85.       RJ=(1.-B)*RJ1+B*RJ2
  86.       ZJ=(1.-B)*ZJ1+B*ZJ2
  87.  
  88.       IF(KIMP.GE.4)WRITE(6,*) ' INTEGRATION IA IB ',IA,IB
  89.  
  90. C>>   LIMITES DE VISIBILITE PROPRE AUX POINTS I ET J
  91. C     ----------------------------------------------
  92.  
  93.       CALL KAVOWN(RI,ZI,RJ,ZJ,DRI,DZI,DRJ,DZJ,KVU,NM,NAL,AL,KIMP)
  94.       IF(KIMP.GE.4) THEN
  95.          WRITE(6,*) ' VISIBILITE PROPRE ',KVU
  96.          IF(KVU.NE.0)CALL UTPRIN(AL,2,NAL)
  97.       ENDIF
  98. C
  99. C>>   TRAITEMENT DES OBSTRUCTEURS
  100. C    ----------------------------
  101.       IF(KVU.NE.0.AND.NOBS.NE.0) THEN
  102.        CALL KAVOTH(RI,ZI,RJ,ZJ,OBS,NOBS,NT,KVU,NAL,NM,AL,BL,KIMP)
  103.       IF(KIMP.GE.4) THEN
  104.          WRITE(6,*) ' OBSTRUCTEURS      ',KVU
  105.          IF(KVU.NE.0)CALL UTPRIN(AL,2,NAL)
  106.       ENDIF
  107.       ENDIF
  108. C
  109. C>>   CALCUL
  110. C    -------
  111.  
  112.       IF(KVU.NE.0) THEN
  113.       CALL KATETA(RI,ZI,RJ,ZJ,DRI,DZI,DRJ,DZJ,NM,NAL,AL,G,KIMP
  114.      &           ,EXTINC,RAD)
  115.       ENDIF
  116.  
  117.        F= F + 4.*G*DB
  118.        IF(KIMP.GE.4)WRITE(6,*) ' IA IB G F ',IA,IB,G,F
  119.  
  120.   30  CONTINUE
  121.  
  122.       FF = FF + F*DA
  123.       IF(KIMP.GE.4) WRITE(6,*) ' IA  FF ',IA,FF
  124.  
  125.   3   CONTINUE
  126.  
  127.       IF(KIMP.GE.4) WRITE(6,*) ' TOTAL FF ',FF
  128.  
  129.       ELSE
  130.  
  131. C>> POINTS DE GAUSS
  132.  
  133.       NA = 1
  134.       NB = 1
  135.       NGA= NG
  136.       NGA2=(NGA+1)/2
  137.       NGB= NG
  138.       NGB2=(NGB+1)/2
  139.       CALL MATG(AG)
  140.  
  141.       IF (AG(1,NGA).LT.1.E-5) THEN
  142.  
  143.       YA(1)=AG(1,NGA)
  144.       HA(1)=AG(2,NGA)
  145.          IF(NGA2.GE.2) THEN
  146.       DO 100 I=1,NGA2-1
  147.       YA(I+1)=AG(2*I+1,NGA)
  148.       YA(NGA2+I)=-YA(I+1)
  149.       HA(I+1)=AG(2*I+2,NGA)
  150.       HA(NGA2+I)=HA(I+1)
  151. 100   CONTINUE
  152.          ENDIF
  153.  
  154.       ELSE
  155.       DO 101 I=1,NGA2
  156.       YA(I)=AG(2*I-1,NGA)
  157.       YA(NGA2+I)=-YA(I)
  158.       HA(I)=AG(2*I,NGA)
  159.       HA(NGA2+I)=HA(I)
  160. 101    CONTINUE
  161.       ENDIF
  162.  
  163.       IF (AG(1,NGB).LT.1.E-5) THEN
  164.  
  165.       YB(1)=AG(1,NGB)
  166.       HB(1)=AG(2,NGB)
  167.          IF(NGB2.GE.2) THEN
  168.       DO 200 I=1,NGB2-1
  169.       YB(I+1)=AG(2*I+1,NGB)
  170.       YB(NGB2+I)=-YA(I+1)
  171.       HB(I+1)=AG(2*I+2,NGB)
  172.       HB(NGB2+I)=HA(I+1)
  173. 200   CONTINUE
  174.          ENDIF
  175.  
  176.       ELSE
  177.       DO 201 I=1,NGB2
  178.       YB(I)=AG(2*I-1,NGB)
  179.       YB(NGA2+I)=-YB(I)
  180.       HB(I)=AG(2*I,NGB)
  181.       HB(NGA2+I)=HB(I)
  182. 201    CONTINUE
  183.       ENDIF
  184.  
  185. C>>   INTEGRATION SUR I : A
  186.  
  187.       FF=0.D0
  188.       DA=1./NA
  189.       DO 1  IA=1,NA
  190.        A = DA/2. + DA*(IA-1)
  191.        DA=1./NA
  192. C bornes
  193.       AL1=A-DA/2.
  194.       AL2=A+DA/2.
  195.  
  196. C>>   GAUSS  SUR I : A
  197.       FA=0.D0
  198.       DO 11 IGA=1,NGA
  199. C YA varie entre -1 et 1.
  200.       ALL= (YA(IGA)+1.)*(AL2-AL1)/2. + AL1
  201.       RI=(1.-ALL)*RI1+ALL*RI2
  202.       ZI=(1.-ALL)*ZI1+ALL*ZI2
  203.  
  204. C>>   INTEGRATION SUR J : B
  205.  
  206.       F=0.D0
  207.       DB=1./NB
  208.       DO 2  IB=1,NB
  209.        B = DB/2. + DB*(IB-1)
  210. C bornes
  211.       BL1=B-DB/2.
  212.       BL2=B+DB/2.
  213.  
  214. C>>   GAUSS  SUR J : B
  215.       FB=0.D0
  216.       DO 21 IGB=1,NGB
  217. C YB varie entre -1 et 1.
  218.       BLL=(YB(IGB)+1.)*(BL2-BL1)/2. + BL1
  219.       RJ=(1.-BLL)*RJ1+BLL*RJ2
  220.       ZJ=(1.-BLL)*ZJ1+BLL*ZJ2
  221.  
  222.       G=0.D0
  223.       IF(KIMP.GE.4)WRITE(6,*) ' INTEGRATION IGA IGB ',IGA,IGB
  224.  
  225. C>>   LIMITES DE VISIBILITE PROPRE AUX POINTS I ET J
  226. C     ----------------------------------------------
  227.  
  228.       CALL KAVOWN(RI,ZI,RJ,ZJ,DRI,DZI,DRJ,DZJ,KVU,NM,NAL,AL,KIMP)
  229.       IF(KIMP.GE.4) THEN
  230.          WRITE(6,*) ' VISIBILITE PROPRE ',KVU
  231.          IF(KVU.NE.0)CALL UTPRIN(AL,2,NAL)
  232.       ENDIF
  233. C
  234. C>>   TRAITEMENT DES OBSTRUCTEURS
  235. C    ----------------------------
  236.       IF(KVU.NE.0.AND.NOBS.NE.0) THEN
  237.        CALL KAVOTH(RI,ZI,RJ,ZJ,OBS,NOBS,NT,KVU,NAL,NM,AL,BL,KIMP)
  238.       IF(KIMP.GE.4) THEN
  239.          WRITE(6,*) ' OBSTRUCTEURS      ',KVU
  240.          IF(KVU.NE.0)CALL UTPRIN(AL,2,NAL)
  241.       ENDIF
  242.       ENDIF
  243. C
  244. C>>   CALCUL
  245. C    -------
  246.  
  247.       IF(KVU.NE.0) THEN
  248.       CALL KATETA(RI,ZI,RJ,ZJ,DRI,DZI,DRJ,DZJ,NM,NAL,AL,G,KIMP
  249.      &           ,EXTINC,RAD)
  250.       ENDIF
  251.  
  252.        FB = FB + 4.*G*HB(IGB)*(BL2-BL1)/2.
  253.  
  254.   21  CONTINUE
  255.        F= F + FB*DB
  256.  
  257.   2   CONTINUE
  258.  
  259.       FA = FA + F*HA(IGA)*(AL2-AL1)/2.
  260.  
  261.   11  CONTINUE
  262.  
  263.       FF = FF + FA*DA
  264.  
  265.   1   CONTINUE
  266.  
  267.       ENDIF
  268.  
  269.       RETURN
  270.       END
  271.  
  272.  
  273.  
  274.  
  275.  
  276.  

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