Cantera  2.1.2
equilibrate.cpp
Go to the documentation of this file.
1 /**
2  * @file equilibrate.cpp Driver routines for the chemical equilibrium solvers.
3  */
4 
5 #include "cantera/equil/equil.h"
9 #include "cantera/base/global.h"
10 
11 namespace Cantera
12 {
13 
14 doublereal equilibrate(MultiPhase& s, const char* XY,
15  doublereal tol, int maxsteps, int maxiter,
16  int loglevel)
17 {
18  if (loglevel > 0) {
19  beginLogGroup("equilibrate",loglevel);
20  addLogEntry("multiphase equilibrate function");
21  beginLogGroup("arguments");
22  addLogEntry("XY",XY);
23  addLogEntry("tol",tol);
24  addLogEntry("maxsteps",maxsteps);
25  addLogEntry("maxiter",maxiter);
26  addLogEntry("loglevel",loglevel);
27  endLogGroup("arguments");
28  }
29  s.init();
30  int ixy = _equilflag(XY);
31  if (ixy == TP || ixy == HP || ixy == SP || ixy == TV) {
32  try {
33  double err = s.equilibrate(ixy, tol, maxsteps, maxiter, loglevel);
34  if (loglevel > 0) {
35  addLogEntry("Success. Error",err);
36  endLogGroup("equilibrate");
37 
38  }
39  return err;
40  } catch (CanteraError& err) {
41  err.save();
42  if (loglevel > 0) {
43  addLogEntry("Failure.",lastErrorMessage());
44  endLogGroup("equilibrate");
45  }
46  throw err;
47  }
48  } else {
49  if (loglevel > 0) {
50  addLogEntry("multiphase equilibrium can be done only for TP, HP, SP, or TV");
51  endLogGroup("equilibrate");
52  }
53  throw CanteraError("equilibrate","unsupported option");
54  return -1.0;
55  }
56  return 0.0;
57 }
58 
59 int equilibrate(thermo_t& s, const char* XY, int solver,
60  doublereal rtol, int maxsteps, int maxiter, int loglevel)
61 {
62  bool redo = true;
63  int retn = -1;
64  int nAttempts = 0;
65  int retnSub = 0;
66 
67  if (loglevel > 0) {
68  beginLogGroup("equilibrate", loglevel);
69  addLogEntry("Single-phase equilibrate function");
70  {
71  beginLogGroup("arguments");
72  addLogEntry("phase",s.id());
73  addLogEntry("XY",XY);
74  addLogEntry("solver",solver);
75  addLogEntry("rtol",rtol);
76  addLogEntry("maxsteps",maxsteps);
77  addLogEntry("maxiter",maxiter);
78  addLogEntry("loglevel",loglevel);
79  endLogGroup("arguments");
80  }
81  }
82  while (redo) {
83  if (solver >= 2) {
84  int printLvlSub = 0;
85  int estimateEquil = 0;
86  try {
87  MultiPhase m;
88  m.addPhase(&s, 1.0);
89  m.init();
90  nAttempts++;
91  vcs_equilibrate(m, XY, estimateEquil, printLvlSub, solver,
92  rtol, maxsteps, maxiter, loglevel-1);
93  redo = false;
94  if (loglevel > 0) {
95  addLogEntry("VCSnonideal solver succeeded.");
96  }
97  retn = nAttempts;
98  } catch (CanteraError& err) {
99  err.save();
100  if (loglevel > 0) {
101  addLogEntry("VCSnonideal solver failed.");
102  }
103  if (nAttempts < 2) {
104  if (loglevel > 0) {
105  addLogEntry("Trying single phase ChemEquil solver.");
106  }
107  solver = -1;
108  } else {
109  if (loglevel > 0) {
110  endLogGroup("equilibrate");
111  }
112  throw err;
113  }
114  }
115  } else if (solver == 1) {
116  try {
117  MultiPhase m;
118  m.addPhase(&s, 1.0);
119  m.init();
120  nAttempts++;
121  equilibrate(m, XY, rtol, maxsteps, maxiter, loglevel-1);
122  redo = false;
123  if (loglevel > 0) {
124  addLogEntry("MultiPhaseEquil solver succeeded.");
125  }
126  retn = nAttempts;
127  } catch (CanteraError& err) {
128  err.save();
129  if (loglevel > 0) {
130  addLogEntry("MultiPhaseEquil solver failed.");
131  }
132  if (nAttempts < 2) {
133  if (loglevel > 0) {
134  addLogEntry("Trying single phase ChemEquil solver.");
135  }
136  solver = -1;
137  } else {
138  if (loglevel > 0) {
139  endLogGroup("equilibrate");
140  }
141  throw err;
142  }
143  }
144  } else { // solver <= 0
145  /*
146  * Call the element potential solver
147  */
148  try {
149  ChemEquil e;
150  e.options.maxIterations = maxsteps;
151  e.options.relTolerance = rtol;
152  nAttempts++;
153  bool useThermoPhaseElementPotentials = true;
154  retnSub = e.equilibrate(s, XY, useThermoPhaseElementPotentials,
155  loglevel-1);
156  if (retnSub < 0) {
157  if (loglevel > 0) {
158  addLogEntry("ChemEquil solver failed.");
159  }
160  if (nAttempts < 2) {
161  if (loglevel > 0) {
162  addLogEntry("Trying MultiPhaseEquil solver.");
163  }
164  solver = 1;
165  } else {
166  throw CanteraError("equilibrate",
167  "Both equilibrium solvers failed");
168  }
169  }
170  retn = nAttempts;
171  s.setElementPotentials(e.elementPotentials());
172  redo = false;
173  if (loglevel > 0) {
174  addLogEntry("ChemEquil solver succeeded.");
175  }
176  }
177 
178  catch (CanteraError& err) {
179  err.save();
180  if (loglevel > 0) {
181  addLogEntry("ChemEquil solver failed.");
182  }
183  // If ChemEquil fails, try the MultiPhase solver
184  if (solver < 0) {
185  if (loglevel > 0) {
186  addLogEntry("Trying MultiPhaseEquil solver.");
187  }
188  solver = 1;
189  } else {
190  redo = false;
191  if (loglevel > 0) {
192  endLogGroup("equilibrate");
193  }
194  throw err;
195  }
196  }
197  }
198  } // while (redo)
199  /*
200  * We are here only for a success
201  */
202  if (loglevel > 0) {
203  endLogGroup("equilibrate");
204  }
205  return retn;
206 }
207 }
Class ChemEquil implements a chemical equilibrium solver for single-phase solutions.
Definition: ChemEquil.h:95
EquilOpt options
Options controlling how the calculation is carried out.
Definition: ChemEquil.h:150
doublereal equilibrate(int XY, doublereal err=1.0e-9, int maxsteps=1000, int maxiter=200, int loglevel=-99)
Set the mixture to a state of chemical equilibrium.
Definition: MultiPhase.cpp:581
doublereal relTolerance
Relative tolerance.
Definition: ChemEquil.h:38
void beginLogGroup(const std::string &title, int loglevel)
Create a new group for log messages.
Definition: global.cpp:90
string lastErrorMessage()
Retrieve the last error message in a string.
Definition: global.cpp:166
This file contains the definition of some high level general equilibration routines.
void addPhase(ThermoPhase *p, doublereal moles)
Add a phase to the mixture.
Definition: MultiPhase.cpp:98
int equilibrate(thermo_t &s, const char *XY, bool useThermoPhaseElementPotentials=false, int loglevel=0)
Definition: ChemEquil.cpp:394
doublereal equilibrate(MultiPhase &s, const char *XY, doublereal tol, int maxsteps, int maxiter, int loglevel)
Equilibrate a MultiPhase object.
Definition: equilibrate.cpp:14
This file contains definitions for utility functions and text for modules, inputfiles, logs, textlogs, HTML_logs (see Input File Handling, Diagnostic Output, Writing messages to the screen and Writing HTML Logfiles).
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:101
A class for multiphase mixtures.
Definition: MultiPhase.h:54
std::string id() const
Return the string id for the phase.
Definition: Phase.cpp:119
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:68
void setElementPotentials(const vector_fp &lambda)
Stores the element potentials in the ThermoPhase object.
void endLogGroup(const std::string &title)
Close the current group of log messages.
Definition: global.cpp:115
Interface class for the vcsnonlinear solver.
void addLogEntry(const std::string &tag, const std::string &value)
Add an entry to an HTML log file.
Definition: global.cpp:95
void init()
Process phases and build atomic composition array.
Definition: MultiPhase.cpp:173
int maxIterations
Maximum number of iterations.
Definition: ChemEquil.h:40
int _equilflag(const char *xy)
map property strings to integers
Definition: ChemEquil.cpp:29
int vcs_equilibrate(thermo_t &s, const char *XY, int estimateEquil, int printLvl, int solver, doublereal rtol, int maxsteps, int maxiter, int loglevel)
Set a single-phase chemical solution to chemical equilibrium.
Chemical equilibrium.
void save()
Function to put this error onto Cantera's error stack.