mt04aa.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\utils\bief\mt04aa.f
00002 !
00073                      SUBROUTINE MT04AA
00074 !                    *****************
00075 !
00076      &( A11 , A12 , A13 ,
00077      &        A22 , A23 ,
00078      &              A33 ,
00079      &  XMUL,SU,SV,U,V,XEL,YEL,SURFAC,IKLE,NELEM,NELMAX)
00080 !
00081 !***********************************************************************
00082 ! BIEF   V6P1                                   21/08/2010
00083 !***********************************************************************
00084 !
00085 !
00086 !
00087 !
00088 !
00089 !
00090 !
00091 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00092 !| A11            |<--| ELEMENTS OF MATRIX
00093 !| A12            |<--| ELEMENTS OF MATRIX
00094 !| A13            |<--| ELEMENTS OF MATRIX
00095 !| A22            |<--| ELEMENTS OF MATRIX
00096 !| A23            |<--| ELEMENTS OF MATRIX
00097 !| A33            |<--| ELEMENTS OF MATRIX
00098 !| IKLE1          |-->| FIRST POINTS OF TRIANGLES
00099 !| IKLE2          |-->| SECOND POINTS OF TRIANGLES
00100 !| IKLE3          |-->| THIRD POINTS OF TRIANGLES
00101 !| NELEM          |-->| NUMBER OF ELEMENTS
00102 !| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
00103 !| SU             |-->| BIEF_OBJ STRUCTURE OF U
00104 !| SURFAC         |-->| AREA OF TRIANGLES
00105 !| SV             |-->| BIEF_OBJ STRUCTURE OF V
00106 !| U              |-->| FUNCTION U USED IN THE FORMULA
00107 !| V              |-->| FUNCTION V USED IN THE FORMULA
00108 !| XEL            |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
00109 !| YEL            |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
00110 !| XMUL           |-->| MULTIPLICATION FACTOR
00111 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00112 !
00113       USE BIEF, EX_MT04AA => MT04AA
00114 !
00115       IMPLICIT NONE
00116       INTEGER LNG,LU
00117       COMMON/INFO/LNG,LU
00118 !
00119 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00120 !
00121       INTEGER, INTENT(IN)             :: NELEM,NELMAX
00122       INTEGER, INTENT(IN)             :: IKLE(NELMAX,*)
00123       DOUBLE PRECISION, INTENT(INOUT) :: A11(*),A12(*),A13(*)
00124       DOUBLE PRECISION, INTENT(INOUT) ::        A22(*),A23(*)
00125       DOUBLE PRECISION, INTENT(INOUT) ::               A33(*)
00126       DOUBLE PRECISION, INTENT(IN)    :: XMUL,U(*),V(*)
00127       TYPE(BIEF_OBJ)  , INTENT(IN)    :: SU,SV
00128       DOUBLE PRECISION, INTENT(IN)    :: XEL(NELMAX,3),YEL(NELMAX,3)
00129       DOUBLE PRECISION, INTENT(IN)    :: SURFAC(NELMAX)
00130 !
00131 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00132 !
00133 !     DECLARATIONS SPECIFIC TO THIS SUBROUTINE
00134 !
00135       INTEGER IELMU,IELMV,IELEM
00136 !
00137       DOUBLE PRECISION SUR48,X2,X3,Y2,Y3,U1,U2,U3,U4,U5,U6
00138       DOUBLE PRECISION V1,V2,V3,V4,V5,V6,AUX
00139       DOUBLE PRECISION AUX1,AUX2,AUX3
00140       DOUBLE PRECISION SUR144,SUR720,U123,V123,ANS1,ANS2,ANS3
00141 !
00142 !-----------------------------------------------------------------------
00143 !
00144       SUR48  = XMUL/48.D0
00145       SUR144 = XMUL/144.D0
00146       SUR720 = XMUL/720.D0
00147 !
00148 !-----------------------------------------------------------------------
00149 !
00150       IELMU = SU%ELM
00151       IELMV = SV%ELM
00152 !
00153 !  CASE WHERE U AND V ARE LINEAR
00154 !
00155       IF(IELMU.EQ.11.AND.IELMV.EQ.11) THEN
00156 !
00157 !   LOOP ON THE ELEMENTS
00158 !
00159       DO IELEM = 1 , NELEM
00160 !
00161       X2  =   XEL(IELEM,2)
00162       X3  =   XEL(IELEM,3)
00163       Y2  =   YEL(IELEM,2)
00164       Y3  =   YEL(IELEM,3)
00165 !
00166       U1   =  U(IKLE(IELEM,1))
00167       U2   =  U(IKLE(IELEM,2))
00168       U3   =  U(IKLE(IELEM,3))
00169 !
00170       V1   =  V(IKLE(IELEM,1))
00171       V2   =  V(IKLE(IELEM,2))
00172       V3   =  V(IKLE(IELEM,3))
00173 !
00174       U123 = U1 + U2 + U3
00175       V123 = V1 + V2 + V3
00176 !
00177       AUX = SUR48 / SURFAC(IELEM)
00178 !
00179       AUX1 = U1*U123+U2**2+U2*U3+U3**2
00180       AUX2 = U1*(V123+V1)+U2*(V123+V2)+U3*(V123+V3)
00181       AUX3 = V1*V123+V2**2+V2*V3+V3**2
00182 !
00183       A12(IELEM) = (  2*Y3*(Y2-Y3)            *AUX1
00184      &                +(X3*(Y3-Y2)+(X3-X2)*Y3)*AUX2
00185      &                +2*X3*(X2-X3)           *AUX3 ) * AUX
00186 !
00187       A13(IELEM) = (  2*Y2*(Y3-Y2)            *AUX1
00188      &                +(2*X2*Y2-X2*Y3-X3*Y2)  *AUX2
00189      &                +2*X2*(X3-X2)           *AUX3 ) * AUX
00190 !
00191       A23(IELEM) = ( -2*Y2*Y3                 *AUX1
00192      &                +(X2*Y3+X3*Y2)          *AUX2
00193      &                -2*X2*X3                *AUX3 ) * AUX
00194 !
00195 !  USES HERE THE 'MAGIC SQUARE' PROPERTIES
00196 !  (SUM OF EACH LINE = 0 FOR EXAMPLE)
00197 !  AND THE SYMMETRICAL PROPERTIES OF THE MATRIX
00198 !
00199       A11(IELEM) = - A12(IELEM) - A13(IELEM)
00200       A22(IELEM) = - A23(IELEM) - A12(IELEM)
00201       A33(IELEM) = - A13(IELEM) - A23(IELEM)
00202 !
00203       ENDDO ! IELEM
00204 !
00205 !-----------------------------------------------------------------------
00206 !  CASE WHERE U AND V ARE QUASI-BUBBLE
00207 !
00208       ELSEIF(IELMU.EQ.12.AND.IELMV.EQ.12) THEN
00209 !
00210 !   LOOP ON THE ELEMENTS
00211 !
00212       DO IELEM = 1 , NELEM
00213 !
00214       X2  =   XEL(IELEM,2)
00215       X3  =   XEL(IELEM,3)
00216       Y2  =   YEL(IELEM,2)
00217       Y3  =   YEL(IELEM,3)
00218 !
00219       U1   =  U(IKLE(IELEM,1))
00220       U2   =  U(IKLE(IELEM,2))
00221       U3   =  U(IKLE(IELEM,3))
00222       U4   =  U(IKLE(IELEM,4))
00223 !
00224       V1   =  V(IKLE(IELEM,1))
00225       V2   =  V(IKLE(IELEM,2))
00226       V3   =  V(IKLE(IELEM,3))
00227       V4   =  V(IKLE(IELEM,4))
00228 !
00229       AUX = SUR144 / SURFAC(IELEM)
00230 !
00231       A12(IELEM) = (X2*((2*((V4+V2)*V3+V3**2+V4**2+V4*V2+V2**2)*
00232      & X3-(2*U3+U4+U2)*V3*Y3-(U3+2*U4+U2)*V4*Y3-(U3+U4+2*U2)*
00233      & V2*Y3)+(2*((V4+V1)*V3+V3**2+V4**2+V4*V1+V1**2)*X3-(2*U3
00234      & +U4+U1)*V3*Y3-(U3+2*U4+U1)*V4*Y3-(U3+U4+2*U1)*V1*Y3)+(
00235      & 2*((V2+V1)*V4+V4**2+V2**2+V2*V1+V1**2)*X3-(2*U4+U2+U1)*
00236      & V4*Y3-(U4+2*U2+U1)*V2*Y3-(U4+U2+2*U1)*V1*Y3))-(2*X3**2
00237      & )*(((V4+V2)*V3+V3**2+V4**2+V4*V2+V2**2)+((V4+V1)*V3+V3**2
00238      & +V4**2+V4*V1+V1**2)+((V2+V1)*V4+V4**2+V2**2+V2*V1+V1**2))
00239      & +X3*(2*Y3-Y2)*(((2*U3+U4+U2)*V3+(U3+2*U4+U2)*V4+(U3+U4
00240      & +2*U2)*V2)+((2*U3+U4+U1)*V3+(U3+2*U4+U1)*V4+(U3+U4+2*
00241      & U1)*V1)+((2*U4+U2+U1)*V4+(U4+2*U2+U1)*V2+(U4+U2+2*U1)*
00242      & V1))+2*Y3*(Y3-Y2)*(-(U4+U2)*U3-(U4+U1)*U3-(U2+U1)*U4-2*
00243      & U3**2-3*U4**2-U4*U2-U4*U1-2*U2**2-U2*U1-2*U1**2))*AUX
00244 !
00245       A13(IELEM) = (-(2*X2**2)*(((V4+V2)*V3+V3**2+V4**2+V4*V2+V2
00246      & **2)+((V4+V1)*V3+V3**2+V4**2+V4*V1+V1**2)+((V2+V1)*V4+V4
00247      & **2+V2**2+V2*V1+V1**2))+X2*((2*((V4+V2)*V3+V3**2+V4**2+
00248      & V4*V2+V2**2)*X3-(2*U3+U4+U2)*(Y3-2*Y2)*V3-(U3+2*U4+U2)
00249      & *(Y3-2*Y2)*V4-(U3+U4+2*U2)*(Y3-2*Y2)*V2)+(2*((V4+V1)*
00250      & V3+V3**2+V4**2+V4*V1+V1**2)*X3-(2*U3+U4+U1)*(Y3-2*Y2)*
00251      & V3-(U3+2*U4+U1)*(Y3-2*Y2)*V4-(U3+U4+2*U1)*(Y3-2*Y2)*
00252      & V1)+(2*((V2+V1)*V4+V4**2+V2**2+V2*V1+V1**2)*X3-(2*U4+U2
00253      & +U1)*(Y3-2*Y2)*V4-(U4+2*U2+U1)*(Y3-2*Y2)*V2-(U4+U2+2*
00254      & U1)*(Y3-2*Y2)*V1))-(X3*Y2)*(((2*U3+U4+U2)*V3+(U3+2*U4+
00255      & U2)*V4+(U3+U4+2*U2)*V2)+((2*U3+U4+U1)*V3+(U3+2*U4+U1)*
00256      & V4+(U3+U4+2*U1)*V1)+((2*U4+U2+U1)*V4+(U4+2*U2+U1)*V2+(
00257      & U4+U2+2*U1)*V1))+2*Y2*(Y3-Y2)*((U4+U2)*U3+(U4+U1)*U3+(
00258      & U2+U1)*U4+2*U3**2+3*U4**2+U4*U2+U4*U1+2*U2**2+U2*U1+2
00259      & *U1**2))*AUX
00260 !
00261       A23(IELEM) = (2*X2*X3*(-2*V3**2-2*V3*V4-V3*V2-V3*V1-3*V4
00262      & **2-2*V4*V2-2*V4*V1-2*V2**2-V2*V1-2*V1**2)+X2*Y3*(4*
00263      & V3*U3+2*V3*U4+V3*U2+V3*U1+2*V4*U3+6*V4*U4+2*V4*U2+2*
00264      & V4*U1+V2*U3+2*V2*U4+4*V2*U2+V2*U1+V1*U3+2*V1*U4+V1*U2+
00265      & 4*V1*U1)+X3*Y2*(4*V3*U3+2*V3*U4+V3*U2+V3*U1+2*V4*U3+6
00266      & *V4*U4+2*V4*U2+2*V4*U1+V2*U3+2*V2*U4+4*V2*U2+V2*U1+V1
00267      & *U3+2*V1*U4+V1*U2+4*V1*U1)+2*Y2*Y3*(-2*U3**2-2*U3*U4
00268      & -U3*U2-U3*U1-3*U4**2-2*U4*U2-2*U4*U1-2*U2**2-U2*U1-2
00269      & *U1**2))*AUX
00270 !
00271 !  USES HERE THE 'MAGIC SQUARE' PROPERTIES
00272 !  (SUM OF EACH LINE = 0 FOR EXAMPLE)
00273 !  AND THE SYMMETRICAL PROPERTIES OF THE MATRIX
00274 !
00275       A11(IELEM) = - A12(IELEM) - A13(IELEM)
00276       A22(IELEM) = - A23(IELEM) - A12(IELEM)
00277       A33(IELEM) = - A13(IELEM) - A23(IELEM)
00278 !
00279       ENDDO ! IELEM
00280 !
00281 !
00282 !-----------------------------------------------------------------------
00283 !  CASE WHERE U AND V ARE P2
00284 !
00285       ELSEIF(IELMU.EQ.13.AND.IELMV.EQ.13) THEN
00286 !
00287 !   LOOP ON THE ELEMENTS
00288 !
00289       DO IELEM = 1 , NELEM
00290 !
00291       X2  =   XEL(IELEM,2)
00292       X3  =   XEL(IELEM,3)
00293       Y2  =   YEL(IELEM,2)
00294       Y3  =   YEL(IELEM,3)
00295 !
00296       U1   =  U(IKLE(IELEM,1))
00297       U2   =  U(IKLE(IELEM,2))
00298       U3   =  U(IKLE(IELEM,3))
00299       U4   =  U(IKLE(IELEM,4))
00300       U5   =  U(IKLE(IELEM,5))
00301       U6   =  U(IKLE(IELEM,6))
00302 !
00303       V1   =  V(IKLE(IELEM,1))
00304       V2   =  V(IKLE(IELEM,2))
00305       V3   =  V(IKLE(IELEM,3))
00306       V4   =  V(IKLE(IELEM,4))
00307       V5   =  V(IKLE(IELEM,5))
00308       V6   =  V(IKLE(IELEM,6))
00309 !
00310       AUX = SUR720 / SURFAC(IELEM)
00311 !
00312       ANS1 = 3.D0*U2**2*Y2**2-4.D0*U4*Y3**2*U3+V2*X2*U3*Y2-
00313      &       4.D0*V1*X3**2*V5+U2*Y3*V3*X3-U2*Y3**2*U3-4.D0*V4*X3**2*V3-
00314      &       V2*X2**2*V3-4.D0*V2*X2**2*V6-V2*X2**2*V1-V2*X3**2*V1-
00315      &       4.D0*U4*Y2**2*U3-V1*X2**2*V3-6.D0*U3**2*Y3*Y2-
00316      &       32.D0*U5**2*Y2*Y3+16.D0*U4*Y3**2*U6+16.D0*V6*X3**2*V5+
00317      &       16.D0*U6*Y3**2*U5+16.D0*V6**2*X3**2+16.D0*U4*Y2**2*U5-
00318      &       4.D0*U1*Y3**2*U5-32.D0*U6**2*Y2*Y3+3.D0*V1**2*X3**2+
00319      &       16.D0*V4*X2**2*V5-4.D0*V1*X2**2*V5+16.D0*V5**2*X3**2-
00320      &       6.D0*U1**2*Y3*Y2-U1*Y3**2*U3-32.D0*V6**2*X3*X2-
00321      &       32.D0*U4**2*Y2*Y3+16.D0*U4*Y2**2*U6+16.D0*V4**2*X2**2-
00322      &       4.D0*U2*Y3**2*U6+16.D0*U5**2*Y3**2+3.D0*V2**2*X3**2-
00323      &       6.D0*V1**2*X2*X3+3.D0*U3**2*Y3**2+8.D0*U2*Y2*U6*Y3+
00324      &       2.D0*V2*X2*V3*X3+2.D0*U2*Y2*U3*Y3-U2*Y2*V1*X3-
00325      &       6.D0*V3**2*X2*X3+16.D0*V4*X3**2*V6-4.D0*U2*Y2*V6*X3+
00326      &       3.D0*V3**2*X3**2+16.D0*U4**2*Y3**2-6.D0*U2**2*Y2*Y3-
00327      &       6.D0*V2**2*X2*X3+16.D0*V4*X2**2*V6-4.D0*V4*X2**2*V3+
00328      &       6.D0*U2*Y2*V2*X3+6.D0*V2*X2*U2*Y3+32.D0*U5*Y2*V5*X3+
00329      &       8.D0*U4*Y3*U3*Y2-6.D0*U1*Y3*V1*X3+4.D0*V4*X3*U3*Y3-
00330      &       16.D0*V4*X3*U6*Y3+U1*Y3*V3*X3-16.D0*U4*Y3*V5*X3
00331       ANS2 = -32.D0*V4*X3*U4*Y3-16.D0*V4*X3*U5*Y3-6.D0*V3*X3*U3*Y3-
00332      &       32.D0*U5*Y3*V5*X3+4.D0*U4*Y3*V3*X3-16.D0*U4*Y3*V6*X3-
00333      &       16.D0*U6*Y3*V5*X3+4.D0*U1*Y3*V5*X3+4.D0*V1*X3*U5*Y3-
00334      &       32.D0*V6*X3*U6*Y3-16.D0*V6*X3*U5*Y3+V1*X3*U3*Y3+
00335      &       4.D0*V2*X3*U6*Y3+4.D0*U2*Y3*V6*X3+U2*Y3*V1*X3+V2*X3*U1*Y3+
00336      &       V2*X3*U3*Y3+3.D0*U1**2*Y3**2+3.D0*U3**2*Y2**2+
00337      &       3.D0*V3**2*X2**2-6.D0*U2*Y3*V2*X3-4.D0*U4*Y3*V3*X2-
00338      &       32.D0*U4*Y3*U6*Y2+32.D0*V6*X3*U6*Y2-32.D0*V5**2*X2*X3+
00339      &       3.D0*U2**2*Y3**2+16.D0*V4**2*X3**2+16.D0*U6*Y3*V5*X2+
00340      &       8.D0*U1*Y3*U5*Y2-4.D0*U1*Y3*V5*X2-32.D0*U6*Y2*U5*Y3-
00341      &       4.D0*U1*Y2*V5*X3+16.D0*U6*Y2*V5*X3-4.D0*V1*X3*U5*Y2+
00342      &       32.D0*V6*X2*U6*Y3+16.D0*V6*X2*U5*Y3-32.D0*V6*X3*V5*X2-
00343      &       4.D0*V1*X2*U5*Y3+8.D0*V1*X2*V5*X3+16.D0*V6*X3*U5*Y2-
00344      &       V1*X3*U3*Y2+6.D0*U1*Y2*V1*X3-V2*X3*U1*Y2-U1*Y2*V3*X3-
00345      &       4.D0*V2*X3*U6*Y2-U2*Y3*V1*X2-4.D0*U2*Y3*V6*X2-V2*X3*U3*Y2-
00346      &       4.D0*U4*Y2*V3*X3+16.D0*U4*Y2*V6*X3-U2*Y3*V3*X2-
00347      &       V1*X2*U3*Y3+8.D0*V2*X2*V6*X3+16.D0*U6**2*Y3**2-
00348      &       V2*X2*U3*Y3-U2*Y2*V3*X3+2.D0*U2*Y2*U1*Y3-V2*X2*U1*Y3+
00349      &       2.D0*V2*X2*V1*X3-4.D0*V2*X2*U6*Y3-32.D0*V4*X2*V5*X3
00350       ANS3 = 16.D0*V4*X2*U5*Y3+32.D0*V4*X2*U4*Y3+2.D0*U1*Y3*U3*Y2-
00351      &       U1*Y3*V3*X2+16.D0*U4*Y3*V6*X2+16.D0*V4*X2*U6*Y3+
00352      &       16.D0*V4*X3*U6*Y2+8.D0*V4*X2*V3*X3-32.D0*V4*X2*V6*X3-
00353      &       4.D0*V4*X2*U3*Y3+2.D0*V1*X2*V3*X3-4.D0*V4*X3*U3*Y2+
00354      &       6.D0*U1*Y3*V1*X2+6.D0*V3*X2*U3*Y3+16.D0*U4*Y3*V5*X2+
00355      &       16.D0*V4*X3*U5*Y2+16.D0*U6*Y2**2*U5-32.D0*U4*Y2*U5*Y3+
00356      &       32.D0*U4*Y2*V4*X3+16.D0*U4*Y2*V5*X3+16.D0*U4*Y3**2*U5-
00357      &       6.D0*U2*Y2*V2*X2+U2*Y2*V1*X2+4.D0*U2*Y2*V6*X2+
00358      &       6.D0*V3*X3*U3*Y2+32.D0*U5*Y3*V5*X2-V1*X3**2*V3
00359      &       -32.D0*U5*Y2*V5*X2+16.D0*V6*X2**2*V5-U1*Y2**2*U3-
00360      &       U2*Y2**2*U1-4.D0*U1*Y2**2*U5+3.D0*V2**2*X2**2-
00361      &       V2*X3**2*V3+16.D0*U6**2*Y2**2+4.D0*U1*Y2*V5*X2-
00362      &       16.D0*U6*Y2*V5*X2-4.D0*U2*Y2**2*U6+3.D0*V1**2*X2**2+
00363      &       16.D0*V6**2*X2**2-32.D0*V4**2*X2*X3+16.D0*V4*X3**2*V5+
00364      &       V1*X2*U3*Y2-6.D0*V1*X2*U1*Y2+4.D0*V4*X2*U3*Y2-
00365      &       16.D0*V4*X2*U6*Y2+U2*Y2*V3*X2+4.D0*V2*X2*U6*Y2+
00366      &       V2*X2*U1*Y2-16.D0*V4*X2*U5*Y2+3.D0*U1**2*Y2**2+
00367      &       U1*Y2*V3*X2+4.D0*U4*Y2*V3*X2-16.D0*U4*Y2*V6*X2-
00368      &       32.D0*V6*X2*U6*Y2-16.D0*V6*X2*U5*Y2+4.D0*V1*X2*U5*Y2
00369       A11(IELEM) = (16.D0*U4**2*Y2**2-U2*Y3**2*U1-4*V2*X3**2*V6-
00370      &             U2*Y2**2*U3+16.D0*V5**2*X2**2-32.D0*U4*Y2*V4*X2-
00371      &             6.D0*V3*X2*U3*Y2-16.D0*U4*Y2*V5*X2+
00372      &             16.D0*U5**2*Y2**2 + ANS1 + ANS2 + ANS3)*2.D0*AUX
00373 !
00374       A22(IELEM) = (-4.D0*U4*Y3**2*U3-4.D0*V1*X3**2*V5+U2*Y3*V3*X3-
00375      &       U2*Y3**2*U3-4.D0*V4*X3**2*V3-V2*X3**2*V1+
00376      &       16.D0*U4*Y3**2*U6+16.D0*V6*X3**2*V5+16.D0*U6*Y3**2*U5+
00377      &       16.D0*V6**2*X3**2-4.D0*U1*Y3**2*U5+3*V1**2*X3**2+
00378      &       16.D0*V5**2*X3**2-U1*Y3**2*U3-4.D0*U2*Y3**2*U6+
00379      &       16.D0*U5**2*Y3**2+3.D0*V2**2*X3**2+3.D0*U3**2*Y3**2+
00380      &       16.D0*V4*X3**2*V6+3.D0*V3**2*X3**2+16.D0*U4**2*Y3**2-
00381      &       6.D0*U1*Y3*V1*X3+4.D0*V4*X3*U3*Y3-16.D0*V4*X3*U6*Y3+
00382      &       U1*Y3*V3*X3-16.D0*U4*Y3*V5*X3-32.D0*V4*X3*U4*Y3-
00383      &       16.D0*V4*X3*U5*Y3-6.D0*V3*X3*U3*Y3-32.D0*U5*Y3*V5*X3+
00384      &       4.D0*U4*Y3*V3*X3-16.D0*U4*Y3*V6*X3-16.D0*U6*Y3*V5*X3+
00385      &       4.D0*U1*Y3*V5*X3+4.D0*V1*X3*U5*Y3-32.D0*V6*X3*U6*Y3-
00386      &       16.D0*V6*X3*U5*Y3+V1*X3*U3*Y3+4.D0*V2*X3*U6*Y3+
00387      &       4.D0*U2*Y3*V6*X3+U2*Y3*V1*X3+V2*X3*U1*Y3+V2*X3*U3*Y3+
00388      &       3.D0*U1**2*Y3**2-6.D0*U2*Y3*V2*X3+3.D0*U2**2*Y3**2+
00389      &       16.D0*V4**2*X3**2+16.D0*U6**2*Y3**2+16.D0*U4*Y3**2*U5-
00390      &       V1*X3**2*V3-V2*X3**2*V3+16.D0*V4*X3**2*V5-U2*Y3**2*U1-
00391      &       4.D0*V2*X3**2*V6)*2.D0*AUX
00392 !
00393       ANS1 = 6.D0*U3**2*Y3*Y2+32.D0*U5**2*Y2*Y3+32.D0*U6**2*Y2*Y3+
00394      &       6.D0*U1**2*Y3*Y2+32.D0*V6**2*X3*X2+32.D0*U4**2*Y2*Y3+
00395      &       6.D0*V1**2*X2*X3-8.D0*U2*Y2*U6*Y3-2.D0*V2*X2*V3*X3-
00396      &       2.D0*U2*Y2*U3*Y3+U2*Y2*V1*X3+6.D0*V3**2*X2*X3+
00397      &       4.D0*U2*Y2*V6*X3+6.D0*U2**2*Y2*Y3+6.D0*V2**2*X2*X3-
00398      &       6.D0*U2*Y2*V2*X3-6.D0*V2*X2*U2*Y3-32.D0*U5*Y2*V5*X3-
00399      &       8.D0*U4*Y3*U3*Y2+4.D0*U4*Y3*V3*X2+32.D0*U4*Y3*U6*Y2-
00400      &       32.D0*V6*X3*U6*Y2+32.D0*V5**2*X2*X3-16.D0*U6*Y3*V5*X2-
00401      &       8.D0*U1*Y3*U5*Y2+4.D0*U1*Y3*V5*X2+32.D0*U6*Y2*U5*Y3+
00402      &       4.D0*U1*Y2*V5*X3-16.D0*U6*Y2*V5*X3+4.D0*V1*X3*U5*Y2-
00403      &       32.D0*V6*X2*U6*Y3-16.D0*V6*X2*U5*Y3+32.D0*V6*X3*V5*X2+
00404      &       4.D0*V1*X2*U5*Y3-8.D0*V1*X2*V5*X3-16.D0*V6*X3*U5*Y2+
00405      &       V1*X3*U3*Y2-6.D0*U1*Y2*V1*X3+V2*X3*U1*Y2+U1*Y2*V3*X3+
00406      &       4.D0*V2*X3*U6*Y2+U2*Y3*V1*X2+4.D0*U2*Y3*V6*X2+V2*X3*U3*Y2+
00407      &       4.D0*U4*Y2*V3*X3-16.D0*U4*Y2*V6*X3+U2*Y3*V3*X2+V1*X2*U3*Y3-
00408      &       8.D0*V2*X2*V6*X3+V2*X2*U3*Y3+U2*Y2*V3*X3-2.D0*U2*Y2*U1*Y3+
00409      &       V2*X2*U1*Y3-2.D0*V2*X2*V1*X3+4.D0*V2*X2*U6*Y3+
00410      &       32.D0*V4*X2*V5*X3-16.D0*V4*X2*U5*Y3-32.D0*V4*X2*U4*Y3-
00411      &       2.D0*U1*Y3*U3*Y2+U1*Y3*V3*X2-16.D0*U4*Y3*V6*X2
00412       A23(IELEM) = -(-16.D0*V4*X2*U6*Y3-16.D0*V4*X3*U6*Y2-
00413      &             8.D0*V4*X2*V3*X3+32.D0*V4*X2*V6*X3+
00414      &             4.D0*V4*X2*U3*Y3-2.D0*V1*X2*V3*X3+
00415      &             4.D0*V4*X3*U3*Y2-6.D0*U1*Y3*V1*X2-
00416      &             6.D0*V3*X2*U3*Y3-16.D0*U4*Y3*V5*X2-
00417      &             16.D0*V4*X3*U5*Y2+32.D0*U4*Y2*U5*Y3-
00418      &             32.D0*U4*Y2*V4*X3-16.D0*U4*Y2*V5*X3-
00419      &             6.D0*V3*X3*U3*Y2-32.D0*U5*Y3*V5*X2+
00420      &             32.D0*V4**2*X2*X3 + ANS1)*AUX
00421 !
00422 !  USES HERE THE 'MAGIC SQUARE' PROPERTIES
00423 !  (SUM OF EACH LINE = 0 FOR EXAMPLE)
00424 !  AND THE SYMMETRICAL PROPERTIES OF THE MATRIX
00425 !
00426       A12(IELEM) = - A22(IELEM) - A23(IELEM)
00427       A13(IELEM) = - A11(IELEM) - A12(IELEM)
00428       A33(IELEM) = - A13(IELEM) - A23(IELEM)
00429 !
00430       ENDDO ! IELEM
00431 !
00432 !
00433 !     OTHER TYPES OF DISCRETISATION FOR U
00434 !
00435 !-----------------------------------------------------------------------
00436 !
00437       ELSE
00438 !
00439         IF(IELMU.EQ.IELMV) THEN
00440         IF (LNG.EQ.1) WRITE(LU,100) IELMU
00441         IF (LNG.EQ.2) WRITE(LU,101) IELMU
00442 100     FORMAT(1X,'MT04AA (BIEF) :',/,
00443      &         1X,'DISCRETISATION DE U ET V : ',1I6,' NON PREVUE')
00444 101     FORMAT(1X,'MT04AA (BIEF) :',/,
00445      &         1X,'DISCRETIZATION OF U AND V : ',1I6,' NOT AVAILABLE')
00446         ELSE
00447         IF (LNG.EQ.1) WRITE(LU,200) IELMU,IELMV
00448         IF (LNG.EQ.2) WRITE(LU,201) IELMU,IELMV
00449 200     FORMAT(1X,'MT04AA (BIEF) :',/,
00450      &         1X,'U ET V DE DISCRETISATIONS DIFFERENTES :',1I6,3X,1I6)
00451 201     FORMAT(1X,'MT04AA (BIEF) :',/,
00452      &         1X,'U AND V OF A DIFFERENT DISCRETISATION:',1I6,3X,1I6)
00453         ENDIF
00454 !
00455         CALL PLANTE(1)
00456         STOP
00457 !
00458       ENDIF
00459 !
00460 !-----------------------------------------------------------------------
00461 !
00462       RETURN
00463       END

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