init_avai.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\sisyphe\init_avai.f
00002 !
00092                      SUBROUTINE INIT_AVAI
00093 !                    ********************
00094 !
00095 !
00096 !***********************************************************************
00097 ! SISYPHE   V6P2                                   21/07/2011
00098 !***********************************************************************
00099 !
00100 !
00101 !
00102 !
00103 !
00104 !
00105 !
00106 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00107 !
00108       USE BIEF
00109       USE DECLARATIONS_TELEMAC
00110       USE DECLARATIONS_SISYPHE
00111       USE INTERFACE_SISYPHE, EX_INIT_AVAI=> INIT_AVAI
00112 !
00113       IMPLICIT NONE
00114 !
00115       INTEGER LNG,LU
00116       COMMON/INFO/LNG,LU
00117 !
00118 !-----------------------------------------------------------------------
00119 !
00120       DOUBLE PRECISION P_DSUM
00121       EXTERNAL         P_DSUM
00122 !
00123       INTEGER I,J,K,NMAXI
00124 !
00125       DOUBLE PRECISION HAUTSED,TEST1
00126 !
00127 !-----------------------------------------------------------------------
00128 !
00129       NMAXI = 0
00130 !
00131 !     THE INITIAL NUMBER OF LAYERS, THEIR THICKNESS AND THEIR COMPOSITION
00132 !     ARE SET BY THE USER
00133 !
00134 !     NOTE: WHEN COMPUTATION CONTINUED INIT_COMPO MUST NOT
00135 !           CHANGE ES AND AVAIL
00136 !
00137       IF(DEBU) THEN
00138 !       TENTATIVE VALUE, THIS IS TO CHECK
00139 !       ADDED BY JMH 30/06/2010
00140         DO J=1,NPOIN
00141           I=1
00142           IF(NOMBLAY.GE.2) THEN
00143             DO K=2,NOMBLAY
00144               IF(ES(J,K).GT.0.D0) I = I + 1
00145             ENDDO
00146           ENDIF
00147           NLAYER%I(J)=I
00148 !         CHECKING ALL LAYERS AND CORRECTING AVAIL
00149 !         DUE TO POSSIBLE SHIFT OF SINGLE PRECISION STORAGE
00150           DO K=1,NLAYER%I(J)
00151             TEST1=0.D0
00152             DO I=1,NSICLA
00153               TEST1=TEST1+AVAIL(J,K,I)
00154             ENDDO
00155             IF(TEST1.GT.1.D-10.AND.(TEST1-1.D0)**2.GT.1.D-10) THEN
00156               DO I=1,NSICLA
00157                 AVAIL(J,K,I)=AVAIL(J,K,I)/TEST1
00158               ENDDO
00159             ENDIF
00160           ENDDO
00161         ENDDO
00162       ELSE
00163 !
00164         CALL INIT_COMPO(IT1%I)
00165 !
00166         DO J=1,NPOIN
00167 !
00168 !       NOMBLAY IS THE MAXIMUM NUMBER OF LAYERS ALLOWED
00169 !
00170         NLAYER%I(J) = IT1%I(J)
00171         IF(NLAYER%I(J).GT.NOMBLAY) THEN
00172           IF(LNG.EQ.1) WRITE(LU,1800) NOMBLAY
00173           IF(LNG.EQ.2) WRITE(LU,1815) NOMBLAY
00174           CALL PLANTE(1)
00175           STOP
00176         ENDIF
00177 !
00178 !       THE HEIGHT OF SEDIMENT (SUM OF ES) MUST NOT BE MORE THAN ZF-ZR
00179 !       IF SO, THE HEIGHT OF THE LAST LAYER IS REDUCED
00180 !       IF THERE ARE LAYERS UNDER ZR, THEY ARE NOT TAKEN INTO ACCOUNT
00181         HAUTSED = 0.D0
00182         DO K=1,IT1%I(J)
00183           IF(HAUTSED + ES(J,K) .GE. ZF%R(J) - ZR%R(J)) THEN
00184             ES(J,K) = ZF%R(J) - ZR%R(J) -  HAUTSED
00185             NLAYER%I(J) = K
00186             HAUTSED = HAUTSED + ES(J,K)
00187             EXIT
00188           ENDIF
00189           HAUTSED = HAUTSED + ES(J,K)
00190         ENDDO
00191 !
00192 !       OTHER LAYERS SET TO 0.D0
00193 !
00194         IF(NLAYER%I(J).LT.NOMBLAY) THEN
00195           DO K=NLAYER%I(J)+1,NOMBLAY
00196             ES(J,K) = 0.D0
00197           ENDDO
00198         ENDIF
00199 !
00200 !       THE THICKNESS OF THE LAST LAYER IS ENLARGED SO THAT
00201 !       THE HEIGHT OF SEDIMENT (SUM OF ES) IS EQUAL TO ZF-ZR
00202 !
00203         IF(HAUTSED.LT.ZF%R(J)-ZR%R(J)) THEN
00204           ES(J,NLAYER%I(J))=ES(J,NLAYER%I(J))+ZF%R(J)-ZR%R(J)-HAUTSED
00205         ENDIF
00206 !
00207         IF(NLAYER%I(J).GT.1) THEN
00208 !         IT IS ASSUMED THAT ELAY0 IS SMALLER THAN THE FIRST STRATUM
00209 !         NEED TO ADD THE CASE WHERE ELAY0 IS LARGER
00210           IF(ELAY0.GT.ES(J,1)) THEN
00211             IF(LNG.EQ.1) THEN
00212               WRITE(LU,*) 'COUCHE ACTIVE TROP GROSSE/STRATIFICATION'
00213             ENDIF
00214             IF(LNG.EQ.2) THEN
00215               WRITE(LU,*) ' ACTIVE LAYER TOO BIG/STRATIFICATION'
00216             ENDIF
00217             CALL PLANTE(1)
00218             STOP
00219           ENDIF
00220         ENDIF
00221 !
00222 !       THE FIRST STRATUM IS SEPARATED BETWEEN ACTIVE LAYER + ACTIVE STRATUM
00223 !
00224         IF(ELAY0.LT.ES(J,1)) THEN
00225           NLAYER%I(J) = NLAYER%I(J) + 1
00226           IF(NLAYER%I(J).GT.NOMBLAY) THEN
00227             IF(LNG.EQ.1) WRITE(LU,1800) NOMBLAY
00228             IF(LNG.EQ.2) WRITE(LU,1815) NOMBLAY
00229             CALL PLANTE(1)
00230             STOP
00231           ENDIF
00232 !         INDICES FOR ES AND AVAIL NEED TO BE OFFSET
00233           IF(NLAYER%I(J).GT.2) THEN
00234             DO K=NLAYER%I(J),3,-1
00235               ES(J,K) = ES(J,K-1)
00236             ENDDO
00237           ENDIF
00238           ES(J,2) = ES(J,1) - ELAY0
00239           ES(J,1) = ELAY0
00240           DO I=1,NSICLA
00241             DO K=NLAYER%I(J),2,-1
00242               AVAIL(J,K,I) = AVAIL(J,K-1,I)
00243             ENDDO
00244           ENDDO
00245         ENDIF
00246 !
00247         ENDDO !J
00248 !
00249       ENDIF
00250 !
00251       NMAXI=0
00252       DO J=1,NPOIN
00253         ELAY%R(J) = ES(J,1)
00254         IF(NLAYER%I(J).GT.1) THEN
00255           ESTRAT%R(J) = ES(J,2)
00256         ENDIF
00257 !
00258 !       UNUSED AVAIL ARE FILLED WITH ZEROS (IS IT USEFUL ???)
00259 !
00260         IF(NLAYER%I(J).LT.NOMBLAY) THEN
00261           DO I = 1, NSICLA
00262             DO K = NLAYER%I(J)+1,NOMBLAY
00263               AVAIL(J,K,I) = 0.D0
00264             ENDDO
00265           ENDDO
00266         ENDIF
00267         IF(NLAYER%I(J).GT.NMAXI) NMAXI = NLAYER%I(J)
00268       ENDDO
00269 !
00270 !     COMPUTES THE TOTAL VOLUME OF SEDIMENTS IN EACH CLASS
00271 !
00272       DO I=1,NSICLA
00273         VOLTOT(I)=0.D0
00274         DO J=1,NPOIN
00275           DO K=1,NLAYER%I(J)
00276             VOLTOT(I) = VOLTOT(I) + ES(J,K)*AVAIL(J,K,I)*VOLU2D%R(J)
00277           ENDDO
00278         ENDDO
00279       ENDDO
00280 !
00281       IF(NCSIZE.GT.1) THEN
00282         DO I=1,NSICLA
00283           VOLTOT(I) = P_DSUM(VOLTOT(I))
00284         ENDDO
00285       ENDIF
00286 !
00287       WRITE(LU,*) 'MAXIMUM INITIAL NUMBER OF LAYERS :',NMAXI
00288       DO I=1,NSICLA
00289         IF(LNG.EQ.1) THEN
00290           WRITE(LU,*)'VOLUME TOTAL DE LA CLASSE ',I ,' :',VOLTOT(I)
00291         ENDIF
00292         IF(LNG.EQ.1) THEN
00293           WRITE(LU,*)'TOTAL VOLUME OF CLASS ',I ,' :',VOLTOT(I)
00294         ENDIF
00295       ENDDO
00296 !
00297 !-----------------------------------------------------------------------
00298 !
00299 1800  FORMAT(1X,'IL Y A PLUS DE ',1I3,' COUCHES DANS LA STRATIFICATION')
00300 1815  FORMAT(1X,'THERE ARE MORE THAN ',1I3,' LAYERS IN STRATIFICATION')
00301 !
00302 !-----------------------------------------------------------------------
00303 !
00304       RETURN
00305       END

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