mt05cc.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\utils\bief\mt05cc.f
00002 !
00062                      SUBROUTINE MT05CC
00063 !                    *****************
00064 !
00065      &( A11 , A12 , A13 , A14 , A15 , A16 ,
00066      &  A21 , A22 , A23 , A24 , A25 , A26 ,
00067      &  A31 , A32 , A33 , A34 , A35 , A36 ,
00068      &  A41 , A42 , A43 , A44 , A45 , A46 ,
00069      &  A51 , A52 , A53 , A54 , A55 , A56 ,
00070      &  A61 , A62 , A63 , A64 , A65 , A66 ,
00071      &  XMUL,SU,SV,U,V,
00072      &  XEL,YEL,IKLE1,IKLE2,IKLE3,IKLE4,IKLE5,IKLE6,
00073      &  NELEM,NELMAX,FORMUL)
00074 !
00075 !***********************************************************************
00076 ! BIEF   V6P1                                   21/08/2010
00077 !***********************************************************************
00078 !
00079 !
00080 !
00081 !
00082 !
00083 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00084 !| A11            |<--| ELEMENTS OF MATRIX
00085 !| ...            |<--| ELEMENTS OF MATRIX
00086 !| A66            |<--| ELEMENTS OF MATRIX
00087 !| FORMUL         |-->| FORMULA DESCRIBING THE RESULTING MATRIX
00088 !| IKLE1          |-->| FIRST POINTS OF TRIANGLES
00089 !| IKLE2          |-->| SECOND POINTS OF TRIANGLES
00090 !| IKLE3          |-->| THIRD POINTS OF TRIANGLES
00091 !| IKLE4          |-->| FOURTH POINTS OF TRIANGLES (QUADRATIC)
00092 !| IKLE5          |-->| FIFTH POINTS OF TRIANGLES (QUADRATIC)
00093 !| IKLE6          |-->| SIXTH POINTS OF TRIANGLES (QUADRATIC)
00094 !| NELEM          |-->| NUMBER OF ELEMENTS
00095 !| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
00096 !| SU             |-->| BIEF_OBJ STRUCTURE OF U
00097 !| SURFAC         |-->| AREA OF TRIANGLES
00098 !| SV             |-->| BIEF_OBJ STRUCTURE OF V
00099 !| U              |-->| FUNCTION U USED IN THE FORMULA
00100 !| V              |-->| FUNCTION V USED IN THE FORMULA
00101 !| XEL            |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
00102 !| YEL            |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
00103 !| XMUL           |-->| MULTIPLICATION FACTOR
00104 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00105 !
00106       USE BIEF!, EX_MT05CC => MT05CC
00107 !
00108       IMPLICIT NONE
00109       INTEGER LNG,LU
00110       COMMON/INFO/LNG,LU
00111 !
00112 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00113 !
00114       INTEGER, INTENT(IN) :: NELEM,NELMAX
00115       INTEGER, INTENT(IN) :: IKLE1(NELMAX),IKLE2(NELMAX)
00116       INTEGER, INTENT(IN) :: IKLE3(NELMAX),IKLE4(NELMAX)
00117       INTEGER, INTENT(IN) :: IKLE5(NELMAX),IKLE6(NELMAX)
00118 !
00119       DOUBLE PRECISION, INTENT(INOUT) :: A11(*),A12(*),A13(*)
00120       DOUBLE PRECISION, INTENT(INOUT) :: A14(*),A15(*),A16(*)
00121       DOUBLE PRECISION, INTENT(INOUT) :: A21(*),A22(*),A23(*)
00122       DOUBLE PRECISION, INTENT(INOUT) :: A24(*),A25(*),A26(*)
00123       DOUBLE PRECISION, INTENT(INOUT) :: A31(*),A32(*),A33(*)
00124       DOUBLE PRECISION, INTENT(INOUT) :: A34(*),A35(*),A36(*)
00125       DOUBLE PRECISION, INTENT(INOUT) :: A41(*),A42(*),A43(*)
00126       DOUBLE PRECISION, INTENT(INOUT) :: A44(*),A45(*),A46(*)
00127       DOUBLE PRECISION, INTENT(INOUT) :: A51(*),A52(*),A53(*)
00128       DOUBLE PRECISION, INTENT(INOUT) :: A54(*),A55(*),A56(*)
00129       DOUBLE PRECISION, INTENT(INOUT) :: A61(*),A62(*),A63(*)
00130       DOUBLE PRECISION, INTENT(INOUT) :: A64(*),A65(*),A66(*)
00131 !
00132       DOUBLE PRECISION, INTENT(IN) :: XMUL
00133       DOUBLE PRECISION, INTENT(IN) :: U(*),V(*)
00134 !
00135 !
00136 !     STRUCTURES OF      U, V
00137       TYPE(BIEF_OBJ), INTENT(IN) :: SU,SV
00138 !
00139       DOUBLE PRECISION, INTENT(IN) :: XEL(NELMAX,3),YEL(NELMAX,3)
00140 !
00141       CHARACTER(LEN=16) :: FORMUL
00142 !
00143 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00144 !
00145 !
00146 !     DECLARATIONS SPECIFIC TO THIS SUBROUTINE
00147 !
00148       INTEGER IELEM,IELMU,IELMV
00149 !
00150       DOUBLE PRECISION X2,X3
00151       DOUBLE PRECISION Y2,Y3
00152       DOUBLE PRECISION U1,U2,U3,U4,U5,U6
00153       DOUBLE PRECISION V1,V2,V3,V4,V5,V6
00154       DOUBLE PRECISION XSU360,XSUR90,XSUR45
00155       DOUBLE PRECISION XSUR2520,XSUR630
00156 !=======================================================================
00157 !
00158 !     EXTRACTS THE TYPE OF ELEMENT FOR VELOCITY
00159 !
00160       IELMU = SU%ELM
00161       IELMV = SV%ELM
00162 !
00163 !-----------------------------------------------------------------------
00164 !
00165       XSU360   = XMUL / 360.D0
00166       XSUR45   = XMUL /  45.D0
00167       XSUR90   = XMUL /  90.D0
00168       XSUR2520 = XMUL /  2520.D0
00169       XSUR630  = XMUL /  630.D0
00170 !
00171 !-----------------------------------------------------------------------
00172 !
00173 !     CASE WHERE U AND V ARE CONSTANT BY ELEMENT
00174 !
00175       IF(IELMU.EQ.10.AND.IELMV.EQ.10) THEN
00176 !
00177 !-----------------------------------------------------------------------
00178 !
00179 !  P1 DISCRETISATION OF THE VELOCITY:
00180 !
00181       DO IELEM = 1 , NELEM
00182 !
00183 !   INITIALISES THE GEOMETRICAL VARIABLES
00184 !
00185         X2  =  XEL(IELEM,2)
00186         X3  =  XEL(IELEM,3)
00187 !
00188         Y2  =  YEL(IELEM,2)
00189         Y3  =  YEL(IELEM,3)
00190 !
00191         U1 = U(IELEM)
00192         U2 = U(IELEM)
00193         U3 = U(IELEM)
00194         V1 = V(IELEM)
00195         V2 = V(IELEM)
00196         V3 = V(IELEM)
00197 !
00198 !  EXTRADIAGONAL TERMS
00199 !
00200         A12(IELEM)=(( 5.D0*V2+V3+V1*6.D0)*X3 +
00201      &             ( -U1*6.D0-U3-U2*5.D0)*Y3) * XSU360
00202         A13(IELEM)=((-V2-V3*5.D0-V1*6.D0)*X2 +
00203      &             (  U1*6.D0+U3*5.D0+U2)*Y2) * XSU360
00204         A14(IELEM)=(( V2*2.D0+V3)*X2+(-V2*2.D0-V3-V1*6.D0)*X3 +
00205      &             ( -U2*2.D0-U3)*Y2+(U1*6.D0+U3+U2*2.D0)*Y3) * XSUR90
00206         A15(IELEM)=((-V3-V2*2.D0)*X2+(V2+V3*2.D0)*X3 +
00207      &             (U2*2.D0+U3)*Y2+(-U2-U3*2.D0)*Y3) * XSUR90
00208         A16(IELEM)=((V2+V3*2.D0+V1*6.D0)*X2+(-V2-V3*2.D0)*X3 +
00209      &             (-U2-U3*2.D0-U1*6.D0)*Y2+(U2+U3*2.D0)*Y3) * XSUR90
00210 !
00211         A21(IELEM)=((V2*6.D0+V3+V1*5.D0)*X2 +
00212      &             (-V2*6.D0-V3-V1*5.D0)*X3 +
00213      &             (-U2*6.D0-U3-U1*5.D0)*Y2 +
00214      &             (U1*5.D0+U3+U2*6.D0)*Y3) * XSU360
00215         A23(IELEM)=((-V2*6.D0-V3*5.D0-V1)*X2 +
00216      &             (  U1+U3*5.D0+U2*6.D0)*Y2) * XSU360
00217         A24(IELEM)=(-V2*X2*6.D0+(V2*6.D0+V3+V1*2.D0)*X3+U2*Y2*6.D0 +
00218      &             (-U2*6.D0-U3-U1*2.D0)*Y3) * XSUR90
00219         A25(IELEM)=(V2*X2*6.D0+(V3*2.D0+V1)*X3-U2*Y2*6.D0 +
00220      &             (-U3*2.D0-U1)*Y3) * XSUR90
00221         A26(IELEM)=((V3-V1)*X2+(-V3*2.D0-V1)*X3+(-U3+U1)*Y2 +
00222      &             (U1+U3*2.D0)*Y3) * XSUR90
00223 !
00224         A31(IELEM)=((V2+V3*6.D0+V1*5.D0)*X2 +
00225      &             (-V2-V3*6.D0-V1*5.D0)*X3 +
00226      &             (-U2-U3*6.D0-U1*5.D0)*Y2 +
00227      &             (U1*5.D0+U3*6.D0+U2)*Y3) * XSU360
00228         A32(IELEM)=((V2*5.D0+V3*6.D0+V1)*X3 +
00229      &             (-U1-U3*6.D0-U2*5.D0)*Y3) * XSU360
00230         A34(IELEM)=((V2*2.D0+V1)*X2+(-V2+V1)*X3+(-U1-U2*2.D0)*Y2+
00231      &             (U2-U1)*Y3)*XSUR90
00232         A35(IELEM)=((-V2*2.D0-V1)*X2-V3*X3*6.D0 +
00233      &              ( U2*2.D0+U1)*Y2+U3*Y3*6.D0) * XSUR90
00234         A36(IELEM)=((-V3*6.D0-V2-V1*2.D0)*X2+V3*X3*6.D0 +
00235      &              ( U2+U3*6.D0+U1*2.D0)*Y2-U3*Y3*6.D0) * XSUR90
00236 !
00237         A41(IELEM)=((-V2*2.D0-V3-V1*6.D0)*X2 +
00238      &              ( V2*2.D0+V3+V1*6.D0)*X3 +
00239      &              ( U1*6.D0+U3+U2*2.D0)*Y2 +
00240      &              (-U1*6.D0-U3-U2*2.D0)*Y3) * XSUR90
00241         A42(IELEM)=((-V2*6.D0-V3-V1*2.D0)*X3 +
00242      &              ( U1*2.D0+U3+U2*6.D0)*Y3) * XSUR90
00243         A43(IELEM)=((-V2*2.D0+V3-V1*2.D0)*X2 +
00244      &              ( U1*2.D0-U3+U2*2.D0)*Y2) * XSUR90
00245         A45(IELEM)=(( 6.D0*V2+2.D0*V3+4.D0*V1)*X2 +
00246      &              (-2.D0*V2-2.D0*V3-2.D0*V1)*X3 +
00247      &              (-6.D0*U2-2.D0*U3-4.D0*U1)*Y2 +
00248      &              (2.D0*U1+2.D0*U3+2.D0*U2)*Y3) * XSUR45
00249         A46(IELEM)=((2.D0*V2+4.D0*V1)*X2 +
00250      &              (2.D0*V3+2.D0*V2+2.D0*V1)*X3 +
00251      &              (-4.D0*U1-2.D0*U2)*Y2 +
00252      &              (-2.D0*U2-2.D0*U3-2.D0*U1)*Y3) * XSUR45
00253 !
00254         A51(IELEM)=((V2*2.D0+V3*2.D0-V1)*X2 +
00255      &             (-V2*2.D0-V3*2.D0+V1)*X3 +
00256      &             (-U2*2.D0-U3*2.D0+U1)*Y2 +
00257      &             (-U1+U3*2.D0+U2*2.D0)*Y3) * XSUR90
00258         A52(IELEM)=((-V2*6.D0-V3*2.D0-V1)*X3 +
00259      &             (U1+U3*2.D0+U2*6.D0)*Y3) * XSUR90
00260         A53(IELEM)=((V2*2.D0+V3*6.D0+V1)*X2 +
00261      &             (-U1-U3*6.D0-U2*2.D0)*Y2) * XSUR90
00262         A54(IELEM)=((-6.D0*V2-4.D0*V3-2.D0*V1)*X2 +
00263      &             (4.D0*V2+2.D0*V3)*X3 +
00264      &             (6.D0*U2+4.D0*U3+2.D0*U1)*Y2 +
00265      &             (-4.D0*U2-2.D0*U3)*Y3) * XSUR45
00266         A56(IELEM)=((-2.D0*V2-4.D0*V3)*X2 +
00267      &             (4.D0*V2+6.D0*V3+2.D0*V1)*X3 +
00268      &             (2.D0*U2+4.D0*U3)*Y2 +
00269      &             (-2.D0*U1-6.D0*U3-4.D0*U2)*Y3) * XSUR45
00270 !
00271         A61(IELEM)=((-V2-V3*2.D0-V1*6.D0)*X2 +
00272      &             (V2+V3*2.D0+V1*6.D0)*X3 +
00273      &             (U2+U3*2.D0+U1*6.D0)*Y2 +
00274      &             (-U2-U3*2.D0-U1*6.D0)*Y3) * XSUR90
00275         A62(IELEM)=((-V2+V3*2.D0+V1*2.D0)*X3 +
00276      &             (-U1*2.D0-U3*2.D0+U2)*Y3) * XSUR90
00277         A63(IELEM)=((V2+V3*6.D0+V1*2.D0)*X2 +
00278      &             (-U1*2.D0-U3*6.D0-U2)*Y2) * XSUR90
00279         A64(IELEM)=((-2.D0*V2-2.D0*V3-2.D0*V1)*X2 +
00280      &             (-4.D0*V1-2.D0*V3)*X3 +
00281      &             (2.D0*U1+2.D0*U3+2.D0*U2)*Y2 +
00282      &             (4.D0*U1+2.D0*U3)*Y3) * XSUR45
00283         A65(IELEM)=((2.D0*V3+2.D0*V2+2.D0*V1)*X2 +
00284      &             (-2.D0*V2-6.D0*V3-4.D0*V1)*X3 +
00285      &             (-2.D0*U2-2.D0*U3-2.D0*U1)*Y2 +
00286      &             (4.D0*U1+6.D0*U3+2.D0*U2)*Y3) * XSUR45
00287 !
00288 !  DIAGONAL TERMS:
00289 !
00290         A11(IELEM) = - A12(IELEM) - A13(IELEM) - A14(IELEM)
00291      &               - A15(IELEM) - A16(IELEM)
00292         A22(IELEM) = - A21(IELEM) - A23(IELEM) - A24(IELEM)
00293      &               - A25(IELEM) - A26(IELEM)
00294         A33(IELEM) = - A31(IELEM) - A32(IELEM) - A34(IELEM)
00295      &               - A35(IELEM) - A36(IELEM)
00296         A44(IELEM) = - A41(IELEM) - A42(IELEM) - A43(IELEM)
00297      &               - A45(IELEM) - A46(IELEM)
00298         A55(IELEM) = - A51(IELEM) - A52(IELEM) - A53(IELEM)
00299      &               - A54(IELEM) - A56(IELEM)
00300         A66(IELEM) = - A61(IELEM) - A62(IELEM) - A63(IELEM)
00301      &               - A64(IELEM) - A65(IELEM)
00302 !
00303       ENDDO ! IELEM
00304 !
00305 !-----------------------------------------------------------------------
00306 !
00307 !     CASE WHERE U AND V ARE LINEAR
00308 !
00309       ELSEIF(IELMU.EQ.11.AND.IELMV.EQ.11) THEN
00310 !
00311 !-----------------------------------------------------------------------
00312 !
00313 !  P1 DISCRETISATION OF THE VELOCITY:
00314 !
00315       DO IELEM = 1 , NELEM
00316 !
00317 !   INITIALISES THE GEOMETRICAL VARIABLES
00318 !
00319         X2  =  XEL(IELEM,2)
00320         X3  =  XEL(IELEM,3)
00321 !
00322         Y2  =  YEL(IELEM,2)
00323         Y3  =  YEL(IELEM,3)
00324 !
00325         U1 = U(IKLE1(IELEM))
00326         U2 = U(IKLE2(IELEM))
00327         U3 = U(IKLE3(IELEM))
00328         V1 = V(IKLE1(IELEM))
00329         V2 = V(IKLE2(IELEM))
00330         V3 = V(IKLE3(IELEM))
00331 !
00332 !  EXTRADIAGONAL TERMS
00333 !
00334         A12(IELEM)=(( 5.D0*V2+V3+V1*6.D0)*X3 +
00335      &             ( -U1*6.D0-U3-U2*5.D0)*Y3) * XSU360
00336         A13(IELEM)=((-V2-V3*5.D0-V1*6.D0)*X2 +
00337      &             (  U1*6.D0+U3*5.D0+U2)*Y2) * XSU360
00338         A14(IELEM)=(( V2*2.D0+V3)*X2+(-V2*2.D0-V3-V1*6.D0)*X3 +
00339      &             ( -U2*2.D0-U3)*Y2+(U1*6.D0+U3+U2*2.D0)*Y3) * XSUR90
00340         A15(IELEM)=((-V3-V2*2.D0)*X2+(V2+V3*2.D0)*X3 +
00341      &             (U2*2.D0+U3)*Y2+(-U2-U3*2.D0)*Y3) * XSUR90
00342         A16(IELEM)=((V2+V3*2.D0+V1*6.D0)*X2+(-V2-V3*2.D0)*X3 +
00343      &             (-U2-U3*2.D0-U1*6.D0)*Y2+(U2+U3*2.D0)*Y3) * XSUR90
00344 !
00345         A21(IELEM)=((V2*6.D0+V3+V1*5.D0)*X2 +
00346      &             (-V2*6.D0-V3-V1*5.D0)*X3 +
00347      &             (-U2*6.D0-U3-U1*5.D0)*Y2 +
00348      &             (U1*5.D0+U3+U2*6.D0)*Y3) * XSU360
00349         A23(IELEM)=((-V2*6.D0-V3*5.D0-V1)*X2 +
00350      &             (  U1+U3*5.D0+U2*6.D0)*Y2) * XSU360
00351         A24(IELEM)=(-V2*X2*6.D0+(V2*6.D0+V3+V1*2.D0)*X3+U2*Y2*6.D0 +
00352      &             (-U2*6.D0-U3-U1*2.D0)*Y3) * XSUR90
00353         A25(IELEM)=(V2*X2*6.D0+(V3*2.D0+V1)*X3-U2*Y2*6.D0 +
00354      &             (-U3*2.D0-U1)*Y3) * XSUR90
00355         A26(IELEM)=((V3-V1)*X2+(-V3*2.D0-V1)*X3+(-U3+U1)*Y2 +
00356      &             (U1+U3*2.D0)*Y3) * XSUR90
00357 !
00358         A31(IELEM)=((V2+V3*6.D0+V1*5.D0)*X2 +
00359      &             (-V2-V3*6.D0-V1*5.D0)*X3 +
00360      &             (-U2-U3*6.D0-U1*5.D0)*Y2 +
00361      &             (U1*5.D0+U3*6.D0+U2)*Y3) * XSU360
00362         A32(IELEM)=((V2*5.D0+V3*6.D0+V1)*X3 +
00363      &             (-U1-U3*6.D0-U2*5.D0)*Y3) * XSU360
00364         A34(IELEM)=((V2*2.D0+V1)*X2+(-V2+V1)*X3+(-U1-U2*2.D0)*Y2+
00365      &             (U2-U1)*Y3)*XSUR90
00366         A35(IELEM)=((-V2*2.D0-V1)*X2-V3*X3*6.D0 +
00367      &              ( U2*2.D0+U1)*Y2+U3*Y3*6.D0) * XSUR90
00368         A36(IELEM)=((-V3*6.D0-V2-V1*2.D0)*X2+V3*X3*6.D0 +
00369      &              ( U2+U3*6.D0+U1*2.D0)*Y2-U3*Y3*6.D0) * XSUR90
00370 !
00371         A41(IELEM)=((-V2*2.D0-V3-V1*6.D0)*X2 +
00372      &              ( V2*2.D0+V3+V1*6.D0)*X3 +
00373      &              ( U1*6.D0+U3+U2*2.D0)*Y2 +
00374      &              (-U1*6.D0-U3-U2*2.D0)*Y3) * XSUR90
00375         A42(IELEM)=((-V2*6.D0-V3-V1*2.D0)*X3 +
00376      &              ( U1*2.D0+U3+U2*6.D0)*Y3) * XSUR90
00377         A43(IELEM)=((-V2*2.D0+V3-V1*2.D0)*X2 +
00378      &              ( U1*2.D0-U3+U2*2.D0)*Y2) * XSUR90
00379         A45(IELEM)=(( 6.D0*V2+2.D0*V3+4.D0*V1)*X2 +
00380      &              (-2.D0*V2-2.D0*V3-2.D0*V1)*X3 +
00381      &              (-6.D0*U2-2.D0*U3-4.D0*U1)*Y2 +
00382      &              (2.D0*U1+2.D0*U3+2.D0*U2)*Y3) * XSUR45
00383         A46(IELEM)=((2.D0*V2+4.D0*V1)*X2 +
00384      &              (2.D0*V3+2.D0*V2+2.D0*V1)*X3 +
00385      &              (-4.D0*U1-2.D0*U2)*Y2 +
00386      &              (-2.D0*U2-2.D0*U3-2.D0*U1)*Y3) * XSUR45
00387 !
00388         A51(IELEM)=((V2*2.D0+V3*2.D0-V1)*X2 +
00389      &             (-V2*2.D0-V3*2.D0+V1)*X3 +
00390      &             (-U2*2.D0-U3*2.D0+U1)*Y2 +
00391      &             (-U1+U3*2.D0+U2*2.D0)*Y3) * XSUR90
00392         A52(IELEM)=((-V2*6.D0-V3*2.D0-V1)*X3 +
00393      &             (U1+U3*2.D0+U2*6.D0)*Y3) * XSUR90
00394         A53(IELEM)=((V2*2.D0+V3*6.D0+V1)*X2 +
00395      &             (-U1-U3*6.D0-U2*2.D0)*Y2) * XSUR90
00396         A54(IELEM)=((-6.D0*V2-4.D0*V3-2.D0*V1)*X2 +
00397      &             (4.D0*V2+2.D0*V3)*X3 +
00398      &             (6.D0*U2+4.D0*U3+2.D0*U1)*Y2 +
00399      &             (-4.D0*U2-2.D0*U3)*Y3) * XSUR45
00400         A56(IELEM)=((-2.D0*V2-4.D0*V3)*X2 +
00401      &             (4.D0*V2+6.D0*V3+2.D0*V1)*X3 +
00402      &             (2.D0*U2+4.D0*U3)*Y2 +
00403      &             (-2.D0*U1-6.D0*U3-4.D0*U2)*Y3) * XSUR45
00404 !
00405         A61(IELEM)=((-V2-V3*2.D0-V1*6.D0)*X2 +
00406      &             (V2+V3*2.D0+V1*6.D0)*X3 +
00407      &             (U2+U3*2.D0+U1*6.D0)*Y2 +
00408      &             (-U2-U3*2.D0-U1*6.D0)*Y3) * XSUR90
00409         A62(IELEM)=((-V2+V3*2.D0+V1*2.D0)*X3 +
00410      &             (-U1*2.D0-U3*2.D0+U2)*Y3) * XSUR90
00411         A63(IELEM)=((V2+V3*6.D0+V1*2.D0)*X2 +
00412      &             (-U1*2.D0-U3*6.D0-U2)*Y2) * XSUR90
00413         A64(IELEM)=((-2.D0*V2-2.D0*V3-2.D0*V1)*X2 +
00414      &             (-4.D0*V1-2.D0*V3)*X3 +
00415      &             (2.D0*U1+2.D0*U3+2.D0*U2)*Y2 +
00416      &             (4.D0*U1+2.D0*U3)*Y3) * XSUR45
00417         A65(IELEM)=((2.D0*V3+2.D0*V2+2.D0*V1)*X2 +
00418      &             (-2.D0*V2-6.D0*V3-4.D0*V1)*X3 +
00419      &             (-2.D0*U2-2.D0*U3-2.D0*U1)*Y2 +
00420      &             (4.D0*U1+6.D0*U3+2.D0*U2)*Y3) * XSUR45
00421 !
00422 !  DIAGONAL TERMS:
00423 !
00424         A11(IELEM) = - A12(IELEM) - A13(IELEM) - A14(IELEM)
00425      &               - A15(IELEM) - A16(IELEM)
00426         A22(IELEM) = - A21(IELEM) - A23(IELEM) - A24(IELEM)
00427      &               - A25(IELEM) - A26(IELEM)
00428         A33(IELEM) = - A31(IELEM) - A32(IELEM) - A34(IELEM)
00429      &               - A35(IELEM) - A36(IELEM)
00430         A44(IELEM) = - A41(IELEM) - A42(IELEM) - A43(IELEM)
00431      &               - A45(IELEM) - A46(IELEM)
00432         A55(IELEM) = - A51(IELEM) - A52(IELEM) - A53(IELEM)
00433      &               - A54(IELEM) - A56(IELEM)
00434         A66(IELEM) = - A61(IELEM) - A62(IELEM) - A63(IELEM)
00435      &               - A64(IELEM) - A65(IELEM)
00436       ENDDO ! IELEM
00437 !
00438 !-----------------------------------------------------------------------
00439 !
00440       ELSEIF(IELMU.EQ.13.AND.IELMV.EQ.13) THEN
00441 !
00442 !-----------------------------------------------------------------------
00443 !
00444 !  P2 DISCRETISATION OF THE VELOCITY:
00445 !
00446       IF(FORMUL(16:16).EQ.'N') THEN
00447 !
00448 !  N SCHEME
00449 !
00450         IF (LNG.EQ.1) WRITE(LU,12)
00451         IF (LNG.EQ.2) WRITE(LU,13)
00452         CALL PLANTE(1)
00453       ELSE
00454 !
00455       DO IELEM = 1 , NELEM
00456 !
00457 !  TRADITIONAL METHOD
00458 !
00459 !
00460 !   INITIALISES THE GEOMETRICAL VARIABLES
00461 !
00462         X2  =  XEL(IELEM,2)
00463         X3  =  XEL(IELEM,3)
00464 !
00465         Y2  =  YEL(IELEM,2)
00466         Y3  =  YEL(IELEM,3)
00467 !
00468         U1 = U(IKLE1(IELEM))
00469         U2 = U(IKLE2(IELEM))
00470         U3 = U(IKLE3(IELEM))
00471         U4 = U(IKLE4(IELEM))
00472         U5 = U(IKLE5(IELEM))
00473         U6 = U(IKLE6(IELEM))
00474         V1 = V(IKLE1(IELEM))
00475         V2 = V(IKLE2(IELEM))
00476         V3 = V(IKLE3(IELEM))
00477         V4 = V(IKLE4(IELEM))
00478         V5 = V(IKLE5(IELEM))
00479         V6 = V(IKLE6(IELEM))
00480 !
00481 !  EXTRADIAGONAL TERMS
00482 !
00483       A12(IELEM)= ((V5*20.D0+4.D0*8.D0*V4-11.D0*V3+16.D0*V6 +
00484      &             V1*18.D0+V2*9.D0)*X3+(-16.D0*U6-U5*20.D0+
00485      &             11.D0*U3-32.D0*U4-U1*18.D0-U2*9.D0)*Y3) * XSUR2520
00486       A13(IELEM)= ((-16.D0*V4-V5*20.D0-4.D0*8.D0*V6-V3*9.D0 -
00487      &             V1*18.D0+11.D0*V2)*X2+(U1*18.D0+U5*20.D0 +
00488      &             U3*9.D0+32.D0*U6+16.D0*U4-11.D0*U2)*Y2) * XSUR2520
00489       A14(IELEM)= ((-V3+8.D0*V4+4.D0*V6-V1*6.D0+12.D0*V5       +
00490      &              2.D0*2.D0*V2)*X2+(-20.D0*V4+V3*5.D0-16.D0*V6-
00491      &              24.D0*V1-8.D0*V5)*X3+(-4.D0*U2-8.D0*U4-12.D0*U5+
00492      &              U3+U1*6.D0-4.D0*U6)*Y2+(16.D0*U6+20.D0*U4+
00493      &              24.D0*U1-U3*5.D0+4.D0*2.D0*U5)*Y3) * XSUR630
00494       A15(IELEM)= ((V3-12.D0*V5-4.D0*V6+V1*6.D0-8.D0*V4-4.D0*V2)*X2 +
00495      &             (4.D0*V4-V1*6.D0-V2+12.D0*V5+4.D0*V3+8.D0*V6)*X3 +
00496      &             (4.D0*U2+12.D0*U5-U1*6.D0+8.D0*U4+4.D0*U6-U3)*Y2 +
00497      &             (-8.D0*U6-4.D0*U4+U1*6.D0+U2-4.D0*U3-12.D0*U5)*Y3)
00498      &             * XSUR630
00499       A16(IELEM)= (( 16.D0*V4+20.D0*V6+24.D0*V1+8.D0*V5-V2*5.D0)*X2 +
00500      &             (-8.D0*V6-4.D0*V4+V1*6.D0-4.D0*V3-12.D0*V5+V2)*X3+
00501      &             (-24.D0*U1-8.D0*U5+U2*5.D0-20.D0*U6-16.D0*U4)*Y2 +
00502      &             (8.D0*U6+4.D0*U4-U1*6.D0+4.D0*U3+12.D0*U5-U2)*Y3)
00503      &             * XSUR630
00504 !
00505       A21(IELEM)= ((-11.D0*V3+2.D0*8.D0*V5+V6*20.D0+V1*9.D0+32.D0*V4+
00506      &             V2*18.D0)*X2+(-4.D0*8.D0*V4-V1*9.D0-V2*18.D0-
00507      &             2.D0*8.D0*V5+11.D0*V3-V6*20.D0)*X3+(-U2*18.D0-
00508      &             2.D0*8.D0*U5-U1*9.D0-32.D0*U4-U6*20.D0+11.D0*U3)*Y2+
00509      &             (U6*20.D0+4.D0*8.D0*U4+U1*9.D0+U2*18.D0-11.D0*U3+
00510      &             2.D0*8.D0*U5)*Y3) * XSUR2520
00511       A23(IELEM)= ((-2.D0*8.D0*V4-4.D0*8.D0*V5-V6*20.D0-V3*9.D0 +
00512      &             11.D0*V1-V2*18.D0)*X2+(-11.D0*U1+4.D0*8.D0*U5 +
00513      &             U3*9.D0+U6*20.D0+2.D0*8.D0*U4+U2*18.D0)*Y2)
00514      &             *XSUR2520
00515       A24(IELEM)=((4.D0*V3-12.D0*V4+4.D0*V1+4.D0*V6-12.D0*V5-
00516      &            V2*30.D0)*X2+(20.D0*V4-V3*5.D0+8.D0*V6+24.D0*V2+
00517      &            16.D0*V5)*X3+(U2*30.D0+12.D0*U4+12.D0*U5-4.D0*U3-
00518      &            4.D0*U1-4.D0*U6)*Y2+(-4.D0*2.D0*U6-20.D0*U4-24.D0*U2+
00519      &            U3*5.D0-16.D0*U5)*Y3)*XSUR630
00520       A25(IELEM)=((-4.D0*V3+12.D0*V5-4.D0*V6-4.D0*V1+12.D0*V4+
00521      &            V2*30.D0)*X2+(4.D0*V4-V1-V2*6.D0+8.D0*V5+4.D0*V3+
00522      &            12.D0*V6)*X3+(-U2*30.D0-12.D0*U5+4.D0*U1-12.D0*U4+
00523      &            4.D0*U6+4.D0*U3)*Y2+(-12.D0*U6-4.D0*U4+U1+U2*6.D0-
00524      &            4.D0*U3-8.D0*U5)*Y3)*XSUR630
00525       A26(IELEM)=((4.D0*V5-4.D0*V4-V1*5.D0+V3*5.D0)*X2+(-12.D0*V6+
00526      &            V2*6.D0-4.D0*V4+V1-4.D0*V3-8.D0*V5)*X3+(-4.D0*U5-
00527      &            U3*5.D0+U1*5.D0+4.D0*U4)*Y2+(12.D0*U6+4.D0*U4-U2*6.D0+
00528      &            8.D0*U5+4.D0*U3-U1)*Y3)*XSUR630
00529 !
00530       A31(IELEM)=((V3*18.D0+16.D0*V5+32.D0*V6+V1*9.D0+V4*20.D0-
00531      &            11.D0*V2)*X2+(-V4*20.D0-V1*9.D0+11.D0*V2-2.D0*8.D0*V5-
00532      &            V3*18.D0-4.D0*8.D0*V6)*X3+(11.D0*U2-2.D0*8.D0*U5-
00533      &            U1*9.D0-U4*20.D0-4.D0*8.D0*U6-U3*18.D0)*Y2+
00534      &            (4.D0*8.D0*U6+U4*20.D0+U1*9.D0-11.D0*U2+U3*18.D0+
00535      &            16.D0*U5)*Y3)*XSUR2520
00536       A32(IELEM)=((32.D0*V5+V4*20.D0+V3*18.D0+16.D0*V6-11.D0*V1+
00537      &            V2*9.D0)*X3+(-16.D0*U6-32.D0*U5-U3*18.D0-U4*20.D0+
00538      &            11.D0*U1-U2*9.D0)*Y3)*XSUR2520
00539       A34(IELEM)=((-V3*6.D0+8.D0*V5+4.D0*V6+12.D0*V4-V1+4.D0*V2)*X2+
00540      &           (-V2*5.D0-4.D0*V5+V1*5.D0+4.D0*V6)*X3+(U1-8.D0*U5+
00541      &           U3*6.D0-4.D0*U6-4.D0*U2-12.D0*U4)*Y2+(-4.D0*U6+U2*5.D0+
00542      &           4.D0*U5-U1*5.D0)*Y3)*XSUR630
00543       A35(IELEM)=((V3*6.D0-8.D0*V5-4.D0*V6+V1-12.D0*V4-4.D0*V2)*X2+
00544      &           (4.D0*V4+4.D0*V1+4.D0*V2-12.D0*V5-V3*30.D0-
00545      &           12.D0*V6)*X3+(4.D0*U2+8.D0*U5-U1+12.D0*U4+4.D0*U6-
00546      &           U3*6.D0)*Y2+(12.D0*U6-4.D0*U4-4.D0*U1-4.D0*U2+U3*30.D0+
00547      &           12.D0*U5)*Y3)*XSUR630
00548       A36(IELEM)=((-24.D0*V3-8.D0*V4-20.D0*V6-16.D0*V5+V2*5.D0)*X2+
00549      &           (-4.D0*V4+V3*30.D0+12.D0*V6-4.D0*V2+12.D0*V5-
00550      &           4.D0*V1)*X3+(-U2*5.D0+16.D0*U5+24.D0*U3+8.D0*U4+
00551      &           20.D0*U6)*Y2+(-12.D0*U6+4.D0*U4+4.D0*U2-U3*30.D0+
00552      &           4.D0*U1-12.D0*U5)*Y3)*XSUR630
00553 !
00554       A41(IELEM)=(-U3*Y2*5.D0-20.D0*V6*X2-V3*X3*5.D0+12.*V1*X3-
00555      &           12.D0*V1*X2-12.D0*U1*Y3+12.D0*U1*Y2-20.D0*U6*Y3+
00556      &           V3*X2*5.D0+U3*Y3*5.D0+20.D0*V6*X3+20.*U6*Y2-
00557      &           40.D0*U4*Y3+40.*V4*X3-4.D0*U5*Y3+4.D0*V5*X3+
00558      &           40.D0*U4*Y2+4.D0*U5*Y2-4.D0*V5*X2-40.*V4*X2-
00559      &           8.D0*V2*X3+8.D0*U2*Y3+8.D0*V2*X2-
00560      &           8.D0*U2*Y2)*XSUR630
00561       A42(IELEM)=((-20.D0*V5-40.D0*V4+5.D0*V3-4.D0*V6+8.D0*V1-
00562      &           12.D0*V2)*X3+(4.D0*U6+20.D0*U5-5.D0*U3+40.D0*U4-
00563      &           8.D0*U1+12.D0*U2)*Y3) * XSUR630
00564       A43(IELEM)=((-24.D0*V4+4.D0*V5+4.D0*V6+3.D0*V3-4.D0*V1-
00565      &           4.D0*V2)*X2+(4.D0*U1-4.D0*U5-3.D0*U3-4.D0*U6+
00566      &           24.D0*U4+4.D0*U2)*Y2) * XSUR630
00567       A45(IELEM)=(12.D0*U3*Y2+32.D0*V6*X2+4.D0*V3*X3+4.D0*V1*X3-
00568      &           8.D0*V1*X2-4.D0*U1*Y3+8.D0*U1*Y2+32.D0*U6*Y3-
00569      &           12.D0*V3*X2-4.D0*U3*Y3-32.D0*V6*X3-32.D0*U6*Y2+
00570      &           32.D0*U4*Y3-32.D0*V4*X3+32.D0*U5*Y3-32.D0*V5*X3-
00571      &           96.D0*U4*Y2-48.D0*U5*Y2+48.D0*V5*X2+96.D0*V4*X2+
00572      &           4.D0*V2*X3-4.D0*U2*Y3+12.D0*V2*X2-12.D0*U2*Y2)
00573      &           * XSUR630
00574       A46(IELEM)=((-8.D0*V3+16.D0*V6+16.D0*V1+64.D0*V4-4.D0*V2)*X2+
00575      &           (-4.D0*V1+32.D0*V4-4.D0*V3+32.D0*V6+32.D0*V5-
00576      &           4.D0*V2)*X3+(8.D0*U3-16.D0*U1+4.D0*U2-16.D0*U6-
00577      &           64.D0*U4)*Y2+(-32.D0*U6+4.D0*U1-32.D0*U5+4.D0*U3-
00578      &           32.D0*U4+4.D0*U2)*Y3)*XSUR630
00579 !
00580       A51(IELEM)=((4.D0*V3+24.D0*V5-4.D0*V6-V1*3.D0-4.D0*V4+
00581      &           4.D0*V2)*X2+(4.D0*V4+V1*3.D0-4.D0*V2-24.D0*V5-
00582      &           4.D0*V3+4.D0*V6)*X3+(-4.D0*U2-24.D0*U5+U1*3.D0+
00583      &           4.D0*U4+4.D0*U6-4.D0*U3)*Y2+(-4.D0*U6-4.D0*U4-U1*3.D0+
00584      &           4.D0*U2+4.D0*U3+24.D0*U5)*Y3)*XSUR630
00585       A52(IELEM)=((-40.D0*V5-20.D0*V4+8.D0*V3-4.D0*V6+V1*5.D0-
00586      &            12.D0*V2)*X3+(4.D0*U6+40.D0*U5-8.D0*U3+
00587      &            20.D0*U4-U1*5.D0+12.D0*U2)*Y3)*XSUR630
00588       A53(IELEM)=((4.D0*V4+40.D0*V5+20.D0*V6+2.D0*6.D0*V3-V1*5.D0-
00589      &           8.D0*V2)*X2+(U1*5.D0-40.D0*U5-12.D0*U3-20.D0*U6-
00590      &           4.D0*U4+8.D0*U2)*Y2)*XSUR630
00591       A54(IELEM)=((8.D0*V3-96.D0*V5-32.D0*V6+12.D0*V1-48.D0*V4-
00592      &           12.D0*V2)*X2+(16.D0*V4-8.D0*V1+16.D0*V2+64.D0*V5-
00593      &           4.D0*V3)*X3+(-12.D0*U1+96.D0*U5-8.D0*U3+12.D0*U2+
00594      &           32.D0*U6+48.D0*U4)*Y2+(8.D0*U1+4.D0*U3-16.D0*U4-
00595      &           64.D0*U5-16.D0*U2)*Y3)*XSUR630
00596       A56(IELEM)=((-16.D0*V3-16.D0*V6+8.D0*V1-64.D0*V5+4.D0*V2)*X2+
00597      &           (32.D0*V4+12.D0*V3-12.D0*V1+48.D0*V6+96.D0*V5-
00598      &           8.D0*V2)*X3+(-8.D0*U1-4.D0*U2+16.D0*U3+16.D0*U6+
00599      &           64.D0*U5)*Y2+(-48.D0*U6+12.D0*U1+8.D0*U2-12.D0*U3-
00600      &           32.D0*U4-96.D0*U5)*Y3)*XSUR630
00601 !
00602       A61(IELEM)=((8.D0*V3-4.D0*V5-40.D0*V6-12.D0*V1-20.D0*V4+
00603      &           V2*5.D0)*X2+(20.D0*V4+12.D0*V1-V2*5.D0+4.D0*V5-
00604      &           8.D0*V3+40.D0*V6)*X3+(-U2*5.D0+4.D0*U5+12.D0*U1+
00605      &           20.D0*U4+40.D0*U6-8.D0*U3)*Y2+(-40.D0*U6-20.D0*U4-
00606      &           12.D0*U1+U2*5.D0+8.D0*U3-4.D0*U5)*Y3)*XSUR630
00607       A62(IELEM)=((-4.D0*V5-4.D0*V4+4.D0*V3+24.D0*V6+4.D0*V1-
00608      &           V2*3.D0)*X3+(-24.D0*U6+4.D0*U5-4.D0*U3+4.D0*U4-
00609      &           4.D0*U1+U2*3.D0)*Y3)*XSUR630
00610       A63(IELEM)=((4.D0*V4+20.D0*V5+40.D0*V6+12.D0*V3-8.D0*V1-
00611      &           V2*5.D0)*X2+(8.D0*U1-20.D0*U5-12.D0*U3-40.D0*U6-
00612      &           4.D0*U4+U2*5.D0)*Y2)*XSUR630
00613       A64(IELEM)=((-32.D0*V4+4.D0*V1+4.D0*V2-32.D0*V5+4.D0*V3-
00614      &           32.D0*V6)*X2+(-16.D0*V4-16.D0*V1+8.D0*V2+4.D0*V3-
00615      &           64.D0*V6)*X3+(32.D0*U6+32.D0*U4-4.D0*U1-4.D0*U2-
00616      &           4.D0*U3+32.D0*U5)*Y2+(64.D0*U6-8.D0*U2+16.D0*U1-
00617      &           4.D0*U3+16.D0*U4)*Y3)*XSUR630
00618       A65(IELEM)=((-4.D0*V1+32.D0*V4-4.D0*V3+32.D0*V6+32.D0*V5-
00619      &           4.D0*V2)*X2+(-32.D0*V4+8.D0*V1+12.D0*V2-48.D0*V5-
00620      &           12.D0*V3-96.D0*V6)*X3+(-32.D0*U6+4.D0*U1-32.D0*U5+
00621      &           4.D0*U3-32.D0*U4+4.D0*U2)*Y2+(96.D0*U6+32.D0*U4-
00622      &           8.D0*U1-12.D0*U2+12.D0*U3+48.D0*U5)*Y3)*XSUR630
00623 !
00624 !  THE DIAGONAL TERMS ARE OBTAINED BY MEANS OF THE 'MAGIC SQUARE':
00625 !
00626         A11(IELEM) = - A12(IELEM) - A13(IELEM) - A14(IELEM)
00627      &               - A15(IELEM) - A16(IELEM)
00628         A22(IELEM) = - A21(IELEM) - A23(IELEM) - A24(IELEM)
00629      &               - A25(IELEM) - A26(IELEM)
00630         A33(IELEM) = - A31(IELEM) - A32(IELEM) - A34(IELEM)
00631      &               - A35(IELEM) - A36(IELEM)
00632         A44(IELEM) = - A41(IELEM) - A42(IELEM) - A43(IELEM)
00633      &               - A45(IELEM) - A46(IELEM)
00634         A55(IELEM) = - A51(IELEM) - A52(IELEM) - A53(IELEM)
00635      &               - A54(IELEM) - A56(IELEM)
00636         A66(IELEM) = - A61(IELEM) - A62(IELEM) - A63(IELEM)
00637      &               - A64(IELEM) - A65(IELEM)
00638 !
00639       ENDDO ! IELEM
00640 !
00641       ENDIF
00642 !
00643 !-----------------------------------------------------------------------
00644 !
00645       ELSE
00646 !
00647         IF (LNG.EQ.1) WRITE(LU,10) IELMU,IELMV
00648         IF (LNG.EQ.2) WRITE(LU,11) IELMU,IELMV
00649 10    FORMAT(1X,'MT05CC (BIEF) : TYPES DE VITESSES NON PREVU : ',2I6)
00650 11    FORMAT(1X,
00651      & 'MT05CC (BIEF) : TYPES OF VELOCITIES NOT AVAILABLE : ',2I6)
00652 12    FORMAT(1X,'MT05CC (BIEF) : SCHEMA N NON PREVU  ')
00653 13    FORMAT(1X,
00654      & 'MT05CC (BIEF) : N SCHEMES NOT AVAILABLE  ')
00655         CALL PLANTE(1)
00656         STOP
00657 !
00658       ENDIF
00659 !
00660 !-----------------------------------------------------------------------
00661 !
00662       RETURN
00663       END

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