Reactions

ReactionsReactions is a small handy freeware for generation of all possible reactions between specified substances (pure minerals and components of solid solutions, oxides and chemical elements) with selection of the arbitrary set of linearly independent reactions. This program may be useful for education, for analysis of diagrams (petrogenetic grids and results of multi-equilibrium geothermobarometry), for balancing of reactions involving minerals with real compositions, for evaluation of solid solution end-members, etc.

 The general approach is a sequence of relatively simple procedures that has been described in many publications (for example, in the chapter 5 of Frank Spear's monography Metamorphic Phase Equilibria and Pressure-Temperature-Time Paths (1993) or in the web-tutorial of Dexter Perkins at Teaching Phase Equilibria) and used in such programs as TWQ (Berman, 1991) and THERMOCALC (Powell & Holland, 1988). Nevertheless, I've decided to give a detailed description of this important algorithm. The description is accompanied by a rather complicated example of the system comprising almandine (Alm), sillimanite (Sill), annite (Ann), ferrosilite (Fs), muscovite (Mu), phlogopite (Phl), pyrope (Prp) and quartz (Qtz).

ScreenShot

 Calculation of reactions is performed in several stages.

 First of all, let's write a compositional matrix for all the substances that may be involved in the reactions:


AlmSillAnnFsMuPhlPrpQtz
K00101100
Fe30320000
Mg00000330
Al22103120
Si31323331
H00202200
O1251261212122

 Then transform it to the so-called Reduced Row-Echelon Form (RREF) by the Gauss-Jordan elimination method. However, let's modify a little the commonly used procedure to produce the trapezoidal matrix with the left square identity part (with "1" on the main diagonal and "0" in other cells) and deleted zero rows (which correspond to linearly dependent elements). You can find the Pascal implementation in the documentation.

 All the matrix transformations are based on elementary equivalent operations:

  1. Swapping of rows and columns
  2. Multiplication of a row by a constant
  3. Addition of rows

AlmSillAnnFsMuPhlPrpQtz
Fe1010.6670000Fe/3, moved up
K00101100
Mg00000330
Al02‑1‑1.3333120Al-Fe*2
Si01003331Si-Fe*3
H00202200
O050-21212122O-Fe*12

AlmSillAnnFsMuPhlPrpQtz
Fe1010.6670000
Al01‑0.5‑0.6671.50.510Al/2, moved up
K00101100
Mg00000330
Si000.50.6671.52.521Si-Al
H00202200
O002.51.3334.59.572O-Al*5

AlmSillAnnFsMuPhlPrpQtz
Fe1000.667‑1‑100Fe-K
Al010‑0.6672110Al+K/2
K00101100
Mg00000330
Si0000.6671221Si-K/2
H00000000H-K*2, deleted
O0001.3332772O-K*2.5

AlmSillAnnFsMuPhlPrpQtz
Fe1000‑2‑3‑2‑1Fe-Si/3*2
Al01003331Al+Si/3*2
K00101100
Si00011.5331.5Si*3/2, moved up
Mg00000330
O00000330O-Si/3*4

AlmSillAnnFsPhlMuPrpQtzPhl moved to the left
Fe1000‑3‑2‑2‑1
Al01003331
K00101100
Si000131.531.5
Mg00003030
O00003030

AlmSillAnnFsPhlMuPrpQtz
Fe10000‑21‑1Fe+Mg*3
Al01000301Al-Mg*3
K001001‑10K-Mg
Si000101.501.5Si-Mg*3
Mg00001010Mg/3
O00000000O-Mg*3, deleted

 The resulting trapezoidal matrix consists of only linearly independent rows (elements), its left part is a strict identity matrix and the right rectangular part contains calculated coefficients for some set of linearly independent reactions (with 3 reactions in this case):

AlmSillAnnFsPhlMuPrpQtz
Fe10000‑21‑1
Al01000301
K001001‑10
Si000101.501.5
Mg00001010

 Thus, to get the matrix of coefficients of these independent reactions we should simply transpose the right rectangular part and augment it by the negative identity matrix at the right side:

AlmSillAnnFsPhlMuPrpQtz
‑2311.50‑100
10‑1010‑10
‑1101.5000‑1

 The same matrix with integer coefficients (after multiplication by appropriate constants):

AlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
‑46230‑2001 4Alm+2Mu=6Sill+2Ann+3Fs(Phl,Prp,Qtz)
10‑1010‑102 Ann+Prp=Alm+Phl(Sill,Fs,Mu,Qtz)
‑2203000‑23 2Alm+2Qtz=2Sill+3Fs(Ann,Phl,Mu,Prp)

 All the following calculations use only this matrix with coefficients of linearly independent reactions: all other (dependent) reactions are different combinations of this small set.

 The last column of this table contains lists of substances not involved in reactions: it is a commonly accepted notation for reactions in complex systems (the "ID" of reaction). It should be noted that this notation plays a key role at generation of dependent reactions because calculations are based on reduction of coefficients to zero just for absent substances.

 As we can conclude from inspection of the compositional matrix transformed to RREF, the maximum number of substances M involved in any reaction should be more than the number of independent elements (rows) E by 1, i.e. M=E+1 (M=5+1=6 for this example system). Therefore, the number of absent substances N can not be less than total number of substances Q minus this maximum number of reaction substances M, i.e. N=Q-E-1 (N=8-5-1=2 for this example system).

 The procedure for generation of reactions is based on evaluation of all possible combinations by N substances in the set of Q substances. The total number of such combinations C can be estimated using the formula:
$$C = \frac{Q!}{M! \times N!}$$
Thus, we should evaluate and check 28 combinations (by 2 substances in the set of 8 substances) for our example system. However, all linearly independent reactions have got numbers of absent substances (N) >2 (3 and 4) because of so-called degeneracy: elimination of some substances due to appropriate stoichiometric relationships. Degeneracy leads to correspondence of several combinations to a single reaction, i.e. the number of possible dependent reactions should be (essentially) less than the total number of combinations. At the same time, it's not necessary to calculate reaction coefficients for every generated combination: just check whether its set of absent substances is a subset of the analogous set for one of already calculated reactions. If so, this combination can be skipped. Therefore every reaction record should also contain some data structure and functionality for effective operations with sets (e.g. based on bit masks).

 Calculation of reaction coefficients is performed by transformation of the matrix of linearly independent reactions: placing of columns with absent substances to its left side and reducing to a trapezoidal matrix by the Gauss elimination method. This method differs a bit from the Gauss-Jordan method used before: the left part of the resulting matrix is not a strict identity matrix, i.e. values on the main diagonal are not necessarily "1" and values above the main diagonal are not necessarily "0" (it is just an upper triangular matrix). The first values of the last row in this matrix will be nullified, i.e. the last row will contain coefficients of searched reaction (which you can convert to integers). The first step can be optimized if array with numbers of columns is used (to avoid moving whole columns).

 So, let's calculate all the remaining linearly dependent reactions in our system:


CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
1(Alm,Sill)‑46230‑200
10‑1010‑10
-2203000‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑1‑11.5010‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
00‑1.33320.6670.667‑0.667‑2
 
00‑2311‑1‑342Ann+Prp+3Qtz=3Fs+Phl+Mu(Alm,Sill)

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
2(Alm,Ann)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑1‑11.5010‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑400‑222‑2
 
0‑200‑111‑152Sill+Phl+Qtz=Mu+Prp(Alm,Ann,Fs)
 
3(Alm,Fs) – subset of reaction #5

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
4(Alm,Phl)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑1‑11.5010‑2
 
0‑2‑23020‑462Sill+2Ann+4Qtz=3Fs+2Mu(Alm,Phl,Prp)

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
5(Alm,Mu)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑1‑11.5010‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
02‑2320‑2‑272Ann+2Prp+2Qtz=2Sill+3Fs+2Phl(Alm,Mu)
 
6(Alm,Prp) – subset of reaction #6

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
7(Alm,Qtz)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
01.5‑0.50.751‑0.5‑10
0‑1‑11.5010‑2
 
‑46230‑200
0‑1‑11.5010‑2
01.5‑0.50.751‑0.5‑10
 
06‑234‑2‑4082Ann+2Phl+4Prp=6Sill+3Fs+4Phl(Alm,Qtz)

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
8(Sill,Ann)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
10‑1010‑10
‑0.6670‑0.667200.6670‑2
 
‑46230‑200
10‑1010‑10
‑1.333002‑0.6670.6670.667‑2
 
‑2003‑111‑392Alm+Phl+3Qtz=3Fs+Mu+Prp(Sill,Ann)
 
9(Sill,Fs) – subset of reaction #2

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
10(Sill,Phl)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
10‑1010‑10
‑0.6670‑0.667200.6670‑2
 
‑10‑13010‑310Alm+Ann+3Qtz=3Fs+Mu(Sill,Phl,Prp)
 
11(Sill,Mu) – subset of reaction #2
12(Sill,Prp) – subset of reaction #10
13(Sill,Qtz) – subset of reaction #2
14(Ann,Fs) – subset of reaction #5
15(Ann,Phl) – subset of reaction #3
16(Ann,Mu) – subset of reaction #3
17(Ann,Prp) – subset of reaction #3

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
18(Ann,Qtz)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
‑1301.51‑1‑10
‑2203000‑2
 
‑46230‑200
‑2203000‑2
‑1301.51‑1‑10
 
‑26032‑2‑20112Alm+2Mu+2Prp=6Sill+3Fs+2Phl(Ann,Qtz)

CombinationAbsent substancesAlmSillAnnFsPhlMuPrpQtz#ReactionAbsent substances
19(Fs,Phl)‑46230‑200
10‑1010‑10
‑2203000‑2
 
‑46230‑200
10‑1010‑10
2‑4‑20020‑2
 
1‑2‑10010‑1122Sill+Ann+Qtz=Alm+Mu(Fs,Phl,Prp)
 
20(Fs,Mu) – subset of reaction #2
21(Fs,Prp) – subset of reaction #12
22(Fs,Qtz) – subset of reaction #2
23(Phl,Mu) – subset of reaction #3
24(Phl,Prp) – subset of reaction #10
25(Phl,Qtz) – subset of reaction #1
26(Mu,Prp) – subset of reaction #3
27(Mu,Qtz) – subset of reaction #2
28(Prp,Qtz) – subset of reaction #1

 Thus, only 9 of 28 possible combinations (without taking into account the calculated set of 3 linearly independent reactions) provide new dependent reactions in this system, whereas the remaining 19 combinations correspond to already found (degenerate) reactions.

 Detection of linear dependency of some reaction on the specified set of linearly independent reactions (IR set) is also based on transformation of the coefficients' matrix of all these reactions to a trapezoidal form: if the tested reaction corresponds to combination of reactions from the IR set, the last row will be nullified after application of the Gauss elimination method.