mt02pt.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\utils\bief\mt02pt.f
00002 !
00088                      SUBROUTINE MT02PT
00089 !                    *****************
00090 !
00091      &(T,XM,XMUL,SF,SG,SH,F,G,H,X,Y,Z,IKLE,NELEM,NELMAX,INCHYD)
00092 !
00093 !***********************************************************************
00094 ! BIEF   V6P2                                   21/08/2010
00095 !***********************************************************************
00096 !
00097 !
00098 !
00099 !
00100 !
00101 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00102 !| F              |-->| FUNCTION USED IN THE FORMULA
00103 !| FORMUL         |-->| FORMULA DESCRIBING THE RESULTING MATRIX
00104 !| G              |-->| FUNCTION USED IN THE FORMULA
00105 !| H              |-->| FUNCTION USED IN THE FORMULA
00106 !| IKLE           |-->| CONNECTIVITY TABLE.
00107 !| INCHYD         |-->| IF YES, TREATS HYDROSTATIC INCONSISTENCIES
00108 !| NELEM          |-->| NUMBER OF ELEMENTS
00109 !| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
00110 !| NPLAN          |-->| NUMBER OF PLANES IN THE MESH OF PRISMS
00111 !| SF             |-->| STRUCTURE OF FUNCTIONS F
00112 !| SG             |-->| STRUCTURE OF FUNCTIONS G
00113 !| SH             |-->| STRUCTURE OF FUNCTIONS H
00114 !| SURFAC         |-->| AREA OF 2D ELEMENTS
00115 !| T              |<->| WORK ARRAY FOR ELEMENT BY ELEMENT DIAGONAL
00116 !| X              |-->| ABSCISSAE OF POINTS
00117 !| Y              |-->| ORDINATES OF POINTS
00118 !| Z              |-->| ELEVATIONS OF POINTS
00119 !| XM             |<->| OFF-DIAGONAL TERMS
00120 !| XMUL           |-->| COEFFICIENT FOR MULTIPLICATION
00121 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00122 !
00123       USE BIEF, EX_MT02PT => MT02PT
00124       USE DECLARATIONS_TELEMAC, ONLY : TETRA
00125 !
00126       IMPLICIT NONE
00127       INTEGER LNG,LU
00128       COMMON/INFO/LNG,LU
00129 !
00130 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00131 !
00132       INTEGER, INTENT(IN)             :: NELEM,NELMAX
00133       INTEGER, INTENT(IN)             :: IKLE(NELMAX,6)
00134       DOUBLE PRECISION, INTENT(INOUT) :: T(NELMAX,6),XM(NELMAX,15)
00135       DOUBLE PRECISION, INTENT(IN)    :: XMUL
00136       DOUBLE PRECISION, INTENT(IN)    :: F(*),G(*),H(*)
00137       TYPE(BIEF_OBJ), INTENT(IN)      :: SF,SG,SH
00138       DOUBLE PRECISION, INTENT(IN)    :: X(*),Y(*),Z(*)
00139       LOGICAL, INTENT(IN)             :: INCHYD
00140 !
00141 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00142 !
00143 !     SPECIFIC DECLARATIONS
00144 !
00145       DOUBLE PRECISION X2,Y2,Z2,X3,Y3,Z3,X4,Y4,Z4
00146       DOUBLE PRECISION DIAG1,DIAG2,DIAG3,DIAG4
00147       DOUBLE PRECISION EXTR12,EXTR13,EXTR14,EXTR23,EXTR24,EXTR34
00148       INTEGER IT1,IT2,IT3,IT4,I,I1,I2,I3,I4,I5,I6,NUM1,NUM2,NUM3,NUM4
00149       INTEGER IGLOB(6),SENS(3),STO(6,6),IELEM
00150 !
00151       DOUBLE PRECISION SUR24,HTOT,VTOT,WTOT,COEF
00152       DOUBLE PRECISION T1,T3,T5,T7,T9,T11,T13,T15,T17,T19,T21,T23
00153       DOUBLE PRECISION T35,T49,T28,T42,T51,T54,SURF
00154 !
00155 !     PRISM SPLITTING (SEE EXPLANATIONS IN DECLARATIONS_TELEMAC)
00156 !
00157 !     INTEGER TETRA(2,2,2,3,4)
00158 !     DATA TETRA / 0,1,1,1,1,1,1,0,0,4,4,4,4,4,4,0,0,6,4,5,5,4,6,0,
00159 !    &             0,2,2,2,2,2,2,0,0,6,6,6,6,6,6,0,0,3,1,2,2,1,3,0,
00160 !    &             0,3,3,3,3,3,3,0,0,5,5,5,5,5,5,0,0,2,3,4,1,6,5,0,
00161 !    &             0,4,5,4,6,6,5,0,0,2,3,3,1,2,1,0,0,4,5,3,6,2,1,0 /
00162 !
00163 !     STORAGE CONVENTION FOR THE PRISM EXTRADIAGONAL TERMS
00164 !
00165       DATA STO / 00 , 01 , 02 , 03 , 04 , 05 ,
00166      &           01 , 00 , 06 , 07 , 08 , 09 ,
00167      &           02 , 06 , 00 , 10 , 11 , 12 ,
00168      &           03 , 07 , 10 , 00 , 13 , 14 ,
00169      &           04 , 08 , 11 , 13 , 00 , 15 ,
00170      &           05 , 09 , 12 , 14 , 15 , 00 /
00171 !
00172 !***********************************************************************
00173 !
00174 !     DIFFERENT WAYS TO SPLIT PRISMS (TO ENSURE MATCHING OF TETRAHEDRONS)
00175 !
00176 !     TETRA(2,2,2,3,4)
00177 !
00178 !     FIRST 3 DIMENSIONS : TYPE OF FACE
00179 !                      1 : CUT RECTANGLE BETWEEN  LOW-LEFT AND HIGH-RIGHT
00180 !                      2 : CUT RECTANGLE BETWEEN  HIGH-LEFT AND LOW-RIGHT
00181 !
00182 !     4TH DIMENSION : NUMBER OF TETRAHEDRON
00183 !     5TH DIMENSION : 4 POINTS OF THE TETRAHEDRON (IN LOCAL PRISM NUMBERING)
00184 !
00185 !     1 1 2 CUTTING
00186 !
00187 !     TETRA(1,1,2,1,1)= 1
00188 !     TETRA(1,1,2,1,2)= 2
00189 !     TETRA(1,1,2,1,3)= 3
00190 !     TETRA(1,1,2,1,4)= 6
00191 !
00192 !     TETRA(1,1,2,2,1)= 4
00193 !     TETRA(1,1,2,2,2)= 6
00194 !     TETRA(1,1,2,2,3)= 5
00195 !     TETRA(1,1,2,2,4)= 1
00196 !
00197 !     TETRA(1,1,2,3,1)= 5
00198 !     TETRA(1,1,2,3,2)= 2
00199 !     TETRA(1,1,2,3,3)= 1
00200 !     TETRA(1,1,2,3,4)= 6
00201 !
00202 !     2 1 1 CUTTING
00203 !
00204 !     TETRA(2,1,1,1,1)= 1
00205 !     TETRA(2,1,1,1,2)= 2
00206 !     TETRA(2,1,1,1,3)= 3
00207 !     TETRA(2,1,1,1,4)= 4
00208 !
00209 !     TETRA(2,1,1,2,1)= 4
00210 !     TETRA(2,1,1,2,2)= 6
00211 !     TETRA(2,1,1,2,3)= 5
00212 !     TETRA(2,1,1,2,4)= 2
00213 !
00214 !     TETRA(2,1,1,3,1)= 6
00215 !     TETRA(2,1,1,3,2)= 3
00216 !     TETRA(2,1,1,3,3)= 2
00217 !     TETRA(2,1,1,3,4)= 4
00218 !
00219 !     1 2 1 CUTTING
00220 !
00221 !     TETRA(1,2,1,1,1)= 1
00222 !     TETRA(1,2,1,1,2)= 2
00223 !     TETRA(1,2,1,1,3)= 3
00224 !     TETRA(1,2,1,1,4)= 5
00225 !
00226 !     TETRA(1,2,1,2,1)= 4
00227 !     TETRA(1,2,1,2,2)= 6
00228 !     TETRA(1,2,1,2,3)= 5
00229 !     TETRA(1,2,1,2,4)= 3
00230 !
00231 !     TETRA(1,2,1,3,1)= 4
00232 !     TETRA(1,2,1,3,2)= 1
00233 !     TETRA(1,2,1,3,3)= 3
00234 !     TETRA(1,2,1,3,4)= 5
00235 !
00236 !     2 2 1 CUTTING
00237 !
00238 !     TETRA(2,2,1,1,1)= 1
00239 !     TETRA(2,2,1,1,2)= 2
00240 !     TETRA(2,2,1,1,3)= 3
00241 !     TETRA(2,2,1,1,4)= 4
00242 !
00243 !     TETRA(2,2,1,2,1)= 4
00244 !     TETRA(2,2,1,2,2)= 6
00245 !     TETRA(2,2,1,2,3)= 5
00246 !     TETRA(2,2,1,2,4)= 3
00247 !
00248 !     TETRA(2,2,1,3,1)= 5
00249 !     TETRA(2,2,1,3,2)= 2
00250 !     TETRA(2,2,1,3,3)= 4
00251 !     TETRA(2,2,1,3,4)= 3
00252 !
00253 !     1 2 2 CUTTING
00254 !
00255 !     TETRA(1,2,2,1,1)= 1
00256 !     TETRA(1,2,2,1,2)= 2
00257 !     TETRA(1,2,2,1,3)= 3
00258 !     TETRA(1,2,2,1,4)= 5
00259 !
00260 !     TETRA(1,2,2,2,1)= 4
00261 !     TETRA(1,2,2,2,2)= 6
00262 !     TETRA(1,2,2,2,3)= 5
00263 !     TETRA(1,2,2,2,4)= 1
00264 !
00265 !     TETRA(1,2,2,3,1)= 6
00266 !     TETRA(1,2,2,3,2)= 3
00267 !     TETRA(1,2,2,3,3)= 5
00268 !     TETRA(1,2,2,3,4)= 1
00269 !
00270 !     2 1 2 CUTTING
00271 !
00272 !     TETRA(2,1,2,1,1)= 1
00273 !     TETRA(2,1,2,1,2)= 2
00274 !     TETRA(2,1,2,1,3)= 3
00275 !     TETRA(2,1,2,1,4)= 6
00276 !
00277 !     TETRA(2,1,2,2,1)= 4
00278 !     TETRA(2,1,2,2,2)= 6
00279 !     TETRA(2,1,2,2,3)= 5
00280 !     TETRA(2,1,2,2,4)= 2
00281 !
00282 !     TETRA(2,1,2,3,1)= 4
00283 !     TETRA(2,1,2,3,2)= 1
00284 !     TETRA(2,1,2,3,3)= 6
00285 !     TETRA(2,1,2,3,4)= 2
00286 !
00287 !-----------------------------------------------------------------------
00288 !
00289       SUR24=1.D0/24.D0
00290 !
00291       IF(SF%ELM.EQ.41.AND.SG%ELM.EQ.41.AND.SH%ELM.EQ.41) THEN
00292 !
00293 !-----------------------------------------------------------------------
00294 !
00295 !   LINEAR DISCRETISATION OF DIFFUSION COEFFICIENTS
00296 !
00297 !   LOOP ON THE PRISMS
00298 !
00299       DO IELEM=1,NELEM
00300 !
00301       IGLOB(1)=IKLE(IELEM,1)
00302       IGLOB(2)=IKLE(IELEM,2)
00303       IGLOB(3)=IKLE(IELEM,3)
00304       IGLOB(4)=IKLE(IELEM,4)
00305       IGLOB(5)=IKLE(IELEM,5)
00306       IGLOB(6)=IKLE(IELEM,6)
00307 !
00308       IF(IGLOB(1).GT.IGLOB(2)) THEN
00309         SENS(1)=1
00310       ELSE
00311         SENS(1)=2
00312       ENDIF
00313       IF(IGLOB(2).GT.IGLOB(3)) THEN
00314         SENS(2)=1
00315       ELSE
00316         SENS(2)=2
00317       ENDIF
00318       IF(IGLOB(3).GT.IGLOB(1)) THEN
00319         SENS(3)=1
00320       ELSE
00321         SENS(3)=2
00322       ENDIF
00323 !
00324 !     FOOTPRINT OF THE PRISM
00325 !
00326       SURF=0.5D0*
00327      &    ((X(IGLOB(2))-X(IGLOB(1)))*(Y(IGLOB(3))-Y(IGLOB(1)))
00328      &    -(X(IGLOB(3))-X(IGLOB(1)))*(Y(IGLOB(2))-Y(IGLOB(1))))
00329 !
00330 ! INITIALISES TO 0.D0
00331 !
00332       T(IELEM,1)=0.D0
00333       T(IELEM,2)=0.D0
00334       T(IELEM,3)=0.D0
00335       T(IELEM,4)=0.D0
00336       T(IELEM,5)=0.D0
00337       T(IELEM,6)=0.D0
00338 !
00339       XM(IELEM, 1)= 0.D0
00340       XM(IELEM, 2)= 0.D0
00341       XM(IELEM, 3)= 0.D0
00342       XM(IELEM, 4)= 0.D0
00343       XM(IELEM, 5)= 0.D0
00344 !
00345       XM(IELEM, 6)= 0.D0
00346       XM(IELEM, 7)= 0.D0
00347       XM(IELEM, 8)= 0.D0
00348       XM(IELEM, 9)= 0.D0
00349 !
00350       XM(IELEM,10)= 0.D0
00351       XM(IELEM,11)= 0.D0
00352       XM(IELEM,12)= 0.D0
00353 !
00354       XM(IELEM,13)= 0.D0
00355       XM(IELEM,14)= 0.D0
00356       XM(IELEM,15)= 0.D0
00357 !
00358 !-----------------------------------------------------------------------
00359 !     LOOP OVER  TI
00360 !
00361       DO I=1,3
00362 !
00363 !     TETRAHEDRON POINTS NUMBERS IN THE PRISM NUMBERING
00364 !
00365       NUM1=TETRA(SENS(1),SENS(2),SENS(3),I,1)
00366       NUM2=TETRA(SENS(1),SENS(2),SENS(3),I,2)
00367       NUM3=TETRA(SENS(1),SENS(2),SENS(3),I,3)
00368       NUM4=TETRA(SENS(1),SENS(2),SENS(3),I,4)
00369 !
00370 !     GLOBAL NUMBERS OF THE TETRAHEDRON POINTS
00371 !
00372       IT1=IGLOB(NUM1)
00373       IT2=IGLOB(NUM2)
00374       IT3=IGLOB(NUM3)
00375       IT4=IGLOB(NUM4)
00376 !
00377 ! VISCOSITY ALONG X Y AND Z
00378 !
00379       HTOT=F(IT1)+F(IT2)+F(IT3)+F(IT4)
00380       VTOT=G(IT1)+G(IT2)+G(IT3)+G(IT4)
00381       WTOT=H(IT1)+H(IT2)+H(IT3)+H(IT4)
00382 !
00383       X2=X(IT2)-X(IT1)
00384       Y2=Y(IT2)-Y(IT1)
00385       Z2=Z(IT2)-Z(IT1)
00386       X3=X(IT3)-X(IT1)
00387       Y3=Y(IT3)-Y(IT1)
00388       Z3=Z(IT3)-Z(IT1)
00389       X4=X(IT4)-X(IT1)
00390       Y4=Y(IT4)-Y(IT1)
00391       Z4=Z(IT4)-Z(IT1)
00392 !
00393 !-----------------------------------------------------------------------
00394 !    COEF:  THANKS MAPLE...
00395 !-----------------------------------------------------------------------
00396 !
00397       T1  = X2*Y3
00398       T3  = X2*Y4
00399       T5  = X3*Y2
00400       T7  = X4*Y2
00401       T9  = X3*Z2
00402       T11 = X4*Z2
00403 !     T13 = 4 TIMES THE TETRAHEDRON VOLUME ?
00404       T13 = T1*Z4-T3*Z3-T5*Z4+T7*Z3+T9*Y4-T11*Y3
00405 !
00406       T15 = -Y2*Z3+Y3*Z2
00407       T17 =  X2*Z4-T11
00408       T19 = -Y3*Z4+Y4*Z3
00409       T21 =  X2*Z3-T9
00410       T23 = -Y2*Z4+Y4*Z2
00411       T35 =  X3*Z4-X4*Z3
00412       T49 =  X3*Y4-X4*Y3
00413 !
00414 !     IF WIDTH MORE THAN 0.01 M
00415 !
00416       COEF=XMUL*SUR24/MAX(T13,0.01D0*SURF)
00417 !
00418       T28 = -T19+T23-T15
00419       T42 = -T35+T17-T21
00420       T51 = T3-T7
00421       T54 = T49-T3+T7+T1-T5
00422 !
00423       DIAG1  =COEF*( HTOT*T28**2 +VTOT*T42**2 +WTOT*T54**2)
00424       DIAG2  =COEF*( HTOT*T19**2 +VTOT*T35**2 +WTOT*T49**2)
00425       DIAG3  =COEF*( HTOT*T23**2 +VTOT*T17**2 +WTOT*T51**2)
00426       EXTR12 =COEF*( HTOT*T28*T19+VTOT*T42*T35-WTOT*T54*T49)
00427       EXTR13 =COEF*(-HTOT*T28*T23-VTOT*T42*T17+WTOT*T54*T51)
00428       EXTR23 =COEF*(-HTOT*T19*T23-VTOT*T35*T17-WTOT*T49*T51)
00429 !
00430 !     DEDUCED FROM PROPERTIES OF THE DIFFUSION MATRIX
00431 !
00432       EXTR14 = -(EXTR13+EXTR12+DIAG1)
00433       EXTR24 = -(EXTR23+DIAG2+EXTR12)
00434       EXTR34 = -(DIAG3+EXTR23+EXTR13)
00435       DIAG4  = -(EXTR14+EXTR24+EXTR34)
00436 !
00437 ! ASSEMBLES ON PRISM
00438 !
00439       T(IELEM,NUM1) = T(IELEM,NUM1)+ DIAG1
00440       T(IELEM,NUM2) = T(IELEM,NUM2)+ DIAG2
00441       T(IELEM,NUM3) = T(IELEM,NUM3)+ DIAG3
00442       T(IELEM,NUM4) = T(IELEM,NUM4)+ DIAG4
00443 !
00444       XM(IELEM,STO(NUM1,NUM2))=XM(IELEM,STO(NUM1,NUM2))+EXTR12
00445       XM(IELEM,STO(NUM1,NUM3))=XM(IELEM,STO(NUM1,NUM3))+EXTR13
00446       XM(IELEM,STO(NUM1,NUM4))=XM(IELEM,STO(NUM1,NUM4))+EXTR14
00447       XM(IELEM,STO(NUM2,NUM3))=XM(IELEM,STO(NUM2,NUM3))+EXTR23
00448       XM(IELEM,STO(NUM2,NUM4))=XM(IELEM,STO(NUM2,NUM4))+EXTR24
00449       XM(IELEM,STO(NUM3,NUM4))=XM(IELEM,STO(NUM3,NUM4))+EXTR34
00450 !
00451       ENDDO ! I
00452 !
00453 !---------------------------------------------------------------
00454 !
00455       ENDDO ! IELEM
00456 !
00457       ELSEIF(SF%ELM.EQ.40.AND.SG%ELM.EQ.40.AND.SH%ELM.EQ.40) THEN
00458 !
00459 !
00460 !-----------------------------------------------------------------------
00461 !
00462 !   P0 DISCRETISATION OF DIFFUSION COEFFICIENTS (CONSTANT ON A PRISM)
00463 !
00464 !   LOOP ON THE PRISMS
00465 !
00466       DO IELEM=1,NELEM
00467 !
00468       IGLOB(1)=IKLE(IELEM,1)
00469       IGLOB(2)=IKLE(IELEM,2)
00470       IGLOB(3)=IKLE(IELEM,3)
00471       IGLOB(4)=IKLE(IELEM,4)
00472       IGLOB(5)=IKLE(IELEM,5)
00473       IGLOB(6)=IKLE(IELEM,6)
00474 !
00475       IF(IGLOB(1).GT.IGLOB(2)) THEN
00476         SENS(1)=1
00477       ELSE
00478         SENS(1)=2
00479       ENDIF
00480       IF(IGLOB(2).GT.IGLOB(3)) THEN
00481         SENS(2)=1
00482       ELSE
00483         SENS(2)=2
00484       ENDIF
00485       IF(IGLOB(3).GT.IGLOB(1)) THEN
00486         SENS(3)=1
00487       ELSE
00488         SENS(3)=2
00489       ENDIF
00490 !
00491 !     FOOTPRINT OF THE PRISM
00492 !
00493       SURF=0.5D0*
00494      &    ((X(IGLOB(2))-X(IGLOB(1)))*(Y(IGLOB(3))-Y(IGLOB(1)))
00495      &    -(X(IGLOB(3))-X(IGLOB(1)))*(Y(IGLOB(2))-Y(IGLOB(1))))
00496 !
00497 ! INITIALISES TO 0.D0
00498 !
00499       T(IELEM,1)=0.D0
00500       T(IELEM,2)=0.D0
00501       T(IELEM,3)=0.D0
00502       T(IELEM,4)=0.D0
00503       T(IELEM,5)=0.D0
00504       T(IELEM,6)=0.D0
00505 !
00506       XM(IELEM, 1)= 0.D0
00507       XM(IELEM, 2)= 0.D0
00508       XM(IELEM, 3)= 0.D0
00509       XM(IELEM, 4)= 0.D0
00510       XM(IELEM, 5)= 0.D0
00511 !
00512       XM(IELEM, 6)= 0.D0
00513       XM(IELEM, 7)= 0.D0
00514       XM(IELEM, 8)= 0.D0
00515       XM(IELEM, 9)= 0.D0
00516 !
00517       XM(IELEM,10)= 0.D0
00518       XM(IELEM,11)= 0.D0
00519       XM(IELEM,12)= 0.D0
00520 !
00521       XM(IELEM,13)= 0.D0
00522       XM(IELEM,14)= 0.D0
00523       XM(IELEM,15)= 0.D0
00524 !
00525 !-----------------------------------------------------------------------
00526 !     LOOP OVER  TI
00527 !
00528       DO I=1,3
00529 !
00530 !     TETRAHEDRON POINTS NUMBERS IN THE PRISM NUMBERING
00531 !
00532       NUM1=TETRA(SENS(1),SENS(2),SENS(3),I,1)
00533       NUM2=TETRA(SENS(1),SENS(2),SENS(3),I,2)
00534       NUM3=TETRA(SENS(1),SENS(2),SENS(3),I,3)
00535       NUM4=TETRA(SENS(1),SENS(2),SENS(3),I,4)
00536 !
00537 !     GLOBAL NUMBERS OF THE TETRAHEDRON POINTS
00538 !
00539       IT1=IGLOB(NUM1)
00540       IT2=IGLOB(NUM2)
00541       IT3=IGLOB(NUM3)
00542       IT4=IGLOB(NUM4)
00543 !
00544 ! VISCOSITY ALONG X Y AND Z
00545 !
00546       HTOT=4*F(IELEM)
00547       VTOT=4*G(IELEM)
00548       WTOT=4*H(IELEM)
00549 !
00550       X2=X(IT2)-X(IT1)
00551       Y2=Y(IT2)-Y(IT1)
00552       Z2=Z(IT2)-Z(IT1)
00553       X3=X(IT3)-X(IT1)
00554       Y3=Y(IT3)-Y(IT1)
00555       Z3=Z(IT3)-Z(IT1)
00556       X4=X(IT4)-X(IT1)
00557       Y4=Y(IT4)-Y(IT1)
00558       Z4=Z(IT4)-Z(IT1)
00559 !
00560 !-----------------------------------------------------------------------
00561 !    COEF:  THANKS MAPLE...
00562 !-----------------------------------------------------------------------
00563 !
00564       T1  = X2*Y3
00565       T3  = X2*Y4
00566       T5  = X3*Y2
00567       T7  = X4*Y2
00568       T9  = X3*Z2
00569       T11 = X4*Z2
00570 !     T13 = 4 TIMES THE TETRAHEDRON VOLUME ?
00571       T13 = T1*Z4-T3*Z3-T5*Z4+T7*Z3+T9*Y4-T11*Y3
00572 !
00573       T15 = -Y2*Z3+Y3*Z2
00574       T17 =  X2*Z4-T11
00575       T19 = -Y3*Z4+Y4*Z3
00576       T21 =  X2*Z3-T9
00577       T23 = -Y2*Z4+Y4*Z2
00578       T35 =  X3*Z4-X4*Z3
00579       T49 =  X3*Y4-X4*Y3
00580 !
00581 !     IF WIDTH MORE THAN 0.01 M
00582 !
00583       COEF=XMUL*SUR24/MAX(T13,0.01D0*SURF)
00584 !
00585       T28 = -T19+T23-T15
00586       T42 = -T35+T17-T21
00587       T51 = T3-T7
00588       T54 = T49-T3+T7+T1-T5
00589 !
00590       DIAG1  =COEF*( HTOT*T28**2 +VTOT*T42**2 +WTOT*T54**2)
00591       DIAG2  =COEF*( HTOT*T19**2 +VTOT*T35**2 +WTOT*T49**2)
00592       DIAG3  =COEF*( HTOT*T23**2 +VTOT*T17**2 +WTOT*T51**2)
00593       EXTR12 =COEF*( HTOT*T28*T19+VTOT*T42*T35-WTOT*T54*T49)
00594       EXTR13 =COEF*(-HTOT*T28*T23-VTOT*T42*T17+WTOT*T54*T51)
00595       EXTR23 =COEF*(-HTOT*T19*T23-VTOT*T35*T17-WTOT*T49*T51)
00596 !
00597 !     DEDUCED FROM PROPERTIES OF THE DIFFUSION MATRIX
00598 !
00599       EXTR14 = -(EXTR13+EXTR12+DIAG1)
00600       EXTR24 = -(EXTR23+DIAG2+EXTR12)
00601       EXTR34 = -(DIAG3+EXTR23+EXTR13)
00602       DIAG4  = -(EXTR14+EXTR24+EXTR34)
00603 !
00604 ! ASSEMBLES ON PRISM
00605 !
00606       T(IELEM,NUM1) = T(IELEM,NUM1)+ DIAG1
00607       T(IELEM,NUM2) = T(IELEM,NUM2)+ DIAG2
00608       T(IELEM,NUM3) = T(IELEM,NUM3)+ DIAG3
00609       T(IELEM,NUM4) = T(IELEM,NUM4)+ DIAG4
00610 !
00611       XM(IELEM,STO(NUM1,NUM2))=XM(IELEM,STO(NUM1,NUM2))+EXTR12
00612       XM(IELEM,STO(NUM1,NUM3))=XM(IELEM,STO(NUM1,NUM3))+EXTR13
00613       XM(IELEM,STO(NUM1,NUM4))=XM(IELEM,STO(NUM1,NUM4))+EXTR14
00614       XM(IELEM,STO(NUM2,NUM3))=XM(IELEM,STO(NUM2,NUM3))+EXTR23
00615       XM(IELEM,STO(NUM2,NUM4))=XM(IELEM,STO(NUM2,NUM4))+EXTR24
00616       XM(IELEM,STO(NUM3,NUM4))=XM(IELEM,STO(NUM3,NUM4))+EXTR34
00617 !
00618       ENDDO ! I
00619 !
00620 !---------------------------------------------------------------
00621 !
00622       ENDDO ! IELEM
00623 !
00624 !-----------------------------------------------------------------------
00625 !
00626 !
00627       ELSE
00628 !
00629         IF (LNG.EQ.1) WRITE(LU,1000) SF%ELM,SG%ELM,SH%ELM
00630         IF (LNG.EQ.2) WRITE(LU,1001) SF%ELM,SG%ELM,SH%ELM
00631 1000    FORMAT(1X,'MT02PT (BIEF) : MAUVAIS TYPE DE F,G OU H : ',
00632      &  I6,1X,I6,1X,I6)
00633 1001    FORMAT(1X,'MT02PT (BIEF) : WRONG TYPE OF F,G OR H: ',
00634      &  I6,1X,I6,1X,I6)
00635         CALL PLANTE(1)
00636         STOP
00637 !
00638       ENDIF
00639 !
00640 !-----------------------------------------------------------------------
00641 !
00642 !  TREATMENT OF HYDROSTATIC INCONSISTENCIES
00643 !
00644       IF(INCHYD) THEN
00645 !
00646       DO IELEM=1,NELEM
00647 !
00648         I1=IKLE(IELEM,1)
00649         I2=IKLE(IELEM,2)
00650         I3=IKLE(IELEM,3)
00651         I4=IKLE(IELEM,4)
00652         I5=IKLE(IELEM,5)
00653         I6=IKLE(IELEM,6)
00654 !
00655         IF(MAX(Z(I1),Z(I2),Z(I3)).GT.MIN(Z(I4),Z(I5),Z(I6))) THEN
00656 !
00657           T(IELEM,1)  =0.D0
00658           T(IELEM,2)  =0.D0
00659           T(IELEM,3)  =0.D0
00660           T(IELEM,4)  =0.D0
00661           T(IELEM,5)  =0.D0
00662           T(IELEM,6)  =0.D0
00663           XM(IELEM, 1)=0.D0
00664           XM(IELEM, 2)=0.D0
00665           XM(IELEM, 3)=0.D0
00666           XM(IELEM, 4)=0.D0
00667           XM(IELEM, 5)=0.D0
00668           XM(IELEM, 6)=0.D0
00669           XM(IELEM, 7)=0.D0
00670           XM(IELEM, 8)=0.D0
00671           XM(IELEM, 9)=0.D0
00672           XM(IELEM,10)=0.D0
00673           XM(IELEM,11)=0.D0
00674           XM(IELEM,12)=0.D0
00675           XM(IELEM,13)=0.D0
00676           XM(IELEM,14)=0.D0
00677           XM(IELEM,15)=0.D0
00678 !
00679         ENDIF
00680 !
00681       ENDDO ! IELEM
00682 !
00683       ENDIF
00684 !
00685 !-----------------------------------------------------------------------
00686 !
00687       RETURN
00688       END

Generated on Fri Aug 31 2013 18:12:58 by S.E.Bourban (HRW) using doxygen 1.7.0