Sources Esope | Cast3M

chmslv
C CHMSLV    SOURCE    CHAT      05/01/12    22:00:06     5004
      SUBROUTINE CHMSLV(IDSCHI,SP2,ITMAX,EPS,ICALCLOG,IEM)
C========================================================================
C
C      SP ISSU DE TRIOEF ( TRSOLV)
C
C Modif PhM
C a) Prise en charge du calcul en log de concentration
c    (argument ICALCLOG)
C
C========================================================================
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
      SEGMENT IDSCHI
           REAL*8 GK(NYDIM),AA(NYDIM,NXDIM),FF(NZDIM,NPDIM)
           INTEGER IDX(NXDIM),IDY(NYDIM),IDZ(NZDIM),IDP(NPDIM),NN(6)
           INTEGER IDECY(NYDIM),IONZ(NXDIM)
           CHARACTER*32 NAME(NXDIM),NAMESP(NYDIM)
      ENDSEGMENT
      SEGMENT SP2
           REAL*8 GX(NXDIM),XX(NXDIM),GS(NZDIM),SS(NZDIM)
           REAL*8 TOT(NXDIM),TOTAQ(NXDIM),TOTFIX(NXDIM),GKS(NZDIM)
           REAL*8 YY(NXDIM),ZZ(NXDIM,NXDIM),CC(NYDIM),GC(NYDIM)
      ENDSEGMENT
      NXDIM=IDX(/1)
      NYDIM=IDY(/1)
      NZDIM=IDZ(/1)
      NPDIM=IDP(/1)
      ITER=0
 
 
c      write(6,*)  'ICALCLOG',ICALCLOG
 
 
C
*************************************************************************
C PRISE EN COMPTE DES SOLUTIONS SOLIDES : METHODE 1
      N4S=0
*      IF(NZDIM.NE.0)THEN
*        I3S=NN(1)+NN(2)+NN(3)+1
*        I4S=NN(1)+NN(2)+NN(3)+NN(4)
*        DO 13 I7S=I3S,I4S
*          IDY7=IDY(I7S)
*          CALL CHIADY(IDZ,NZDIM,IDY7,ID7)
*          IF(ID7.NE.0)THEN
*            N4S=N4S+1
*            CALL CHMREX(IDSCHI,LGKMOD,LGKTMP,IDY7,4,5)
*            CALL CHMREX(IDSCHI,LGKMOD,LGKTMP,IDY7,5,4)
*          ENDIF
*   13 CONTINUE
*      ENDIF
**************************************************************************
 
      NC=NN(1)+NN(2)
      NX=NXDIM-NN(3)-NN(4)+N4S
C      write(6,*)'chmslv nc nx gk',nc,nx
C      write(6,120)(idy(j),gk(j),j=1,nc)
      DO 1 J=1,NX
      YY(J)=0.D0
 1    CONTINUE
 
 
 
1000  CONTINUE
  120 FORMAT(6(1X,I6,1PD15.6))
C  COMPLEXES
      DO 2 I=1,NC
      V=GK(I)
*      write(6,*)'chmslv idy',idy(i),'gk',gk(i)
      DO 3 J=1,NX
      V=V+AA(I,J)*GX(J)
*      write(6,*)'chmslv idx',idx(j)
*      write(6,*)'chmslv aa',aa(i,j),'gx',gx(j)
3     CONTINUE
      GC(I)=V
C
      CC(I)=10.D0**GC(I)
*      write(6,*)'chmslv idy',IDY(I),'c',CC(I),'logc',GC(I)
2     CONTINUE
 
 
C MOLE BALANCE
 
C      ----------------- A REVOIR----------
C     IF(IADSORB.GT.0) CALL TRSURT(SP1,SP2,SURFACE)
C      --------------------------------------
      DO 201 J=1,NX
      V=-TOT(J)
      DO 200 I=1,NC
      V=V+AA(I,J)*CC(I)
200   CONTINUE
      YY(J)=V
C
201   CONTINUE
 
C CALCUL Z (Jacobienne pour les iterations de newton)
      DO 301 I=1,NX
      DO 300 J=1,NX
      ZZ(I,J)=0.D0
300   CONTINUE
301   CONTINUE
      DO 410 J=1,NX
      DO 405 I=1,NC
      DO 400 K=1,NX
c modif PhM
         IF (ICALCLOG.EQ.1) THEN
            ZZ(J,K)=ZZ(J,K)+AA(I,J)*AA(I,K)*CC(I)*LOG(10.)
         ELSE
            ZZ(J,K)=ZZ(J,K)+AA(I,J)*AA(I,K)*CC(I)/XX(K)
         END IF
c modif PhM
400   CONTINUE
405   CONTINUE
C      write(6,*)'chmslv iter zz',iter,'IDX ',IDX(J)
c      write(6,120)(idx(k),ZZ(J,K),K=1,NX)
410   CONTINUE
C
C
C      ----------------- A REVOIR----------
C     IF (IADSORB.GT.0) CALL TRSURZ(SP1,SP2,SURFACE)
C      ----------------------------------------------
 
C CONVERGENCE TEST
 
      DO 800 J=1,NX
      TJ = TOT(J)
      VMAX=ABS(TJ)
      DO 810 I=1,NC
      IF(ABS(AA(I,J)*CC(I)).GE.VMAX) THEN
              VMAX=ABS(AA(I,J)*CC(I))
      ENDIF
810   CONTINUE
      YJ = YY(J)
C                   pour test et impression
C     idximp=idx(j)
C      WRITE(6,*)' CHMSLV J  YJ VMAX EPS ',J,YJ,VMAX,EPS
c       write(6,*) '    J  YJ/VMAX ',J,ABS(YJ)/VMAX
      IF (ABS(YJ)/VMAX.GT.EPS) GOTO 840
800   CONTINUE
 
 
c       write(6,*) 'convergence en ', iter, ' etapes'
c       write(6,120) (idx(j),gx(j),j=1,nx)
 
      RETURN
 
840   CONTINUE
 
       ITER=ITER+1
      IF(ITER.GT.ITMAX) THEN
C     WRITE(6,*)' CHMSLV YJ VMAX IDX TJ',YJ,VMAX,IDXIMP,TJ
C       WRITE(6,*) ' PROBLEME DE CONVERGENCE EN CHIMIE AU NOEUD ',NI
C       WRITE(6,*) ' ITER: ',ITER,'> ITMAX: ',ITMAX
C       WRITE(6,*) ' L ETAT DU PROBLEME EST LE SUIVANT: '
C       WRITE(6,*) 'EXECUTION CHIMIE TERMINEE ** ERREUR 7'
        IEM = 1
        RETURN
      ENDIF
 
C     WRITE(6,*) 'CONVERGENCE TEST OK'
C
C ITERATE
 
C Resolution du systeme lineaire
      CALL CHMSMQ(ZZ,YY,NX,NXDIM,IEM)
C      write(6,*)'chmslv  xx',iter
c       write(6,120)(idx(j),yy(j),j=1,nx)
 
 
 
      IF(IEM.NE.0)RETURN
C
 
C
C Calcul de la concentration a l'etape l
C
      DO 500 J=1,NX
c modif PhM
c Deux moyen de calculer la concentration :
c C^l = 10^{LC^{l-1} - dLC^{l}} ou
c C^l = C^{l-1} * 10^{-dLC^{l}}
c            XX(J) = 10.D0**(GX(J))
         IF (ICALCLOG.EQ.1) THEN
            GX(J) = GX(J)-YY(J)
            XX(J) = XX(J) *  10.D0**(-YY(J))
         ELSE
            XX(J)=XX(J)-YY(J)
            IF (XX(J).LE.0.D0) THEN
               XX(J)=(XX(J)+YY(J))/10.D0
            ENDIF
c modif PhM
c            IF (ABS(XX(J)).LT.1.D-60) XX(J)=1.D-60
            IF (ABS(XX(J)).LT.1.D-100) XX(J)=1.D-100
c modif PhM
C      WRITE(6,510) XX(J)
            GX(J)=LOG10(XX(J))
         END IF
C modif PhM
 
*      write(6,*)'chmslv idx',idx(j),'gx',gx(j)
 
 500     CONTINUE
 510     FORMAT ('XX(J): ', 1PE19.10)
      GO TO 1000
      END