* Program SupcrtRxn * Reads reaction stoichiometries from Soltherm and writes Supcrt * reaction file(s) for Supcrt92 and SupcrtHP. * Originally written in 1996 for Soltherm.dat * Updated in 2008 to read Soltherm.xpt * Further updates to add component species * SupcrtRxn uses two name index tables to find equivalent names for * Sprons and Soltherm: * - Mins_ptx.tbl Minerals and Gases * - Aqss_pec.tbl Aqueous species * Copyright 1996, 2008, 2016, 2026 James Palandri * James Palandri hereby disclaims all copyright interest in the program * “SupcrtRxn” (which reformats text data) written by James Palandri * Signature of James Palandri, 24 January 2024 * James Palandri * This file is part of SupcrtRxn. * SupcrtRxn is free software: you can redistribute it and/or modify it under * the terms of the GNU General Public License as published by the * Free Software Foundation, either version 3 of the License, or (at your * option) any later version. SupcrtRxn is distributed in the hope that it * will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General * Public License for more details. You should have received a copy of the GNU * General Public License along with SupcrtRxn. If not, see * . Implicit none integer icompo parameter (icompo = 66) integer ispchg(icompo) integer numspc, iread, istoch(10) integer flsize, mintot, iloop, jloop integer nmin, ngas, nwater, naqs, minflg integer inp1, inp2, iout1, iout2, iout3 REAL*8 RADIUS(ICOMPO), RSTOCH(10), DRADIU, DCHARG REAL*8 SPCWT(ICOMPO), H2OCOE, MOLWT CHARACTER SKPREC, SKPRXN CHARACTER*3 CON_FLAG, GAS_FLAG CHARACTER*8 SPCNAM(ICOMPO), BLANK8, CHECK, CHRFLG, ref CHARACTER*20 MINNAM, SOLNAM CHARACTER*12 HPABBR, BLNK12, RXNFIL CHARACTER*20 BLNK20, HPNAME, SPECIE(ICOMPO) CHARACTER*32 DUMMY32 CHARACTER*80 DUMMY INP1 = 20 !soltherm INP2 = 21 !*.tbl IOUT1 = 10 !scratch IOUT2 = 11 !rxn file iout3 = 12 !missing minerals FLSIZE = 800 DATA BLANK8 /' '/ DATA BLNK12 /' '/ DATA BLNK20 /' '/ DATA SOLNAM /'********************'/ DATA SPECIE(1) /'H+ '/ DATA SPECIE(2) /'H2O '/ DATA SPECIE(3) /'Cl- '/ DATA SPECIE(4) /'SO4-2 '/ DATA SPECIE(5) /'HCO3- '/ DATA SPECIE(6) /'HS- '/ DATA SPECIE(7) /'H4SiO4,aq '/ DATA SPECIE(8) /'Al+3 '/ DATA SPECIE(9) /'Ca+2 '/ DATA SPECIE(10) /'Mg+2 '/ DATA SPECIE(11) /'Fe+2 '/ DATA SPECIE(12) /'K+ '/ DATA SPECIE(13) /'Na+ '/ DATA SPECIE(14) /'Mn+2 '/ DATA SPECIE(15) /'Zn+2 '/ DATA SPECIE(16) /'Cu+ '/ DATA SPECIE(17) /'Pb+2 '/ DATA SPECIE(18) /'Ag+ '/ DATA SPECIE(19) /'Au+ '/ !** now Au+ or Au+3, rewrite DATA SPECIE(20) /'Hg+2 '/ !** now Hg++, rewrite DATA SPECIE(21) /'Sr+2 '/ DATA SPECIE(22) /'Ba+2 '/ DATA SPECIE(23) /'F- '/ DATA SPECIE(24) /'Sb(OH)3,AQ '/ !-- not in sprons DATA SPECIE(25) /'H2AsO3- '/ !** now H2AsO3-, but no As mins DATA SPECIE(26) /'Ni+2 '/ DATA SPECIE(27) /'HPO4-2 '/ DATA SPECIE(28) /'Co+2 '/ DATA SPECIE(29) /'MoO4-2 '/ DATA SPECIE(30) /'UO2+2 '/ !-- not in sprons DATA SPECIE(31) /'O2,AQ '/ DATA SPECIE(32) /'NH4+ '/ DATA SPECIE(33) /'B(OH)3,AQ '/ DATA SPECIE(34) /'Ti(OH)4,aq '/ DATA SPECIE(35) /'ACETATE,AQ '/ DATA SPECIE(36) /'OXALATE,AQ '/ DATA SPECIE(37) /'SUCCINATE,AQ '/ DATA SPECIE(38) /'MALONATE,AQ '/ DATA SPECIE(39) /'Pd+2 '/ DATA SPECIE(40) /'Pt+2 '/ DATA SPECIE(41) /'HTe- '/ DATA SPECIE(42) /'Sn+2 '/ DATA SPECIE(43) /'HSe- '/ DATA SPECIE(44) /'WO4-- '/ DATA SPECIE(45) /'METHANE,AQ '/ DATA SPECIE(46) /'H2,AQ '/ DATA SPECIE(47) /'N2,AQ '/ DATA SPECIE(48) /'NH3,AQ '/ DATA SPECIE(49) /'Sc+3 '/ DATA SPECIE(50) /'Y+3 '/ DATA SPECIE(51) /'La+3 '/ DATA SPECIE(52) /'Ce+3 '/ DATA SPECIE(53) /'Pr+3 '/ DATA SPECIE(54) /'Nd+3 '/ DATA SPECIE(55) /'Pm+3 '/ DATA SPECIE(56) /'Sm+3 '/ DATA SPECIE(57) /'Eu+3 '/ DATA SPECIE(58) /'Gd+3 '/ DATA SPECIE(59) /'Tb+3 '/ DATA SPECIE(60) /'Dy+3 '/ DATA SPECIE(61) /'Ho+3 '/ DATA SPECIE(62) /'Er+3 '/ DATA SPECIE(63) /'Tm+3 '/ DATA SPECIE(64) /'Yb+3 '/ DATA SPECIE(65) /'Lu+3 '/ DATA SPECIE(66) /'Li+ '/ OPEN (INP1, FILE = 'SOLTHERM.XPT', STATUS = 'OLD') 2 WRITE (*,'(//,1X,A30,A21,/)') 'Select the data type for which', + ' to create reactions: ' WRITE (*,'(3x,A,A)') ' Minerals & gases, name index table =', + ' Mins_ptx.tbl, (Enter 1)' WRITE (*,'(3x,A,A/)') ' Aqueous species, name index table =', + ' Aqss_ptx.tbl, (Enter 2)' READ (*,'(I1)') MINFLG IF (MINFLG.EQ.1) THEN rxnfil(1:12) = 'SupcMins.rxn' open (inp2, file = 'Mins_ptx.tbl', status = 'old') Open (iout3, file = 'Miss-Srx-min.lst', status = 'unknown') write (iout3,'(/,1x,A,/,A,/)') + 'These Soltherm full P-T grid items have no record in SpronsHP', + ' and so also no record in the name index table Mins_ptx.tbl' ELSE IF (MINFLG.EQ.2) THEN rxnfil(1:12) = 'SupcAqss.rxn' open (inp2, file = 'Aqss_ptx.tbl', status = 'old') Open (iout3, file = 'Miss-Srx-aqs.lst', status = 'unknown') write (iout3,'(/,1x,A,/,A,/)') + 'These Soltherm full P-T grid items have no record in SpronsHP', + ' and so also no record in the name index table Mins_ptx.tbl' ELSE goto 2 ENDIF **** Read soltherm.dat until component species data is found 10 READ (INP1,'(A32)') DUMMY32 IF (DUMMY32.NE.'H+ +1 308 ') GOTO 10 BACKSPACE (INP1) **** Read component aqueous species and properties DO 15 ILOOP = 1,ICOMPO READ (INP1,810) SPCNAM(ILOOP),ISPCHG(ILOOP), + RADIUS(ILOOP),SPCWT(ILOOP) 810 FORMAT (A8,T21,I2,T28,F4.2,T38,F10.5) 15 CONTINUE IF (MINFLG.EQ.2) THEN * Do aqueous species here. Start two lines after components. READ (INP1,'(A1)') SKPREC READ (INP1,'(A1)') SKPREC 20 OPEN (IOUT1, STATUS = 'SCRATCH') MINTOT = 0 DO 100 ILOOP = 1, FLSIZE IF (SOLNAM.NE.BLNK20) THEN **** Rewind name table and find start of data. 35 REWIND INP2 23 READ (INP2,'(A1)') SKPREC IF (SKPREC.EQ.'*') GOTO 23 BACKSPACE INP2 **** Read aqueous species reactions from soltherm.dat 30 READ (INP1,'(A20,T91,A3)') SOLNAM, CON_FLAG **** Test if done IF (SOLNAM.EQ.BLNK20.AND.ILOOP.EQ.1) GOTO 700 IF (SOLNAM.NE.BLNK20) THEN IF (CON_FLAG.NE.'All') THEN READ (INP1,'(4/)') GOTO 30 ENDIF READ (INP1,'(T74,A8)') ref READ (INP1,821) NUMSPC, + (RSTOCH(IREAD),ISTOCH(IREAD), IREAD = 1, NUMSPC) 821 FORMAT (I2,1X,12(F10.3,I2)) READ (INP1,'(37/)') **** Read lookup table file, find HPNAME 40 READ (INP2,'(A20,12X,A20)') MINNAM, HPNAME **** If match not found proceed to the next aqs species, line 35 IF (MINNAM.NE.SOLNAM) THEN IF (HPNAME.EQ.BLNK20) then write (iout3,'(4x,a20,2x,a8)') solnam, ref GOTO 35 endif GOTO 40 ENDIF 45 WRITE (IOUT1,'(A1,A20)') ' ', MINNAM NMIN = 0 NGAS = 0 NWATER = 0 **** Plus one to soltherm numspec to account for the species itself. NAQS = NUMSPC + 1 **** Check if water present, Set nwater to 1 if yes, 0 if no. **** Save coefficient, decrement NAQS. DO 50 JLOOP = 1, NUMSPC IF (ISTOCH(JLOOP).EQ.2) THEN NWATER = 1 H2OCOE = RSTOCH(JLOOP) NAQS = NAQS - 1 ENDIF 50 CONTINUE **** Write nm, na, ng, nw and stoichiometry to scratch file. WRITE (IOUT1,'(4I4)') NMIN,NAQS,NGAS,NWATER WRITE (IOUT1,950) -1.000,HPNAME DO 60 JLOOP = 1, NUMSPC IF (ISTOCH(JLOOP).EQ.2) GOTO 60 WRITE(IOUT1,950) RSTOCH(JLOOP),SPECIE(ISTOCH(JLOOP)) 60 CONTINUE IF (NWATER.EQ.1) WRITE (IOUT1,950) H2OCOE, SPECIE(2) H2OCOE = 0 WRITE (IOUT1,'(A8)') BLANK8 MINTOT = MINTOT + 1 ENDIF ENDIF 100 CONTINUE **** Write eof marker, write contents to output file, close scratch file. WRITE (IOUT1,900) '********' REWIND IOUT1 CALL RXNWRT(MINTOT, CHRFLG, IOUT1, IOUT2, RXNFIL) CLOSE (IOUT1) IF (SOLNAM.EQ.BLNK20) GOTO 700 GOTO 20 ELSE *====================================================================== **** Gases and minerals here. **** Find start of gas and mineral data in soltherm. 110 READ (INP1,800) CHRFLG IF (CHRFLG.NE.'C123 COE') GOTO 110 READ (INP1,'(/)') *** Do gases first. NMIN = 0 NGAS = 1 **** Read reactions from soltherm and write to scratch file. 115 OPEN (IOUT1, STATUS = 'scratch') MINTOT = 0 DO 200 ILOOP = 1, FLSIZE IF (SOLNAM.NE.BLNK20) THEN **** Rewind name table and find start of data. 130 REWIND INP2 140 READ (INP2,'(A1)') SKPREC IF (SKPREC.EQ.'*') GOTO 140 BACKSPACE INP2 **** Reset species numbers and coefficients to zeros DO 145 JLOOP = 1, 10 RSTOCH(JLOOP) = 0.0 ISTOCH(JLOOP) = 0 145 CONTINUE **** Read mineral or gas name from SOLTHERM 300 READ (INP1,'(A20,T81,A3,T91,A3)')SOLNAM,GAS_FLAG,CON_FLAG **** Change ngas and nmin to minerals upon reaching gas buffers IF (SOLNAM.EQ.'fO2-0.7 ') THEN NGAS = 0 NMIN = 1 ENDIF **** Test if done IF (SOLNAM.EQ.BLNK20.AND.ILOOP.EQ.1) GOTO 700 IF (SOLNAM.NE.BLNK20) THEN **** Skip L-V minerals IF (CON_FLAG.NE.'All') THEN READ (INP1,*) !skip wt, vol and dens data READ (INP1,*) !skip stochiometry **** Skip gas fugacity coeffs if present IF (GAS_FLAG.EQ.'GAS') THEN READ (INP1,'(A8)') CHECK IF (.NOT.(CHECK(1:3).EQ.'B,C')) BACKSPACE INP1 ENDIF READ (INP1,*) !skip log K's READ (INP1,*) !skip blank line GOTO 300 **** Ready to read stochiomitry ENDIF READ (INP1,'(T74,A8)') ref READ (INP1,821) NUMSPC, + (RSTOCH(IREAD),ISTOCH(IREAD), IREAD = 1, NUMSPC) IF (GAS_FLAG.EQ.'GAS') THEN READ (INP1,'(A8)') CHECK IF (.NOT.(CHECK(1:3).EQ.'B,C')) BACKSPACE INP1 ENDIF READ (INP1,'(37/)') **** Read lookup table file, find match between minnam (supcrt) and solnam (soltherm) 150 READ (INP2,'(A20,12X,A20)') MINNAM,HPNAME **** If match not found proceed to next gas or min, line 130 IF (MINNAM.NE.SOLNAM) THEN IF (HPNAME.EQ.BLNK20) then write (iout3,'(4x,a20,2x,a8)') solnam, ref GOTO 130 endif GOTO 150 ENDIF WRITE (IOUT1,'(A1,A20)') ' ', MINNAM NWATER = 0 NAQS = NUMSPC **** Check if water present, Set nwater to 1 if yes. **** Save coefficient, decrement NAQS. DO 160 JLOOP = 1, NUMSPC IF (ISTOCH(JLOOP).EQ.2) THEN NWATER = 1 H2OCOE = RSTOCH(JLOOP) NAQS = NAQS - 1 ENDIF 160 CONTINUE **** Write nm, na, ng, nw and stoichiometry in the scratch file. WRITE (IOUT1,'(4I4)') NMIN,NAQS,NGAS,NWATER **** Need to write species/phases in proper order: minerals, aqs species, gases, H2O IF (NMIN.GT.0) WRITE (IOUT1,950) -1.000,HPNAME DO 180 JLOOP = 1, NUMSPC IF (ISTOCH(JLOOP).EQ.2) GOTO 180 WRITE(IOUT1,950) RSTOCH(JLOOP),SPECIE(ISTOCH(JLOOP)) 180 CONTINUE IF (NGAS.GT.0) WRITE (IOUT1,950) -1.000,HPNAME IF (NWATER.EQ.1) WRITE (IOUT1,950) H2OCOE, SPECIE(2) H2OCOE = 0 WRITE (IOUT1,'(A8)') BLANK8 MINTOT = MINTOT + 1 ENDIF ENDIF 200 CONTINUE **** Write eof marker, write contents to output file, close scratch file. WRITE (IOUT1,900) '********' REWIND IOUT1 CALL RXNWRT(MINTOT, CHRFLG, IOUT1, IOUT2, RXNFIL) CLOSE (IOUT1) IF (SOLNAM.EQ.BLNK20) GOTO 700 GOTO 115 ENDIF 700 WRITE(*,'(1X,4X,A20)') 'EXECUTION COMPLETED.' *===================================================================== **** Input formats for SOLTHERM.DAT 800 FORMAT (A8) 803 FORMAT (A1) **** Output formats for SUP.RXN 900 FORMAT (A8) 902 FORMAT (A8,F8.2,I2,1X,10(F7.3,I2)) 930 FORMAT (A8,I2,F6.1,F10.5) 950 FORMAT (1X,F9.3,2X,A20,2X,A30) 999 END *===================================================================== * Subroutine RXNWRT * Write SupcrtRxn with the format header and the total number of * reactions followed by the contents of the scratch file. *===================================================================== SUBROUTINE RXNWRT (MINTOT, CHRFLG, IOUT1, IOUT2, RXNFIL) INTEGER MINTOT, IOUT1, IOUT2, ITEMP CHARACTER*8 CHRFLG CHARACTER*12 RXNFIL CHARACTER*80 LINE IOUT1 = 10 IOUT2 = 11 WRITE (*,'(A21,A12)') 'WRITING RXNFIL: ', RXNFIL 10 OPEN (IOUT2, FILE = RXNFIL, STATUS = 'UNKNOWN') WRITE (IOUT2,'(A10,A50)')' Line 1: ', +'nreac, iwet (free format) ' WRITE (IOUT2,'(A10,A50)')' Line 2: ', +'[blank] (free format) ' WRITE (IOUT2,'(A10,A50)')' Line 3: ', +'descriptive title (a80) ' WRITE (IOUT2,'(A10,A50)')' Line 4: ', +'nm, na, ng, nw (free format) ' WRITE (IOUT2,'(A10,A50)')' nm Lines:', +' coeff mname mform (1X,F9.3,2X,A20,2X,A30) ' WRITE (IOUT2,'(A10,A50)')' ng Lines:', +' coeff aname aform (1X,F9.3,2X,A20,2X,A30) ' WRITE (IOUT2,'(A10,A50)')' na Lines:', +' coeff gname gform (1X,F9.3,2X,A20,2X,A30) ' WRITE (IOUT2,'(A10,A50)')' [1 Line: ', +' coeff H2O H2O] (1X,F9.3,2X,A20,2X,A30) ' WRITE (IOUT2,'(A60)')' ' WRITE (IOUT2,'(A37)')'*** each of the nreac reaction blocks' WRITE (IOUT2,'(A37)')'*** contains 3+nm+ng+na+nw lines ' WRITE (IOUT2,'(A60)')' ' WRITE (IOUT2,'(A10,A50)')'**********', +'********************************************* ' WRITE (IOUT2,'(A60)')' ' WRITE (IOUT2,'(I3,I3)') MINTOT,1 WRITE (IOUT2,'(A60)')' ' 1000 READ (IOUT1,'(A8)') CHRFLG IF (CHRFLG.NE.'********') THEN BACKSPACE IOUT1 READ (IOUT1,'(A80)') LINE WRITE (IOUT2,'(A80)') LINE GOTO 1000 ELSE CLOSE (IOUT2) ENDIF END *=====================================================================