{VERSION 5 0 "IBM INTEL NT" "5.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }1 1 0 0 0 0 1 0 1 0 2 2 0 1 }{PSTYLE "Headi ng 1" -1 3 1 {CSTYLE "" -1 -1 "Times" 1 18 0 0 0 1 2 1 2 2 2 2 1 1 1 1 }1 1 0 0 8 4 1 0 1 0 2 2 0 1 }} {SECT 0 {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}{PARA 0 " " 0 "" {TEXT -1 75 "Formulation of ordered compound (speciation) model for omphacitic pyroxene " }}{PARA 0 "" 0 "" {TEXT -1 82 "that can be \+ described as a mixture of two components (c=2) Jadeite (NaAlSi2O6) and " }}{PARA 0 "" 0 "" {TEXT -1 87 "diopside (CaMgSi2O6) as a mixture of \+ three species with the following site populations:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 61 " \+ Site:" }}{PARA 0 "" 0 "" {TEXT -1 78 " \+ 1 2 3 \+ 4" }}{PARA 0 "" 0 "" {TEXT -1 68 " \+ M2a M2b M1a M1b" }}{PARA 0 "" 0 "" {TEXT -1 66 " 1 Diopside Ca Ca Mg Mg" } }{PARA 0 "" 0 "" {TEXT -1 60 "Species: 2 Omphacite Na C a Al Mg" }}{PARA 0 "" 0 "" {TEXT -1 69 " 3 Jadeite Na Na Al Al" }}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 97 "each site has a multipl icity of 1/2 per formula unit. The use of such models for solid soluti ons " }}{PARA 0 "" 0 "" {TEXT -1 90 "has been championed in petrology \+ largely by Powell & Holland, but such models are standard" }}{PARA 0 " " 0 "" {TEXT -1 88 "practice in metallurgy and materials science. Powe ll & Holland formulate the speciation " }}{PARA 0 "" 0 "" {TEXT -1 93 "problem in terms of species activities/fugacities, here we first deve lop a different approach" }}{PARA 0 "" 0 "" {TEXT -1 97 "where the pro blem is formulated from the total free energy of the solution, and the n for didactic" }}{PARA 0 "" 0 "" {TEXT -1 52 "purposes we develop the activity/fugacity approach. " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 94 "The first step is to formulate the soluti on model in terms of the species fractions y[1]..y[3]" }}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 33 "Formulation of the solution model" }} {EXCHG {PARA 0 "" 0 "" {TEXT -1 103 "in general the model is formulate d as the sum of three terms, i.e., G_bar = Gmech - T Sconf + Gexcess " }}{PARA 0 "" 0 "" {TEXT -1 108 "each of which is expected to be a fu nction of the properties of the s species used to describe the solutio n." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 103 "He re we designate the free energy of the solution as G_bar to remind us \+ that the free energy is defined" }}{PARA 0 "" 0 "" {TEXT -1 99 "here p er mole of solution species, i.e., it is the molar free energy, not th e absolute free energy." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 32 "Gmech := sum(G0[i]*y[i],i=1..s);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 104 "If we are only concerned with internal p hase equilibria we are free to set two (i.e., c) of the species " }} {PARA 0 "" 0 "" {TEXT -1 85 "energies arbitrarily, so here we set G0[1 ]=G0[3]=0, then the Delta G of the reaction " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 22 "1/2 Di + 1/2 Jd = Omph" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 103 "is G0[2] , here we assume this is not a funciton of pressure or temperature, i. e., it is the enthalpy of" }}{PARA 0 "" 0 "" {TEXT -1 84 "formation of ordered omphacite species from diopside and jadeite, a model paramete r." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 97 "Mor e generally, in a c-component solution it must be possible to describe the composition of any " }}{PARA 0 "" 0 "" {TEXT -1 96 "species as a \+ mixture of c other species (provided these have independent compositio ns), thus if " }}{PARA 0 "" 0 "" {TEXT -1 101 "a model involves s spec ies it is always possible to write s-c linearly independent reactions, and the" }}{PARA 0 "" 0 "" {TEXT -1 91 "thermodynamic properties of t hese reactions are in effect parameters of the solution model." } {MPLTEXT 1 0 0 "" }{TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "s := 3; G0[1] := 0; G0[3] := 0; G0[2] := Delta; " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 42 "The confi gurational entropy on any site is" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "Smix := -R*sum(m[j]*sum(z[j,i]*ln(z[j,i]),i=1..n[j]),j=1..q);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 114 "where q is the number of distinct sites on which mixing occurs (4, i.e., m2a, m2b...) and m[j] is the m ultiplicity" }}{PARA 0 "" 0 "" {TEXT -1 115 "of each site per formula \+ unit (1/2) and n[i] is the number of kinds of atoms that mix on site j (2 for each site in" }}{PARA 0 "" 0 "" {TEXT -1 10 "our Cpx). " }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 110 "The next step in the model formulation is to express the atomic fractions z[j, i] in terms of species fractions" }}{PARA 0 "" 0 "" {TEXT -1 108 "y[k] using the site and compositional notation from the table above, e.g., z[1,1] is the atomic fraction of " }}{PARA 0 "" 0 "" {TEXT -1 15 "Ca \+ on site M2a:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 56 "q := 4: m[1 ] := 1/2:m[2] := 1/2:m[3] := 1/2:m[4] := 1/2:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "n[1] := 2:n[2] := 2:n[3] := 2:n[4] := 2:" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 14 "z[1,1] := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "z[1,2] := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 " z[2,1] := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "z[2,2] := ...;" } }{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "z[3,1] := ...;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 14 "z[3,2] := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "z[4,1] := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "z[4,2] := ... ;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "Smix;" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 95 "The excess term is the final component of the model, he re we will assume that the excess energy" }}{PARA 0 "" 0 "" {TEXT -1 60 "involves only pair-wise interactions of the form W*y[1]*y[2]" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 58 "Gex := W[1,2]*y[1]*y[2]+W[2,3]*y[3] *y[2]+W[1,3]*y[1]*y[3];" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 46 "The to tal free energy of the solution is then:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "G_bar := Gmech - T*Smix + Gex;" }}}}{SECT 0 {PARA 3 " " 0 "" {TEXT -1 48 "Speciation as a function of chemical composition" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 98 "the stable speciation for any gi ven bulk composition must be that which minimizes the free energy " }} {PARA 0 "" 0 "" {TEXT -1 81 "of the solution, and the necessary condit ion for a minimimum is diff(f(x),x) = 0." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 101 "to obtain the speciation as a fun ction of the bulk composition of the solution we now need to convert" }}{PARA 0 "" 0 "" {TEXT -1 96 "the above expression into an expression with a single unknown that can be solved for a parameter" }}{PARA 0 " " 0 "" {TEXT -1 94 "that determines the speciation (i.e., the degree o f order in our cpx). To this end we may use:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "eq1 := 1 = y[1] + y[2] + y[3];" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 89 "and if the bulk composition of the pyroxene is \+ defined by X, the mole fraction of Jadeite" }}{PARA 0 "" 0 "" {TEXT -1 32 "component (NOT Jadeite species)." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "eq2 := X = ...;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 31 "then we have for y[1] and y[3]:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 59 "Y1 := solve(...,...); Y3 := solve(subs(y[1]=Y1,...),y [3]); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 88 "So we now have now have a function for G that we can minimize at any specified values of" }} {PARA 0 "" 0 "" {TEXT -1 96 "X,T,W[1,2],W[1,3],W[2,3],Delta for y[2], \+ we can the back substitute this value into the previous" }}{PARA 0 "" 0 "" {TEXT -1 34 "equations to obtain y[1] and y[3]:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "dG := diff(subs(...,G_bar),...);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 98 "test if it works with an ideal sol ution (W's=0), are all values of Delta (i.e. G0[omph] feasible)?" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 91 "fsolve(subs(R=8.314,X=0.5,T= 1000,W[1,2]=0e3,W[1,3]=0e3,W[2,3]=0e3,Delta=...,dG),y[2]=0..1);" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 110 "now make a loop to compute and pl ot the speciation as a function of X, and, anticipating that we will w ant to " }}{PARA 0 "" 0 "" {TEXT -1 71 "compute phase relations for th e Cpx solution also calculate the molar G" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 89 "to simplify things assign the solu tion model parameters dummy values outside of the loop:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 49 "w12 := ...; w13 := ...; w23 := ...; delta := ...;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "pi := 0; dX := 0. 025: XX := dX: R := 8.314: T:= 1000:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "while XX < 1 do:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 " pi := \+ pi + 1:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 90 " y2[pi] := fsolve(subs (X=XX,W[1,2]=w12,W[1,3]=w13,W[2,3]=w23,Delta=delta,dG),y[2]=0..1):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " y1[pi] := ...:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " y3[pi] := ...:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 " x[pi] := XX:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " g[pi ] := ...:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " XX := XX + dX:" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "end do:" }}{PARA 0 "" 0 "" {TEXT -1 69 "put the points into arrays for plotting and free energy minimizati on:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "y1_points := [seq([x[i], y1 [i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "y2_points := [s eq([x[i], y2[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "y3_p oints := [seq([x[i], y3[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "g_x_points := [seq([x[i], g[i]],i=1..pi)]:" }}{PARA 0 "" 0 "" {TEXT -1 54 "plot the speciation as a function of bulk composition:" } }{PARA 0 "> " 0 "" {MPLTEXT 1 0 189 "plot([y1_points,y1_points,y2_poin ts,y2_points,y3_points,y3_points], x=0..1, style=[point,line],symbol=[ circle,cross],color=[blue,red,yellow],axes=boxed,thickness=3,labels=[ \"X(JD)\",\"y[i]\"]);" }}{PARA 0 "" 0 "" {TEXT -1 83 "do free energy m inimization to determine the (internal) phase relations of the Cpx:" } }{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "with(simplex):" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 44 "stable_points := convexhull( g_x_points ):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 231 "plot([g_x_points,stable_points,sta ble_points], x=0..1, style=[point,point,line],symbol=[circle,cross],co lor=[blue,red,yellow],axes=boxed,thickness=3,legend=[\"Compounds\",\"S table compounds\",\"Minimum G surface\"],labels=[\"X(Jd)\",\"G\"]);" } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 84 " Find solution model parameters that would allow omphacitic pyroxene to coexist with " }}{PARA 0 "" 0 "" {TEXT -1 64 "diopside-rich or jadeit e-rich pyroxene as is sometimes the case." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 85 "Plot the equilirbium composition of the Cpx in a ternary compos ition diagram, to this" }}{PARA 0 "" 0 "" {TEXT -1 85 "end it may be h elpful to note that two orthogonal cordinates y[2]-y[3] transform to a " }}{PARA 0 "" 0 "" {TEXT -1 34 "equilateral triangle as: " } }{PARA 0 "" 0 "" {TEXT -1 11 " " }}{PARA 0 "" 0 "" {TEXT -1 43 " y[1]' = y[1] + y[2]/2" }}{PARA 0 "" 0 "" {TEXT -1 46 " y[2]' = y[2] * sin(Pi/3)" }}} {EXCHG }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 39 "Derivation of species activities from G" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 92 "Fu gacities and activities are an artificial thermodynamic concept introd uced by noting that " }}{PARA 0 "" 0 "" {TEXT -1 90 "the total free en ergy of a solution may always be expressed as a linear combination of \+ the" }}{PARA 0 "" 0 "" {TEXT -1 76 "partial molar free energies of it' s constituent species g[j] = diff(G,n[j]):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "eq3 := Gtotal = sum(g[i]*y[i],i=1..ss);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 96 "where the activity of the species j relat ive to the pure species at the pressure and temperature" }}{PARA 0 "" 0 "" {TEXT -1 28 "of interest is defined from:" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 58 "eq4 := g[i] = G0[i] + RT*ln(aideal[i]) + RT*ln (gamma[i]) ;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 94 "where activity is normally defined as the produce aideal[i]*gamma[i] and aideal[i] is r elated " }}{PARA 0 "" 0 "" {TEXT -1 95 "to the partial molar entropy o f species i, without going into detail the ideal activity can be " }} {PARA 0 "" 0 "" {TEXT -1 98 "written as the product the of the site fr actions of the atoms present in the pure species, raised " }}{PARA 0 " " 0 "" {TEXT -1 102 "to a power corresponding to the site multiplicity . if this product is not unity for the pure species, " }}{PARA 0 "" 0 "" {TEXT -1 100 "the product must be normalized by the value of the pr oduct for the pure species (i.e., when y[i]=1, " }}{PARA 0 "" 0 "" {TEXT -1 26 "aideal[i] must = 1), e.g.," }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 95 "aideal[di] = z[m2a,Ca]^(1/2)*z[m2b,C a]^(1/2)*z[m1a,Mg]^(1/2)*z[m2a,Mg]/Normalization_factor[di]" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 64 "and substutitin g the site populations for pure diopside we have:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 96 "1 = 1^(1/2)*1^(1/2)*1^(1/ 2)*1/Normalization_factor[di], i.e., the normalization factor is unity ." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 5 "thus: " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 97 "aideal[1] := z[1,...]^( 1/2)*z[2,...]^(1/2)*z[3,...]^(1/2)*z[4,...]^(1/2)/Normalization_factor [1];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 31 "Normalization_facto r[1] := ...;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 97 "aideal[2] : = z[1,...]^(1/2)*z[2,...]^(1/2)*z[3,...]^(1/2)*z[4,...]^(1/2)/Normaliz ation_factor[2];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "Normali zation_factor[2] := 1;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 97 "a ideal[3] := z[1,...]^(1/2)*z[2,...]^(1/2)*z[3,...]^(1/2)*z[4,...]^(1/2 )/Normalization_factor[3];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "Normalization_factor[3] := 1;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 94 "comparison of out previous expression for the free energy of a sol ution with eq3 and eq4 shows" }}{PARA 0 "" 0 "" {TEXT -1 25 "that (RTl ng[i] is R*T*ln:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "eq5 := \+ GGex/RT = sum(lng[i]*y[i],i=1..ss);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 97 "but because we know that the partial molar energies for any arb itrary speciation substituted into" }}{PARA 0 "" 0 "" {TEXT -1 95 "eq5 define a plane that extrapolates to RTlng[i] when y[i]=1, eq5 can be \+ rearranged to obtain a" }}{PARA 0 "" 0 "" {TEXT -1 66 "general express ion for the activity of each species in a solution:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 62 "eq6 := lng[j] = (GGex/RT + sum((K[i,j]-y[ i])*lng[i],i=1..ss));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 73 "where K[ i,j] is the Kronecker delta (i.e., = 0 when i=j and 1 otherwise)." }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "K := array(1..3,1..3,[[0,1,1],[1,0, 1],[1,1,0]]):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 65 "lnaj := ln(aideal[ j]) + subs(lng[i]=Diff(GGex,y[i])/RT,rhs(eq6));" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 85 "then substituting our function for G of Cpx as a fun ction of its species composition:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "lna[1] := subs(j=1,ss=3,lnaj):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 67 "lna[1] := collect(eval(subs(GGex=Gex,Diff=diff,lna[1] )),\{RT,y[1]\});" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "lna[2] \+ := subs(j=2,ss=3,lnaj):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 67 "lna[2] : = collect(eval(subs(GGex=Gex,Diff=diff,lna[2])),\{RT,y[2]\});" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "lna[3] := subs(j=3,ss=3,lnaj ):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 67 "lna[3] := collect(eval(subs(G Gex=Gex,Diff=diff,lna[3])),\{RT,y[3]\});" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 70 "Speciation \+ problems in terms of an activity/fugacity formulation for G" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 84 "To solve the speciation problem we note that for any equilibrium speciation th e free" }}{PARA 0 "" 0 "" {TEXT -1 88 "energy change for every indepen dent reaction that can be written among the species must " }}{PARA 0 " " 0 "" {TEXT -1 63 "be zero. For our Cpx there is only one such reacti on for which:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "eq7 := 0 = Delta + RT*lnK;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "lnK := \+ simplify(.....);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 92 "And substitut ing our former mass balance constraints for y[1] and y[3] as functions of y[2]" }}{PARA 0 "" 0 "" {TEXT -1 88 "and bulk composition X, eq7 \+ becomes a univariate function in y[2] that may be solved for" }}{PARA 0 "" 0 "" {TEXT -1 47 "a given X and set of solution model parameters: " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 99 "if we remember that that partial molar energies (basically the activities) \+ are derivatives of G_bar" }}{PARA 0 "" 0 "" {TEXT -1 98 "it is apparen t that basically we have again solved for the condition diff(G_bar,y[2 ]) = 0, but in " }}{PARA 0 "" 0 "" {TEXT -1 91 "a somewhat less transp arent way than previously. In fact by the previous method we evaluate " }}{PARA 0 "" 0 "" {TEXT -1 101 "the derivative once, whereas here to obtain the activities we must do nine differentiations of G_bar." }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 55 "eq8 := simplify(subs(RT=R*T,y[1]=Y1,y[3]=Y3,rhs(eq7)));" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 17 "test if it works:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 92 "fsolve(subs(R=8.314,X=0.5,T=1000,W[ 1,2]=0e3,W[1,3]=0e3,W[2,3]=0e3,Delta=...,eq8),y[2]=0..1);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 67 "make a loop to calculate speciation and a ctivity as a function of X" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 51 "w12 := 0e3; w13 := 0e3; w23 := 0e3; delta := -20e3;" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 53 "pi := 0; dX := 0.025: XX := dX: R := 8.314: T: = 1000:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "while XX < 1 do:" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 " pi := pi + 1:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 91 " y2[pi] := fsolve(subs(X=XX,W[1,2]=w12,W[1,3]=w1 3,W[2,3]=w23,Delta=delta,eq8),y[2]=0..1):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 27 " y1[pi] := XX-1/2*y2[pi]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 29 " y3[pi] := 1-XX-1/2*y2[pi]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 " x[pi] := XX:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 112 " a1[pi] := subs(y[1]=y1[pi],y[2]=y2[pi],y[3]=y3[pi],W[1,2]=w12, W[1,3]=w13,W[2,3]=w23,Delta=delta,exp(lna[1])):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 112 " a2[pi] := s ubs(y[1]=y1[pi],y[2]=y2[pi],y[3]=y3[pi],W[1,2]=w12,W[1,3]=w13,W[2,3]=w 23,Delta=delta,exp(lna[2])):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 112 " \+ a3[pi] := subs(y[1]=y1[pi],y[2]=y2[pi],y[3]=y3[pi],W[1,2]=w12,W[1,3]= w13,W[2,3]=w23,Delta=delta,exp(lna[3])):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " XX := XX + dX:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "end do:" }}{PARA 0 "" 0 "" {TEXT -1 69 "put the points into arrays \+ for plotting and free energy minimization:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y1_points := [seq([x[i], y1[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y2_points := [seq([x[i], y2[i]],i=1..pi)] :" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y3_points := [seq([x[i], y3[i] ],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "a1_points := [seq([ x[i], a1[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "a2_point s := [seq([x[i], a2[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "a3_points := [seq([x[i], a3[i]],i=1..pi)]:" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 213 "plot([y1_points,y2_points,y3_points,a1_points ,a2_points,a3_points], x=0..1, style=[line,line,line,point,point,point ],color=[red,blue,yellow,red,blue,yellow],axes=boxed,thickness=3,label s=[\"X(JD)\",\"y[i] or a[i]\"]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 58 "COH-fluid s peciation as a function of chemical composition" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 94 "The only difference in fluid and s olid speciation calculations arises due to the convention by" }}{PARA 0 "" 0 "" {TEXT -1 105 "which the free energy of a fluid is described, i.e., ordinarily in terms of fugacities f[i] which account" }}{PARA 0 "" 0 "" {TEXT -1 106 "for the difference in energy from a pure speci es i ideal gas at the reference pressure and the temperature" }}{PARA 0 "" 0 "" {TEXT -1 107 "of interest, the fugacity may also be expresse d as p*y[i]*phi[i], where phi[i] is the fugacity coefficient," }} {PARA 0 "" 0 "" {TEXT -1 76 "which is generally a function of pressure , temperature, and speciation. Thus" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 54 "eq1 := Gtotal = sum(y[i]*(G0[i]+RT*ln(f[i])),i=1..ss) ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 63 "eq2 := Gtotal = sum(y[ i]*(G0[i]+RT*ln(p*y[i]*phi[i])),i=1..ss);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 97 "In the present case we consider a graphite saturated COH \+ fluid composed of five dominant species:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 6 "1 - H2" }}{PARA 0 "" 0 "" {TEXT -1 7 "2 - CH4" }}{PARA 0 "" 0 "" {TEXT -1 7 "3 - H2O" }}{PARA 0 "" 0 "" {TEXT -1 6 "4 - CO" }}{PARA 0 "" 0 "" {TEXT -1 7 "5 - CO2" }}{PARA 0 " " 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 96 "for notational cla rity we equate named and indexed species fractions and fugacities as \+ follows:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "s := 5;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 67 "yh2 := y[1]; ych4 := y[2]; y h2o := y[3]; yco := y[4]; yco2 := y[5];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 67 "fh2 := f[1]; fch4 := f[2]; fh2o := f[3]; fco := f[4]; fco2 := f[5];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 114 "f[1] := \+ p*y[1]*phi[1]; f[2] := p*y[2]*phi[2]; f[3] := p*y[3]*phi[3]; f[4] := p *y[4]*phi[4]; f[5] := p*y[5]*phi[5];" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 97 "because we are concerned with g raphite saturated conditions we may reduce the problem to one with" }} {PARA 0 "" 0 "" {TEXT -1 83 "two components, O and H, by making the Le gendre transformation Omega = G - n_C*G_C." }}{PARA 0 "" 0 "" {TEXT -1 99 "For a system with c=2 and s=5 there are three (s-c) linearly in dependent reactions possible between" }}{PARA 0 "" 0 "" {TEXT -1 50 "t he species (and carbon), which we may specify as " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 60 "1) C + 2 H2 = CH4 \+ => Delta_1" }}{PARA 0 "" 0 "" {TEXT -1 59 "2) 2 CO = C + CO2 => Delta_2" }}{PARA 0 " " 0 "" {TEXT -1 52 "3) C + 2 H2O = CH4 + CO2 => Delta _3" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 98 "the free energy change of these reactions (which is identical in G or Ome ga, provided the chemical" }}{PARA 0 "" 0 "" {TEXT -1 98 "potential of carbon is determined by saturation), gives us three equations in term s of the unknown" }}{PARA 0 "" 0 "" {TEXT -1 17 "speciation, i.e.:" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "eq3 := 0 = Delta_1 + RT*ln(fch4/fh2^2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "eq4 := 0 = Delta_2 + RT*ln(fco2/fco^2);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "eq5 := 0 = Delta_3 + RT*ln(f co2*fch4/fh2o^2);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 74 "we obtain tw o additional equations by the general mass balance constraint:" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 28 "eq6 := 1 = sum(y[i],i=1..s); " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 89 "and the constraint on the bul k composition of the system specified here by X = nO/(nO+nH)" }}{PARA 0 "" 0 "" {TEXT -1 77 "where nO and nH can be expressed in terms of th e fractions of the species as:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "nO := ...;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "nH := ...;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "eq7 := X = nO/(nO+ nH);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 92 "so we now have five equations, eq3-eq7, in five unknowns (y[1-5 ]) for any specified T, P, X," }}{PARA 0 "" 0 "" {TEXT -1 97 "and ther modynamic parameters phi[1-5] and Delta_[1-3]. In the general case tha t phi is a function" }}{PARA 0 "" 0 "" {TEXT -1 101 "of speciation we \+ would solve these iteratively, but here we simplify the problem by ass uming ideality" }}{PARA 0 "" 0 "" {TEXT -1 50 "in which case the fugac ity coefficients are unity:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 95 "in maple it is possible to group equations into sets and then solve the entire set of equations" }}{PARA 0 "" 0 "" {TEXT -1 92 "simultaneously with fsolve (numeric) or msolve (analytic) , this is more compact so we employ" }}{PARA 0 "" 0 "" {TEXT -1 17 "th is method here:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 29 "eqs := \{eq3,eq4,eq5,eq6,eq7\}:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "unkowns := \{y[1],y[2],y[3],y[4],y[ 5]\}" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 121 "neqs := subs(phi[1 ]=1,phi[2]=1,phi[3]=1,phi[4]=1,phi[5]=1,Delta_1=...,Delta_2=...,Delta_ 3=...,P=1000,RT=8314,p=1000,eqs);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 35 "fsolve(subs(X=0.33,neqs),unknowns);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 69 "find plausible values for the free energy chang es Delta_1 ... Delta_4" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 87 "Given a reaction CO + 1/2 O2 = CO2 for wh ich the free energy change is Delta_4 compute " }}{PARA 0 "" 0 "" {TEXT -1 90 "the fugacity of oxygen from the above speciation calculat ion. What value of Delta_4 would " }}{PARA 0 "" 0 "" {TEXT -1 69 "just ify our implicit assumption that O2 is not a significant species?" }} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 101 "and lastly we would like to be able to plot the speciati on, or perhaps fugacities, as a function of X" }}{PARA 0 "" 0 "" {TEXT -1 27 "with a do-loop of the form:" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 33 "pi := 0; dX := 0.09: XX := 0.01: " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "while XX < 1 do:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 " pi := pi + 1:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 " ans := fsolve(subs(X=XX,neqs),unknowns);" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " y1[pi] := rhs(ans[1]):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " y2[pi] := rhs(ans[2]):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " y3[pi] := rhs(ans[3]):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " y4[pi] := rhs(ans[4]):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " y5[pi] := rhs(ans[4]):" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " x[pi] := XX:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 " XX := XX + dX:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "end do;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y1_points := [seq ([x[i], y1[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y2_poi nts := [seq([x[i], y2[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y3_points := [seq([x[i], y3[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y4_points := [seq([x[i], y5[i]],i=1..pi)]:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "y5_points := [seq([x[i], y5[i]],i=1..pi)] :" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 179 "plot([y1_points,y2_po ints,y3_points,y4_points,y5_points], x=0..1, style=[line],symbol=[circ le],color=[blue,red,yellow,green,orange],axes=boxed,thickness=3,labels =[\"X(O)\",\"y[i]\"]);" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }} {PARA 0 "" 0 "" {TEXT -1 94 "It turns out that this maple strategy doe s not work for plotting because the order that fsolve" }}{PARA 0 "" 0 "" {TEXT -1 90 "outputs the species fractions changes between solution s, probably there is some simple way" }}{PARA 0 "" 0 "" {TEXT -1 95 "a round this problem, but an alternative which i suggest as an excercise here, is to explicitly " }}{PARA 0 "" 0 "" {TEXT -1 98 "solve the equ ations to obtain a single function of X and one species fraction, once this is solved" }}{PARA 0 "" 0 "" {TEXT -1 24 "the remaining fraction s " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 46 "COH-fluid speciation as a function of fugacity" }}}}{MARK "5 32 5 0" 24 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }