1 : C CHMSLV SOURCE CHAT 05/01/12 22:00:06 5004 2 : SUBROUTINE CHMSLV(IDSCHI,SP2,ITMAX,EPS,ICALCLOG,IEM) 3 : C======================================================================== 4 : C 5 : C SP ISSU DE TRIOEF ( TRSOLV) 6 : C 7 : C Modif PhM 8 : C a) Prise en charge du calcul en log de concentration 9 : c (argument ICALCLOG) 10 : C 11 : C======================================================================== 12 : IMPLICIT INTEGER(I-N) 13 : IMPLICIT REAL*8(A-H,O-Z) 14 : SEGMENT IDSCHI 15 : REAL*8 GK(NYDIM),AA(NYDIM,NXDIM),FF(NZDIM,NPDIM) 16 : INTEGER IDX(NXDIM),IDY(NYDIM),IDZ(NZDIM),IDP(NPDIM),NN(6) 17 : INTEGER IDECY(NYDIM),IONZ(NXDIM) 18 : CHARACTER*32 NAME(NXDIM),NAMESP(NYDIM) 19 : ENDSEGMENT 20 : SEGMENT SP2 21 : REAL*8 GX(NXDIM),XX(NXDIM),GS(NZDIM),SS(NZDIM) 22 : REAL*8 TOT(NXDIM),TOTAQ(NXDIM),TOTFIX(NXDIM),GKS(NZDIM) 23 : REAL*8 YY(NXDIM),ZZ(NXDIM,NXDIM),CC(NYDIM),GC(NYDIM) 24 : ENDSEGMENT 25 : NXDIM=IDX(/1) 26 : NYDIM=IDY(/1) 27 : NZDIM=IDZ(/1) 28 : NPDIM=IDP(/1) 29 : ITER=0 30 : 31 : 32 : c write(6,*) 'ICALCLOG',ICALCLOG 33 : 34 : 35 : C 36 : ************************************************************************* 37 : C PRISE EN COMPTE DES SOLUTIONS SOLIDES : METHODE 1 38 : N4S=0 39 : * IF(NZDIM.NE.0)THEN 40 : * I3S=NN(1)+NN(2)+NN(3)+1 41 : * I4S=NN(1)+NN(2)+NN(3)+NN(4) 42 : * DO 13 I7S=I3S,I4S 43 : * IDY7=IDY(I7S) 43 : * CALL CHIADY(IDZ,NZDIM,IDY7,ID7) 44 : * IF(ID7.NE.0)THEN 45 : * N4S=N4S+1 45 : * CALL CHMREX(IDSCHI,LGKMOD,LGKTMP,IDY7,4,5) 45 : * CALL CHMREX(IDSCHI,LGKMOD,LGKTMP,IDY7,5,4) 46 : * ENDIF 47 : * 13 CONTINUE 48 : * ENDIF 49 : ************************************************************************** 50 : 51 : NC=NN(1)+NN(2) 52 : NX=NXDIM-NN(3)-NN(4)+N4S 53 : C write(6,*)'chmslv nc nx gk',nc,nx 54 : C write(6,120)(idy(j),gk(j),j=1,nc) 55 : DO 1 J=1,NX 56 : YY(J)=0.D0 57 : 1 CONTINUE 58 : 59 : 60 : 61 : 1000 CONTINUE 62 : 120 FORMAT(6(1X,I6,1PD15.6)) 63 : C COMPLEXES 64 : DO 2 I=1,NC 65 : V=GK(I) 66 : * write(6,*)'chmslv idy',idy(i),'gk',gk(i) 67 : DO 3 J=1,NX 68 : V=V+AA(I,J)*GX(J) 69 : * write(6,*)'chmslv idx',idx(j) 70 : * write(6,*)'chmslv aa',aa(i,j),'gx',gx(j) 71 : 3 CONTINUE 72 : GC(I)=V 73 : C 74 : CC(I)=10.D0**GC(I) 75 : * write(6,*)'chmslv idy',IDY(I),'c',CC(I),'logc',GC(I) 76 : 2 CONTINUE 77 : 78 : 79 : C MOLE BALANCE 80 : 81 : C ----------------- A REVOIR---------- 81 : C IF(IADSORB.GT.0) CALL TRSURT(SP1,SP2,SURFACE) 82 : C -------------------------------------- 83 : DO 201 J=1,NX 84 : V=-TOT(J) 85 : DO 200 I=1,NC 86 : V=V+AA(I,J)*CC(I) 87 : 200 CONTINUE 88 : YY(J)=V 89 : C 90 : 201 CONTINUE 91 : 92 : C CALCUL Z (Jacobienne pour les iterations de newton) 93 : DO 301 I=1,NX 94 : DO 300 J=1,NX 95 : ZZ(I,J)=0.D0 96 : 300 CONTINUE 97 : 301 CONTINUE 98 : DO 410 J=1,NX 99 : DO 405 I=1,NC 100 : DO 400 K=1,NX 101 : c modif PhM 102 : IF (ICALCLOG.EQ.1) THEN 103 : ZZ(J,K)=ZZ(J,K)+AA(I,J)*AA(I,K)*CC(I)*LOG(10.) 104 : ELSE 105 : ZZ(J,K)=ZZ(J,K)+AA(I,J)*AA(I,K)*CC(I)/XX(K) 106 : END IF 107 : c modif PhM 108 : 400 CONTINUE 109 : 405 CONTINUE 110 : C write(6,*)'chmslv iter zz',iter,'IDX ',IDX(J) 111 : c write(6,120)(idx(k),ZZ(J,K),K=1,NX) 112 : 410 CONTINUE 113 : C 114 : C 115 : C ----------------- A REVOIR---------- 115 : C IF (IADSORB.GT.0) CALL TRSURZ(SP1,SP2,SURFACE) 116 : C ---------------------------------------------- 117 : 118 : C CONVERGENCE TEST 119 : 120 : DO 800 J=1,NX 121 : TJ = TOT(J) 122 : VMAX=ABS(TJ) 123 : DO 810 I=1,NC 124 : IF(ABS(AA(I,J)*CC(I)).GE.VMAX) THEN 125 : VMAX=ABS(AA(I,J)*CC(I)) 126 : ENDIF 127 : 810 CONTINUE 128 : YJ = YY(J) 129 : C pour test et impression 130 : C idximp=idx(j) 131 : C WRITE(6,*)' CHMSLV J YJ VMAX EPS ',J,YJ,VMAX,EPS 132 : c write(6,*) ' J YJ/VMAX ',J,ABS(YJ)/VMAX 133 : IF (ABS(YJ)/VMAX.GT.EPS) GOTO 840 134 : 800 CONTINUE 135 : 136 : 137 : c write(6,*) 'convergence en ', iter, ' etapes' 138 : c write(6,120) (idx(j),gx(j),j=1,nx) 139 : 140 : RETURN 141 : 142 : 840 CONTINUE 143 : 144 : ITER=ITER+1 145 : IF(ITER.GT.ITMAX) THEN 146 : C WRITE(6,*)' CHMSLV YJ VMAX IDX TJ',YJ,VMAX,IDXIMP,TJ 147 : C WRITE(6,*) ' PROBLEME DE CONVERGENCE EN CHIMIE AU NOEUD ',NI 148 : C WRITE(6,*) ' ITER: ',ITER,'> ITMAX: ',ITMAX 149 : C WRITE(6,*) ' L ETAT DU PROBLEME EST LE SUIVANT: ' 150 : C WRITE(6,*) 'EXECUTION CHIMIE TERMINEE ** ERREUR 7' 151 : IEM = 1 152 : RETURN 153 : ENDIF 154 : 155 : C WRITE(6,*) 'CONVERGENCE TEST OK' 156 : C 157 : C ITERATE 158 : 159 : C Resolution du systeme lineaire 159 : CALL CHMSMQ(ZZ,YY,NX,NXDIM,IEM) 160 : C write(6,*)'chmslv xx',iter 161 : c write(6,120)(idx(j),yy(j),j=1,nx) 162 : 163 : 164 : 165 : IF(IEM.NE.0)RETURN 166 : C 167 : 168 : C 169 : C Calcul de la concentration a l'etape l 170 : C 171 : DO 500 J=1,NX 172 : c modif PhM 173 : c Deux moyen de calculer la concentration : 174 : c C^l = 10^{LC^{l-1} - dLC^{l}} ou 175 : c C^l = C^{l-1} * 10^{-dLC^{l}} 176 : c XX(J) = 10.D0**(GX(J)) 177 : IF (ICALCLOG.EQ.1) THEN 178 : GX(J) = GX(J)-YY(J) 179 : XX(J) = XX(J) * 10.D0**(-YY(J)) 180 : ELSE 181 : XX(J)=XX(J)-YY(J) 182 : IF (XX(J).LE.0.D0) THEN 183 : XX(J)=(XX(J)+YY(J))/10.D0 184 : ENDIF 185 : c modif PhM 186 : c IF (ABS(XX(J)).LT.1.D-60) XX(J)=1.D-60 187 : IF (ABS(XX(J)).LT.1.D-100) XX(J)=1.D-100 188 : c modif PhM 189 : C WRITE(6,510) XX(J) 190 : GX(J)=LOG10(XX(J)) 191 : END IF 192 : C modif PhM 193 : 194 : * write(6,*)'chmslv idx',idx(j),'gx',gx(j) 195 : 196 : 500 CONTINUE 197 : 510 FORMAT ('XX(J): ', 1PE19.10) 198 : GO TO 1000 199 : END 200 : 201 : 202 : 203 : 204 : 205 : 206 : 207 : 208 : 209 : 210 : 211 : 212 :
© Cast3M 2003 - All rights reserved.
Disclaimer