Télécharger la plaquette Cast3M


fr en

Télécharger profch.eso

Retour à la liste

Numérotation des lignes : 

profch

  1. C PROFCH    SOURCE    PV        16/11/17    22:01:15     9180           
  2.       SUBROUTINE PROFCH(MATRAK,IMPR,KPIMP,NIMP)
  3. C************************************************************************
  4. C
  5. C     OBJET :
  6. C        CE SOUS PROGRAMME PREPARE L'INVERSION DE LA MATRICE DE PRESSION
  7. C
  8. C        KTYPI= 1 OU 5, METHODE DIRECTE
  9. C
  10. C                       1/ IOPTI = 0 PAS DE RENUMEROTATION
  11. C                          IOPTI = 1 RENUMEROTATION DE LA MATRICE
  12. C                                    DE PRESSION
  13. C                       2/ CALCUL DU PROFIL ET DE LA TAILLE
  14. C                          DE LA DEMI-MATRICE
  15. C
  16. C        MATRAK : objet de type MATRAK
  17. C                 Si KIZCL=0 On cree IZL
  18. C                 Sinon rien
  19. C        IMPR 1 ou 2 niveau d'impressions
  20. C        KPIMP = 0 pas de pression imposee
  21. C        KPIMP = 1 pression imposée en un point (le dernier)
  22. C        KPIMP = 2 pression moyenne imposée
  23. C        NIMP      numero de la pression imposée
  24. C
  25. C************************************************************************
  26.       IMPLICIT INTEGER(I-N)
  27.       IMPLICIT REAL*8 (A-H,O-Z)
  28.  
  29. C***
  30. C-INC SMMATRAKANC
  31. C*************************************************************************
  32. C
  33. C    REPERAGE ET STOKAGE DES MATRICES ELEMENTAIRES  puis assemblees
  34. C
  35.  
  36. * LGEOC SPG de la pression et/ou des multiplicateurs de Lagrange
  37. * (points CENTRE ) pour chaque operateur de contrainte
  38. * KGEOC SPG pour la totalite des points CENTRE.
  39. * KGEOS SPG pour la totalite des points SOMMET (Diagonale vitesse)
  40. * KLEMC Connectivites de l'ensemble des contraintes
  41. * LIZAFM(NBSOUS) contient les pointeurs IZAFM des sous-zones
  42.  
  43.       SEGMENT MATRAK
  44.       INTEGER LGEOC(NBOP),IDEBS(NBOP),IFINS(NBOP)
  45.       INTEGER LIZAFM(NBSOUS)
  46.       INTEGER IKAM0 (NBSOUS)
  47.       INTEGER IMEM  (NBELC)
  48.       INTEGER KLEMC,KGEOS,KGEOC,KDIAG,KCAC,KIZCL,KIZGC
  49.       ENDSEGMENT
  50.  
  51.       SEGMENT IZAFM
  52.       REAL*8  AM(NNELP,NP,IESP),RPGI(NELAX)
  53.       ENDSEGMENT
  54.  
  55.       POINTEUR IPMJ.IZAFM,IPMK.IZAFM
  56.  
  57. C*******************************************************************
  59. -INC SMCOORD
  60. -INC SMLENTI
  61.       POINTEUR IZIPAD.MLENTI
  62. -INC SMELEME
  63.       POINTEUR IZK.MELEME
  64. -INC SMTABLE
  65.       POINTEUR MTABP.MTABLE
  66.       PARAMETER (IBLSIZ=50000)
  67.  
  68.       SEGMENT IZTRAV
  69.          INTEGER ITAB(NJA)
  70.       ENDSEGMENT
  71. C
  72. C IVK = Tableau de travail
  73.       PARAMETER (MAXIVK=100)
  74.       INTEGER IVK(MAXIVK)
  75. C***
  76.       CHARACTER*8 TYPE
  77.       CHARACTER*4 LMOT(1)
  78.       DATA LMOT/'IMPR'/
  79.       DATA IOPTI/1/
  80. C***
  81.  
  82.       SEGACT MATRAK*MOD
  83.       IF(KIZCL.NE.0)THEN
  84.       SEGDES MATRAK
  85.       RETURN
  86.       ENDIF
  87.  
  88.       KGEO=KGEOC
  89.       CALL KNBEL(KGEO,NBELC)
  90.       MELEME=KLEMC
  91.  
  92. C
  93. C--- CALNOE renvoie dans IZTKB un tableau contenant pour chaque noeud
  94. C    du domaine la liste des {l{ments poss{dant ce noeud. Les noeuds
  95. C    sont rep{r{s dans la num{rotation locale(domaniale), les {l{ments
  96. C    sont pris dans l'ordre fourni par la num{rotation naturelle.
  97. C    Il est possible, partant de ce tableau, de construire pour chaque
  98. C    {l{ment la liste des {l{ments voisins. Ceci est tr}s utile pour
  99. C    le calcul de la matrice de pression ou Akl.ne.0 si k et l sont
  100. C    deux {l{ments voisins.
  101. C
  102. C    Cette technique n'est pas (encore) utilis{e si KTYPI = 1 ou 5.
  103. C    Les maillages trait{s dans ce cas {tant "petits", ce n'est pas
  104. C    trop grave. Par contre lorsque KTYPI=2,3,4, c'est @ dire pour les
  105. C    tr}s gros maillages, son emploi est obligatoire.
  106. C
  107. C    Cet appel ne doit etre fait que si on stocke la matrice de
  108. C    pression (KTPI.LE.5).
  109. C
  110.       CALL KALNO0(MELEME,IZK,MLENTI,IZIPAD,IRET)
  111.       IF(IRET.EQ.0)GO TO 90
  112.  
  113. C---  Activation de tous les MELEMEs du domaine
  114.  
  115.       SEGACT MELEME,IZK,MLENTI
  116.       NBSOUS=LISOUS(/1)
  117.       IF(NBSOUS.EQ.0)NBSOUS=1
  118.       DO 400 NS=1,NBSOUS
  119.       IF(NBSOUS.EQ.1)IPT1=MELEME
  120.       IF(NBSOUS.NE.1)IPT1=LISOUS(NS)
  121.       SEGACT IPT1
  122.  400  CONTINUE
  123.  
  124. C
  125. C************************************************************************
  126. C     C H O L E S K I
  127. C************************************************************************
  128. C
  129.       NJA=NBELC
  130.       NJAN=NBELC
  131.       NJAB=NBELC
  132.  
  133.       SEGINI IZL,IZTRAV
  134.       KIZCL=IZL
  135.  
  136.       LONG=0
  137.  
  138.       CALL INITI(NUAN,NBELC,0)
  139.       CALL INITI(NUNA,NBELC,0)
  140.       NUAN(1)=1
  141.       NUNA(1)=1
  142.  
  143. C********* ON STOKE LES POINTEURS DES MELEME
  144.  
  145.       KI=1
  146.       DO 40 NS=1,NBSOUS
  147.       IF(NBSOUS.EQ.1)IPT1=MELEME
  148.       IF(NBSOUS.NE.1)IPT1=LISOUS(NS)
  149.  
  150.       NEL=IPT1.NUM(/2)
  151.       IPOINT=IPT1
  152.       CALL INITI(IMEL(KI),NEL,IPOINT)
  153.  
  154.       DO 41 I=1,NEL
  155.       IMJ(KI+I-1)=I
  156.  41   CONTINUE
  157.  
  158.       KI=KI+NEL
  159.  
  160.  40   CONTINUE
  161.  
  162. C*********
  163.  
  164.       IF(IOPTI.EQ.1)GO TO 33
  165.  
  166. C********* ON NE RENUMEROTE PAS
  167.  
  168.       IF(IMPR.GE.2)WRITE(6,*)' PAS DE RENUMEROTATION DE LA PRESSION '
  169.       DO 32 I=1,NBELC
  170.       NUAN(I)=I
  171.       NUNA(I)=I
  172.  32   CONTINUE
  173.       GO TO 34
  174.  
  175. C********* OPTIMISATION DE LA NUMEROTAION DE LA MATRICE DE PRESSION
  176.  33   CONTINUE
  177.       IF(IMPR.GE.2)WRITE(6,*)' ON RENUMEROTE LA MATRICE DE PRESSION '
  178. C*********
  179.  
  180. C     write(6,*)' SUB PROFCH : IPADL :'
  181. C     nbp=ipadl(/1)
  182. C     write(6,1001)(IPADL(II),II=1,NBP)
  183. C     im2=lect(/1)
  184. C     im1=izk.num(/1)
  185. C     write(6,*)' IM1=',IM1 ,' Taille IZK.NUM'
  186. C     do 1873 i=1,im2
  187. C     nbb=lect(i)
  188. C     write(6,*)' I=->',i,' nb=',lect(i)
  189. C     write(6,1001)(izk.num(ii,i),ii=1,nbb)
  190. C1873 continue
  191.  
  192.       KC=1
  193.  
  194.       DO 30 KK=1,NBELC
  195.  
  196.       K=NUNA(KK)
  197.  
  198.       IF(K.NE.0)GO TO 310
  199. C   ON A PLUS D'ELEMENT NOUVELLEMENT NUMEROTE
  200.       WRITE(6,*)' EST CE LA DISCONTINUITE ',K,KC,KK
  201.       KC=KC+1
  202.       DO 311 KKA=1,NBELC
  203.       IF(NUAN(KKA).EQ.0)GO TO 312
  204. C SI NUAN(KKA)=0 C'EST QUE L'ELEMENT KKA N'A PAS ETE NUMEROTE
  205.  311  CONTINUE
  206.       WRITE(6,*)' IL YA UN PROBLEME DE NUMEROTATION'
  207.       CALL ARRET(0)
  208.  312  CONTINUE
  209.       NUAN(KKA)=KC
  210.       NUNA(KC)=KKA
  211.       K=NUNA(KC)
  212.  310  CONTINUE
  213.  
  214.       IPT1=IMEL(K)
  215.       NP=IPT1.NUM(/1)
  216.       DO 31 I=1,NP
  217.       NU=IPT1.NUM(I,IMJ(K))
  218. C     NU=IPT1.NUM(I,K)
  219.       NU=IZIPAD.LECT(NU)
  220.       NELV=LECT(NU)
  221.       DO 31 NN=1,NELV
  222.       KN=IZK.NUM(NN,NU)
  223.  
  224.       IF(NUAN(KN).NE.0)GO TO 31
  225.       KC=KC+1
  226.       NUAN(KN)=KC
  227.       NUNA(KC)=KN
  228.  31   CONTINUE
  229.  
  230.  30   CONTINUE
  231.  
  232.       DO 35 I=1,NBELC
  233.       ITAB(I)=IMJ(NUNA(I))
  234.  35   CONTINUE
  235.       CALL RSETI(IMJ,ITAB,NBELC)
  236.  
  237.       DO 36 I=1,NBELC
  238.       ITAB(I)=IMEL(NUNA(I))
  239.  36   CONTINUE
  240.       CALL RSETI(IMEL,ITAB,NBELC)
  241.  
  242.       SEGSUP IZTRAV
  243. C             **********************************************
  244.  34   CONTINUE
  245.  
  246.       IF(IMPR.GE.2)THEN
  247.       WRITE(6,1929)
  248.  1929 FORMAT(/10X,' SUB PROFIM : RENUMEROTAION ZONE PRESSION ')
  249.       WRITE(6,*)' CORRESPONDANCE NUMER. ANC. --> NUMER. NOUVELLE '
  250.       WRITE(6,1930)(N,NUAN(N),N=1,NBELC)
  251.       WRITE(6,*)' CORRESPONDANCE NUMER. NOUVELLE --> NUMER. ANC. '
  252.       WRITE(6,1930)(N,NUNA(N),N=1,NBELC)
  253.  1930 FORMAT(7(1X,'**',I5,3X,I5,1X,'**'))
  254.       ENDIF
  255. C
  256.       DO 10 KK=1,NBELC
  257.       K=IMJ(KK)
  258.       IPT1=IMEL(KK)
  259.       NP=IPT1.NUM(/1)
  260.       K0=KK
  261.       DO 11 I=1,NP
  262.       NU=IPT1.NUM(I,K)
  263.       NU=IZIPAD.LECT(NU)
  264.       NELV=LECT(NU)
  265.       DO 11 NN=1,NELV
  266.       KN=IZK.NUM(NN,NU)
  267.       KN=NUAN(KN)
  268.       IF(KN.LT.K0)K0=KN
  269.  11   CONTINUE
  270.       IF(KK.NE.NIMP.OR.KPIMP.NE.2)THEN
  271.       KZA(KK)=KK-K0+1
  272.       ELSE
  273.       KZA(KK)=NBELC-1
  274.       ENDIF
  275.       LONG=LONG+KZA(KK)
  276.  10   CONTINUE
  277.  
  278. C
  279. C La matrice est stockée par blocs de plusieurs lignes. En effet, il faut
  280. C l'équivalent de 900 opérations pour activer/désactiver un segment. Si
  281. C une ligne de la matrice comprend 50 éléments, on fait 100 opérations
  282. C utiles et 900 inutiles lors de la résolution.
  283. C La taille des blocs est fixée par le paramètre IBLSIZ. La diagonale est
  284. C stockée dans un segment à part.
  285. C
  286. C Il nous faut constituer un certain nombre de tableaux. Le premier segment
  287. C décrit la matrice. Il contient un tableau donnant les numéros des premières
  288. C lignes de chaque bloc et un tableau de pointeurs sur des segments
  289. C descripteurs de bloc.
  290. C
  291. C Il faut d'abord déterminer le nombre de blocs.
  292. C
  293.  
  294.       NBLK=1
  295.       IF(LONG-NBELC.GT.IBLSIZ) THEN
  296.       ILBL=0
  297.       DO 12 KK=2,NBELC
  298.          LA=KZA(KK)-1
  299.          IF(LA.GT.IBLSIZ) THEN
  300.             NBLK=NBLK+2
  301.             ILBL=0
  302.             GOTO 12
  303.          ENDIF
  304.          ILBL=ILBL+LA
  305.          IF(ILBL.GT.IBLSIZ) THEN
  306.             NBLK=NBLK+1
  307.             ILBL=LA
  308.             GOTO 12
  309.          ENDIF
  310.  12   CONTINUE
  311.       ENDIF
  312. C
  313. C On peut allouer le segment descripteur de la matrice. On va pouvoir
  314. C remplir des numéros de ligne des débuts de bloc.
  315. C
  316.       SEGINI IDMAT
  317.       IF(NBLK.EQ.1) THEN
  318.          NLDBLK(1)=2
  319.          NLDBLK(2)=NBELC+1
  320.       ELSE
  321.          IBLK=1
  322.          ILBL=0
  323.          NLDBLK(1)=2
  324.          NLDBLK(NBLK+1)=NBELC+1
  325.          DO 13 KK=2,NBELC
  326.             LA=KZA(KK)-1
  327.             IF(LA.GT.IBLSIZ) THEN
  328.               NLDBLK(IBLK+1)=KK
  329.               NLDBLK(IBLK+2)=KK+1
  330.               IBLK=IBLK+2
  331.               ILBL=0
  332.               GOTO 13
  333.             ENDIF
  334.             ILBL=ILBL+LA
  335.             IF(ILBL.GT.IBLSIZ) THEN
  336.               NLDBLK(IBLK+1)=KK
  337.               ILBL=LA
  338.               IBLK=IBLK+1
  339.               GOTO 13
  340.             ENDIF
  341.  13      CONTINUE
  342.       ENDIF
  343. C
  344. C On peut maintenant remplir le segment descripteur de chaque bloc.
  345. C
  346.       DO 15 IBLK=1,NBLK
  347.       NLBLK=NLDBLK(IBLK+1)-NLDBLK(IBLK)
  348.       SEGINI IDBLK
  349.       IDESCR(IBLK)=IDBLK
  350.       ILON=0
  351.       IDEBLK(1)=1
  352.       DO 16 KK=NLDBLK(IBLK),NLDBLK(IBLK+1)-1
  353.          LA=KZA(KK)-1
  354.          ILON=ILON+LA
  355.          IDEBLK(KK-NLDBLK(IBLK)+2)=IDEBLK(KK-NLDBLK(IBLK)+1)+LA
  356.  16   CONTINUE
  357.       SEGDES IDBLK
  358.  15   CONTINUE
  359.       SEGDES IDMAT
  360. C
  361. C On stocke dans KZA1 la référence au segment IDMAT à partir duquel tout se
  362. C retrouve.
  363. C
  364. C     write(6,*)' sub PROFCH on met IDMAT -> KZA1 ',IDMAT
  365.       KZA1=IDMAT
  366.  14   CONTINUE
  367. C
  368.       IF(IMPR.GE.2)WRITE(6,*)(KZA(M),M=1,NBELC),LONG
  369. C
  370.       IF(IMPR.GE.1)WRITE(6,200)LONG
  371.       SEGDES IZL
  372.       SEGSUP IZK,IZIPAD,MLENTI
  373.  
  374. C---  D{sactivation de tous les MELEMEs
  375.  
  376.       DO 51 NS=1,NBSOUS
  377.       IF(NBSOUS.EQ.1)IPT1=MELEME
  378.       IF(NBSOUS.NE.1)IPT1=LISOUS(NS)
  379.       SEGDES IPT1
  380.  51   CONTINUE
  381.       SEGDES MELEME
  382.  
  383.       RETURN
  384.  90   CONTINUE
  385.       WRITE(6,*)' Retour anormal de PROFCH'
  386.       RETURN
  387.  
  388. C************************************************************************
  389.  100  FORMAT(40X,1H*,9X,'.... TAILLE DE LA MATRICE DE PRESSION ',I7,
  390.      *' MOTS',13X,1H*/40X,1H*,72X,1H*)
  391.  200  FORMAT(10X,1H*,9X,'.... TAILLE DE LA MATRICE DE PRESSION ',I7,
  392.      *' MOTS',13X,1H*)
  393.  1001 FORMAT(20(1X,I5))
  394.       END
  395.  
  396.  
  397.  
  398.  
  399.  
  400.  

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