cost_function.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\telemac2d\cost_function.f
00002 !
00110                      SUBROUTINE COST_FUNCTION
00111 !                    ************************
00112 !
00113      &(JCOUT,OPTION,MODE)
00114 !
00115 !***********************************************************************
00116 ! TELEMAC2D   V6P1                                   21/08/2010
00117 !***********************************************************************
00118 !
00119 !
00120 !
00121 !
00122 !
00123 !
00124 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00125 !| JCOUT          |<->| COST FUNCTION
00126 !| MODE           |-->| FCT: COST FUNCTION
00127 !|                |   | GRD: GRADIENT COST FUNCTION
00128 !|                |   | RHS: RIGHT HAND SIDE
00129 !| OPTION         |-->| 1: COST FUNCTION COMPUTED WITH DEPTH
00130 !|                |   | 2: COST FUNCTION COMPUTED WITH CELERITY
00131 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00132 !
00133       USE BIEF
00134       USE DECLARATIONS_TELEMAC2D
00135       USE INTERFACE_TELEMAC2D, EX_COST_FUNCTION => COST_FUNCTION
00136 !
00137       IMPLICIT NONE
00138 !
00139 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00140 !
00141       DOUBLE PRECISION , INTENT(INOUT) :: JCOUT
00142       INTEGER , INTENT(IN)             :: OPTION
00143       CHARACTER(LEN=3) , INTENT(IN)    :: MODE
00144 !
00145 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00146 !
00147       COMMON/INFO/LNG,LU
00148       INTEGER LNG,LU
00149 !
00150       INTEGER I,J
00151       DOUBLE PRECISION C
00152 !
00153 !=======================================================================
00154 !
00155       IF(MODE.EQ.'FCT') THEN
00156 !
00157 !       HERE U,V AND H ARE GIVEN BY A CALL TO LITENR
00158 !
00159         IF(OPTION.EQ.1) THEN
00160 !
00161           DO I=1,NPOIN
00162 !
00163             JCOUT = JCOUT + ALPHA1%R(I) * (H%R(I)-HD%R(I))**2
00164      &                    + ALPHA2%R(I) * (U%R(I)-UD%R(I))**2
00165      &                    + ALPHA3%R(I) * (V%R(I)-VD%R(I))**2
00166 !
00167           ENDDO
00168 !
00169         ELSEIF(OPTION.EQ.2) THEN
00170 !
00171           DO I=1,NPOIN
00172 !
00173             JCOUT = JCOUT + ALPHA1%R(I) * GRAV *
00174      &              (SQRT(MAX(H%R(I),0.D0))-SQRT(MAX(HD%R(I),0.D0)))**2
00175      &                    + ALPHA2%R(I) * (U%R(I)-UD%R(I))**2
00176      &                    + ALPHA3%R(I) * (V%R(I)-VD%R(I))**2
00177 !
00178           ENDDO
00179 !
00180         ELSE
00181 !
00182           IF(LNG.EQ.1) THEN
00183             WRITE(LU,*) 'COST_FUNCTION : OPTION NON PREVUE : ',MODE
00184           ENDIF
00185           IF(LNG.EQ.2) THEN
00186             WRITE(LU,*) 'COST_FUNCTION: UNEXPECTED OPTION : ',MODE
00187           ENDIF
00188           CALL PLANTE(1)
00189           STOP
00190 !
00191 !       TEST ON OPTION
00192         ENDIF
00193 !
00194 !=======================================================================
00195 !
00196       ELSEIF(MODE.EQ.'GRD') THEN
00197 !
00198       IF(    INCLU2(ESTIME,'FROTTEMENT')
00199      &   .OR.INCLU2(ESTIME,'FRICTION'  )  ) THEN
00200 !
00201 !     IDENTIFICATION OF FRICTION
00202 !
00203       IF (KFROT.EQ.3.OR.KFROT.EQ.2) THEN
00204 !
00205       CALL SLOPES(TE3,ZF,MESH)
00206       CALL VECTOR(T1,'=','MASBAS          ',U%ELM,1.D0,
00207      &            T3,T3,T3,T3,T3,T3,MESH,.TRUE.,TE3)
00208 !
00209       CALL FRICTI(T3,T4,T2,T2,UN,VN,HN,CF,MESH,T2,T5,VERTIC,UNSV2D,MSK,
00210      &            MASKEL,HFROT)
00211 !
00212       CALL OS( 'X=XY    ' , T3 , T1 , T1 , C )
00213       CALL OS( 'X=XY    ' , T4 , T1 , T1 , C )
00214 !
00215       CALL OS( 'X=CXYZ  ' , T3 , QQ , UU , -2.D0 )
00216       CALL OS( 'X=CXYZ  ' , T4 , RR , VV , -2.D0 )
00217 !
00218       CALL OS( 'X=X+Y   ' , T3 , T4 , T4 , C )
00219       CALL OS( 'X=Y/Z   ' , T3 , T3 , CHESTR , C )
00220 !
00221       ELSE
00222 !
00223         IF(LNG.EQ.1) THEN
00224           WRITE(LU,*) 'COST_FUNCTION : FROTTEMENT NON TRAITE : ',KFROT
00225         ENDIF
00226         IF(LNG.EQ.2) THEN
00227           WRITE(LU,*) 'COST_FUNCTION: UNEXPECTED FRICTION LAW: ',KFROT
00228         ENDIF
00229         CALL PLANTE(1)
00230         STOP
00231 !
00232       ENDIF
00233 !
00234       ELSE
00235 !
00236         IF(LNG.EQ.1) THEN
00237           WRITE(LU,*) 'COST_FUNCTION : PARAMETRE NON PREVU :'
00238           WRITE(LU,*) ESTIME
00239           WRITE(LU,*) 'VERIFIER LE MOT CLEF : ESTIMATION DE PARAMETRES'
00240         ENDIF
00241         IF(LNG.EQ.2) THEN
00242           WRITE(LU,*) 'COST_FUNCTION: UNEXPECTED PARAMETER :'
00243           WRITE(LU,*) ESTIME
00244           WRITE(LU,*) 'CHECK THE KEY-WORD: PARAMETER ESTIMATION'
00245         ENDIF
00246         CALL PLANTE(1)
00247         STOP
00248 !
00249       ENDIF
00250 !
00251 ! COMPUTATION OF A GRADIENT FOR EVERY ZONE, AFTER BUILDING
00252 ! A GRADIENT VECTOR FOR EVERY POINT (IN T3)
00253 !
00254       IF(NZONE.GT.0) THEN
00255 !       IF ONE IDENTIFIES ONLY ONE PARAMETER S FOR A SET OF NODES
00256 !       GRADIENT DJ/DS IS THE SUM OF THE GRADIENTS OF EACH
00257 !       NODE OF THE SET
00258         DO J=1,NZONE
00259           DO I=1,NPOIN
00260             IF(ZONE%I(I).EQ.J) GRADJ%R(J)=GRADJ%R(J)+T3%R(I)
00261           ENDDO
00262         ENDDO
00263       ELSE
00264 !       NOTE JMH: HERE IT IS SUPPOSED THAT NPARAM = NPOIN
00265         CALL OS('X=X+Y   ',GRADJ,T3,T3,0.D0)
00266       ENDIF
00267 !
00268 !=======================================================================
00269 !
00270       ELSEIF(MODE.EQ.'RHS') THEN
00271 !
00272 !           IT    IT    IT
00273 !  TERMS 2 W   ( X   - M   ) OR EQUIVALENT DEPENDING ON THE OPTION
00274 !           IP    IP    IP
00275 !
00276       IF(OPTION.EQ.1) THEN
00277 !
00278         CALL OS( 'X=Y-Z   ', CV1 , HH , HD , C )
00279         CALL OS( 'X=Y-Z   ', CV2 , UU , UD , C )
00280         CALL OS( 'X=Y-Z   ', CV3 , VV , VD , C )
00281 !
00282         CALL OS( 'X=CXY   ', CV1 , ALPHA1 , ALPHA1 , 2.D0 )
00283         CALL OS( 'X=CXY   ', CV2 , ALPHA2 , ALPHA2 , 2.D0 )
00284         CALL OS( 'X=CXY   ', CV3 , ALPHA3 , ALPHA3 , 2.D0 )
00285 !
00286       ELSEIF(OPTION.EQ.2) THEN
00287 !
00288 !       HERE COST FUNCTION COMPUTED WITH CELERITY INSTEAD OF DEPTH
00289         CALL OS( 'X=SQR(Y)', T1  , HH , T1 , C   )
00290         CALL OS( 'X=SQR(Y)', T2  , HD , T2 , C   )
00291         CALL OS( 'X=Y-Z   ', T3  , T1 , T2 , C )
00292         CALL OS( 'X=Y/Z   ', CV1 , T3 , T1 , C )
00293         CALL OS( 'X=Y-Z   ', CV2 , UU , UD , C )
00294         CALL OS( 'X=Y-Z   ', CV3 , VV , VD , C )
00295 !
00296         CALL OS( 'X=CXY   ', CV1 , ALPHA1 , ALPHA1 , GRAV )
00297         CALL OS( 'X=CXY   ', CV2 , ALPHA2 , ALPHA2 , 2.D0 )
00298         CALL OS( 'X=CXY   ', CV3 , ALPHA3 , ALPHA3 , 2.D0 )
00299 !
00300       ENDIF
00301 !
00302 !=======================================================================
00303 !
00304       ELSE
00305 !
00306 !=======================================================================
00307 !
00308         IF(LNG.EQ.1) WRITE(LU,*) 'COST_FUNCTION : MODE NON PREVU : ',
00309      &                           MODE
00310         IF(LNG.EQ.2) WRITE(LU,*) 'COST_FUNCTION: UNEXPECTED MODE : ',
00311      &                           MODE
00312         CALL PLANTE(1)
00313         STOP
00314 !
00315 !=======================================================================
00316 !
00317 !     TEST ON MODE
00318       ENDIF
00319 !
00320 !-----------------------------------------------------------------------
00321 !
00322       RETURN
00323       END

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