@@@@@ @@ @ @ @@@ @ @ @@ @@@ @ @ @ @ @ @ @@ @@@@@ @ @@ @ @ @@@ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @@ @@ @ @ @ @ @ @ @ @ @ @@ @ @ @ @ @ @@ @ @ @ @ @ @ @ @@ @@@@ @ @@ @ @ @ @ @ @ @ @@@@ @ @ @ @ @ @ @ @@@ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @@ @ @ @ @ @ @@ @ @@ @@ @@@ @ @ @@ @@@ @ @ @ @@ @@@@ @ @ @ @ @@ @ @ @ @ @@ @ @ TOUGH3 IS A PROGRAM FOR MULTIPHASE MULTICOMPONENT FLOW IN PERMEABLE MEDIA, INCLUDING HEAT FLOW. IT IS A MEMBER OF THE MULKOM FAMILY OF CODES, DEVELOPED AT LAWRENCE BERKELEY LABORATORY BY KARSTEN PRUESS ET AL. ******************************************************************************** ********************** TOUGH3 VER. 1.0.0 (January 2018) ********************* ******************************************************************************** SUMMARY OF DISK FILES FILE *MESH* EXISTS --- OPEN AS AN OLD FILE WARNING: FILES *MESHA* AND *MESHB* EXIST --- information in file MESH or ELEME and CONNE blocks in input file will be ignored FILE *INCON* EXISTS --- OPEN AS AN OLD FILE FILE *GENER* EXISTS --- OPEN AS AN OLD FILE FILE *SAVE* EXISTS --- OPEN AS AN OLD FILE FILE *TABLE* EXISTS --- OPEN AS AN OLD FILE =================================================================================================================================== PROBLEM TITLE: Lorain_CCS (Reading MULTI block) (Reading INDEX block) (Reading START block) HAVE READ UNKNOWN BLOCK LABEL " " --- IGNORE THIS, AND CONTINUE READING INPUT DATA HAVE READ UNKNOWN BLOCK LABEL " " --- IGNORE THIS, AND CONTINUE READING INPUT DATA (Reading PARAM block) (Reading ROCKS block) DOMAIN NO. 1 MATERIAL NAME -- ROME DOMAIN NO. 2 MATERIAL NAME -- MT SI (Reading SELEC block) (Reading TIMES block) HAVE READ UNKNOWN BLOCK LABEL " " --- IGNORE THIS, AND CONTINUE READING INPUT DATA (Reading GENER block) WRITE FILE *GENER* FROM INPUT DATA (Reading ELEME block) WRITE FILE *MESH* FROM INPUT DATA (Reading CONNE block) (Reading INCON block) WRITE FILE *INCON* FROM INPUT DATA (Reading ENDCY block) ************************************************************************************ * EVALUATE FLOATING POINT ARITHMETIC * ************************************************************************************ * * * FLOATING POINT PROCESSOR HAS APPROXIMATELY 15 SIGNIFICANT DIGITS * * * * DEFAULT VALUE OF INCREMENT FACTOR FOR NUMERICAL DERIVATIVES IS DFAC = 0.1490E-07 * * DEFAULT VALUE FOR DFAC WILL BE USED * * * ************************************************************************************ PARAMETERS FOR FLEXIBLE DIMENSIONING OF MAJOR ARRAYS (MAIN PROGRAM) ARE AS FOLLOWS MNEL = 200000 MNCON = 190000 MNEQ = 4 MNK = 3 MNPH = 4 MNB = 6 MNOGN = 10 MGTAB = 0 =================================================================================================================================== MAXIMUM NUMBER OF VOLUME ELEMENTS (GRID BLOCKS): MNEL = 200000 MAXIMUM NUMBER OF CONNECTIONS (INTERFACES): MNCON = 190000 MAXIMUM LENGTH OF PRIMARY VARIABLE ARRAYS: MPRIM = 800000 MAXIMUM NUMBER OF GENERATION ITEMS (SINKS/SOURCES): MNOGN = 10 MAXIMUM NUMBER OF TABULAR (TIME-DEPENDENT) GENERATION DATA: MGTAB = 0 LENGTH OF SECONDARY PARAMETER ARRAY: MNSEC = 38000000 MESH HAS 200000 ELEMENTS AND 190000 CONNECTIONS (INTERFACES) BETWEEN THEM GENER HAS 10 SINKS/SOURCES END OF PART 1 INPUT JOB --- ELAPSED TIME = 0.0000 SECONDS *********************************************************************************************************************************** * ARRAY *MOP* ALLOWS TO GENERATE MORE PRINTOUT IN VARIOUS SUBROUTINES, AND TO MAKE SOME CALCULATIONAL CHOICES. * *********************************************************************************************************************************** MOP(1) = 0 *** ALLOWS TO GENERATE A SHORT PRINTOUT FOR EACH NEWTON-RAPHSON ITERATION = 0, 1, OR 2: GENERATE 0, 1, OR 2 LINES OF PRINTOUT MORE PRINTOUT IS GENERATED FOR MOP(I) > 0 IN THE FOLLOWING SUBROUTINES (THE LARGER MOP IS, THE MORE WILL BE PRINTED). MOP(2) = 0 *** CYCIT MOP(3) = 0 *** MULTI MOP(4) = 1 *** QU MOP(5) = 0 *** EOS MOP(6) = 0 *** LINEQ MOP(7) = 0 *** IF UNEQUAL ZERO, WILL GENERATE A PRINTOUT OF INPUT DATA CALCULATIONAL CHOICES OFFERED BY MOP ARE AS FOLLOWS: MOP(8) = 0 *** IF ISOT IS NEGATIVE, CHOOSES OPTION FOR REDUCING FRACTURE-MATRIX INTERFACE AREA. MOP(9) = 1 *** CHOOSES FLUID COMPOSITION ON WITHDRAWAL (PRODUCTION). = 0: ACCORDING TO RELATIVE MOBILITIES. = 1: ACCORDING TO COMPOSITION IN PRODUCING ELEMENT. MOP(10) = 0 *** CHOOSES INTERPOLATION FORMULA FOR DEPENDENCE OF THERMAL CONDUCTIVITY ON LIQUID SATURATION (SL). = 0: K = KDRY + SQRT(SL)*(KWET-KDRY) = 1: K = KDRY + SL*(KWET-KDRY) = 2: K = C0+C1*T+C2*Sw+C3*POR MOP(11) = 0 *** CHOOSES EVALUATION OF MOBILITY AND ABSOLUTE PERMEABILITY AT INTERFACES. = 0: MOBILITIES ARE UPSTREAM WEIGHTED WITH WUP. (DEFAULT IS WUP = 1.0). PERMEABILITY IS UPSTREAM WEIGHTED. = 1: MOBILITIES ARE AVERAGED BETWEEN ADJACENT ELEMENTS. PERMEABILITY IS UPSTREAM WEIGHTED. = 2: MOBILITIES ARE UPSTREAM WEIGHTED WITH WUP. (DEFAULT IS WUP = 1.0). PERMEABILITY IS HARMONIC WEIGHTED. = 3: MOBILITIES ARE AVERAGED BETWEEN ADJACENT ELEMENTS. PERMEABILITY IS HARMONIC WEIGHTED. = 4: MOBILITY * PERMEABILITY PRODUCT IS HARMONIC WEIGHTED. MOP(12) = 0 *** CHOOSES PROCEDURE FOR INTERPOLATING GENERATION RATES FROM A TIME TABLE. = 0: TRIPLE LINEAR INTERPOLATION. = 1: "STEP FUNCTION" OPTION. = 2: RIGOROUS STEP RATE OPTION. MOP(13) = 0 *** DEFINES CONTENT OF INCON AND SAVE FILE. = 0: STANDARD CONTENT. = 2: READS PARAMETERS OF HYSTERESIS MODEL FROM FILE INCON. MOP(15) = 0 *** ALLOWS TO SELECT A SEMI-ANALYTICAL HEAT EXCHANGE CALCULATION WITH CONFINING BEDS. = 0: NO SEMI-ANALYTICAL HEAT EXCHANGE. = 1: SEMI-ANALYTICAL LINEAR HEAT EXCHANGE ENGAGED. INITIAL TEMPERATURE OF CONFINING BEDS IS UNIFORM. = 2: SEMI-ANALYTICAL LINEAR HEAT EXCHANGE ENGAGED. INITIAL TEMPERATURE OF CONFINING LAYERS IS NOT UNIFORM. = 5: SEMI-ANALYTICAL RADIAL HEAT EXCHANGE WITH PROPERTIES GIVEN IN MATERIAL QLOSS = 6: SEMI-ANALYTICAL RADIAL HEAT EXCHANGE WITH DEPTH-DEPENDENT PROPERTIES (DEPTH, RADIUS, TEMPERATURE, CONDUCTIVITY, DENSITY, CAPACITY) PROVIDED ON FILE radqloss.dat MOP(16) = 3 *** PERMITS TO CHOOSE TIME STEP SELECTION OPTION = 0: AUTOMATIC TIME STEPPING BASED ON MAXIMUM CHANGE IN SATURATION. = 1: AUTOMATIC TIME STEPPING BASED ON NUMBER OF ITERATIONS NEEDED FOR CONVERGENCE. > 1: INCREASE TIME STEP BY AT LEAST A FACTOR 2, IF CONVERGENCE OCCURS IN .LE. MOP(16) ITERATIONS. MOP(17) = 0 *** HANDLES TIME STEPPING AFTER LINEAR EQUATION SOLVER FAILURE. = 0: NO TIME STEP REDUCTION DESPITE LINEAR EQUATION SOLUTION FAILURE. = 9: REDUCE TIME STEP AFTER LINEAR EQUATION SOLUTION FAILURE. MOP(18) = 1 *** ALLOWS TO SELECT HANDLING OF INTERFACE DENSITY. = 0: PERFORM UPSTREAM WEIGHTING FOR INTERFACE DENSITY. > 0: COMPUTE INTERFACE DENSITY AS AVERAGE OF THE TWO GRID BLOCK DENSITIES. HOWEVER, WHEN ONE OF THE TWO PHASE SATURATIONS IS ZERO, DO UPSTREAM WEIGHTING. MOP(21) = 3 *** PERMITS TO SELECT LINEAR EQUATION SOLVER FROM PACKAGE < 2: DEFAULTS TO MOP(21) = 3 = 2: DSLUBC: BI-CONJUGATE GRADIENT SOLVER; PRECONDITIONER: INCOMPLETE LU FACTORIZATION = 3: DSLUCS: BI-CONJUGATE GRADIENT SOLVER - LANCZOS TYPE; PRECONDITIONER: INCOMPLETE LU FACTORIZATION = 4: DSLUGM: GENERALIZED MINIMUM RESIDUAL CONJUGATE GRADIENTS; PRECONDITIONER: INCOMPLETE LU FACTORIZATION = 5: DLUSTB: STABILIZED BI-CONJUGATE GRADIENT SOLVER; PRECONDITIONER: INCOMPLETE LU FACTORIZATION = 6: LUBAND: DIRECT SOLVER USING LU DECOMPOSITION = 7: AZTEC: PARALLEL ITERATIVE SOLVER = 8: PETSc: PARALLEL ITERATIVE SOLVER MOP(24) = 0 *** PERMITS TO SELECT HANDLING OF MULTIPHASE DIFFUSIVE FLUXES AT INTERFACES = 0: HARMONIC WEIGHTING OF FULLY-COUPLED EFFECTIVE MULTIPHASE DIFFUSIVITY = 1: SEPARATE HARMONIC WEIGHTING FOR EACH PHASE FLUX *********************************************************************************************************************************** *********************************************************************************************************************************** &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& Summary of capabilities for random permeability modification &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& Modification of absolute permeability on a grid block-by-grid block basis will be made when a domain "SEED " is present in data block "ROCKS", as follows. k ---> k' = k*m Here, k is the absolute permeability specified for the reservoir domain to which the grid block belongs. Parameter m is a "permeability modifier" which can be internally generated or externally prescribed by the user on a block-by-block basis. When permeability modification is in effect, the strength of capillary pressure will, following Leverett (1941), automatically be scaled as Pcap ---> Pcap' = Pcap/SQRT(m). User-supplied permeability modifiers have to be entered as parameter "PMX" in columns 41-50 of an ELEMEnt record. Permeability modification options are selected through parameters in data block "ROCKS". &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& Summary of available permeability modification options (with s - random number between 0 and 1; PMX - user-supplied modifiers in data block "ELEME"): (1) externally supplied: m = PMX - PER(2) (2) "linear" (DROK.ne.0): m = PER(1) * s - PER(2) (3) "logarithmic" (DROK.eq.0): m = exp(- PER(1) * s) - PER(2) &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& &&&& if a domain "SEED " is present, permeability modification will be made &&&& if no domain "SEED " is present, no permeability modification will be made >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< >>>>>>>>>>>>>>>>>>>>>>>>>>> domain = "SEED " is not present, no permeability modification will be made <<<<<<<<<<<<<<<<<<<<<<<<<< >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< Data provided in domain "SEED " are used to select the following options. DROK = *** random number seed for internal generation of "linear" permeability modifiers. = 0: (default) no internal generation of "linear" permeability modifiers. > 0: perform "linear" permeability modification; random modifiers are generated internally with DROK as seed. POR = *** random number seed for internal generation of "logarithmic" permeability modifiers, = 0: (default) no internal generation of "logarithmic" permeability modifiers. > 0: perform "logarithmic" permeability modification; random modifiers are generated internally with POR as seed. &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& &&&&& note: if both DROK and POR are specified as non-zero, DROK takes precedence &&&&& &&&&& if both DROK and POR are zero, permeability modifiers as supplied through "ELEME" data will be used &&&&& &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& PER(1) = *** scale factor (optional) for internally generated permeability modifiers. = 0: (defaults to PER(1) = 1): permeability modifiers are generated as random numbers in the interval (0, 1). > 0: permeability modifiers are generated as random numbers in the interval (0, PER(1)). PER(2) = *** shift (optional) for internal or external permeability modifiers. = 0: (default) no shift is applied to permeability modifiers. > 0: permeability modifiers are shifted according to m' = m - PER(2). All m' < 0 are set equal to zero. FILE *CO2TAB* EXISTS --- OPEN AS AN OLD FILE ************************************************* * CO2 PROPERTIES READ FROM FILE CO2TAB * * NPK=127; NTK= 51 * * P(1) = 0.100000E+06; P(127) = 0.600000E+08 * * T(1) = 0.304000E+01; T( 51) = 0.103040E+03 * ************************************************* *********************************************************************************************************************************** * ECO2M: EQUATION OF STATE FOR MIXTURES OF H2O/NaCl/CO2 * *********************************************************************************************************************************** OPTIONS SELECTED ARE: (NK,NEQ,NPH,NB) = (3,4,4, 6) NK = 3 - NUMBER OF FLUID COMPONENTS NEQ = 4 - NUMBER OF EQUATIONS PER GRID BLOCK NPH = 4 - NUMBER OF PHASES THAT CAN BE PRESENT NB = 6 - NUMBER OF SECONDARY PARAMETERS (OTHER THAN COMPONENT MASS FRACTIONS) FOR NB = 6, DIFFUSION IS "OFF", FOR NB = 8, DIFFUSION IS "ON" AVAILABLE OPTIONS: (NK,NEQ,NPH,NB)=(3,4,4,6 OR 8) -WATER, NACL, NCG; NON-ISOTH.; VAR. (P, XS OR SS+10, X3 OR SAQ OR SG, T OR SG) (3,3,4,6 OR 8) -WATER, NACL, NCG; ISOTHERMAL; VAR. (P, XS OR SS+10, X3 OR SAQ OR SG, T) DEFAULT OPTIONS ARE (3,4,4,6) - NON-ISOTHERMAL, DIFFUSION "OFF" *********************************************************************************************************************************** THE PRIMARY VARIABLES ARE P - PRESSURE; XS - SALT MASS FRACTION IN LIQUID; SS SOLID SATURATION; X3 - NCG MASS FRACTION SAQ - AQ. PHASE SATURATION; SG - GAS PHASE SATURATION T - TEMPERATURE ****************************** ********************************************************* * COMPONENTS * * FLUID PHASE CONDITION PRIMARY VARIABLES * ****************************** ********************************************************* * * * * * # 1 - WATER * * SINGLE-PHASE P, XSM, X3, T * * * * * * # 2 - NACL * * TWO-PHASE W/AQUEOUS P, XSM, SAQ, T * * * * * * # 3 - CO2 * * TWO-PHASE W/O AQUEOUS P, XSM, SG, T * * * * * * # 4 - HEAT * * THREE FLUID PHASES P, XSM, SAQ, SG * * * * * ****************************** ********************************************************* *********************************************************************************************************************************** * ARRAY *IE* ALLOWS TO MAKE CHOICES AMONG DIFFERENT EOS OPTIONS * *********************************************************************************************************************************** IE(1) = 1: *** NUMBER OF ADDITIONAL RECORDS READ IN DATA BLOCK SELEC. IE(3) = 0: *** ALLOWS CHOICE OF OPTIONS FOR CALCULATION OF BRINE VISCOSITY. = 0: PHILLIPS ET AL. (1981) (DEFAULT). = 1: PALLISER AND MCKIBBIN (1998). = 2: MAO AND SUN (2006) = 3: POTTER (1976). IE(4) = 1: *** ALLOWS TO CHOOSE AMONG DIFFERENT CORRELATIONS FOR THE DENSITY OF COMPRESSED BRINE. = 1: ANDERSEN ET AL. (1992). (DEFAULT=0) = 2: PRITCHETT (1993). = 3: BRINE COMPRESSIBILITY EQUAL TO WATER COMPRESSIBILITY AT THE SAME REDUCED TEMPERATURE. = 4: BRINE COMPRESSIBILITY EQUAL TO WATER COMPRESSIBILITY AT THE SAME TEMPERATURE = 5: BATZLE AND WANG (1992). = 6: DRIESNER (2007). IE(8) = 0: *** ALLOWS CHOICE OF THE MANNER IN WHICH PHASE ASSIGNMENTS ARE HANDLED FOR SUPERCRITICAL TEMPERATURES. RECALL THAT FOR T > TCRIT, POINTS WITH P.GE.PCRIT ARE ARBITRARILY ASSIGNED AS LIQUID, POINTS WITH P.LT.PCRIT AS GAS. = 0: A-L <==> A-G TRANSITIONS WILL ALWAYS BE MADE, BASED ON COMPARING P WITH PCRIT. .NE. 0: A-L <==> A-G TRANSITIONS FOR SUPERCRITICAL T WILL ONLY BE MADE WHEN ITER.LE.IE(8). IE(9) = 0: *** ALLOWS TO CHOOSE A FINITE WINDOW FOR CERTAIN PHASE TRANSITIONS. A "HAIRTRIGGER" CRITERION WILL BE APPLIED ONLY FOR ITER < IE(9) FOR ITER.GE.IE(9), LIQ ==> L-G, A-L ==> A-L-G, AND A-G ==> A-L-G TRANSITIONS WILL ONLY BE MADE WHEN PCO2 IS OUTSIDE A FINITE WINDOW. THE LIQ ==> L-G AND A-L ==> A-L-G TRANSITIONS WILL ONLY BE MADE WHEN PCO2 < (1-FE(3))*PSAT,CO2. THE A-G ==> A-L-G TRANSITION WILL ONLY BE MADE WHEN PCO2 > (1+FE(3))*PSAT,CO2. The L ==> A-L transition will only be made when X3 < YCO2*(1-FE(3)). THE A ==> A-L, A-G TRANSITION WILL ONLY BE MADE WHEN X3 > XCO2AQ*(1.+FE(4)) IE(11) = 3: *** ALLOWS CHOICE OF DEPENDENCE OF PERMEABILITY ON ACTIVE FLOW POROSITY. = 0: PERMEABILITY DOES NOT VARY WITH FLOW POROSITY. = 1: PERMEABILITY VARIES AS PHIF**FE(1). = 2: SERIES FRACTURE MODEL. = 3: SERIES TUBE MODEL. IE(12) = 0: *** ALLOWS CHOICE OF MODEL FOR WATER SOLUBILITY IN CO2. = 0: AFTER SPYCHER AND PRUESS (2005). = 1: EVAPORATION MODEL (WATER DENSITY IN GAS BASED ON DENSITY OF VAPOR IN EQUILIBRIUM WITH BRINE). IE(13) = 0: *** ALLOWS CHOICE OF DEPENDENCE OF BRINE DENSITY ON DISSOLVED CO2. = 0: BRINE DENSITY VARIES WITH DISSOLVED CO2, ACCORDING TO THE CORRELATION FOR CO2 MOLAR VOLUME IN SOLUTION FROM J. GARCIA (LBNL-49023). = 1: BRINE DENSITY IS INDEPENDENT OF CO2 CONCENTRATION. IE(14) = 0: *** ALLOWS CHOICE OF TREATMENT OF THERMOPHYSICAL PROPERTIES DEPENDENCE ON SALT CONTENT. = 0: FULL DEPENDENCE (DEFAULT). > 0: SUPPRESSION OF ALL THERMOPHYSICAL PROPERTIES DEPENDENCE OTHER THAN SALT SOLUBILITY. IE(15) = 1: *** ALLOWS CHOICE OF CORRELATION FOR VAPOR SATURATED BRINE ENTHALPY. = 1: MICHAELIDES (1981). = 2: MILLER (1978). (OBSOLETE) = 3: PHILLIPS ET AL. (1981). = 4: LORENZ, MARIC AND RIRSHL (2000) (DEFAULT=0). = 5: DRIESNER (2007). fe(3) =0.000000E+00: *** if set to a small non-zero number, will specify a finite window for some phase changes [ie(9)] fe(4) =0.000000E+00 (winfac1): *** set to non-zero number, for finite window for a ==> a-l and a ==> a-g [ie(9)] *********************************************************************************************************************************** +++++++++ CANNOT FIND PARAMETERS AT ELEMENT *00001* XX(M) = 0.000000E+00 0.000000E+00 0.000000E+00 0.000000E+00 }}}}}}}}}}} ERRONEOUS DATA INITIALIZATION }}}}}}}}}} STOP EXECUTION---------- excessive aqueous phase saturation change at element *00001*; initial Saq = 0.000000E+00