1 : C RLENCT SOURCE CHAT 05/01/13 03:01:32 5004 2 : SUBROUTINE RLENCT(MELFL,MELSOM,MLEPOI,MLECOE,MLEPOF,MLECOF) 3 : C************************************************************************ 4 : C 5 : C PROJET : CASTEM 2000 6 : C 7 : C NOM : RLENCT 8 : C 9 : C DESCRIPTION : Appelle par GRADI2 10 : C 11 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI) 12 : C 13 : C AUTEUR : A. BECCANTINI 14 : C 15 : C************************************************************************ 16 : C 17 : C 18 : C This subroutine computes the coefficients to compute the gradient 19 : C at intefaces with respect to the values on its neighbors. 20 : C The neighbors are 'CENTRE' points or 'boundary condition' points. 21 : C 22 : C**** Inputs: 23 : C 24 : C MELFL = 'FACEL' meleme 25 : C 26 : C MELSOM = 'SOMMET' meleme 27 : C 28 : C MLEPOI = list of integers. 29 : C MLEPOI.LECT(i) points to the list neighbors of 30 : C MELSOM.NUM(1,I). Neighbors are 'CENTRE' points or 31 : C 'boundary condition' points 32 : C MLECOE = list of integers. 33 : C MLECOE.LECT(i) points to the list of real of coeffients 34 : C to compute the vertex values 35 : C 36 : C MLEPOF = list of integers. 37 : C MLEPOI.LECT(i) points to the list neighbors of 38 : C MELFL.NUM(2,i). Neighbors are 'CENTRE' points or 39 : C 'SOMMET' points. 40 : C MLECOF = list of integers. 41 : C MLECOE.LECT(i) points to the matrix of coeffients to 42 : C compute the gradient with respect the neighbors value 43 : C 44 : C**** Output: 45 : C 46 : C MLEPOF = list of integers. 47 : C MLEPOI.LECT(i) points to the list neighbors of 48 : C MELFL.NUM(2,i). Neighbors are 'CENTRE' points or 49 : C 'boundary condition' points. The first one is the 50 : C 'FACE' point itself. 51 : C MLECOE = list of integers. 52 : C MLECOE.LECT(i) points to the matrix of coeffients to 53 : C compute the gradient 54 : C 55 : C**** Variables de COOPTIO 56 : C 57 : C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR 58 : C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ 59 : C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES 60 : C & ,IECHO, IIMPI, IOSPI 61 : C & ,IDIM 62 : CC & ,MCOORD 63 : C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE 64 : C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU 65 : C & ,NORINC,NORVAL,NORIND,NORVAD 66 : C & ,NUCROU, IPSAUV, IFICLE, IPREFI 67 : CC 68 : IMPLICIT INTEGER(I-N) 69 : -INC CCOPTIO 70 : -INC SMCOORD 71 : -INC SMLENTI 72 : -INC SMLREEL 73 : -INC SMELEME 74 : C 75 : INTEGER NTP, NFAC, IFAC, NVOIF, IVOIF, NGV, NLS, NLV 76 : & ,NVOIS,IVOIS,NGVS, LAST, LAST0, NSOMM, IPOS 77 : & ,I1,ICELL,NGF 78 : REAL*8 CELL 79 : INTEGER JG 80 : INTEGER N1,N2 81 : SEGMENT MATRIX 82 : REAL*8 MAT(N1,N2) 83 : ENDSEGMENT 84 : POINTEUR MELSOM.MELEME, MLEPOI.MLENTI,MLECOE.MLENTI, MELFL.MELEME 85 : & ,MLEPOF.MLENTI,MLECOF.MLENTI 86 : & ,MATCOE.MATRIX,MATCO1.MATRIX,MLRCOE.MLREEL,MLREST.MLENTI 87 : & ,MLESOM.MLENTI,MLPOSI.MLENTI 88 : C 89 : C 90 : SEGACT MELFL 91 : C 92 : NTP=MCOORD.XCOOR(/1)/(IDIM+1) 93 : C 94 : C**** Chaining list 95 : C LAST 96 : C MLREST(NTP) 97 : C 98 : JG=NTP 99 : SEGINI MLREST 100 : LAST=-1 101 : C 102 : C**** Position of a point in the list of neighbors 103 : C MLREST is used to clean it at the end 104 : C 105 : JG=NTP 106 : SEGINI MLPOSI 107 : C 108 : C**** We create the MLENTI for the sommets 109 : C 109 : CALL KRIPAD(MELSOM,MLESOM) 110 : IF(IERR .NE. 0)GOTO 9999 111 : C En KRIPAD 112 : C SEGACT MELSOM 113 : C SEGINI MLESOM 114 : C 115 : SEGACT MLEPOF*MOD 116 : NFAC=MLEPOF.LECT(/1) 117 : SEGACT MLECOF*MOD 118 : C 119 : SEGACT MLEPOI 120 : SEGACT MLECOE 121 : C 122 : NSOMM=MLEPOI.LECT(/1) 123 : DO I1=1,NSOMM,1 124 : MLENTI=MLEPOI.LECT(I1) 125 : SEGACT MLENTI 126 : MLREEL=MLECOE.LECT(I1) 127 : SEGACT MLREEL 128 : ENDDO 129 : C 130 : DO IFAC=1,NFAC,1 131 : NGF=MELFL.NUM(2,IFAC) 132 : MLENT1=MLEPOF.LECT(IFAC) 133 : SEGACT MLENT1 134 : NVOIF=MLENT1.LECT(/1) 135 : C 136 : C******* We fill MLREST, MLPOSI 137 : C 138 : LAST=-1 139 : IPOS=1 140 : MLREST.LECT(NGF)=LAST 141 : LAST=NGF 142 : DO IVOIF=1,NVOIF,1 143 : NGV=MLENT1.LECT(IVOIF) 144 : C 145 : C********** First of all, we have to check if this is a 146 : C 'SOMMET' point. In that case we have to replace 147 : C it by its neighbors. 148 : C 149 : NLS=MLESOM.LECT(NGV) 150 : C 151 : IF(NLS .GT. 0)THEN 152 : C 'SOMMET' 153 : MLENT2=MLEPOI.LECT(NLS) 154 : NVOIS=MLENT2.LECT(/1) 155 : DO IVOIS=1,NVOIS,1 156 : NGVS=MLENT2.LECT(IVOIS) 157 : NLV=MLREST.LECT(NGVS) 158 : IF(NLV .EQ. 0)THEN 159 : C 160 : C**************** New point 161 : C 162 : IPOS=IPOS+1 163 : MLREST.LECT(NGVS)=LAST 164 : LAST=NGVS 165 : ENDIF 166 : ENDDO 167 : ELSE 168 : C 'CENTRE' 169 : NLV=MLREST.LECT(NGV) 170 : IF(NLV .EQ. 0)THEN 171 : C 172 : C************* New point 173 : C 174 : IPOS=IPOS+1 175 : MLREST.LECT(NGV)=LAST 176 : LAST=NGV 177 : ENDIF 178 : ENDIF 179 : ENDDO 180 : C 181 : C********** We create the new list of neighbors 182 : C 183 : JG=IPOS 184 : SEGINI MLENTI 185 : MLEPOF.LECT(IFAC)=MLENTI 186 : LAST0=LAST 187 : DO IVOIF=1,IPOS,1 188 : I1=(IPOS-IVOIF)+1 189 : MLENTI.LECT(I1)=LAST 190 : MLPOSI.LECT(LAST)=I1 191 : LAST=MLREST.LECT(LAST) 192 : ENDDO 193 : IF(LAST .NE. -1)THEN 194 : WRITE(IOIMP,*) 'subroutine rlenct.eso' 194 : CALL ERREUR(5) 195 : ENDIF 196 : LAST=LAST0 197 : C 198 : C******* Summarizing 199 : C 200 : C MLENTI: list of the new 'FACE' neighbors 201 : C MLENT1: list of the old 'FACE' neighbors 202 : C MLENT2: is free. It has been used and it will be used 203 : C for the 'SOMMET' neighbors. 204 : C MLPOSI: position of the new neighbors into MLENTI 205 : C MLREST + LAST : chaining list, used to clean MLPOSI 206 : C 206 : C******* Let us call 207 : C MATCOE: matrix of the 'FACE' coeff (IDIM+1,*) 208 : C MATCO1: matrix of the old 'FACE' coeff. (IDIM+1,*) 209 : C MLRCOE: list of the 'SOMMET coeff 210 : C 211 : C 212 : MATCO1=MLECOF.LECT(IFAC) 213 : SEGACT MATCO1 214 : N1=IDIM+1 215 : N2=MLENTI.LECT(/1) 216 : SEGINI MATCOE 217 : MLECOF.LECT(IFAC)=MATCOE 218 : C 219 : C******* Loop on the old 'FACE' neighbors 220 : C 221 : NVOIF=MLENT1.LECT(/1) 222 : DO IVOIF=1,NVOIF,1 223 : NGV=MLENT1.LECT(IVOIF) 224 : NLS=MLESOM.LECT(NGV) 225 : C 226 : IF(NLS .GT. 0)THEN 227 : C 'SOMMET' 228 : MLENT2=MLEPOI.LECT(NLS) 229 : MLRCOE=MLECOE.LECT(NLS) 230 : NVOIS=MLENT2.LECT(/1) 231 : DO IVOIS=1,NVOIS,1 232 : NGVS=MLENT2.LECT(IVOIS) 233 : IPOS=MLPOSI.LECT(NGVS) 234 : IF(IPOS .EQ. 0)THEN 235 : WRITE(IOIMP,*) 'subroutine rlenct.eso' 235 : CALL ERREUR(5) 236 : ENDIF 237 : DO I1=1,IDIM+1,1 238 : MATCOE.MAT(I1,IPOS)=MATCOE.MAT(I1,IPOS)+ 239 : & (MATCO1.MAT(I1,IVOIF)*MLRCOE.PROG(IVOIS)) 240 : ENDDO 241 : ENDDO 242 : ELSE 243 : C 'CENTRE' 244 : IPOS=MLPOSI.LECT(NGV) 245 : DO I1=1,IDIM+1,1 246 : MATCOE.MAT(I1,IPOS)=MATCOE.MAT(I1,IPOS)+ 247 : & MATCO1.MAT(I1,IVOIF) 248 : ENDDO 249 : ENDIF 250 : ENDDO 251 : C 252 : CC 253 : CC******* Test 254 : CC 255 : C ipos=mlenti.lect(/1) 256 : C write(*,*) 'ngf=',melfl.num(2,ifac) 257 : C write(*,*) 'ntvois=',ipos 258 : C write(*,*) 'nvois=',(mlenti.lect(ivoif),ivoif=1,ipos,1) 259 : C write(*,*) 'Position=', 260 : C & (mlposi.lect(mlenti.lect(ivoif)),ivoif=1,ipos,1) 261 : C write(*,*) 'coeff(1) =',(matcoe.mat(1,ivoif),ivoif=1,ipos,1) 262 : C write(*,*) 'coeff(2) =',(matcoe.mat(2,ivoif),ivoif=1,ipos,1) 263 : C write(*,*) 'coeff(3) =',(matcoe.mat(3,ivoif),ivoif=1,ipos,1) 264 : C if(idim.eq.3) write(*,*) 'coeff(4)=', 265 : C & (matcoe.mat(4,ivoif),ivoif=1,ipos,1) 266 : C cell=0.0D0 267 : C do ivoif=1,ipos,1 268 : C cell=cell+matcoe.mat(1,ivoif) 269 : C enddo 270 : C write(*,*) 'sum=',cell 271 : C if(abs(cell-1.0d0) .gt. 1.0d-10)then 272 : CC It must be true if I just consider Dirichlet B.C. 272 : C call erreur(5) 273 : C goto 9999 274 : C endif 275 : C 276 : C******* We clean MLPOSI and MLREST 277 : C 278 : NVOIS=MLENTI.LECT(/1) 279 : DO IVOIF=1,NVOIS,1 280 : MLPOSI.LECT(LAST)=0 281 : ICELL=LAST 282 : LAST=MLREST.LECT(ICELL) 283 : MLREST.LECT(ICELL)=0 284 : ENDDO 285 : IF(LAST .NE. -1)THEN 286 : WRITE(IOIMP,*) 'subroutine rlenct.eso' 286 : CALL ERREUR(5) 287 : ENDIF 288 : C 289 : SEGSUP MATCO1 290 : SEGSUP MLENT1 291 : SEGDES MATCOE 292 : SEGDES MLENTI 293 : C 294 : ENDDO 295 : CC 296 : CC******* Test 297 : CC 298 : C do ifac=1,nfac,1 299 : C mlenti=mlepof.lect(ifac) 300 : C matcoe=mlecof.lect(ifac) 301 : C segact mlenti 302 : C segact matcoe 303 : C ipos=mlenti.lect(/1) 304 : C write(*,*) 'ngf=',melfl.num(2,ifac) 305 : C write(*,*) 'ntvois=',ipos 306 : C write(*,*) 'nvois=',(mlenti.lect(ivoif),ivoif=1,ipos,1) 307 : C write(*,*) 'coeff(1) =',(matcoe.mat(1,ivoif),ivoif=1,ipos,1) 308 : C write(*,*) 'coeff(2) =',(matcoe.mat(2,ivoif),ivoif=1,ipos,1) 309 : C write(*,*) 'coeff(3) =',(matcoe.mat(3,ivoif),ivoif=1,ipos,1) 310 : C if(idim.eq.3) write(*,*) 'coeff(4)=', 311 : C & (matcoe.mat(4,ivoif),ivoif=1,ipos,1) 312 : C cell=0.0D0 313 : C do ivoif=1,ipos,1 314 : C cell=cell+matcoe.mat(1,ivoif) 315 : C enddo 316 : C write(*,*) 'sum=',cell 317 : C if(abs(cell-1.0d0) .gt. 1.0d-10)then 318 : CC It must be true if I just consider Dirichlet B.C. 318 : C call erreur(5) 319 : C goto 9999 320 : C endif 321 : C segdes mlenti 322 : C segdes matcoe 323 : C enddo 324 : C 325 : SEGDES MELFL 326 : C 327 : SEGSUP MLREST 328 : SEGSUP MLPOSI 329 : C 330 : SEGDES MELSOM 331 : SEGSUP MLESOM 332 : C 333 : SEGDES MLEPOF 334 : SEGDES MLECOF 335 : C 336 : NSOMM=MLEPOI.LECT(/1) 337 : DO I1=1,NSOMM,1 338 : MLENTI=MLEPOI.LECT(I1) 339 : SEGSUP MLENTI 340 : MLREEL=MLECOE.LECT(I1) 341 : SEGSUP MLREEL 342 : ENDDO 343 : SEGSUP MLEPOI 344 : SEGSUP MLECOE 345 : C 346 : 9999 CONTINUE 347 : RETURN 348 : END 349 : 350 : 351 :
© Cast3M 2003 - All rights reserved.
Disclaimer