-
|
If I have a pseudobinary system that's PbCl2 - ZnCl2 with
I can't get any equilibrium calculations to converge and get two different behaviors depending on my starting points. I think the root cause for both cases are the same - there's one independent chemical potential. Are there any workaround and/or potential solutions for this? cc: @richardotis Example code(requires from pycalphad import Database, equilibrium, Model, variables as v
# Assumes ModelMQMQA is defined, requires the cs_dat_support_linear branch
dbf = Database('CuZnFeCl-Viitala (1).dat')
phases = ['LIQUIDSOLN', 'PBCL2(S)', 'ZNCL2(S)']
print(equilibrium(dbf, ['PB', 'ZN', 'CL'], phases, {v.N: 1, v.P: 101325, v.T: 500, v.X('PB'): 0.0666666666, v.X('CL'): 0.6666666666}, model={'LIQUIDSOLN': ModelMQMQA(dbf, ['PB', 'ZN', 'CL'], 'LIQUIDSOLN')}, verbose=True, calc_opts={'pdens': 10_000}))Stoichiometric phases as starting pointWith stoichiometric phases as the starting point, the only variables are the phase fractions for each stoichiometric phase, so not too surprising that there aren't enough degrees of freedom. Liquid starting pointIf I have the liquid as a starting point, it also gets stuck. |
Beta Was this translation helpful? Give feedback.
Replies: 2 comments
-
|
One possible-but-not-very-good workaround is to "nudge" the composition off the psuedobinary and add another phase that could be used to mass balance and define an chemical potential. Here I choose CL gas: from pycalphad import Database, equilibrium, Model, variables as v
# Assumes ModelMQMQA is defined, requires the cs_dat_support_linear branch
dbf = Database('CuZnFeCl-Viitala (1).dat')
phases = ['LIQUIDSOLN', 'PBCL2(S)', 'ZNCL2(S)', 'CL2(G)']
print(equilibrium(dbf, ['PB', 'ZN', 'CL'], phases, {v.N: 1, v.P: 101325, v.T: 600, v.X('PB'): 0.0666666666, v.X('CL'): 0.66667}, model={'LIQUIDSOLN': ModelMQMQA(dbf, ['PB', 'ZN', 'CL'], 'LIQUIDSOLN')}, verbose=True, calc_opts={'pdens': 10_000}))verbose solver output: This might work, depending on the calculation needs. The phase fraction of Cl gas could be made vanishingly small (by approaching the pseudobinary composition line) and neglected if phase fractions are of interest. The total energy and phase compositions and internal degrees of freedom (if any) should be quite close. On the other hand, the pure element chemical potentials will not be meaningful (similar to the pure element chemical potentials if only one stoichiometric compound is stable in a binary system). This makes me wonder if it would be effective to try treating the pseudobinary case in a similar way that we do stoichiometric compounds. |
Beta Was this translation helpful? Give feedback.
-
|
It seems like the new minimizer fixes this. from pycalphad import Database, calculate, equilibrium, variables as v
import pycalphad
print(pycalphad.__version__)
TDB_CA_MG_O = """
ELEMENT MG ALPHA 0 0 0 !
ELEMENT CA BETA 0 0 0 !
ELEMENT O GAMMA 0 0 0 !
PHASE LIQUID % 2 1 1 !
CONSTITUENT LIQUID : CA MG : O : !
PHASE HALITE % 2 1 1 !
CONSTITUENT HALITE : CA MG : O : !
PARAMETER G(LIQUID,CA:O;0) 1 0; 10000 N !
PARAMETER G(LIQUID,MG:O;0) 1 0; 10000 N !
PARAMETER G(HALITE,CA:O;0) 1 12000 - 10*T; 10000 N !
PARAMETER G(HALITE,MG:O;0) 1 8000 - 10*T; 10000 N !
"""
dbf = Database(TDB_CA_MG_O)
comps = ["CA", "MG", "O"]
phases = ["LIQUID", "HALITE"]
conds = {v.P: 101325, v.T: 300, v.N: 1, v.X('O'): 0.5, v.X('MG'): 0.25}
eq = equilibrium(dbf, comps, phases, conds)
#print(eq)
print("Phase ", eq.Phase.values.squeeze())
print("NP ", eq.NP.values.squeeze())
print("X ", eq.X.values.squeeze())
print("Y ", eq.Y.values.squeeze())
print("MU ", eq.MU.values.squeeze())
TDB_AL_CR_O = """
ELEMENT AL BETA 0 0 0 !
ELEMENT CR ALPHA 0 0 0 !
ELEMENT O GAMMA 0 0 0 !
PHASE LIQUID % 2 2 3 !
CONSTITUENT LIQUID : AL CR : O : !
PARAMETER G(LIQUID,AL:O;0) 1 0; 10000 N !
PARAMETER G(LIQUID,CR:O;0) 1 0; 10000 N !
PHASE CORUND % 2 2 3 !
CONSTITUENT CORUND : AL CR : O : !
PARAMETER G(CORUND,AL:O;0) 1 12000 - 10*T; 10000 N !
PARAMETER G(CORUND,CR:O;0) 1 8000 - 10*T; 10000 N !
"""
dbf = Database(TDB_AL_CR_O)
comps = ["AL", "CR", "O"]
phases = ["LIQUID", "CORUND"]
conds = {v.P: 101325, v.T: 300, v.N: 1, v.X('O'): 0.6, v.X('CR'): 0.25}
eq = equilibrium(dbf, comps, phases, conds)
#print(eq)
print("Phase ", eq.Phase.values.squeeze())
print("NP ", eq.NP.values.squeeze())
print("X ", eq.X.values.squeeze())
print("Y ", eq.Y.values.squeeze())
print("MU ", eq.MU.values.squeeze())0.8.5
|
Beta Was this translation helpful? Give feedback.
It seems like the new minimizer fixes this.