mt12ab.f

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

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