Berkeley Lab
Bringing Science Solutions to the World

TOUGH Bugs & Fixes

The TOUGH codes are continually updated in response to new scientific challenges and user needs. While the TOUGH codes are thoroughly tested before release, exhaustive testing is not possible. Consequently, developers and users occasionally detect programming bugs in a released code version. These bugs are fixed as they become known. However, such bug fixes cannot immediately be distributed to all licensed users due to the intricate review, approval, licensing, and distribution process at Berkeley Lab.

This page contains a living list of known bugs. Certain bugs may be ignored, as they only affect (often rarely used) program options or combination of features; others may be fixed with only a few code changes; yet others may require more subtantial recoding.

Please consult the list to see whether the identified issues affect your work. If so, please follow the instructions or request an updated version of the code. If you suspect there is a bug in TOUGH, please contact the developers and provide them with a detailed description of the error as well as the conditions and computing environment under which the error occurs.

Please also consult the following pages:

#
Date		 Title/
Description		 Affected
Module/File		 Proposed Fix

14

06/21

 Phase condition bug in EWASG TOUGH3/eoswasg.f
Action: search for “(Palliser & McKibbin, 1998)”

Go back 10 lines replace the section from “CALL HALITE(TX,XEQ,XeqG)” to right before calling SATB subroutine by the following (changes indicated by red):

CALL HALITE(TX,XEQ,XeqG)
IF (XSM.GT.1.5D0) THEN
SS=XSM-10.D0
XS=XEQ
XSL=XEQ
ELSE
SS=0.D0
XSL=XSM
XS=XSM
END IF
IF(NK.EQ.3) THEN
CALL HENRY(TX,XX(3),PA,XSL,RKH)
IF (IGOOD.NE.0) GO TO 200
ELSE
PA=0.D0
END IF
click to enlarge
13
06/21		 GOFT printout “TOUGH3/Input_Output.f90”

There are two places in the file for this fix.
First place (between line 1450~1550):
Changed text:

!      IELL=LEN_TRIM(LINE(:14))-5  ! yq, replaced by the next two lines. Works for element names 5-8 characters
IELL=5
IF (MOP2(2).NE.0) IELL=MOP2(2)
BACKSPACE(NUINPUT)
!      CFI3='(Ax,A5,yX,I5)’  ! yq remove S_GOFT, to be backward compatible. replaced by the next line
CFI3='(Ax,yX,5X,I5)’
WRITE(CFI3(3:3),'(I1)’)       IELL
!      WRITE(CFI3(8:8),'(I1)’)       10-IELL      ! yq replaced by the next line
WRITE(CFI3(5:5),'(I1)’)       10-IELL
do 181 i=1,IGOFT
!         READ(NUINPUT,CFI3) egoft(i),S_GOFT(I),IGOFTF(I)   ! yq for backward compatible, replaced by next line

Reference figure to help locate the place:

click to enlarge
Second place need to be changed (line 3700~3750)

find the blocks in local numbering for the time dependent print out (GOFT)

N0=0
ALLOCATE(NGOFT(IGOFT),STAT=iI)
do  iI=1,igoft
do n=1,NOGN
IF(ELEG(N).EQ.egoft(iI)) then
!                IF(TRIM(S_GOFT(iI)).EQ.”.OR.S_GOFT(iI)(:5).EQ.SOURCE(N)) THEN   ! yq, S_GOFT is never used anywhere
N0=N0+1
ngoft(N0)=n
igoftf(N0)=igoftf(iI)
GOTO 7883
!                ENDIF                      ! yq, comment out

Reference figure to help locate the place:

click to enlarge
12
06/21		 Handling PSAT greater than PL during iteration “TOUGH3/thrmb.f90” around line 309.

The solution is to remove the IF THEN clause as indicated below.

In file ” thrmb.f90″ around line 309, remove the IF THEN clause as indicated below

! IF (PL.GT.PSAT) THEN
CALL COWAT(TX,PSATW+DP,DW,UW)
IF (IGOOD.NE.0) RETURN
! ENDIF
DW0=DW
11

01/21

 COFT output when ICOFTF = -1

Bug in printing flow for phases and components using COFT block

 TOUGH3/Input_Output.f90

In Input_Output.f90, subroutine fgtab, near the end of the section
“tabulate CONNE data”

Replace:
if(flo(nnp+jj).ge.0.) then
xfm(nk*(jj-1)+j)=flo(nnp+jj)*xm2(nk*(j-1)+j)
else
xfm(nk*(jj-1)+j)=flo(nnp+jj)*xm1(nk*(j-1)+j)
endif

By:
if(flo(nnp+jj).ge.0.) then
xfm(nk*(jj-1)+j)=flo(nnp+jj)*xm2(nk*(jj-1)+j)
else
xfm(nk*(jj-1)+j)=flo(nnp+jj)*xm1(nk*(jj-1)+j)
endif
10

07/20

 Windows 10 and Ubuntu Terminal Ubuntu Terminal

If your operating system is Windows 10, a better option than Cygwin is Ubuntu Terminal. Updated instruction for TOUGH3 installation on Ubuntu Terminal can be found here (Installation File, READ ME File v2)

(previous version/superseded July 2020: READ ME File v1)
9

06/20

 Element-by-element heterogeneity

Bug in identifying connection for permeability modification when MOP2(20)>0.

 TOUGH3, MULTI.f90 (Option 1) In Multi.f90 around line 590, make the following change.

Replace:
do i = 1,2
if (iran .ne. 0) then
perl(i) = perl(i)*pm(nexidx(i))
end if
! Assign anisotropic permeability for each grid block
if (mop2(20) .gt. 0) then
USRX=USERX(ISO,N)
CALL HETERO(nexidx(i),nmat(i),nloc2(i),iso,USRX,perl(i),ph(i))
endif
! ysw for anisotropic fracrtures 6/30/95
if(ifaniso(nmat(i)).eq.1) perl(i)=perl(i)+PERF(iso,nmat(i))
end do

By:
do i = 1,2
! Assign anisotropic permeability for each grid block
! change from conne. index N to connected elem. index nexidx
USRX=USERX(ISO,nexidx(i))
CALL HETERO(nexidx(i),nmat(i),nloc2(i),iso,USRX,perl(i),ph(i))
! ysw for anisotropic fractures 6/30/95
if(ifaniso(nmat(i)).eq.1) perl(i)=perl(i)+PERF(iso,nmat(i))
end do

(Option 2) Request& the updated version of the code
8

11/17

 CO2 solubility model in ECO2N V2.0

In TOUGH2/ECO2N V2.0 there is a bug wherein variables TEMP1 and TEMP2 are not defined. These variables control the transition from a low-temperature CO2 solubility model to a high-temperature model. When they are not defined, on most computers they will be assumed to be zero, which implies that the high-temperature solubility model will be used for all temperatures, which will lead to small inaccuracies for low-temperature systems (T below 99 deg C). Also, the option IE(16)=1 (invoke the ECO2N V1 CO2 solubility model) will not work.

 TOUGH2 V2.1, CO2Proper_new.f (Option 1) In the section “Block Data Thermodata” add the line

DATA TEMP1,TEMP2/99.0D0,109.0D0/

(Option 2) Request& the updated version of the executable
7

10/15

 Phase change in ECO2N V2.0

Bug in specifying DELX when phase change to single-phase gas. Note this bug may only affect the convergence behavior in cases where no NaCl is present and a grid blcok dries our completely.

 TOUGH2 V2.1, eco2n_v20.f (Option 1) Make the following changes to SUBROUTINE EOS, file eco2n.f, Lines 1551 and 1553:

case (4)
if(SX > 1.0-dfac-SS) then
DELX(NLOC+3)=-DFAC
else
DELX(NLOC+3)=DFAC
endif
SX=SX+DELX(NLOC+3)
end select

(Option 2) Download the updated version of eco2.f and recompile

(Option 3) Request& the updated version of the executable

6

05/15

 Hysteresis Restart

When doing a restart using the hysteresis module with MOP(13)=2, history parameters for hysteresis (read in the INCON file) are not assigned correctly if the keyword START is present in the input file

 TOUGH2 V2.1, t2fm.f (only if IRP=ICP=12)

iTOUGH2/hysteresis module, t2f.f (only if IRP=ICP=12)

 (Option 1) Make the following changes to SUBROUTINE RFILE, file t2f.f (for iTOUGH2) of t2fm.f (for TOUGH2 V2.1):

DO 2005 J=1,NEL
...unchanged material not shown...

CFiHysteresisB
IF (IHYSTER.LT.0) THEN
SOLDH(J)=SLD
HOLD(J)=HLD
ICURV(J)=ICRV
SOROLD(J)=SORE
STRN1OLD(J)=STR1
STRN2OLD(J)=STR2
STRN3OLD(J)=STR3
HTRN1OLD(J)=HTR1
HTRN2OLD(J)=HTR2
HTRN3OLD(J)=HTR3
SCMP1OLD(J)=SCM1
SCMP2OLD(J)=SCM2
BB2OLD(J)=BBB2
HMINOLD(J)=HMIN2
SHMINOLD(J)=SHMIN2
HMINPOLD(J)=HMINP2
ENDIF
CFiHysteresisE
GOTO 2002
2005 CONTINUE

(Option 2) Download the updated version of t2fm.f and recompile (for TOUGH2 V2.1 only)

(Option 3) Request& the updated version of the executable

5

06/15

 AZTEC Option Broadcast

AZTEC options were not being broadcast to all ranks

 TOUGH2.01-MP

File Paral_Subs.f90

 Download file Paral_Subs.f90 and recompile

4

06/15

 MOLES/SEC vs. KG/S

The header preceding the printout of diffusive fluxes shows units of MOLES/SEC. This is only correct for EOS modules that are based on molar fractions (e.g., TMVOC); for modules that are based on mass fractions, the correct units are KG/S. This is a reporting error; the calculations are correct.

 TOUGH2 V2.1

All EOS modules (except TMVOC)

File t2fm.f

Subroutine OUTDF

 (Option 1) Since this is a reporting error rather than a computational issue, it may be sufficient for the user to be aware of the correct units

(Option 2) Make following changes to SUBROUTINE OUTDF, file t2fm.f:

C-----PRINT A SHORT HEADER.
C
IF(EOSN(1).EQ.'TMVOC1') THEN
PRINT 5000,H1,TITLE,KCYC,ITER,SUMTIM,'MOLES'
ELSE
PRINT 5000,H1,TITLE,KCYC,ITER,SUMTIM,'KG'
ENDIF
C

...

5000 FORMAT(A1/10X,A80//80X,'KCYC =',I5,' - ITER =',I5,
& ' - TIME =',E12.6/
&' MASS FLOW RATES (',A,'/S) FROM DIFFUSION'/)

(Option 3) Download the updated version of t2fm.f and recompile (for TOUGH2 V2.1 only)

(Option 4) Request& the updated version of the executable

3

06/15

 T2WELL-ECO2N Heat Flow Calculation

(1) Heat flow between the wellbore and surrounding formation calculated using Ramey’s solution was not correctly scaled into the energy balance equation.

(2) Accounting for potential energy required all wellbore elements to be connected monotonically downward, leading to errors if U-shaped boreholes are defined.

(3) See Notes for additional changes and capability enhancements.

 T2WELL-ECO2N Request& the updated version of T2WELL-ECO2N

2

02/15

 TOUGH2-MP

Minor clean-up, bug fixes, and updates; see Notes; implemented in TOUGH2.01-MP (includes Fix #1)

 TOUGH2-MP TOUGH2-MP has been replaced by TOUGH3 and is no longer supported.

1

02/15

 AZTEC Memory Operation

A faulty memory operation in the AZTEC library caused premature run termination

 AZTEC (external library to TOUGH2-MP)

File az_zort.c

 (Option 1) In the Aztec2.1 source, edit file az_sort.c (Aztec2.1/lib/az_sort.c).

Change line 1433:
memcpy(ptr1, ptr2, thelength);

to:
memmove(ptr1, ptr2, thelength);

(Option 2) Download the updated file az_sort.c and recompile