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

removed include files that were suppose to have been moved.

Browse source
This commit is contained in:
Harry Moffat 2012-04-05 19:17:01 +00:00
parent 25ba149aab
commit 458fbf579d
8 changed files with 3 additions and 4252 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

127
src/equil/vcs_DoubleStarStar.h
View file

@ -1,127 +0,0 @@
/**
* @file vcs_DoubleStarStar.h
*
* Header file for class DoubleStarStar
*/
#ifndef VCS_DOUBLESTARSTAR_H
#define VCS_DOUBLESTARSTAR_H
#include <vector>
namespace VCSnonideal
{
using std::size_t;
//! A class for 2D double arrays stored in column-major
//! (Fortran-compatible) form.
/*!
* In this form, the data entry for an n row, m col
* matrix is
* index = i + (n-1) * j
* where
* Matrix[j][i]
* i = row
* j = column
* The way this is instantiated is via the constructor:
* DoubleStarStar Dmatrix(mcol, mrow);
*
* The way this is referenced is via the notation:
* Dmatrix[icol][irow]
*/
class DoubleStarStar
{
public:
//! Default constructor. Create an empty array.
DoubleStarStar();
//! Constructor.
/*!
* Create an \c nrow by \c mcol double array, and initialize
* all elements to \c v.
*
* @param mcol Number of columns
* @param nrow Number of rows
*/
DoubleStarStar(size_t mcol, size_t nrow, double v = 0.0);
//! copy constructor
/*!
* @param y object to be copied
*/
DoubleStarStar(const DoubleStarStar& y);
/// assignment operator
/*!
* @param y object to be copied
*/
DoubleStarStar& operator=(const DoubleStarStar& y);
//! Resize the array, and fill the new entries with 'v'
/*!
* @param mrow This is the number of columns in the new matrix
* @param ncol This is the number of rows
* @param v Default fill value -> defaults to zero.
*/
void resize(size_t mcol, size_t nrow, double v = 0.0);
//! Pointer to the top of the column
/*!
* @param jcol This is the jth column
*
* @return returns the pointer to the top of the jth column
*/
double* operator[](size_t jcol);
//! Returns a const Pointer to the top of the jth column
/*!
* @param jcol This is the jth column
*
* @return returns the pointer to the top of the jth column
*/
const double* operator[](size_t jcol) const;
//! Returns a double ** pointer to the base address
/*!
* This is the second way to get to the data
* This returns a double ** which can later be used in
* Dmatrix[icol][irow] notation to get to the data
*/
double* const* baseDataAddr();
//! Returns a const double ** pointer to the base address
/*!
* This is the second way to get to the data
* This returns a double ** which can later be used in
* Dmatrix[icol][irow] notation to get to the data
*/
double const* const* constBaseDataAddr() const;
//! Number of rows
size_t nRows() const;
//! Number of columns
size_t nColumns() const;
private:
//! Storage area
std::vector<double> m_data;
//! Vector of addresses for the top of the columns
/*!
* Length = mcol
*/
std::vector<double*> m_colAddr;
//! number of rows
size_t m_nrows;
//! number of columns
size_t m_ncols;
};
}
#endif

4
src/equil/vcs_IntStarStar.cpp
View file

@ -3,12 +3,12 @@
*
* Header file for class IntStarStar
*/
#include "vcs_IntStarStar.h"
#include "cantera/equil/vcs_IntStarStar.h"
namespace VCSnonideal
{
//!Default constructor. Create an empty array.
//Default constructor. Create an empty array.
IntStarStar::IntStarStar() :
m_nrows(0),
m_ncols(0)

101
src/equil/vcs_IntStarStar.h
View file

@ -1,101 +0,0 @@
/**
* @file vcs_IntStarStar.h
*
* Header file for class IntStarStar
*/
#ifndef VCS_INTSTARSTAR_H
#define VCS_INTSTARSTAR_H
#include <vector>
namespace VCSnonideal
{
using std::size_t;
//! A class for 2D int arrays stored in column-major
//! (Fortran-compatible) form.
/*!
* In this form, the data entry for an n row, m col
* matrix is
* index = i + (n-1) * j
* where
* Matrix[j][i]
* i = row
* j = column
*/
class IntStarStar
{
public:
//! Default constructor. Create an empty array.
IntStarStar();
//! Constructor.
/*!
* Create an \c nrow by \c mcol int array, and initialize
* all elements to \c v.
*
* @param mcol Number of columns
* @param nrow Number of rows
*/
IntStarStar(size_t mcol, size_t nrow, int v = 0);
//! copy constructor
IntStarStar(const IntStarStar& y);
/// assignment operator
IntStarStar& operator=(const IntStarStar& y);
//! Resize the array, and fill the new entries with 'v'
/*!
* @param mcol This is the number of columns in the new matrix
* @param nrow This is the number of rows
* @param v Default fill value -> defaults to zero.
*/
void resize(size_t mcol, size_t nrow, int v = 0);
//! Pointer to the top of the column
/*!
* @param jcol Pointer to the top of the jth column
*/
int* operator[](size_t jcol);
//! Pointer to the top of the column
/*!
* @param j Pointer to the top of the jth column
*/
const int* operator[](size_t jcol) const;
//! Returns a int ** pointer to the base address
/*!
* This is the second way to get to the data
* This returns a int ** which can later be used in
* Imatrix[icol][irow] notation to get to the data
*/
int* const* baseDataAddr();
//! Number of rows
size_t nRows() const;
//! Number of columns
size_t nColumns() const;
private:
//! Storage area
std::vector<int> m_data;
std::vector<int*> m_colAddr;
//! number of rows
size_t m_nrows;
//! number of columns
size_t m_ncols;
};
}
#endif

1038
src/equil/vcs_VolPhase.h
View file

File diff suppressed because it is too large Load diff

528
src/equil/vcs_internal.h
View file

@ -1,528 +0,0 @@
/**
* @file vcs_internal.h
* Internal declarations for the VCSnonideal package
*/
/*
* Copyright (2005) Sandia Corporation. Under the terms of
* Contract DE-AC04-94AL85000 with Sandia Corporation, the
* U.S. Government retains certain rights in this software.
*/
#ifndef _VCS_INTERNAL_H
#define _VCS_INTERNAL_H
#include <cstring>
#include "cantera/equil/vcs_defs.h"
#include "vcs_DoubleStarStar.h"
#include "vcs_Exception.h"
#include "cantera/base/global.h"
namespace VCSnonideal
{
using Cantera::npos;
//! Points to the data in a std::vector<> object
#define VCS_DATA_PTR(vvv) (&(vvv[0]))
//! define this Cantera function to replace printf
/*!
* We can replace this with printf easily
*/
#define plogf Cantera::writelogf
//! define this Cantera function to replace cout << endl;
/*!
* We use this to place an endl in the log file, and
* ensure that the IO buffers are flushed.
*/
#define plogendl() Cantera::writelogendl()
//! Global hook for turning on and off time printing.
/*!
* Default is to allow printing. But, you can assign this to zero
* globally to turn off all time printing.
* This is helpful for test suite purposes where you are interested
* in differences in text files.
*/
extern int vcs_timing_print_lvl;
/*
* Forward references
*/
class VCS_SPECIES_THERMO;
class VCS_PROB;
//! Amount of extra printing that is done while in debug mode.
/*!
* 0 -> none
* 1 -> some
* 2 -> alot (default)
* 3 -> everything
*/
//! Class to keep track of time and iterations
/*!
* class keeps all of the counters together.
*/
class VCS_COUNTERS
{
public:
//! Total number of iterations in the main loop
//! of vcs_TP() to solve for thermo equilibrium
int T_Its;
//! Current number of iterations in the main loop
//! of vcs_TP() to solve for thermo equilibrium
int Its;
//! Total number of optimizations of the
//! components basis set done
int T_Basis_Opts;
//! number of optimizations of the components basis set done
int Basis_Opts;
//! Current number of times the initial thermo
//! equilibrium estimator has been called
int T_Calls_Inest;
//! Current number of calls to vcs_TP
int T_Calls_vcs_TP;
//! Current time spent in vcs_TP
double T_Time_vcs_TP;
//! Current time spent in vcs_TP
double Time_vcs_TP;
//! Total Time spent in basopt
double T_Time_basopt;
//! Current Time spent in basopt
double Time_basopt;
//! Time spent in initial estimator
double T_Time_inest;
//! Time spent in the vcs suite of programs
double T_Time_vcs;
};
//! Returns the value of the gas constant in the units specified by parameter
/*!
* @param mu_units Specifies the units.
* - VCS_UNITS_KCALMOL: kcal gmol-1 K-1
* - VCS_UNITS_UNITLESS: 1.0 K-1
* - VCS_UNITS_KJMOL: kJ gmol-1 K-1
* - VCS_UNITS_KELVIN: 1.0 K-1
* - VCS_UNITS_MKS: joules kmol-1 K-1 = kg m2 s-2 kmol-1 K-1
*/
double vcsUtil_gasConstant(int mu_units);
//! Invert an n x n matrix and solve m rhs's
/*!
* Solve a square matrix with multiple right hand sides
*
* \f[
* C X + B = 0;
* \f]
*
* This routine uses Gauss elimination and is optimized for the solution
* of lots of rhs's. A crude form of row pivoting is used here.
* The matrix C is destroyed during the solve.
*
* @return The solution x[] is returned in the matrix <I>B</I>.
* Routine returns an integer representing success:
* - 1 : Matrix is singular
* - 0 : solution is OK
*
*
* @param c Matrix to be inverted. c is in fortran format, i.e., rows
* are the inner loop. Row numbers equal to idem.
* c[i+j*idem] = c_i_j = Matrix to be inverted:
* - i = row number
* - j = column number
*
* @param idem number of row dimensions in c
* @param n Number of rows and columns in c
* @param b Multiple RHS. Note, b is actually the negative of
* most formulations. Row numbers equal to idem.
* b[i+j*idem] = b_i_j = vectors of rhs's:
* - i = row number
* - j = column number
* (each column is a new rhs)
* @param m number of rhs's
*/
int vcsUtil_mlequ(double* c, size_t idem, size_t n, double* b, size_t m);
//! Invert an n x n matrix and solve m rhs's
/*!
* Solve a square matrix with multiple right hand sides
*
* \f[
* C X + B = 0;
* \f]
*
* This routine uses Gauss-Jordan elimination and is optimized for the solution
* of lots of rhs's. Full row and column pivoting is used here. It's been
* shown to be necessary in at least one case.
* The matrix C is destroyed during the solve.
*
* @return The solution x[] is returned in the matrix <I>B</I>.
* Routine returns an integer representing success:
* - 1 : Matrix is singular
* - 0 : solution is OK
*
* @param c Matrix to be inverted. c is in fortran format, i.e., rows
* are the inner loop. Row numbers equal to idem.
* c[i+j*idem] = c_i_j = Matrix to be inverted:
* - i = row number
* - j = column number
*
* @param idem number of row dimensions in c
* @param n Number of rows and columns in c
* @param b Multiple RHS. Note, b is actually the negative of
* most formulations. Row numbers equal to idem.
* b[i+j*idem] = b_i_j = vectors of rhs's:
* - i = row number
* - j = column number
* (each column is a new rhs)
* @param m number of rhs's
*/
int vcsUtil_gaussj(double* c, size_t idem, size_t n, double* b, size_t m);
//! Definition of the function pointer for the root finder
/*!
* see vcsUtil_root1d for a definition of how to use this.
*/
typedef double(*VCS_FUNC_PTR)(double xval, double Vtarget,
int varID, void* fptrPassthrough,
int* err);
//! One dimensional root finder
/*!
*
* This root finder will find the root of a one dimensional
* equation
*
* \f[
* f(x) = 0
* \f]
* where x is a bounded quantity: \f$ x_{min} < x < x_max \f$
*
* The functional to be minimized must have the following call
* structure:
*
* @verbatim
typedef double (*VCS_FUNC_PTR)(double xval, double Vtarget,
int varID, void *fptrPassthrough,
int *err); @endverbatim
*
* xval is the current value of the x variable. Vtarget is the
* requested value of f(x), usually 0. varID is an integer
* that is passed through. fptrPassthrough is a void pointer
* that is passed through. err is a return error indicator.
* err = 0 is the norm. anything else is considered a fatal
* error.
* The return value of the function is the current value of
* f(xval).
*
* @param xmin Minimum permissible value of the x variable
* @param xmax Maximum permissible value of the x parameter
* @param itmax Maximum number of iterations
* @param func function pointer, pointing to the function to be
* minimized
* @param fptrPassthrough Pointer to void that gets passed through
* the rootfinder, unchanged, to the func.
* @param FuncTargVal Target value of the function. This is usually set
* to zero.
* @param varID Variable ID. This is usually set to zero.
* @param xbest Pointer to the initial value of x on input. On output
* This contains the root value.
* @param printLvl Print level of the routine.
*
*
* Following is a nontrial example for vcs_root1d() in which the position of a
* cylinder floating on the water is calculated.
*
* @verbatim
#include <cmath>
#include <cstdlib>
#include "equil/vcs_internal.h"
const double g_cgs = 980.;
const double mass_cyl = 0.066;
const double diam_cyl = 0.048;
const double rad_cyl = diam_cyl / 2.0;
const double len_cyl = 5.46;
const double vol_cyl = Pi * diam_cyl * diam_cyl / 4 * len_cyl;
const double rho_cyl = mass_cyl / vol_cyl;
const double rho_gas = 0.0;
const double rho_liq = 1.0;
const double sigma = 72.88;
// Contact angle in radians
const double alpha1 = 40.0 / 180. * Pi;
double func_vert(double theta1, double h_2, double rho_c) {
double f_grav = - Pi * rad_cyl * rad_cyl * rho_c * g_cgs;
double tmp = rad_cyl * rad_cyl * g_cgs;
double tmp1 = theta1 + sin(theta1) * cos(theta1) - 2.0 * h_2 / rad_cyl * sin(theta1);
double f_buoy = tmp * (Pi * rho_gas + (rho_liq - rho_gas) * tmp1);
double f_sten = 2 * sigma * sin(theta1 + alpha1 - Pi);
double f_net = f_grav + f_buoy + f_sten;
return f_net;
}
double calc_h2_farfield(double theta1) {
double rhs = sigma * (1.0 + cos(alpha1 + theta1));
rhs *= 2.0;
rhs = rhs / (rho_liq - rho_gas) / g_cgs;
double sign = -1.0;
if (alpha1 + theta1 < Pi) sign = 1.0;
double res = sign * sqrt(rhs);
double h2 = res + rad_cyl * cos(theta1);
return h2;
}
double funcZero(double xval, double Vtarget, int varID, void *fptrPassthrough, int *err) {
double theta = xval;
double h2 = calc_h2_farfield(theta);
double fv = func_vert(theta, h2, rho_cyl);
return fv;
}
int main () {
double thetamax = Pi;
double thetamin = 0.0;
int maxit = 1000;
int iconv;
double thetaR = Pi/2.0;
int printLvl = 4;
iconv = VCSnonideal::vcsUtil_root1d(thetamin, thetamax, maxit,
funcZero,
(void *) 0, 0.0, 0,
&thetaR, printLvl);
printf("theta = %g\n", thetaR);
double h2Final = calc_h2_farfield(thetaR);
printf("h2Final = %g\n", h2Final);
return 0;
} @endverbatim
*
*/
int vcsUtil_root1d(double xmin, double xmax, size_t itmax, VCS_FUNC_PTR func,
void* fptrPassthrough,
double FuncTargVal, int varID, double* xbest,
int printLvl = 0);
//! Returns the system wall clock time in seconds
/*!
* @return time in seconds.
*/
double vcs_second();
//! This define turns on using memset and memcpy. I have not run into
//! any systems where this is a problem. It's the fastest way to do
//! low lvl operations where applicable. There are alternative routines
//! available if this ever fails.
#define USE_MEMSET
#ifdef USE_MEMSET
//! Zero a double vector
/*!
* @param vec_to vector of doubles
* @param length length of the vector to zero.
*/
inline void vcs_dzero(double* const vec_to, const size_t length)
{
(void) memset((void*) vec_to, 0, length * sizeof(double));
}
//! Zero an int vector
/*!
* @param vec_to vector of ints
* @param length length of the vector to zero.
*/
inline void vcs_izero(int* const vec_to, const size_t length)
{
(void) memset((void*) vec_to, 0, length * sizeof(int));
}
//! Copy a double vector
/*!
* @param vec_to Vector to copy into. This vector must be dimensioned
* at least as large as the vec_from vector.
* @param vec_from Vector to copy from
* @param length Number of doubles to copy.
*/
inline void vcs_dcopy(double* const vec_to,
const double* const vec_from, const size_t length)
{
(void) memcpy((void*) vec_to, (const void*) vec_from,
(length) * sizeof(double));
}
//! Copy an int vector
/*!
* @param vec_to Vector to copy into. This vector must be dimensioned
* at least as large as the vec_from vector.
* @param vec_from Vector to copy from
* @param length Number of int to copy.
*/
inline void vcs_icopy(int* const vec_to,
const int* const vec_from, const size_t length)
{
(void) memcpy((void*) vec_to, (const void*) vec_from,
(length) * sizeof(int));
}
//! Zero a std double vector
/*!
* @param vec_to vector of doubles
* @param length length of the vector to zero.
*/
inline void vcs_vdzero(std::vector<double> &vec_to, const size_t length)
{
(void) memset((void*)VCS_DATA_PTR(vec_to), 0, (length) * sizeof(double));
}
//! Zero a std int vector
/*!
* @param vec_to vector of ints
* @param length length of the vector to zero.
*/
inline void vcs_vizero(std::vector<int> &vec_to, const size_t length)
{
(void) memset((void*)VCS_DATA_PTR(vec_to), 0, (length) * sizeof(int));
}
//! Copy one std double vector into another
/*!
* This is an inlined function that uses memcpy. memcpy is probably
* the fastest way to do this. This routine requires the vectors to be
* previously dimensioned appropriately. No error checking is done.
*
* @param vec_to Vector to copy into. This vector must be dimensioned
* at least as large as the vec_from vector.
* @param vec_from Vector to copy from
* @param length Number of doubles to copy.
*/
inline void vcs_vdcopy(std::vector<double> & vec_to,
const std::vector<double> & vec_from, size_t length)
{
(void) memcpy((void*)&(vec_to[0]), (const void*) &(vec_from[0]),
(length) * sizeof(double));
}
//! Copy one std integer vector into another
/*!
* This is an inlined function that uses memcpy. memcpy is probably
* the fastest way to do this. This routine requires the
*
* @param vec_to Vector to copy into. This vector must be dimensioned
* at least as large as the vec_from vector.
* @param vec_from Vector to copy from
* @param length Number of integers to copy.
*/
inline void vcs_vicopy(std::vector<int> & vec_to,
const std::vector<int> & vec_from, const int length)
{
(void) memcpy((void*)&(vec_to[0]), (const void*) &(vec_from[0]),
(length) * sizeof(int));
}
#else
extern void vcs_dzero(double* const, const int);
extern void vcs_izero(int* const , const int);
extern void vcs_dcopy(double* const, const double* const, const int);
extern void vcs_icopy(int* const, const int* const, const int);
extern void vcs_vdzero(std::vector<double> &vvv, const int len = -1);
extern void vcs_vizero(std::vector<double> &vvv, const int len = -1);
void vcs_vdcopy(std::vector<double> &vec_to,
const std::vector<double> vec_from, const int len = -1);
void vcs_vicopy(std::vector<int> &vec_to,
const std::vector<int> vec_from, const int len = -1);
#endif
//! determine the l2 norm of a vector of doubles
/*!
* @param vec vector of doubles
*
* @return Returns the l2 norm of the vector
*/
double vcs_l2norm(const std::vector<double> vec);
//! Finds the location of the maximum component in a double vector
/*!
* @param x pointer to a vector of doubles
* @param xSize pointer to a vector of doubles used as a multiplier
* to x[]
* @param j lowest index to search from
* @param n highest index to search from
* @return Return index of the greatest value on X(i) searched
* j <= i < n
*/
size_t vcs_optMax(const double* x, const double* xSize, size_t j, size_t n);
//! Returns the maximum integer in a list
/*!
* @param vector pointer to a vector of ints
* @param length length of the integer vector
*
* @return returns the max integer value in the list
*/
int vcs_max_int(const int* vector, int length);
//! Prints a line consisting of multiple occurrences of the same string
/*!
* This prints a string num times, and then terminate with a
* end of line character
*
* @param str C string that is null terminated
* @param num number of times the string is to be printed
*/
void vcs_print_line(const char* str, int num);
//! Returns a const char string representing the type of the
//! species given by the first argument
/*!
* @param speciesStatus Species status integer representing the type
* of the species.
* @param length Maximum length of the string to be returned.
* Shorter values will yield abbreviated strings.
* Defaults to a value of 100.
*/
const char* vcs_speciesType_string(int speciesStatus, int length = 100);
//! Print a string within a given space limit
/*!
* This routine limits the amount of the string that will be printed to a
* maximum of "space" characters. Printing is done to
* to Cantera's writelog() function.
*
* @param str String, which must be null terminated.
* @param space space limit for the printing.
* @param alignment Alignment of string within the space:
* - 0 centered
* - 1 right aligned
* - 2 left aligned
*/
void vcs_print_stringTrunc(const char* str, size_t space, int alignment);
//! Simple routine to check whether two doubles are equal up to
//! roundoff error
/*!
* Currently it's set to check for 10 digits of
* relative accuracy.
*
* @param d1 first double
* @param d2 second double
*
* @return returns true if the doubles are "equal" and false otherwise
*/
bool vcs_doubleEqual(double d1, double d2);
}
#endif

372
src/equil/vcs_prob.h
View file

@ -1,372 +0,0 @@
/**
* @file vcs_prob.h
* Header for the Interface class for the vcs thermo equilibrium solver package,
*/
/*
* Copyright (2005) Sandia Corporation. Under the terms of
* Contract DE-AC04-94AL85000 with Sandia Corporation, the
* U.S. Government retains certain rights in this software.
*/
#ifndef _VCS_PROB_H
#define _VCS_PROB_H
#include "vcs_DoubleStarStar.h"
#include "vcs_IntStarStar.h"
#include "cantera/equil/vcs_defs.h"
#include <vector>
#include <string>
namespace VCSnonideal
{
class vcs_VolPhase;
class VCS_SPECIES_THERMO;
//! Interface class for the vcs thermo equilibrium solver package,
//! which generally describes the problem to be solved.
/*!
* HKM add:
* HaveEstimate -> 0 no estimate, or estimate that doesn' satisfy elem
* abundances
* 1 have an estimate that satisfies elem_abund.
* 2 Have an estimate that minimizes a subproblem
* and satisfies elem abund.
* solnFound -> True, soln to current problem found and included here
* False, soln has not been found.
*/
class VCS_PROB
{
public:
//! Problem type. I.e., the identity of what is held constant.
/*!
* Currently, T and P are held constant, and this input
* is ignored
*/
int prob_type;
//! Total number of species in the problems
size_t nspecies;
//! Species number used to malloc data structures
size_t NSPECIES0;
//! Number of element constraints in the equilibrium problem
size_t ne;
//! Number of element constraints used to malloc data structures
//! involving elements
size_t NE0;
//! Number of phases in the problem
size_t NPhase;
//! Number of phases used to malloc data structures
size_t NPHASE0;
//! Vector of chemical potentials of the species
/*!
* This is a calculated output quantity
* length = number of species
* units = m_VCS_UnitsFormat;
*/
std::vector<double> m_gibbsSpecies;
//! Total number of moles of the kth species.
/*!
* This is both an input and an output variable.
* On input, this is an estimate of the mole numbers.
* The actual element abundance vector contains the problem specification.
*
* On output, this contains the solution for the total number of moles
* of the kth species.
*
* units = m_VCS_UnitsFormat
*/
std::vector<double> w;
//! Mole fraction vector
/*!
* This is a calculated vector, calculated from w[]
* length number of species.
* -> Take out? -> No, useful for storage of a quantity often needed
*/
std::vector<double> mf;
//! Element abundances for jth element
/*!
* This is input from the input file and is considered a constant from
* thereon within the vcs_solve_TP().
* units = m_VCS_UnitsFormat
*/
std::vector<double> gai;
//! Formula Matrix for the problem
/*!
* FormulaMatrix[j][kspec] = Number of elements, j, in the kspec
* species
*/
DoubleStarStar FormulaMatrix;
//! Specifies the species unknown type
/*!
* There are two types. One is the straightforward
* species, with the mole number w[k], as the
* unknown. The second is the an interfacial
* voltage where w[k] refers to the interfacial
* voltage in volts.
* These species types correspond to metalic
* electrons corresponding to electrodes.
* The voltage and other interfacial conditions
* sets up an interfacial current, which is
* set to zero in this initial treatment.
* Later we may have non-zero interfacial currents.
*/
std::vector<int> SpeciesUnknownType;
//! Temperature (Kelvin)
/*!
* Specification of the temperature for the equilibrium problem
*/
double T;
//! Pressure
/*!
* units given by m_VCS_UnitsFormat
* -> are now PA
*/
double PresPA;
//! Volume of the entire system
/*!
* units given by m_VCS_UnitsFormat
* Note, this is an output variable atm
*/
double Vol;
//! Partial Molar Volumes of species
/*!
* This is a calculated vector, calculated from w[]
* length number of species.
* -> Take out? -> No, useful for storage of a quantity often needed
*/
std::vector<double> VolPM;
//! Units for the chemical potential data, pressure data, volume,
//! and species amounts
/*!
* All internally stored quantities will have these units. Also, printed
* quantities will display in these units.
*
* Chem_Pot Pres vol moles
* ----------------------------------------------------------------------
* -1 VCS_UNITS_KCALMOL = kcal/mol atm cm**3 gmol
* 0 VCS_UNITS_UNITLESS = MU / RT -> no units atm cm**3 gmol
* 1 VCS_UNITS_KJMOL = kJ / mol atm cm**3 gmol
* 2 VCS_UNITS_KELVIN = KELVIN -> MU / R atm cm**3 gmol
* 3 VCS_UNITS_MKS = Joules / Kmol (Cantera) Pa m**3 kmol
* ----------------------------------------------------------------------
*
* see vcs_defs.h for more information
*/
int m_VCS_UnitsFormat;
//! Specification of the initial estimate method
/*!
* iest = Initial estimate: 0 user estimate
* 1 user estimate if satisifies elements
* -1 machine estimate
*/
int iest;
//! Tolerance requirement for major species
double tolmaj;
//! Tolerance requirement for minor species
double tolmin;
//! Mapping between the species and the phases
std::vector<size_t> PhaseID;
//! Vector of strings containing the species names
std::vector<std::string> SpName;
//! vector of strings containing the element names
std::vector<std::string> ElName;
//! vector of Element types
std::vector<int> m_elType;
//! Specifies whether an element constraint is active
/*!
* The default is true
* Length = nelements
*/
std::vector<int> ElActive;
//! Molecular weight of species
/*!
* WtSpecies[k] = molecular weight of species in gm/mol
*/
std::vector<double> WtSpecies;
//! Charge of each species
std::vector<double> Charge;
//! Array of phase structures
std::vector<vcs_VolPhase*> VPhaseList;
// String containing the title of the run
std::string Title;
//! Vector of pointers to thermo structures which identify the model
//! and parameters for evaluating the thermodynamic
//! functions for that particular species
std::vector<VCS_SPECIES_THERMO*> SpeciesThermo;
//! Number of iterations
/*!
* This is an output variable
*/
int m_Iterations;
//! Number of basis optimizations used
/*!
* This is an output variable
*/
int m_NumBasisOptimizations;
//! Print level for print routines
int m_printLvl;
//! Debug print lvl
int vcs_debug_print_lvl;
//! Constructor
/*!
* This constructor initializes the sizes within the object
* to parameter values.
*
* @param nsp number of species
* @param nel number of elements
* @param nph number of phases
*/
VCS_PROB(size_t nsp, size_t nel, size_t nph);
//! Destructor
~VCS_PROB();
//! Resizes all of the phase lists within the structure
/*!
* Note, this doesn't change the number of phases in the problem.
* It will change NPHASE0 if nsp is greater than NPHASE0.
*
* @param nPhase size to dimension all the phase lists to
* @param force If true, this will dimension the size to be equal to nPhase
* even if nPhase is less than the current value of NPHASE0
*/
void resizePhase(size_t nPhase, int force);
//! Resizes all of the species lists within the structure
/*!
* Note, this doesn't change the number of species in the problem.
* It will change NSPECIES0 if nsp is greater than NSPECIES0.
*
* @param nsp size to dimension all the species lists to
* @param force If true, this will dimension the size to be equal to nsp
* even if nsp is less than the current value of NSPECIES0
*/
void resizeSpecies(size_t nsp, int force);
//! Resizes all of the element lists within the structure
/*!
* Note, this doesn't change the number of element constraints in the problem.
* It will change NE0 if nel is greater than NE0.
*
* @param nel size to dimension all the elements lists
* @param force If true, this will dimension the size to be equal to nel
* even if nel is less than the current value of NEL0
*/
void resizeElements(size_t nel, int force);
//! Calculate the element abundance vector
/*!
* Calculates the element abundance vectors from the mole
* numbers
*/
void set_gai();
//! Print out the problem specification in all generality
//! as it currently exists in the VCS_PROB object
/*!
* @param print_lvl Parameter lvl for printing
* 0 - no printing
* 1 - all printing
*/
void prob_report(int print_lvl);
//! Add elements to the local element list
/*!
* This routine sorts through the elements defined in the
* vcs_VolPhase object. It then adds the new elements to
* the VCS_PROB object, and creates a global map, which is
* stored in the vcs_VolPhase object.
* Id and matching of elements is done strictly via the element name,
* with case not mattering.
*
* The routine also fills in the position of the element
* in the vcs_VolPhase object's ElGlobalIndex field.
*
* @param volPhase Object containing the phase to be added.
* The elements in this phase are parsed for
* addition to the global element list
*/
void addPhaseElements(vcs_VolPhase* volPhase);
//! This routine resizes the number of elements in the VCS_PROB object by
//! adding a new element to the end of the element list
/*!
* The element name is added. Formula vector entries ang element
* abundances for the new element are set to zero.
*
* Returns the index number of the new element.
*
* @param elNameNew New name of the element
* @param elType Type of the element
* @param elactive boolean indicating whether the element is active
*
* @return returns the index number of the new element
*/
size_t addElement(const char* elNameNew, int elType, int elactive);
//! This routines adds entries for the formula matrix for one species
/*!
* This routines adds entries for the formula matrix for this object
* for one species
*
* This object also fills in the index filed, IndSpecies, within
* the volPhase object.
*
* @param volPhase object containing the species
* @param k Species number within the volPhase k
* @param kT global Species number within this object
*
*/
size_t addOnePhaseSpecies(vcs_VolPhase* volPhase, size_t k, size_t kT);
void reportCSV(const std::string& reportFile);
//! Set the debug level
/*!
* @param vcs_debug_print_lvl input debug level
*/
void setDebugPrintLvl(int vcs_debug_print_lvl);
};
}
#endif

2083
src/equil/vcs_solve.h
View file

File diff suppressed because it is too large Load diff

2
test_problems/VCSnonideal/NaCl_equil/nacl_equil.cpp
View file

@ -4,7 +4,7 @@
#include "cantera/equilibrium.h"
#include "cantera/equil/vcs_MultiPhaseEquil.h"
#include "equil/vcs_internal.h"
#include "cantera/equil/vcs_internal.h"
#include "cantera/thermo/ThermoFactory.h"
#include "cantera/thermo/IdealGasPhase.h"