Esope source files

   1 : C KSPRJS    SOURCE    MAGN      10/05/18    21:16:07     6675           
   2 :       SUBROUTINE KSPRJS
   3 :      &(MELEME,IPM1,IPM2,IPM3,IAXI,IKAS,INEFMD,KPRE,IZTGG1,IPAD,IK1)
   4 : 
   5 :       IMPLICIT INTEGER(I-N)
   6 :       IMPLICIT REAL*8 (A-H,O-Z)
   7 : C*****************************************************************************
   8 : C
   9 : C     Ce SP calcule les matrices elementaires de divergence alias C
  10 : C
  11 : C  IKAS=1 KMCT calcul de Ct  (Div U)
  12 : C  IKAS=2 KMAC calcul de C uniquement  (Grad p)
  13 : C  IKAS=3 KCCT calcul de C assemblage pour C et Ct
  14 : C
  15 : C*****************************************************************************
  16 :       CHARACTER*8 NOM0
  17 : 
  18 : -INC CCOPTIO
  19 : -INC CCGEOME
  20 : -INC SMCOORD
  21 : -INC SIZFFB
  22 :       POINTEUR IZF1.IZFFM
  23 : -INC SMMATRIK
  24 : -INC SMELEME
  25 : -INC SMLENTI
  26 : -INC SMCHPOI
  27 :       POINTEUR IZTGG1.MPOVAL
  28 : -INC CCREEL
  29 :       DIMENSION KIPM(3),XYZ1(24),UA(3),UB(3)
  30 : 
  31 : 
  32 : C OPERATEUR PRESSION
  33 : C
  34 :       MLENTI=IPAD
  35 :       DEUPI=1.D0
  36 :       IF(IAXI.NE.0)DEUPI=2.D0*XPI
  37 :       IF(IKAS.EQ.3)DEUPI=-DEUPI
  38 : 
  39 :       IF(IDIM.EQ.2)IPM3=IPM1
  40 :       KIPM(1)=IPM1
  41 :       KIPM(2)=IPM2
  42 :       KIPM(3)=IPM3
  43 :       SEGACT MELEME,IPM1*MOD,IPM2*MOD,IPM3*MOD
  44 : 
  45 :       NP=NUM(/1)
  46 :       NEL=NUM(/2)
  47 : 
  48 : C      write(6,*)' INEFMD=',inefmd,'KPRE=',kpre
  49 :       IF(KPRE.NE.2)THEN
  50 :        IF(INEFMD.EQ.3)THEN
  51 :        IF(KPRE.EQ.3)NOM0=NOMS(ITYPEL)//'PRP0'
  52 :        IF(KPRE.EQ.4)NOM0=NOMS(ITYPEL)//'PRP1'
  53 :        IF(KPRE.EQ.5)NOM0=NOMS(ITYPEL)//'PRF1'
  54 :        ELSEIF(INEFMD.EQ.2)THEN
  55 :        IF(KPRE.EQ.3)NOM0=NOMS(ITYPEL)//'MCP0'
  56 :        IF(KPRE.EQ.4)NOM0=NOMS(ITYPEL)//'MCP1'
  57 :        IF(KPRE.EQ.5)NOM0=NOMS(ITYPEL)//'MCF1'
  58 :        ELSEIF(INEFMD.EQ.4)THEN
  59 :        NOM0=NOMS(ITYPEL)//'    '
  60 :        ENDIF
  61 : C      write(6,*)' NOM0=',nom0,' ikas=',ikas,' nel=',nel
  61 :        CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  62 : 
  63 :       SEGACT IZFFM*MOD
  64 :       IZHR=KZHR(1)
  65 :       SEGACT IZHR*MOD
  66 :       NES=GR(/1)
  67 :       NPG=GR(/3)
  68 :       IZF1=KTP(1)
  69 :       SEGACT IZF1*MOD
  70 :       MP1=FN(/1)
  71 : 
  72 :       DO 30 KE=1,NEL
  73 : 
  74 :       DO 12 I=1,NP
  75 :       J=NUM(I,KE)
  76 :       DO 12 N=1,IDIM
  77 :       XYZ(N,I)=XCOOR((J-1)*(IDIM+1)    +N)
  78 :  12   CONTINUE
  79 : 
  79 :       CALL CALJBR
  80 :      &(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE,AJ,SGN)
  81 : C     write(6,*)' Retour caljbr ',mp1,np,npg
  82 : 
  83 :       DO 324 K=1,IDIM
  84 :       IPM4=KIPM(K)
  85 :       KG=(IK1-4)*(K-1)+1
  86 : 
  87 :       DO 323 M=1,MP1
  88 : 
  89 :       DO 323 I=1,NP
  90 :       U=0.D0
  91 :       DO 333 L=1,NPG
  92 : 
  92 :       CALL INITD(UA,3,0.D0)
  92 :       CALL INITD(UB,3,0.D0)
  93 : 
  94 :       DO 533 J=1,NP
  95 :       J1=LECT(NUM(J,KE))
  96 :       UA(KG)=UA(KG)+FN(J,L)*IZTGG1.VPOCHA(J1,KG)
  97 :       UB(K)=UB(K)+HR(K,J,L)*IZTGG1.VPOCHA(J1,KG)
  98 :  533  CONTINUE
  99 : 
 100 :       U=U+FN(M,L)*
 101 :      &(HR(K,I,L)*UA(KG) + FN(I,L)*UB(KG))*PGSQ(L)*DEUPI*RPG(L)
 102 :       IF(IAXI.NE.0.AND.K.EQ.1)THEN
 103 :       U=U+FN(M,L)*FN(I,L)*UA(KG)*PGSQ(L)*DEUPI
 104 :       ENDIF
 105 : 
 106 :  333  CONTINUE
 107 : 
 108 : 
 109 :       if(ikas.ne.2)then
 110 :       IPM4.AM(KE,I,M)=IPM4.AM(KE,I,M)+U
 111 :       else
 112 :       IPM4.AM(KE,M,I)=IPM4.AM(KE,M,I)+U
 113 :       endif
 114 : 
 115 :  323  CONTINUE
 116 :  324  CONTINUE
 117 : 
 118 :  37   CONTINUE
 119 :  30   CONTINUE
 120 : 
 121 :       SEGSUP IZHR,IZFFM
 122 : 
 123 : 
 124 : C CAS MACRO CENTRE
 125 : 
 126 :       ELSEIF(KPRE.EQ.2)THEN
 127 :         NOM0=NOMS(ITYPEL)//'    '
 127 :         CALL KALPBG(NOM0,'FONFORM ',IZFFM)
 128 : 
 129 :         SEGACT IZFFM*MOD
 130 :         IZHR=KZHR(1)
 131 :         SEGACT IZHR*MOD
 132 :         NES=GR(/1)
 133 :         NPG=GR(/3)
 134 :         IZF1=KTP(1)
 135 :         SEGACT IZF1*MOD
 136 :         MPG=IZF1.FN(/2)
 137 : c modif tc initialisation de M à 1 ??? (izf1.fn(/1))
 138 :         M=1
 139 : 
 140 :         NP=GR(/2)
 141 : 
 142 :         DO 40 KE=1,NEL
 143 : 
 144 :         IX=0
 145 :         DO 42 I=1,NP
 146 :         J=NUM(I,KE)
 147 :         DO 42 N=1,IDIM
 148 :         IX=IX+1
 149 :         XYZ1(IX)=XCOOR((J-1)*(IDIM+1)    +N)
 150 :  42     CONTINUE
 151 : 
 151 :       CALL CALJBR(FN,GR,PG,XYZ1,HR,PGSQ,RPG,NES,
 152 :      & IDIM,NP,NPG,IAXI,AIRE,AJ,SGN)
 153 : 
 154 :       DO 424 K=1,IDIM
 155 :       IPM4=KIPM(K)
 156 :       KG=(IK1-4)*(K-1)+1
 157 : 
 158 :       DO 423 I=1,NP
 159 :       U=0.D0
 160 :       DO 433 L=1,NPG
 161 : 
 161 :       CALL INITD(UA,3,0.D0)
 161 :       CALL INITD(UB,3,0.D0)
 162 :       DO 633 J=1,NP
 163 :       J1=LECT(NUM(J,KE))
 164 :       UA(KG)=UA(KG)+FN(J,L)*IZTGG1.VPOCHA(J1,KG)
 165 :       UB(K)=UB(K)+HR(K,J,L)*IZTGG1.VPOCHA(J1,KG)
 166 :  633  CONTINUE
 167 : 
 168 :       U=U+IZF1.FN(M,L)*
 169 :      &(HR(K,I,L)*UA(KG) + FN(I,L)*UB(KG))*PGSQ(L)*DEUPI*RPG(L)
 170 :       IF(IAXI.NE.0.AND.K.EQ.1)THEN
 171 :       U=U+IZF1.FN(M,L)*FN(I,L)*UA(KG)*PGSQ(L)*DEUPI
 172 :       ENDIF
 173 : 
 174 :  433  CONTINUE
 175 : 
 176 :       if(ikas.ne.2)then
 177 :       IPM4.AM(KE,I,1)=IPM4.AM(KE,I,1)+U
 178 :       else
 179 :       IPM4.AM(KE,1,I)=IPM4.AM(KE,1,I)+U
 180 :       endif
 181 : 
 182 :  423  CONTINUE
 183 :  424  CONTINUE
 184 : 
 185 :  40   CONTINUE
 186 : 
 187 :       SEGSUP IZHR,IZFFM
 188 : 
 189 :       ENDIF
 190 : 
 191 :       RETURN
 192 :  1002 FORMAT(10(1X,1PE11.4))
 193 :  1040 FORMAT(1X,'CALCUL MATRICE AM ',I4/10(1X,1PE11.4))
 194 :       END
 195 :  
 196 :  
 197 :  
 198 :  
 199 :  
 200 :  
 201 :  
 202 :  
 203 :  
 204 : 
 205 :  
 206 :  
 207 :  
 208 :  
 209 :  
 210 :