research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Split functions in a file into 2 files, in order to clarify roles.

Browse source
This commit is contained in:
Harry Moffat 2008-01-20 20:29:48 +00:00
parent 4e56b6f7e6
commit c2b8b1d7c8
6 changed files with 375 additions and 174 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

2
Cantera/src/equil/Makefile.in
View file

@ -69,7 +69,7 @@ VCSNONIDEAL_OBJ = vcs_solve_TP.o vcs_VolPhase.o vcs_solve.o vcs_prob.o \
vcs_nasa_poly.o vcs_nondim.o vcs_Exception.o \
vcs_funcVtot.o vcs_inest.o vcs_rearrange.o \
vcs_root1d.o vcs_rxnadj.o vcs_timer_generic.o \
vcs_SpeciesProperties.o \
vcs_SpeciesProperties.o vcs_equilibrate.o \
vcs_prep.o vcs_species_thermo.o vcs_Gibbs.o \
$(DALT_OBJ)

4
Cantera/src/equil/equilibrate.cpp
View file

@ -133,8 +133,8 @@ namespace Cantera {
m->addPhase(&s, 1.0);
m->init();
nAttempts++;
(void) vcs_equilibrate(*m, XY, estimateEquil, printLvlSub,
rtol, maxsteps, maxiter, loglevel-1);
vcs_equilibrate(*m, XY, estimateEquil, printLvlSub, solver,
rtol, maxsteps, maxiter, loglevel-1);
redo = false;
if (loglevel > 0)
addLogEntry("VCSnonideal solver succeeded.");

168
Cantera/src/equil/vcs_MultiPhaseEquil.cpp
View file

@ -1,9 +1,8 @@
/**
* @file vcs_MultiPhaseEquil.cpp
* Driver routines for equilibrium solvers
* Driver routine for the VCSnonideal equilibrium solver package
*/
/*
*
* $Id$
*/
/*
@ -11,6 +10,7 @@
* Contract DE-AC04-94AL85000 with Sandia Corporation, the
* U.S. Government retains certain rights in this software.
*/
#include "vcs_MultiPhaseEquil.h"
#include "vcs_prob.h"
#include "vcs_internal.h"
@ -1618,171 +1618,7 @@ namespace Cantera {
return VCS_SUCCESS;
}
/*
* Set a single-phase chemical solution to chemical equilibrium.
* This is a convenience function that uses one or the other of
* the two chemical equilibrium solvers.
*
* @param s The object to set to an equilibrium state
*
* @param XY An integer specifying the two properties to be held
* constant.
*
* @param estimateEquil Boolean indicating whether the solver
* should estimate its own initial condition.
* If false, the initial mole fraction vector
* in the %ThermoPhase object is used as the
* initial condition.
*
* @param printLvl Determines the amount of printing that
* gets sent to stdout from the vcs package
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param solver The equilibrium solver to use. If solver = 0,
* the ChemEquil solver will be used, and if
* solver = 1, the vcs_MultiPhaseEquil solver will
* be used (slower than ChemEquil,
* but more stable). If solver < 0 (default, then
* ChemEquil will be tried first, and if it fails
* vcs_MultiPhaseEquil will be tried.
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or
* pressure iterations to take when T and/or P is
* not held fixed.
*
* @param loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose
* messages are written as loglevel increases. The
* messages are written to a file in HTML format for viewing
* in a web browser. @see HTML_logs
*/
int vcs_equilibrate(thermo_t& s, const char* XY,
bool estimateEquil, int printLvl,
int solver,
doublereal rtol, int maxsteps, int maxiter,
int loglevel) {
MultiPhase* m = 0;
bool redo = true;
int retn = 1;
beginLogGroup("equilibrate", loglevel);
addLogEntry("Single-phase equilibrate function");
{
beginLogGroup("arguments");
addLogEntry("phase",s.id());
addLogEntry("XY",XY);
addLogEntry("solver",solver);
addLogEntry("rtol",rtol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
}
if (solver > 0) {
m = new MultiPhase;
try {
/*
* Set the kmoles of the phase to 1.0, arbitrarily.
* It actually doesn't matter.
*/
m->addPhase(&s, 1.0);
m->init();
retn = vcs_equilibrate(*m, XY, estimateEquil, printLvl,
rtol, maxsteps, maxiter, loglevel);
redo = false;
addLogEntry("MultiPhaseEquil solver succeeded.");
delete m;
}
catch (CanteraError err) {
addLogEntry("MultiPhaseEquil solver failed.");
delete m;
throw err;
}
} else {
throw CanteraError("vcs_equilibrate",
"ChemEquil not implemented in this interface yet");
}
/*
* We are here only for a success
*/
endLogGroup("equilibrate");
return retn;
}
int vcs_equilibrate(MultiPhase& s, const char* XY,
bool estimateEquil, int printLvl,
doublereal tol, int maxsteps, int maxiter,
int loglevel) {
int ixy = _equilflag(XY);
int retn = vcs_equilibrate_1(s, ixy, estimateEquil, printLvl,
tol, maxsteps, maxiter, loglevel);
return retn;
};
/*
* Set a multiphase mixture to a state of chemical equilibrium.
* This is the top-level driver for multiphase equilibrium. It
* doesn't do much more than call the equilibrate method of class
* MultiPhase, except that it adds some messages to the logfile,
* if loglevel is set > 0.
*
* @return Returns the value of 1 if successful.
*
* @ingroup equil
*/
int vcs_equilibrate_1(MultiPhase& s, int ixy,
bool estimateEquil, int printLvl,
doublereal tol, int maxsteps, int maxiter, int loglevel) {
static int counter = 0;
int retn = 1;
beginLogGroup("equilibrate",loglevel);
addLogEntry("multiphase equilibrate function");
beginLogGroup("arguments");
addLogEntry("XY",ixy);
addLogEntry("tol",tol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
int printLvlSub = MAX(0, printLvl-1);
s.init();
try {
vcs_MultiPhaseEquil *eqsolve = new vcs_MultiPhaseEquil(&s, printLvlSub);
int err = eqsolve->equilibrate(ixy, estimateEquil, printLvlSub,
tol, maxsteps, maxiter);
if (err != 0) {
retn = 0;
}
addLogEntry("Success. Error", err);
endLogGroup("equilibrate");
// hard code a csv output file.
if (printLvl > 0) {
string reportFile = "vcs_equilibrate_res.csv";
if (counter > 0) {
reportFile = "vcs_equilibrate_res_" + int2str(counter) + ".csv";
}
eqsolve->reportCSV(reportFile);
counter++;
}
delete eqsolve;
}
catch (CanteraError e) {
addLogEntry("Failure.", lastErrorMessage());
endLogGroup("equilibrate");
throw e;
}
return retn;
}
}

16
Cantera/src/equil/vcs_MultiPhaseEquil.h
View file

@ -111,6 +111,10 @@ namespace Cantera {
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param solver Determines which solver is used.
* - 1 MultiPhaseEquil solver
* - 2 VCSnonideal Solver (default)
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
@ -129,6 +133,7 @@ namespace Cantera {
*/
int vcs_equilibrate(MultiPhase& s, const char* XY,
bool estimateEquil = false, int printLvl = 0,
int solver = 2,
doublereal rtol = 1.0e-9, int maxsteps = 1000,
int maxiter = 100, int loglevel = -99);
@ -157,6 +162,10 @@ namespace Cantera {
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param solver Determines which solver is used.
* - 1 MultiPhaseEquil solver
* - 2 VCSnonideal Solver (default)
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
@ -174,9 +183,10 @@ namespace Cantera {
* @ingroup equilfunctions
*/
int vcs_equilibrate_1(MultiPhase& s, int ixy,
bool estimateEquil = false, int printLvl = 0,
doublereal rtol = 1.0e-9, int maxsteps = 1000,
int maxiter = 100, int loglevel = -99);
bool estimateEquil = false, int printLvl = 0,
int solver = 2,
doublereal rtol = 1.0e-9, int maxsteps = 1000,
int maxiter = 100, int loglevel = -99);
//! Cantera's Interface to the Multiphase chemical equilibrium solver.
/*!

357
Cantera/src/equil/vcs_equilibrate.cpp Normal file
View file

@ -0,0 +1,357 @@
/**
* @file vcs_equilibrate.cpp
* Driver routines for equilibrium solvers
*/
/*
*
* $Id$
*/
/*
* Copywrite (2006) Sandia Corporation. Under the terms of
* Contract DE-AC04-94AL85000 with Sandia Corporation, the
* U.S. Government retains certain rights in this software.
*/
#include "vcs_MultiPhaseEquil.h"
#include "vcs_prob.h"
#include "vcs_internal.h"
#include "vcs_VolPhase.h"
#include "vcs_species_thermo.h"
#include "vcs_SpeciesProperties.h"
#include "vcs_VolPhase.h"
#include "vcs_nasa_poly.h"
#include "vcs_solve.h"
#include "equil.h"
#include "ct_defs.h"
#include "mix_defs.h"
#include "speciesThermoTypes.h"
#ifdef WITH_IDEAL_SOLUTIONS
#include "IdealSolidSolnPhase.h"
#endif
#ifdef WITH_ELECTROLYTES
#include "IdealMolalSoln.h"
#endif
#include "ChemEquil.h"
#include <string>
#include <vector>
using namespace Cantera;
using namespace std;
namespace Cantera {
/*
* Set a single-phase chemical solution to chemical equilibrium.
* This is a convenience function that uses one or the other of
* the two chemical equilibrium solvers.
*
* @param s The object to set to an equilibrium state
*
* @param XY An integer specifying the two properties to be held
* constant.
*
* @param estimateEquil Boolean indicating whether the solver
* should estimate its own initial condition.
* If false, the initial mole fraction vector
* in the %ThermoPhase object is used as the
* initial condition.
*
* @param printLvl Determines the amount of printing that
* gets sent to stdout from the vcs package
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param solver The equilibrium solver to use. If solver = 0,
* the ChemEquil solver will be used, and if
* solver = 1, the vcs_MultiPhaseEquil solver will
* be used (slower than ChemEquil,
* but more stable). If solver < 0 (default, then
* ChemEquil will be tried first, and if it fails
* vcs_MultiPhaseEquil will be tried.
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or
* pressure iterations to take when T and/or P is
* not held fixed.
*
* @param loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose
* messages are written as loglevel increases. The
* messages are written to a file in HTML format for viewing
* in a web browser. @see HTML_logs
*/
int vcs_equilibrate(thermo_t& s, const char* XY,
bool estimateEquil, int printLvl,
int solver,
doublereal rtol, int maxsteps, int maxiter,
int loglevel) {
MultiPhase* m = 0;
int retn = 1;
int retnSub = 0;
beginLogGroup("equilibrate", loglevel);
// retry:
addLogEntry("Single-phase equilibrate function");
{
beginLogGroup("arguments");
addLogEntry("phase",s.id());
addLogEntry("XY",XY);
addLogEntry("solver",solver);
addLogEntry("rtol",rtol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
}
if (solver == 2) {
m = new MultiPhase;
try {
/*
* Set the kmoles of the phase to 1.0, arbitrarily.
* It actually doesn't matter.
*/
m->addPhase(&s, 1.0);
m->init();
retn = vcs_equilibrate(*m, XY, estimateEquil, printLvl, solver,
rtol, maxsteps, maxiter, loglevel);
addLogEntry("MultiPhaseEquil solver succeeded.");
delete m;
}
catch (CanteraError &err) {
addLogEntry("MultiPhaseEquil solver failed.");
delete m;
throw err;
}
} else if (solver == 1) {
m = new MultiPhase;
try {
m->addPhase(&s, 1.0);
m->init();
(void) equilibrate(*m, XY, rtol, maxsteps, maxiter, loglevel-1);
if (loglevel > 0)
addLogEntry("MultiPhaseEquil solver succeeded.");
delete m;
retn = 1;
}
catch (CanteraError &err) {
if (loglevel > 0) {
addLogEntry("MultiPhaseEquil solver failed.");
}
delete m;
throw err;
}
} else if (solver == 0) {
ChemEquil *e = new ChemEquil;
try {
e->options.maxIterations = maxsteps;
e->options.relTolerance = rtol;
retnSub = e->equilibrate(s,XY,loglevel-1);
if (retnSub < 0) {
if (loglevel > 0) {
addLogEntry("ChemEquil solver failed.");
}
delete e;
throw CanteraError("equilibrate",
"ChemEquil equilibrium solver failed");
}
retn = 1;
s.setElementPotentials(e->elementPotentials());
delete e;
if (loglevel > 0) {
addLogEntry("ChemEquil solver succeeded.");
}
}
catch (CanteraError &err) {
if (loglevel > 0) {
addLogEntry("ChemEquil solver failed.");
}
delete e;
throw err;
}
} else {
throw CanteraError("vcs_equilibrate",
"unknown solver");
}
/*
* We are here only for a success
*/
endLogGroup("equilibrate");
return retn;
}
// Set a multi-phase chemical solution to chemical equilibrium.
/*
* This function uses the vcs_MultiPhaseEquil interface to the
* vcs solver.
* The function uses the element abundance vector that is
* currently consistent with the composition within the phases
* themselves. Two other thermodynamic quantities, determined by the
* XY string, are held constant during the equilibration.
*
* @param s The object to set to an equilibrium state
*
* @param XY A character string specifying the two properties to
* be held constant
*
* @param estimateEquil Boolean indicating whether the solver
* should estimate its own initial condition.
* If false, the initial mole fraction vector
* in the %ThermoPhase object is used as the
* initial condition.
*
* @param printLvl Determines the amount of printing that
* gets sent to stdout from the vcs package
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or
* pressure iterations to take when T and/or P is
* not held fixed.
*
* @param loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose
* messages are written as loglevel increases. The
* messages are written to a file in HTML format for viewing
* in a web browser. @see HTML_logs
*
* @ingroup equilfunctions
*/
int vcs_equilibrate(MultiPhase& s, const char* XY,
bool estimateEquil, int printLvl, int solver,
doublereal tol, int maxsteps, int maxiter,
int loglevel) {
int ixy = _equilflag(XY);
int retn = vcs_equilibrate_1(s, ixy, estimateEquil, printLvl, solver,
tol, maxsteps, maxiter, loglevel);
return retn;
};
// Set a multi-phase chemical solution to chemical equilibrium.
/*
* This function uses the vcs_MultiPhaseEquil interface to the
* vcs solver.
* The function uses the element abundance vector that is
* currently consistent with the composition within the phases
* themselves. Two other thermodynamic quantities, determined by the
* XY string, are held constant during the equilibration.
*
* @param s The object to set to an equilibrium state
*
* @param XY An integer specifying the two properties to be held
* constant.
*
* @param estimateEquil Boolean indicating whether the solver
* should estimate its own initial condition.
* If false, the initial mole fraction vector
* in the %ThermoPhase object is used as the
* initial condition.
*
* @param printLvl Determines the amount of printing that
* gets sent to stdout from the vcs package
* (Note, you may have to compile with debug
* flags to get some printing).
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or
* pressure iterations to take when T and/or P is
* not held fixed.
*
* @param loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose
* messages are written as loglevel increases. The
* messages are written to a file in HTML format for viewing
* in a web browser. @see HTML_logs
*
* @ingroup equilfunctions
*/
int vcs_equilibrate_1(MultiPhase& s, int ixy,
bool estimateEquil, int printLvl, int solver,
doublereal tol, int maxsteps, int maxiter, int loglevel) {
static int counter = 0;
int retn = 1;
beginLogGroup("equilibrate",loglevel);
addLogEntry("multiphase equilibrate function");
beginLogGroup("arguments");
addLogEntry("XY",ixy);
addLogEntry("tol",tol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
int printLvlSub = MAX(0, printLvl-1);
s.init();
if (solver == 2) {
try {
vcs_MultiPhaseEquil *eqsolve = new vcs_MultiPhaseEquil(&s, printLvlSub);
int err = eqsolve->equilibrate(ixy, estimateEquil, printLvlSub,
tol, maxsteps, maxiter);
if (err != 0) {
retn = 0;
}
addLogEntry("Success. Error", err);
endLogGroup("equilibrate");
// hard code a csv output file.
if (printLvl > 0) {
string reportFile = "vcs_equilibrate_res.csv";
if (counter > 0) {
reportFile = "vcs_equilibrate_res_" + int2str(counter) + ".csv";
}
eqsolve->reportCSV(reportFile);
counter++;
}
delete eqsolve;
}
catch (CanteraError &e) {
addLogEntry("Failure.", lastErrorMessage());
endLogGroup("equilibrate");
throw e;
}
} else if (solver == 1) {
if (ixy == TP || ixy == HP || ixy == SP || ixy == TV) {
try {
double err = s.equilibrate(ixy, tol, maxsteps, maxiter, loglevel);
if (loglevel > 0) {
addLogEntry("Success. Error",err);
endLogGroup("equilibrate");
}
return 0;
}
catch (CanteraError &e) {
if (loglevel > 0) {
addLogEntry("Failure.",lastErrorMessage());
endLogGroup("equilibrate");
}
throw e;
}
}
else {
if (loglevel > 0) {
addLogEntry("multiphase equilibrium can be done only for TP, HP, SP, or TV");
endLogGroup("equilibrate");
}
throw CanteraError("equilibrate","unsupported option");
//return -1.0;
}
} else {
throw CanteraError("vcs_equilibrate_1", "unknown solver");
}
return retn;
}
}

2
Cantera/src/equil/vcs_solve.h
View file

@ -779,9 +779,7 @@ public:
VCS_COUNTERS *m_VCount;
#ifdef DEBUG_MODE
int vcs_debug_print_lvl;
#endif
//! Units for the chemical potential data:
/*!