Iniciar sesiónCerrar
Correo:
Contraseña:
Entrar
Recordar sesión en este navegador
Iniciar sesiónIniciar sesiónCrear cuentaCrear cuenta
LWP » Foros » Fortran » CUAL ES EL ERRO EN MI BUCLE

Fortran - CUAL ES EL ERRO EN MI BUCLE

Volver
Nuevo Tema
<<>>
 
Vista:

CUAL ES EL ERRO EN MI BUCLE

Publicado por jlito (1 intervención) el 28/07/2014 04:18:50
Estaría muy agradecido si me ayudan, no veo el error en mi bucle de incremento de temperatura, estoy usando el FORCE 2.0

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
EXTERNAL DERIVS
      DIMENSION X0(1:3),X1(1:3),F(1:3),XEND(1:4),XWRK(1:4,1:3),X(1:3)
      CHARACTER*64 filename
      REAL ZETA,WF,WPS,S,T0,TF,PUNTOS,MEK0,INCREMT,T,T1,M1,M2,M3,M4,
     +SUPER,F1,F2,F3,F4,F5,DECIS
      REAL PASOS,H,M
C DEFINO EL VECTOR X CON LOS VALORES DE X0 COMO 0
50    X0(1)=0.
      X0(2)=0.
      X0(3)=0.
      Q=0.
      M=0.
80    PRINT*,'INTRODUCE EL PESO DEL CATALIZADOR (g)'
      READ(*,*)WP
      IF(WP.LE.0.)THEN
      PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
      GO TO 80
      ELSE
      ENDIF
      WRITE(*,*) ' '
85    PRINT*,'INTRODUCE LA VELOCIDAD ESPACIAL (h-1)'
      READ(*,*)S
      IF (S.LE.0.) THEN
      PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
      GO TO 85
      ELSE
      ENDIF
      WRITE(*,*) ' '
90    PRINT*,'INTRODUCE LA CONCENTRACION INICIAL DE MEK (ppm)'
      READ(*,*)MEK0
      IF(MEK0.LE.0.)THEN
      PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
      GOTO 90
      ELSE
      ENDIF
      WRITE(*,*) ' '
92    PRINT*,'INTRODUCE LA TEMPERATURA INICIAL (§C)'
      READ(*,*)TO
      IF(.0.GT.TO)THEN
       PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
      GOTO 92
 
      ENDIF
 
94    PRINT*,'INTRODUCE LA TEMPERATURA FINAL (§C)'
      READ(*,*)TF
      IF((.0.GT.TF).OR.(TO.GT.TF))THEN
       PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
       GOTO 92
       IF(.0.GT.T0)THEN
        PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
        GOTO 94
       ENDIF
 
      ENDIF
      WRITE(*,*) ' '
96    PRINT*,'INTRODUCE EL NUMERO DE PUNTOS DE LA CURVA DE LIGHT-OFF'
      READ(*,*) PUNTOS
      IF (PUNTOS.LE.2.) THEN
      PRINT*,'ME PARECE QUE TE HAS EQUIVOCADO'
      GOTO 96
      ENDIF
      INCREMT=(TF-T0)/PUNTOS
      PRINT*,'INTRODUCE EL NUMERO DE PASOS POR PUNTO QUE QUIERES'
      READ(*,*)PASOS
      H=WP/PASOS
      PRINT*,'LA LONGITUD DEL PASO SERA',H,'g'
      WRITE(*,*) ' '
 
 
97    WRITE (*,*)'TU REACTOR ES MENBRANA O LECHO FIJO'
      WRITE(*,*) ' '
      WRITE (*,*)'1-MENBRANA'
      WRITE (*,*)'2-LECHO FIJO'
      READ(*,*)DECIS
      IF (DECIS.EQ.2.) THEN
      M1=0.344561752
      M2=-0.582216338
      M3=0.469444001
      M4=-0.764382122
      PRINT*,'¨CUAL ES SU SUPERFICIE  ESPECIFICA?'
      READ(*,*) SUPER
      F5=SUPER/5.55
      F1=F5**M1
      F2=F5**M2
      F3=F5**M3
      F4=F5**M4
 
      ELSE
      ENDIF
      IF(DECIS.EQ.1.) THEN
      M1=0.339163953
      M2=-0.687032826
      M3=0.179069907
      M4=-3.59974550
      PRINT*,'¨CUAL  ES SU SUPERFICIE ESPECIFICA?'
      READ(*,*) SUPER
      F5=SUPER/5.55
      F1=F5**M1
      F2=F5**M2
      F3=F5**M3
      F4=F5**M4
 
      ELSE
      ENDIF
      IF((DECIS.NE.1.).AND.(DECIS.NE.2.))THEN
      GOTO 97
 
      ELSE
      END IF
 
      N=3.
 
      WRITE(*,'(A)')'INTRODUCE EL NOMBRE DEL ARCHIVO, POR FAVOR'
      READ (*,'(A)') filename
      OPEN (1, FILE=filename, STATUS='unknown')
      WRITE (1,100)
100   FORMAT (15X,'RUNGE-KUTTA DE CUARTO ORDEN PARA UN SISTEMA DE')
      WRITE (1,110)
110   FORMAT (25X,'ECUACIONES DIFERENCIALES ORDINARIAS',/,5X,74(1H*),//)

      WRITE (1,111)
111   FORMAT (10X,'VELOCIDAD ESPACIAL',5X,'MEK INICIAL',5X,'W FE2O3',/,
     +5X,74(1H-),//)

      WRITE (1,112)S,MEK0,WP
112   FORMAT (15X,F8.2,10X,F8.2,10X,F8.2,/,5X,74(1H-),//)

      WRITE (1,113)
113   FORMAT (10X,'TEMPERATURA (§C)',10X,'MEK',10X,'CO',10X,'CO2',/,
     +5X,74(1H*),//)

      WRITE (1,140)T,(X0(I),I=1,N)
140   FORMAT (15X,F8.2,10X,3(F11.3,5X))
 
C INICIO BUCLE INCREMENTO DE TEMPERATURA
 
      T=T0
 
98    WF=WP
 
      X0(1)=MEK0
      X0(2)=0.
      X0(3)=0.
      Q=S*WP
 
      W0=0.
 
      XEND(1)=0.
      XEND(2)=0.
      XEND(3)=0.
 
      X(1)=MEK0
      X(2)=0.
      X(3)=0.
 
 
 
 
 
C     ESCRIBE LOS RESULTADOS EN PANTALLA
      WRITE(*,170) T,X0(1),X0(2),X0(3)
170   FORMAT (1X,'T=',F10.2,'MEK=',F14.3,'CO=',F14.3,'CO2=',F14.3)
      WRITE (1,150) T,(X0(I),I=1,N)
150   FORMAT(14X,F8.2,7X,3(F11.3,4X))
      IF (T.LT.TF) THEN
      T1=T+INCREMT
      T=T1
          M=M+1
      IF (M.GE.20.) THEN
       M=0.
      ELSE
      ENDIF
      GOTO 98
      ELSE
      ENDIF
 
      WRITE (1,160)
160   FORMAT (5X,74(1H-))
 
      DO K=1,4
        DO J=1,3
		XWRK(K,J)=0.
        END DO
      END DO
 
5     CALL RKSIST (DERIVS,W0,H,X0,XEND,XWRK,F,B,T,Q,F1,F2,F3,F4)
      W0=W0+H
      DO I=1,N
      X0(I)=XEND(I)
      IF (X0(I).LE.1.) THEN
      TEND=T-INCREMT
      PRINT*,'SE HA ALCANZADO COMBUSTION COMPLETA A',TEND,'§C'
      GOTO 200
      ENDIF
      END DO
 
      IF (W0.LT.WF) THEN
      GOTO 5
      ELSE
      ENDIF
 
200   WRITE (1,*) CHAR (12)
      CLOSE (UNIT=1)
      WRITE (*,*)''
      WRITE (*,*)'QUIERES VOLVER A INTENTARLO?'
      WRITE (*,*)''
      WRITE (*,*)'1-SI'
      WRITE (*,*)'2-NO'
      READ (*,*)ZETA
      IF (ZETA.EQ.1) THEN
      GOTO 50
 
      ELSE
      ENDIF
      IF((ZETA.NE.1).AND.(ZETA.NE.2)) THEN
      GOTO 200
 
      ELSE
      ENDIF
 
      STOP
      END
C ========================================================
C   CON ESTE SUBPROGRAMA RESOLVEMOS LAS ECUACIONES
C ========================================================
      SUBROUTINE RKSIST (DERIVS,W0,H,X0,XEND,XWRK,F,N,T,Q,F1,F2,F3,F4)
      DIMENSION X0(1:N),XEND(1:N),XWRK(1:4,1:N),F(1:N)
C     CALCULA LA PRIMERA DERIVADA
      CALL DERIVS (X0,W0,F,N,T,Q,F1,F2,F3,F4)
      DO I=1,N
        XWRK(1,I)=H*F(I)
        XEND(I)=X0(I)+XWRK(1,I)/2
      IF (XEND(I).LT.0.) THEN
           XEND(I)=0
      ENDIF
      END DO
C     SEGUNDA ESTIMADA
      CALL DERIVS (XEND,W0+H/2,F,N,T,Q,F1,F2,F3,F4)
      DO I=1,N
      XWRK(2,I)=H*F(I)
      XEND(I)=X0(I)+XWRK(2,I)/2
      IF (XEND(I).LT.0.) THEN
      XEND(I)=0
      ENDIF
      END DO
C     TERCERA ESTIMADA
      CALL DERIVS (XEND,W0+H/2,F,N,T,Q,F1,F2,F3,F4)
      DO I=1,N
      XWRK(3,I)=H*F(I)
      XEND(I)=X0(I)+XWRK(3,I)
      IF (XEND(I).LT.0.) THEN
      XEND(I)=0
      ENDIF
      END DO
C     CUARTA ESTIMADA
      CALL DERIVS (XEND,W0+H,F,N,T,Q,F1,F2,F3,F4)
      DO I=1,N
        XWRK(4,I)=H*F(I)
      END DO
C     CALCULA EL VALOR DE X AL FINAL DEL INTERVALO
      DO I=1,N
      XEND(I)=X0(I)+(XWRK(1,I)+2*XWRK(2,I)+2*XWRK(3,I)+XWRK(4,I))/6
      IF (XEND(I).LT.0.) THEN
      XEND(I)=0
      ENDIF
      END DO
 
      RETURN
      END
C
C	ESTE PROGRAMA DEFINE LA FUNCION EN TERMINOS F y X C
C
      SUBROUTINE DERIVS (X,W,F,N,T,Q,F1,F2,F3,F4)
      IMPLICIT REAL (K)
      DIMENSION X(1:N),F(1:N)
 
      REAL RI,R2,R3
 
C
C     INTRODUCIMOS LAS DECLARACIONES DE VARIABLES DEL PROBLEMA*
C     Aqu¡ se colocan las constantes del modelo cinético calculadas a la temperatura de     referencia (en este ejemplo aplicadas a Fe2O3)
      R=8.314
 
      KA0=0.00111344075
      EA=299666.9361
      KB0=26.3408816
      EB=1918.31555
      KC0=2249.78113
      EC=-40977.9019
      KC20=1.09651168E-11
      EC2=647663.9403
      K10=0.00137229485
      E1=-36119.5068
      K20=0.000488923772
      E2=36158.3548
 
C CALCULAMOS LAS CONSTANTES  DE REACCiÓN Y DE ADSORCIÓN
 
      KA=KA0*EXP((-EA/R)*((1/(T+273.15))-( 1/473. 15)))
      KB=KB0*EXP((-EB/R)*((1/(T+273.15))-(1/473.15)))
      KC=KC0*EXP((-EC/R)*((1/(T+273.15))-(1/473.15)))
      KC2=KC20*EXP((-EC2/R)*((1/(T+273. 15))-(1/473.15)))
      K1=K10*EXP((-E1/R)*((1/(T+273.15))-(1/473.15)))
      K2=K20*EXP((-E2/R)*((1/(T+273.15))-(1/473.15)))
 
 
 
C
C     INICIO DE ESCRITURA DE LAS ECUACIONES DEL SISTEMA C
C
C     CALCULO DE VELOCIDAD DE MEK A CO2
      R1=F1*(KA*X(1))/((1+K1*X(1)+K2*X(2))**2)
 
C    CALCULO DE VELOCIDAD DE MEK A CO
      R2=F2*(KB*X(1))/((1+K1*X(1)+K2*X(2))**2)
 
C     CALCULO DE VELOCIDAD DE CO A C02
      R3=F3*KC*X(2)/((I+K1*X(1)+K2*X(2))**2)-F4*KC2*X(3)
 
C    CALCULO DE dMEK/dW
      F(1)=(-R1-R2)/(Q)
 
C	CALCULO DE dCO/dW F(2)=«(R2*4.)-RJ)/(Q)
      F(2)=((R2*4.)-R3)/(Q)
 
C	CALCULO DE dC02/dW F(3)=(RI*4.+RJ)/(Q)
 
      F(3)=(R1*4.+R3)/(Q)
 
      RETURN
      END
Valora esta pregunta
Me gusta: Está pregunta es útil y esta claraNo me gusta: Está pregunta no esta clara o no es útil
0
Responder

CUAL ES EL ERRO EN MI BUCLE

Publicado por Sebastian (10 intervenciones) el 28/07/2014 08:23:25
Hola jlito

Quizás si das más detalles acerca del error, ¿cuál es el mensaje al compilar? etc, sería más fácil ayudarte.
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar

CUAL ES EL ERRO EN MI BUCLE

Publicado por Sebastian (10 intervenciones) el 30/07/2014 07:43:08
En la sentencia

1
5 CALL RKSIST (DERIVS,W0,H,X0,XEND,XWRK,F,B,T,Q,F1,F2,F3,F4)

parece ser que la variable B no es del tipo adecuado. al sustituir B --> N el programa corre perfectamente.

Saludos
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar

CUAL ES EL ERRO EN MI BUCLE

Publicado por Capitan Kirk (19 intervenciones) el 25/08/2014 08:12:58
Profundizando en lo que te ha comentado Sebastian, ocurre que en la llamada a RKSIST estás pasando una variable B, de la que no he visto declaración ni inicialización alguna. En la rutina RKSIST estás utilizando esa variable (que corresponde al parámetro N) como dimensión para una serie de matrices. Tenemos dos posibles causas del error:

1. Si mal no recuerdo, salvo que haya una declaración explícita, las variables cuyos nombres comienzan por I, J, K, L, M, N son, por defecto, de tipo entero, y las que empiezan por cualquier otra letra son de tipo real. En tu caso, en la llamada a RKSIST estarías pasando un valor de tipo real en el lugar que corresponde a un valor de tipo entero.

2. La variable B no tiene declaración ni inicialización alguna, por tanto cuando la pasas como parámetro en la llamada a RKSIST es probable que le esté asignando a B un valor de cero o un valor no válido para el uso que se hace de esa variable dentro de RKSIST.

Y, como se te ha comentado, da más detalles del error; esencialmente, si ha sido un error de compilación o de ejecución, y el mensaje en concreto.

Saludos,
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar