cvtrvf.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\utils\bief\cvtrvf.f
00002 !
00084                      SUBROUTINE CVTRVF
00085 !                    *****************
00086 !
00087      &(F,FN,FSCEXP,DIFT,CONV,H,HN,HPROP,UCONV,VCONV,DM1,ZCONV,SOLSYS,
00088      & VISC,VISC_S,SM,SMH,YASMH,SMI,YASMI,FBOR,MASKTR,MESH,
00089      & T1,T2,T3,T4,T5,T6,T7,T8,HNT,HT,AGGLOH,TE1,DT,ENTET,BILAN,
00090      & OPDTRA,MSK,MASKEL,S,MASSOU,OPTSOU,LIMTRA,KDIR,KDDL,NPTFR,FLBOR,
00091      & YAFLBOR,VOLU2D,V2DPAR,UNSV2D,IOPT,FLBORTRA,MASKPT,
00092      & RAIN,PLUIE,TRAIN,OPTADV)
00093 !
00094 !***********************************************************************
00095 ! BIEF   V7P0                                   21/08/2010
00096 !***********************************************************************
00097 !
00098 !
00099 !
00100 !
00101 !
00102 !
00103 !
00104 !
00105 !
00106 !
00107 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00108 !| AGGLOH         |-->| MASS-LUMPING IN CONTINUITY EQUATION
00109 !| BILAN          |-->| LOGICAL TRIGGERING A MASS BALANCE INFORMATION
00110 !| CONV           |-->| LOGICAL, IF YES THERE IS ADVECTION OF F
00111 !| DIFT           |-->| LOGICAL, IF YES THERE IS DIFFUSION OF F
00112 !| DM1            |-->| THE PIECE-WISE CONSTANT PART OF ADVECTION FIELD
00113 !|                |   | IS DM1*GRAD(ZCONV), SEE SOLSYS.
00114 !| DT             |-->| TIME STEP
00115 !| ENTET          |-->| LOGICAL, IF YES INFORMATION IS GIVEN ON MASS
00116 !|                |   | CONSERVATION.
00117 !| F              |<--| F AT TIME T(N+1)
00118 !| FBOR           |-->| DIRICHLET CONDITIONS ON F.
00119 !| FLBOR          |-->| FLUXES AT BOUNDARIES
00120 !| FLBORTRA       |<->| TRACER FLUXES AT BOUNDARIES
00121 !| FN             |-->| F AT TIME T(N)
00122 !| FSCEXP         |-->| EXPLICIT PART OF THE SOURCE TERM
00123 !|                |   | EQUAL TO ZERO EVERYWHERE BUT ON SOURCES
00124 !|                |   | WHERE THERE IS FSCE - (1-TETAT) FN
00125 !|                |   | SEE DIFSOU
00126 !| HNT,HT         |<--| WORK ARRAYS (MODIFIED DEPTHS TO TAKE MASS-LUMPING
00127 !|                |   | INTO ACCOUNT)
00128 !| HPROP          |-->| PROPAGATION DEPTH (DONE IN CVDFTR).
00129 !| IOPT           |-->| OPTIONS FOR COMPUTATION (NUMBER BETWEEN 0 AND 13)
00130 !|                |   | THE TENS (IOPT2, I.E. 0 OR 1):
00131 !|                |   | 0: UCONV OBEYS THE CONTINUITY EQUATION
00132 !|                |   | 1: UCONV DOES NOT OBEY THE CONTINUITY EQUATION
00133 !|                |   | THE UNITS (IOPT1, I.E. 0 TO 3): VARIANT FOR FLUXES
00134 !|                |   | 0: CONSTANT PER ELEMENT = 0
00135 !|                |   | 1: CHI-TUAN PHAM'S CONSTANT
00136 !|                |   | 2: N SCHEME
00137 !|                |   | 3: PSI SCHEME
00138 !| KDDL           |-->| CONVENTION FOR DEGREE OF FREEDOM
00139 !| KDIR           |-->| CONVENTION FOR DIRICHLET POINT
00140 !| LIMTRA         |-->| BOUNDARY CONDITIONS ON BOOUNDARY POINTS
00141 !| MASKEL         |-->| MASKING OF ELEMENTS
00142 !|                |   | =1. : NORMAL   =0. : MASKED ELEMENT
00143 !| MASKPT         |-->| MASKING PER POINT.
00144 !| MASSOU         |-->| MASS OF TRACER ADDED BY SOURCE TERM
00145 !|                |   | SEE DIFSOU
00146 !| MESH           |-->| MESH STRUCTURE
00147 !| MSK            |-->| IF YES, THERE IS MASKED ELEMENTS.
00148 !| NPTFR          |-->| NUMBER OF BOUNDARY POINTS
00149 !| OPDADV         |-->| SCHEME OPTION FOR THE ADVECTION OF TRACERS
00150 !|                |   | WITH N SCHEME:
00151 !|                |   |  1: EXPLICIT
00152 !|                |   |  2: IMPLICIT
00153 !|                |   |  3: PREDICTOR-CORRECTOR 2ND ORDER IN TIME (MONOTONICITY NOT PROVED)
00154 !|                |   |  4: IMPLICIT PREDICTOR EXPLICIT CORRECTOR 2ND ORDER IN TIME
00155 !|                |   | WITH PSI SCHEME:
00156 !|                |   |  1: EXPLICIT
00157 !|                |   |  2: PREDICTOR-CORRECTOR 1ST ORDER IN TIME
00158 !|                |   |  3: PREDICTOR-CORRECTOR 2ND ORDER IN TIME (MONOTONICITY NOT PROVED)
00159 !|                |   |  4: IMPLICIT PREDICTOR EXPLICIT CORRECTOR 2ND ORDER IN TIME
00160 !| OPDTRA         |-->| OPTION FOR THE DIFFUSION OF TRACERS
00161 !| OPTSOU         |-->| TYPE OF SOURCES
00162 !|                |   | 1: NORMAL
00163 !|                |   | 2: DIRAC
00164 !| S              |-->| VOID STRUCTURE
00165 !| SM             |-->| SOURCE TERMS.
00166 !| SMH            |-->| SOURCE TERM IN CONTINUITY EQUATION
00167 !| SMI            |-->| IMPLICIT SOURCE TERM
00168 !| SOLSYS         |-->| 1 OR 2. IF 2 ADVECTION FIELD IS UCONV + DM1*GRAD(ZCONV)
00169 !| T1             |<->| WORK BIEF_OBJ STRUCTURE
00170 !| T2             |<->| WORK BIEF_OBJ STRUCTURE
00171 !| T3             |<->| WORK BIEF_OBJ STRUCTURE
00172 !| T4             |<->| WORK BIEF_OBJ STRUCTURE
00173 !| T5             |<->| WORK BIEF_OBJ STRUCTURE
00174 !| T6             |<->| WORK BIEF_OBJ STRUCTURE
00175 !| T7             |<->| WORK BIEF_OBJ STRUCTURE
00176 !| T8             |<->| WORK BIEF_OBJ STRUCTURE
00177 !| TE1            |<->| WORK BIEF_OBJ STRUCTURE FOR ELEMENTS
00178 !| TRAIN          |-->| VALUE OF TRACER IN THE RAIN
00179 !| UCONV,VCONV    |-->| ADVECTION VELOCITY FIELD
00180 !| UNSV2D         |-->| INVERSE OF INTEGRALS OF TEST FUNCTIONS
00181 !| VOLU2D         |-->| INTEGRAL OF TEST FUNCTIONS, NOT ASSEMBLED IN PARALLEL
00182 !| V2DPAR         |-->| INTEGRAL OF TEST FUNCTIONS, ASSEMBLED IN PARALLEL
00183 !| VISC           |-->| VISCOSITY COEFFICIENTS ALONG X,Y AND Z .
00184 !|                |   | IF P0 : PER ELEMENT
00185 !|                |   | IF P1 : PERR POINT
00186 !| VISC_S         |<->| WORK ARRAY FOR SAVING VISC
00187 !| YAFLBOR        |-->| IF YES FLBOR IS GIVEN
00188 !| YASMH          |-->| IF YES, SMH MUST BE TAKEN INTO ACCOUNT
00189 !| YASMI          |-->| IF YES, SMI MUST BE TAKEN INTO ACCOUNT
00190 !| ZCONV          |-->| THE PIECE-WISE CONSTANT PART OF ADVECTION FIELD
00191 !|                |   | IS DM1*GRAD(ZCONV), SEE SOLSYS.
00192 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00193 !
00194       USE BIEF, EX_CVTRVF => CVTRVF
00195 !
00196       IMPLICIT NONE
00197       INTEGER LNG,LU
00198       COMMON/INFO/LNG,LU
00199 !
00200 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00201 !
00202       INTEGER, INTENT(IN)             :: OPDTRA,OPTSOU,KDIR,NPTFR,SOLSYS
00203       INTEGER, INTENT(IN)             :: KDDL,IOPT,OPTADV
00204       INTEGER, INTENT(INOUT)          :: LIMTRA(NPTFR)
00205       DOUBLE PRECISION, INTENT(IN)    :: DT,AGGLOH,TRAIN
00206       DOUBLE PRECISION, INTENT(INOUT) :: MASSOU
00207       LOGICAL, INTENT(IN)             :: BILAN,CONV,YASMH,YAFLBOR
00208       LOGICAL, INTENT(IN)             :: DIFT,MSK,ENTET,YASMI,RAIN
00209       TYPE(BIEF_OBJ), INTENT(IN)      :: MASKEL,H,HN,DM1,ZCONV,MASKPT
00210       TYPE(BIEF_OBJ), INTENT(IN)      :: VOLU2D,V2DPAR,UNSV2D,HPROP
00211       TYPE(BIEF_OBJ), INTENT(INOUT)   :: F,SM,HNT,HT
00212       TYPE(BIEF_OBJ), INTENT(IN)      :: UCONV,VCONV,FN,SMI,SMH
00213       TYPE(BIEF_OBJ), INTENT(INOUT)   :: FBOR
00214       TYPE(BIEF_OBJ), INTENT(INOUT)   :: TE1,FLBORTRA
00215       TYPE(BIEF_OBJ), INTENT(INOUT)   :: T1,T2,T4,T5,T6,T7,T8
00216       TYPE(BIEF_OBJ), INTENT(IN)      :: FSCEXP,S,MASKTR
00217       TYPE(BIEF_OBJ), INTENT(IN)      :: VISC_S,VISC,PLUIE
00218       TYPE(BIEF_MESH) :: MESH
00219       TYPE(BIEF_OBJ), INTENT(IN),    TARGET :: FLBOR
00220       TYPE(BIEF_OBJ), INTENT(INOUT), TARGET :: T3
00221 !
00222 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00223 !
00224       INTEGER IELMF,I,IOPT1,IOPT2,N,VARIANT,IELEM,I1,I2,I3,ITER
00225 !
00226 !-----------------------------------------------------------------------
00227 !
00228       DOUBLE PRECISION DT_REMAIN,DDT,TDT,SECU,COEMIN,COESOU
00229 !
00230       CHARACTER(LEN=16) FORMUL
00231       DOUBLE PRECISION, POINTER, DIMENSION(:) :: SAVE_HT,SAVE_HNT
00232       DOUBLE PRECISION, POINTER, DIMENSION(:) :: FXMAT,FXMATPAR
00233       DOUBLE PRECISION, POINTER, DIMENSION(:) :: FMIN,FMAX
00234       TYPE(BIEF_OBJ), POINTER                 :: FLBOUND
00235 !
00236       DOUBLE PRECISION P_DMIN,C
00237       EXTERNAL         P_DMIN
00238 !
00239       INTEGER NITMAX,NIT
00240       DATA NITMAX/200/
00241 !
00242 !-----------------------------------------------------------------------
00243 !
00244       SAVE_HT =>HT%R
00245       SAVE_HNT=>HNT%R
00246       FXMAT=>MESH%MSEG%X%R(1:MESH%NSEG)
00247 !     IN PARALLEL MODE, ASSEMBLED AND NON ASSEMBLED VERSIONS ARE DIFFERENT
00248       IF(NCSIZE.GT.1) THEN
00249         FXMATPAR=>MESH%MSEG%X%R(MESH%NSEG+1:2*MESH%NSEG)
00250       ELSE
00251         FXMATPAR=>MESH%MSEG%X%R(1:MESH%NSEG)
00252       ENDIF
00253 !
00254 !-----------------------------------------------------------------------
00255 !
00256 !     EXTRACTS THE OPTIONS
00257 !
00258       IOPT2=IOPT/10
00259       IOPT1=IOPT-10*IOPT2
00260 !
00261 !-----------------------------------------------------------------------
00262 !
00263 !     IELMF = F%ELM
00264 !     FORCED TO LINEAR
00265       IELMF=11
00266 !
00267 !     TAKES MASS-LUMPING INTO ACCOUNT IN THE CONTINUITY EQUATION
00268 !
00269       IF(ABS(1.D0-AGGLOH).GT.1.D-8) THEN
00270         CALL VECTOR(HT ,'=','MASVEC          ',IELMF,
00271      &              1.D0-AGGLOH,H ,S,S,S,S,S,MESH,MSK,MASKEL)
00272         CALL VECTOR(HNT,'=','MASVEC          ',IELMF,
00273      &              1.D0-AGGLOH,HN,S,S,S,S,S,MESH,MSK,MASKEL)
00274         IF(NCSIZE.GT.1) THEN
00275           CALL PARCOM(HT ,2,MESH)
00276           CALL PARCOM(HNT,2,MESH)
00277         ENDIF
00278         CALL OS('X=YZ    ',X=HT ,Y=HT ,Z=UNSV2D)
00279         CALL OS('X=YZ    ',X=HNT,Y=HNT,Z=UNSV2D)
00280         CALL OS('X=X+CY  ',X=HT ,Y=H  ,C=AGGLOH)
00281         CALL OS('X=X+CY  ',X=HNT,Y=HN ,C=AGGLOH)
00282       ELSE
00283 !       CALL OS('X=Y     ',X=HT ,Y=H )
00284 !       CALL OS('X=Y     ',X=HNT,Y=HN)
00285         HT%R =>H%R
00286         HNT%R=>HN%R
00287       ENDIF
00288 !
00289 !     IF NO FLBOR IS GIVEN, IT IS COMPUTED HERE IN T3
00290 !
00291       IF(YAFLBOR) THEN
00292         FLBOUND => FLBOR
00293       ELSE
00294 !       MASK=5 FOR NON NEUMANN BOUNDARIES IN DIFFIN
00295         CALL VECTOR(T3,'=','FLUBDF          ',1,1.D0,HPROP,HPROP,HPROP,
00296      &              UCONV,VCONV,VCONV,MESH,.TRUE.,MASKTR%ADR(5)%P)
00297         FLBOUND => T3
00298       ENDIF
00299 !
00300 !-----------------------------------------------------------------------
00301 !
00302 !     CORRECTION OF THE BOUNDARY CONDITIONS
00303 !
00304 !-----------------------------------------------------------------------
00305 !
00306 !     A SIMILAR CORRECTION IS DONE IN DIFFIN, BUT IT MAY BE INCOMPATIBLE
00307 !     AS U*N MAY NOT BE OF THE SAME SIGN AS FLBOR HERE (RARE BUT ALREADY
00308 !     SEEN).
00309 !
00310       CALL CPSTVC(HN,T7)
00311 !     JUST IN CASE INTERNAL POINTS HAVE NON INITIALISED VALUES
00312       CALL OS('X=0     ',X=T7)
00313       DO I=1,MESH%NPTFR
00314         N=MESH%NBOR%I(I)
00315         T7%R(N)=FLBOUND%R(I)
00316       ENDDO
00317       IF(NCSIZE.GT.1) CALL PARCOM(T7,2,MESH)
00318       DO I=1,MESH%NPTFR
00319         N=MESH%NBOR%I(I)
00320         IF(LIMTRA(I).EQ.KDIR.AND.T7%R(N).GT.0.D0) THEN
00321           LIMTRA(I)=KDDL
00322         ELSEIF(LIMTRA(I).EQ.KDDL.AND.T7%R(N).LT.0.D0) THEN
00323           LIMTRA(I)=KDIR
00324           FBOR%R(I)=FN%R(N)
00325         ENDIF
00326       ENDDO
00327 !
00328 !     T7 MAY NOW BE REUSED
00329 !
00330 !-----------------------------------------------------------------------
00331 !
00332 !     INITIALISES THE TRACER FLUX AT THE BOUNDARY
00333 !
00334 
00335       DO I=1,MESH%NPTFR
00336         IF(LIMTRA(I).EQ.KDIR) THEN
00337 !         FLBOR IS NOT ASSEMBLED IN PARALLEL MODE
00338           FLBORTRA%R(I)=FLBOUND%R(I)*FBOR%R(I)
00339         ELSE
00340 !         FOR KDDL, WILL BE DONE IN TVF OR TVF_2 OR TVF_IMP
00341           FLBORTRA%R(I)=0.D0
00342         ENDIF
00343       ENDDO
00344 !
00345 !     COMPUTES THE FLUXES PHIIJ = FXMAT
00346 !
00347       FORMUL='HUGRADP         '
00348       IF(SOLSYS.EQ.2) FORMUL(8:8)='2'
00349       CALL VECTOR(T2,'=',FORMUL,IELMF,-1.D0,
00350      &            HPROP,DM1,ZCONV,UCONV,VCONV,VCONV,MESH,MSK,MASKEL)
00351 !                 T2 AS HUGRADP IS NOT USED AS AN ASSEMBLED VECTOR
00352 !                 BUT TO GET THE NON ASSEMBLED FORM MESH%W
00353       NIT=0
00354       DT_REMAIN=DT
00355       TDT=0.D0
00356       CALL CPSTVC(H ,T5)
00357       CALL CPSTVC(H ,T4)
00358       CALL CPSTVC(F,T8)
00359 !
00360 !     T4 WILL BE F PROGRESSIVELY UPDATED
00361 !     T5 WILL BE THE DEPTH AT THE END OF THE SUB-TIMESTEP
00362 !     (INITIALISED HERE TO CALL  CFLVF)
00363 !
00364       DO I=1,HN%DIM1
00365         T4%R(I)=FN%R(I)
00366         T5%R(I)=HNT%R(I)
00367       ENDDO
00368 !
00369 !     T1 WILL BE THE DEPTH ACCORDING TO THE CONTINUITY EQUATION
00370 !
00371       IF(IOPT2.EQ.1) THEN
00372         DO I=1,HN%DIM1
00373           T1%R(I)=HNT%R(I)
00374         ENDDO
00375       ENDIF
00376 !
00377 100   CONTINUE
00378       NIT=NIT+1
00379 !
00380 !---------------------------------------
00381 ! VARIOUS OPTIONS TO COMPUTE THE FLUXES
00382 !---------------------------------------
00383 !
00384       IF(NIT.EQ.1.OR.IOPT1.EQ.3) THEN
00385         CALL FLUX_EF_VF(FXMAT,MESH%W%R,MESH%NSEG,MESH%NELEM,
00386      &                  MESH%ELTSEG%I,MESH%ORISEG%I,
00387      &                  MESH%IKLE%I,.TRUE.,2    ,T4)
00388 !                                          IOPT1 HERE FORCED TO N SCHEME
00389 !       CANCELS FLUXES TO AND FROM MASKED POINTS
00390         IF(MSK) THEN
00391           CALL FLUX_MASK(FXMAT,MESH%NSEG,
00392      &                   MESH%GLOSEG%I,MESH%GLOSEG%DIM1,MASKPT%R)
00393         ENDIF
00394 !       ASSEMBLES THE FLUXES AT INTERFACES IN PARALLEL MODE, THIS
00395 !       IS FOR UPWINDING (STORED IN SECOND DIMENSION OF MESH%MSEG)
00396         IF(NCSIZE.GT.1) THEN
00397           CALL OV('X=Y     ',FXMATPAR,FXMAT,FXMAT,0.D0,MESH%NSEG)
00398           CALL PARCOM2_SEG(FXMATPAR,FXMATPAR,FXMATPAR,
00399      &                     MESH%NSEG,1,2,1,MESH,1,11)
00400         ENDIF
00401       ENDIF
00402 !
00403 !--------------------------------------------
00404 ! DETERMINES THE LARGEST ADMISSIBLE TIMESTEP
00405 !--------------------------------------------
00406 !
00407 !     THIS COULD BE PUT OUTSIDE THE LOOP, BUT T7 USED LATER IN THE LOOP...
00408 !
00409 !     IN CFLVF, T7 WILL BE FLBOR WITH A DIMENSION NPOIN
00410       CALL OS('X=0     ',X=T7)
00411       DO I=1,MESH%NPTFR
00412         N=MESH%NBOR%I(I)
00413         T7%R(N)=FLBOUND%R(I)
00414       ENDDO
00415       IF(NCSIZE.GT.1) CALL PARCOM(T7,2,MESH)
00416 !
00417 !     MASKS FLBOR IF(MSK)
00418 !
00419       IF(MSK) CALL OS('X=XY    ',X=T7,Y=MASKPT)
00420 !
00421 !     COMPUTES THE MAXIMUM TIMESTEP ENSURING MONOTONICITY
00422 !     ACCORDING TO THEORY
00423 !
00424       COESOU=0.D0
00425       COEMIN=0.D0
00426       SECU=0.99D0
00427 !
00428       IF(IOPT1.EQ.3.AND.OPTADV.EQ.2) THEN
00429         SECU=1.D0
00430         COEMIN=-1.D0
00431         COESOU=COEMIN
00432       ENDIF
00433 !
00434       CALL CFLVF(DDT,T5%R,HT%R,FXMAT,FXMATPAR,
00435 !                               FLBOR%R(NPOIN)
00436      &           V2DPAR%R,DT_REMAIN,T7%R   ,SMH%R,
00437      &           YASMH,T8,MESH%NSEG,MESH%NPOIN,MESH%NPTFR,
00438      &           MESH%GLOSEG%I,MESH%GLOSEG%DIM1,MESH,MSK,MASKPT,
00439      &           RAIN,PLUIE%R,T4%R,MESH%NELEM,MESH%IKLE%I,
00440      &           LIMTRA,KDIR,FBOR%R,FSCEXP%R,TRAIN,MESH%NBOR%I,
00441      &           T2,T6,SECU,COEMIN,COESOU)
00442 !
00443 !     NOW RECOMPUTING THE PSI FLUXES (THE N FLUXES HAVE BEEN
00444 !     USED FOR THE STABILITY CRITERION).
00445 !
00446       IF(IOPT1.EQ.3) THEN
00447         CALL FLUX_EF_VF(FXMAT,MESH%W%R,MESH%NSEG,MESH%NELEM,
00448      &                  MESH%ELTSEG%I,MESH%ORISEG%I,
00449      &                  MESH%IKLE%I,.TRUE.,IOPT1,T4)
00450 !       CANCELS FLUXES TO AND FROM MASKED POINTS
00451         IF(MSK) THEN
00452           CALL FLUX_MASK(FXMAT,MESH%NSEG,
00453      &                   MESH%GLOSEG%I,MESH%GLOSEG%DIM1,MASKPT%R)
00454         ENDIF
00455 !       ASSEMBLES THE FLUXES AT INTERFACES IN PARALLEL MODE, THIS
00456 !       IS FOR UPWINDING (STORED IN SECOND DIMENSION OF MESH%MSEG)
00457         IF(NCSIZE.GT.1) THEN
00458           CALL OV('X=Y     ',FXMATPAR,FXMAT,FXMAT,0.D0,MESH%NSEG)
00459           CALL PARCOM2_SEG(FXMATPAR,FXMATPAR,FXMATPAR,
00460      &                     MESH%NSEG,1,2,1,MESH,1,11)
00461         ENDIF
00462       ENDIF
00463 !
00464 !
00465       IF(NCSIZE.GT.1) DDT=P_DMIN(DDT)
00466       DDT=MIN(DDT,DT_REMAIN)
00467 !
00468 !     T2 WILL TAKE THE SUCCESSIVE VALUES OF H
00469 !     AT THE BEGINNING OF THE SUB-TIMESTEP
00470 !     WARNING: T2 ALSO USED WITH IOPT2=1, BUT SO FAR PREDICTOR-CORRECTOR
00471 !              NOT USED WITH THIS OPTION
00472 !
00473       IF(IOPT1.EQ.3.AND.OPTADV.EQ.2) THEN
00474         DO I=1,HN%DIM1
00475           T2%R(I)=HNT%R(I)+TDT*(HT%R(I)-HNT%R(I))/DT
00476         ENDDO
00477       ENDIF
00478 !
00479       TDT=TDT+DDT
00480 !
00481 !     T5 WILL TAKE THE SUCCESSIVE VALUES OF H
00482 !     AT THE END OF THE SUB-TIMESTEP
00483 !
00484       DO I=1,HN%DIM1
00485         T5%R(I)=HNT%R(I)+TDT*(HT%R(I)-HNT%R(I))/DT
00486       ENDDO
00487 !
00488 !     IN TVF FACTOR HT/HLIN MAY TRIGGER DIVERGENCE FOR DRY POINTS
00489 !
00490       IF(MSK) THEN
00491         DO I=1,HN%DIM1
00492           IF(MASKPT%R(I).LT.0.5D0) T5%R(I)=HT%R(I)
00493         ENDDO
00494       ENDIF
00495 !
00496 !------------------------------------
00497 !  FINAL RESOLUTION OR PREDICTOR STEP
00498 !------------------------------------
00499 !
00500       CALL TRACVF(F,FN,FSCEXP,HT,HNT,FXMAT,FXMATPAR,V2DPAR,UNSV2D,
00501      &            DDT,FLBOUND,FBOR,SMH,YASMH,T1,T2,T4,T5,T6,T7,T8,
00502      &            MESH,LIMTRA,KDIR,KDDL,OPTSOU,IOPT2,FLBORTRA,MSK,
00503      &            DT,RAIN,PLUIE,TRAIN)
00504 !
00505 !--------------------------------
00506 !  CORRECTOR STEP FOR PSI SCHEME
00507 !--------------------------------
00508 !
00509       IF(IOPT1.EQ.3.AND.OPTADV.EQ.2) THEN
00510 !
00511         CALL FLUX_EF_VF_2(MESH%W%R,MESH%NELEM,
00512      &                    MESH%IKLE%I,IOPT1,MESH%NPOIN,T4,
00513 !    &                    FI_I,FSTAR,
00514      &                    T8%R,F%R,T5%R,MESH%SURFAC%R,DDT)
00515         IF(NCSIZE.GT.1) CALL PARCOM(T8,2,MESH)
00516 !
00517 !       COMPUTE THE SECOND ORDER CORRECTION FORM
00518 !
00519         CALL TVF_2(F%R,FN%R,T4%R,UNSV2D%R,DDT,
00520      &             FLBOUND%R,T7%R,FBOR%R,SMH%R,YASMH,FSCEXP%R,
00521      &             MESH%NPOIN,MESH%NPTFR,MESH%NBOR%I,
00522      &             LIMTRA,KDIR,KDDL,
00523      &             OPTSOU,T5%R,IOPT2,FLBORTRA%R,DDT/DT,RAIN,
00524      &             PLUIE%R,TRAIN,T8%R)
00525 !                                FI_I
00526 !
00527       ENDIF
00528 !
00529 !-----------------
00530 ! END CORRECTOR STEP
00531 !-----------------
00532 !
00533       DO I=1,HN%DIM1
00534 !       T4 IS F(N+1)
00535         T4%R(I)=F%R(I)
00536       ENDDO
00537       IF(IOPT2.EQ.1) THEN
00538         DO I=1,HN%DIM1
00539           T1%R(I)=T2%R(I)
00540         ENDDO
00541       ENDIF
00542 !
00543       DT_REMAIN=DT_REMAIN-DDT
00544 !
00545       IF(DT_REMAIN.NE.0.D0.AND.NIT.LT.NITMAX) GO TO 100
00546 !
00547       IF(NIT.GE.NITMAX) THEN
00548         IF(LNG.EQ.1) WRITE(LU,900) NIT
00549         IF(LNG.EQ.2) WRITE(LU,901) NIT
00550 900     FORMAT(1X,'CVTRVF : ',1I6,' SOUS-ITERATIONS DEMANDEES POUR LE'
00551      &   ,/,1X,   '         SCHEMA VF. DIMINUER LE PAS DE TEMPS')
00552 901     FORMAT(1X,'CVTRVF: ',1I6,' SUB-ITERATIONS REQUIRED FOR THE'
00553      &   ,/,1X,   '         VF SCHEME. DECREASE THE TIME-STEP')
00554         CALL PLANTE(1)
00555         STOP
00556       ELSEIF(ENTET) THEN
00557         IF(LNG.EQ.1) WRITE(LU,902) NIT
00558         IF(LNG.EQ.2) WRITE(LU,903) NIT
00559 902     FORMAT(1X,'CVTRVF (BIEF) : ',1I6,' SOUS-ITERATIONS')
00560 903     FORMAT(1X,'CVTRVF (BIEF): ',1I6,' SUB-ITERATIONS')
00561       ENDIF
00562 !
00563 !-----------------------------------------------------------------------
00564 !
00565 !     EXPLICIT SOURCE TERM
00566 !
00567       DO I = 1,MESH%NPOIN
00568         F%R(I) = F%R(I)+DT*SM%R(I)
00569       ENDDO
00570 !
00571 !     IMPLICIT SOURCE TERM
00572 !
00573       IF(YASMI) THEN
00574         DO I = 1,MESH%NPOIN
00575           F%R(I) = F%R(I)/(1.D0-DT*SMI%R(I)/MAX(H%R(I),1.D-15))
00576         ENDDO
00577       ENDIF
00578 !
00579 !-----------------------------------------------------------------------
00580 !
00581 !     LOCAL MASS BALANCE (FOR CHECKING PURPOSES)
00582 !
00583 !     CALL OS('X=Y-Z   ',X=T7,Y=T5,Z=HT)
00584 !     PRINT*,'DIFFERENCE ENTRE H RECALCULE ET H : ',DOTS(T7,T7)
00585 !     CHECKS THE TRACER EQUATION
00586 !     CALL CPSTVC(FBOR,T4)
00587 !     T4 : F AT THE BOUNDARIES AS TAKEN FOR THE BOUNDARY FLUXES
00588 !     DO I=1,NPTFR
00589 !       IF(LIMTRA(I).EQ.KDIR) THEN
00590 !         T4%R(I)=FBOR%R(I)
00591 !       ELSE
00592 !         T4%R(I)=FN%R(MESH%NBOR%I(I))
00593 !       ENDIF
00594 !     ENDDO
00595 !     CALL OS('X=YZ    ',X=T6,Y=FN,Z=HNT)
00596 !     CALL OS('X=YZ    ',X=T7,Y=F ,Z=HT )
00597 !     MASSETN=P_DOTS(V2DPAR,T6,MESH)
00598 !     MASSET =P_DOTS(V2DPAR,T7,MESH)
00599 !     FXT2   =P_DOTS(FLBOR,T4,MESH)
00600 !     PRINT*,'MASSE INIT: ',MASSETN,' MASSE FINALE: ',MASSET
00601 !     PRINT*,'FLUX: ',FXT2
00602 !     MASSETN = MASSETN - FXT2*DT
00603 !     TSOU=0.D0
00604 !     IF(YASMH) THEN
00605 !       IF(OPTSOU.EQ.1) THEN
00606 !         DO I=1,MESH%NPOIN
00607 !           MASSETN=MASSETN
00608 !    *             +DT*V2DPAR%R(I)*SMH%R(I)*(FSCEXP%R(I)+FN%R(I))
00609 !           TSOU=TSOU+DT*V2DPAR%R(I)*SMH%R(I)*(FSCEXP%R(I)+FN%R(I))
00610 !         ENDDO
00611 !       ELSEIF(OPTSOU.EQ.2) THEN
00612 !         DO I=1,MESH%NPOIN
00613 !           MASSETN=MASSETN
00614 !    *             +DT*SMH%R(I)*(FSCEXP%R(I)+FN%R(I))
00615 !           TSOU=TSOU+DT*SMH%R(I)*(FSCEXP%R(I)+FN%R(I))
00616 !         ENDDO
00617 !       ENDIF
00618 !     ENDIF
00619 !     PRINT*,'CREATION PAR SOURCE : ',TSOU
00620 !     PRINT*,'ERREUR DE MASSE DE TRACEUR VF : ',MASSETN-MASSET
00621 !     CHECKS THE CONTINUITY EQUATION
00622 !     DO I = 1,MESH%NPOIN
00623 !       T5%R(I)=V2DPAR%R(I)*(HT%R(I)-HNT%R(I))
00624 !     ENDDO
00625 !     DO I = 1,MESH%NSEG
00626 !       T5%R(MESH%GLOSEG%I(I)) =
00627 !    *  T5%R(MESH%GLOSEG%I(I)) + DT*MESH%MSEG%X%R(I)
00628 !       T5%R(MESH%GLOSEG%I(I+MESH%NSEG)) =
00629 !    *  T5%R(MESH%GLOSEG%I(I+MESH%NSEG)) - DT*MESH%MSEG%X%R(I)
00630 !     ENDDO
00631 !     DO I = 1,MESH%NPTFR
00632 !       T5%R(MESH%NBOR%I(I))=T5%R(MESH%NBOR%I(I))+DT*FLBOR%R(I)
00633 !     ENDDO
00634 !     IF(YASMH) THEN
00635 !       IF(OPTSOU.EQ.1) THEN
00636 !         DO I = 1,MESH%NPOIN
00637 !           T5%R(I)=T5%R(I)-DT*V2DPAR%R(I)*SMH%R(I)
00638 !         ENDDO
00639 !       ELSEIF(OPTSOU.EQ.2) THEN
00640 !         DO I = 1,MESH%NPOIN
00641 !           T5%R(I)=T5%R(I)-DT*SMH%R(I)
00642 !         ENDDO
00643 !       ENDIF
00644 !     ENDIF
00645 !     MASSET=0.D0
00646 !     MASSETN = 0.D0
00647 !     DO I = 1,MESH%NPOIN
00648 !       MASSET=MASSET+T5%R(I)
00649 !       MASSETN=MAX(MASSETN,ABS(T5%R(I)))
00650 !     ENDDO
00651 !     PRINT*,'ERREUR DE MASSE GLOBALE : ',MASSET,' LOCALE : ',MASSETN
00652 !
00653 !-----------------------------------------------------------------------
00654 !
00655 !     RETURNS POINTERS HT AND HNT
00656 !
00657       HT%R =>SAVE_HT
00658       HNT%R=>SAVE_HNT
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