flux_ef_vf_2.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\utils\bief\flux_ef_vf_2.f
00002 !
00051                      SUBROUTINE FLUX_EF_VF_2
00052 !                    ***********************
00053 !
00054      &(PHIEL,NELEM,IKLE,IOPT,NPOIN,FN,FI_I,FSTAR,H,SU,DDT)
00055 !
00056 !***********************************************************************
00057 ! BIEF   V7P0                                   21/08/2010
00058 !***********************************************************************
00059 !
00060 !
00061 !
00062 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00063 !| FLOW           |<--| FLUXES PER SEGMENTS
00064 !| FN             |-->| OPTIONAL ARGUMENT FOR PSI SCHEME
00065 !| IKLE           |-->| CONNECTIVITY TABLE
00066 !| IOPT           |-->| OPTION FOR THE CONSTANT PER ELEMENT
00067 !| NELEM          |-->| NUMBER OF ELEMENTS
00068 !| NSEG           |-->| NUMBER OF SEGMENTS
00069 !| PHIEL          |-->| PER ELEMENT, FLUXES LEAVING POINTS
00070 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00071 !
00072       USE BIEF, EX_FLUX_EF_VF_2 => FLUX_EF_VF_2
00073 !
00074       IMPLICIT NONE
00075       INTEGER LNG,LU
00076       COMMON/INFO/LNG,LU
00077 !
00078 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00079 !
00080       INTEGER, INTENT(IN)             :: IOPT,NELEM
00081       INTEGER, INTENT(IN)             :: IKLE(NELEM,3)
00082       INTEGER, INTENT(IN)             :: NPOIN
00083       DOUBLE PRECISION, INTENT(IN)    :: PHIEL(NELEM,3)
00084       TYPE(BIEF_OBJ), INTENT(IN)      :: FN
00085       DOUBLE PRECISION, INTENT(INOUT) :: FI_I(NPOIN)
00086       DOUBLE PRECISION, INTENT(IN)    :: FSTAR(NPOIN),H(NPOIN)
00087       DOUBLE PRECISION, INTENT(IN)    :: SU(NELEM),DDT
00088 !
00089 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00090 !
00091       INTEGER IELEM,I1,I2,I3,I
00092       DOUBLE PRECISION K1,K2,K3,F12,F23,F31
00093       DOUBLE PRECISION FN1,FN2,FN3,F21,F32,F13
00094       DOUBLE PRECISION BETA1,BETA2,BETA3,BETA1PSI,BETA2PSI,BETA3PSI
00095       DOUBLE PRECISION FINCORR1,FINCORR2,FINCORR3
00096       DOUBLE PRECISION FI,BETA1FI,BETA2FI,BETA3FI
00097       DOUBLE PRECISION PHITCOR,SUMAX,H1,H2,H3,COEF
00098 !
00099       INTRINSIC MIN,MAX,ABS
00100 !
00101 !-----------------------------------------------------------------------
00102 !
00103       IF(IOPT.EQ.3) THEN
00104 !
00105 !-----------------------------------------------------------------------
00106 !
00107 !       PSI-SCHEME, CORRECTOR
00108 !
00109         DO I=1,NPOIN
00110           FI_I(I)=0.D0
00111         ENDDO
00112 !
00113         DO IELEM = 1,NELEM
00114 !
00115           K1 = -PHIEL(IELEM,1)
00116           K2 = -PHIEL(IELEM,2)
00117           K3 = -PHIEL(IELEM,3)
00118 !
00119           I1=IKLE(IELEM,1)
00120           I2=IKLE(IELEM,2)
00121           I3=IKLE(IELEM,3)
00122 !
00123 !         STARTS WITH N-SCHEME (EQUIVALENT TO LEO POSTMA'S IMPLEMENTATION)
00124 !         FIJ HERE LIKE LAMBDA(I,J) IN BOOK
00125 !
00126           F12=MAX(MIN(K1,-K2),0.D0)
00127           F23=MAX(MIN(K2,-K3),0.D0)
00128           F31=MAX(MIN(K3,-K1),0.D0)
00129           F21=MAX(MIN(K2,-K1),0.D0)
00130           F32=MAX(MIN(K3,-K2),0.D0)
00131           F13=MAX(MIN(K1,-K3),0.D0)
00132 !
00133           FN1=FN%R(I1)
00134           FN2=FN%R(I2)
00135           FN3=FN%R(I3)
00136 !
00137 !         NOT IN MARIO RICCHIUTO'S THEORY !!!!
00138 !
00139 !         FLUXES MODIFIED ACCORDING TO CLASSICAL PSI SCHEME
00140 !
00141           BETA1FI=F12*(FN1-FN2)+F13*(FN1-FN3)
00142           BETA2FI=F21*(FN2-FN1)+F23*(FN2-FN3)
00143           BETA3FI=F31*(FN3-FN1)+F32*(FN3-FN2)
00144           FI=BETA1FI+BETA2FI+BETA3FI
00145           IF(FI.GT.0.D0) THEN
00146             IF(BETA1FI.GT.FI) THEN
00147               F12=F12*FI/BETA1FI
00148               F13=F13*FI/BETA1FI
00149             ELSEIF(BETA1FI.LT.0.D0) THEN
00150               F12=0.D0
00151               F13=0.D0
00152             ENDIF
00153             IF(BETA2FI.GT.FI) THEN
00154               F21=F21*FI/BETA2FI
00155               F23=F23*FI/BETA2FI
00156             ELSEIF(BETA2FI.LT.0.D0) THEN
00157               F21=0.D0
00158               F23=0.D0
00159             ENDIF
00160             IF(BETA3FI.GT.FI) THEN
00161               F31=F31*FI/BETA3FI
00162               F32=F32*FI/BETA3FI
00163             ELSEIF(BETA3FI.LT.0.D0) THEN
00164               F31=0.D0
00165               F32=0.D0
00166             ENDIF
00167           ELSEIF(FI.LT.0.D0) THEN
00168             IF(BETA1FI.LT.FI) THEN
00169               F12=F12*FI/BETA1FI
00170               F13=F13*FI/BETA1FI
00171             ELSEIF(BETA1FI.GT.0.D0) THEN
00172               F12=0.D0
00173               F13=0.D0
00174             ENDIF
00175             IF(BETA2FI.LT.FI) THEN
00176               F21=F21*FI/BETA2FI
00177               F23=F23*FI/BETA2FI
00178             ELSEIF(BETA2FI.GT.0.D0) THEN
00179               F21=0.D0
00180               F23=0.D0
00181             ENDIF
00182             IF(BETA3FI.LT.FI) THEN
00183               F31=F31*FI/BETA3FI
00184               F32=F32*FI/BETA3FI
00185             ELSEIF(BETA3FI.GT.0.D0) THEN
00186               F31=0.D0
00187               F32=0.D0
00188             ENDIF
00189           ELSE
00190             F12=0.D0
00191             F23=0.D0
00192             F31=0.D0
00193             F21=0.D0
00194             F32=0.D0
00195             F13=0.D0
00196           ENDIF
00197 !
00198 !         END OF "NOT IN MARIO RICCHIUTO'S THEORY" !!!!!!
00199 !
00200 !         COMPUTE THE NEW RESIDUAL AND NEW DISTRIBUTION
00201 !
00202           COEF=SU(IELEM)/(DDT*3.D0)
00203           H1=H(I1)*COEF
00204           H2=H(I2)*COEF
00205           H3=H(I3)*COEF
00206 !
00207 !         AS CLASSICAL N SCHEME, BUT WITH DERIVATIVE IN TIME ADDED
00208 !
00209           FINCORR1=H1*(FSTAR(I1)-FN1)+F12*(FN1-FN2)+F13*(FN1-FN3)
00210           FINCORR2=H2*(FSTAR(I2)-FN2)+F21*(FN2-FN1)+F23*(FN2-FN3)
00211           FINCORR3=H3*(FSTAR(I3)-FN3)+F31*(FN3-FN1)+F32*(FN3-FN2)
00212 !
00213 !         PHITCOR IS THE NEW TOTAL RESIDUAL,
00214 !
00215           PHITCOR=FINCORR1+FINCORR2+FINCORR3
00216 !
00217           IF(ABS(PHITCOR).GT.1.D-25) THEN
00218             BETA1=FINCORR1/PHITCOR
00219             BETA2=FINCORR2/PHITCOR
00220             BETA3=FINCORR3/PHITCOR
00221             SUMAX=MAX(0.D0,BETA1)+MAX(0.D0,BETA2)+MAX(0.D0,BETA3)
00222             IF(SUMAX.NE.0.D0) THEN
00223               BETA1PSI=MAX(0.D0,BETA1)/SUMAX
00224               BETA2PSI=MAX(0.D0,BETA2)/SUMAX
00225               BETA3PSI=MAX(0.D0,BETA3)/SUMAX
00226               FI_I(I1)=FI_I(I1)+BETA1PSI*PHITCOR
00227               FI_I(I2)=FI_I(I2)+BETA2PSI*PHITCOR
00228               FI_I(I3)=FI_I(I3)+BETA3PSI*PHITCOR
00229             ENDIF
00230           ENDIF
00231 !
00232         ENDDO
00233 !
00234 !-----------------------------------------------------------------------
00235 !
00236       ELSE
00237 !
00238         IF(LNG.EQ.1) THEN
00239           WRITE(LU,*) 'FLUX_EF_VF_2 :'
00240           WRITE(LU,*) 'OPTION INCONNUE : ',IOPT
00241         ENDIF
00242         IF(LNG.EQ.2) THEN
00243           WRITE(LU,*) 'FLUX_EF_VF_2:'
00244           WRITE(LU,*) 'UNKNOWN OPTION: ',IOPT
00245         ENDIF
00246         CALL PLANTE(1)
00247         STOP
00248 !
00249       ENDIF
00250 !
00251 !-----------------------------------------------------------------------
00252 !
00253       RETURN
00254       END

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