Error termination. Backtrace:

#0 0x7f2c54b2aad0 in ???

#1 0x7f2c54b2b649 in ???

#2 0x7f2c54b2c27f in ???

#3 0x7f2c54d7b1f6 in ???

#4 0x7f2c54d7b4ed in ???

#5 0x55f2418a7fef in ???

#6 0x55f2418a875e in ???

#7 0x7f2c54908d8f in __libc_start_call_main

at ../sysdeps/nptl/libc_start_call_main.h:58

#8 0x7f2c54908e3f in __libc_start_main_impl

at ../csu/libc-start.c:392

#9 0x55f2418a7114 in ???

#10 0xffffffffffffffff in ???

IMPLICIT NONE

INTEGER :: N

PARAMETER (N = 51)

INTEGER :: NM1, NM2, NP1, NP2, INTR, N1

REAL :: DY, DT, CC, TT, AA, PR, SC, TIMEI, STABR, TIMMAX, TIM

REAL, DIMENSION(N) :: U, T, C, Y, A1, B1, C1, D1, XX

INTEGER :: I, J, NINT

OPEN (UNIT=6, FILE='rmigri2024.out', STATUS='NEW')

! Initialization of constants

SC = 0.5

PR = 0.71

TT = 8.0

CC = 2.0

AA = 100.0

TIMEI = 0.5

STABR = 1.0

TIMMAX = 4.0

WRITE (6,10)

10 FORMAT (' N CC TT ')

WRITE (6,20) N, CC, TT, PR

20 FORMAT (I4, 3F8.2)

NM1 = N - 1

DY = 1.0 / REAL(NM1)

DT = STABR * DY * DY

INTR = INT(TIMEI / DT)

NINT = 0

N1 = 0

TIM = 0.0

Y(1) = 0.0

DO I = 1, NM1

Y(I+1) = Y(I) + DY

END DO

DO I = 1, N

T(I) = 0.0

U(I) = 0.0

C(I) = 0.0

END DO

C(1) = 0.5

C(N) = -0.5

U(N) = 0.0

U(1) = 1.0

DO I = 1, NM1/10

N1 = N1 + 1

NINT = NINT + 1

TIM = TIM + DT

CALL FLTEMP(NM1, DY, DT, PR, T, A1, B1, C1, D1, XX)

CALL FLCON(NM1, DY, DT, SC, CC, TT, AA, T, C, A1, B1, C1, D1, XX)

CALL FLVEL(NM1, DY, DT, -0.8, U, T, C, A1, B1, C1, D1, XX)

IF (NINT == INTR) THEN

WRITE (6,30) TIM

30 FORMAT ('Time = ', F8.2)

DO J = 1, N, NM1/10

WRITE (6,40) Y(J), U(J), T(J), C(J)

40 FORMAT (F6.2, 3E20.7)

END DO

CALL HTC(NM1, DY, U, C)

NINT = 0

END IF

IF (TIM >= TIMMAX) EXIT

END DO

CLOSE (UNIT=6)

STOP

END PROGRAM ABCD

SUBROUTINE FLTEMP(NM1, DY, DT, PR, T, A1, B1, C1, D1, XX)

IMPLICIT NONE

INTEGER :: NM1

REAL :: DY, DT, PR, DTYP

REAL, DIMENSION(*) :: T, A1, B1, C1, D1, XX

INTEGER :: I, J

DTYP = DT / (PR * DY * DY)

DO I = 2, NM1

A1(I-1) = DTYP

B1(I-1) = -1.0 - 2.0 * DTYP

C1(I-1) = DTYP

D1(I-1) = -T(I)

END DO

D1(1) = D1(1) - A1(1) * T(1)

D1(NM1) = D1(NM1) - C1(NM1) * T(NM1+1)

CALL THOMAS(NM1, A1, B1, C1, D1, XX)

DO J = 2, NM1

T(J) = XX(J-1)

END DO

RETURN

END SUBROUTINE FLTEMP

SUBROUTINE FLCON(NM1, DY, DT, SC, CC, TT, AA, T, C, A1, B1, C1, D1, XX)

IMPLICIT NONE

INTEGER :: NM1

REAL :: DY, DT, SC, CC, TT, AA, DTYP, CC1

REAL, DIMENSION(*) :: T, C, A1, B1, C1, D1, XX

INTEGER :: I, J

DTYP = DT / (SC * DY * DY)

CC1 = 1.0 / (DY * DY)

DO I = 2, NM1

A1(I-1) = DTYP

B1(I-1) = -1.0 - 2.0 * DTYP

C1(I-1) = DTYP

D1(I-1) = -C(I) - AA * DT * (CC1 * ((C(I) + CC) / (T(I) + TT) * (T(I+1) - 2.0 * T(I) + T(I-1)) &

+ 1.0 / 4.0 * ((C(I+1) - C(I-1)) * (T(I+1) - T(I-1)) / (T(I) + TT)) &

- ((C(I) + CC) * (T(I+1) - T(I-1))**2 / (4.0 * (T(I) + TT)**2))))

END DO

D1(1) = D1(1) - A1(1) * C(1)

D1(NM1) = D1(NM1) - C1(NM1) * C(NM1+1)

CALL THOMAS(NM1, A1, B1, C1, D1, XX)

DO J = 2, NM1

C(J) = XX(J-1)

END DO

RETURN

END SUBROUTINE FLCON

SUBROUTINE FLVEL(NM1, DY, DT, AB, U, T, C, A1, B1, C1, D1, XX)

IMPLICIT NONE

INTEGER :: NM1

REAL :: DY, DT, AB, DTYP

REAL, DIMENSION(*) :: U, T, C, A1, B1, C1, D1, XX

INTEGER :: I, J, N

DTYP = DT / DY**2

DO I = 2, NM1

A1(I-1) = DTYP

B1(I-1) = -1.0 - 2.0 * DTYP

C1(I-1) = DTYP

D1(I-1) = -U(I) - DT * (T(I) + (AB * C(I)))

END DO

D1(1) = D1(1) - A1(1) * U(1)

D1(NM1) = D1(NM1) - C1(NM1) * U(N)

CALL THOMAS(NM1, A1, B1, C1, D1, XX)

DO J = 2, NM1

U(J) = XX(J-1)

END DO

RETURN

END SUBROUTINE FLVEL

SUBROUTINE HTC(NM1, DY, U, C)

IMPLICIT NONE

INTEGER :: NM1

REAL :: DY, DY2, SKIN0, SKIN1, SHK0, SHK1

REAL, DIMENSION(*) :: U, C

DY2 = 2.0 * DY

SKIN0 = (U(2) - U(1)) / DY

SKIN1 = (U(NM1+1) - U(NM1)) / DY

SHK0 = (C(2) - C(1)) / DY

SHK1 = (C(NM1+1) - C(NM1)) / DY

WRITE (6,50) SKIN0, SHK0

50 FORMAT ('Stanton number at y = 0 is', F10.5, ' ', F10.5)

WRITE (6,60) SKIN1, SHK1

60 FORMAT ('Stanton number at y = 1 is', F10.5, ' ', F10.5)

RETURN

END SUBROUTINE HTC

SUBROUTINE THOMAS(NM1, A, B, C, D, X)

IMPLICIT NONE

INTEGER :: NM1, I

REAL, DIMENSION(*) :: A, B, C, D, X

REAL, DIMENSION(NM1) :: CP, DP

CP(1) = C(1) / B(1)

DP(1) = D(1) / B(1)

DO I = 2, NM1

CP(I) = C(I) / (B(I) - A(I) * CP(I-1))

DP(I) = (D(I) - A(I) * DP(I-1)) / (B(I) - A(I) * CP(I-1))

END DO

X(NM1) = DP(NM1)

DO I = NM1-1, 1, -1

X(I) = DP(I) - CP(I) * X(I+1)

END DO

RETURN

END SUBROUTINE THOMAS