MultiPhase.h

Go to the documentation of this file.
00001 /**
00002  * @file MultiPhase.h
00003  * Headers for the \link Cantera::MultiPhase MultiPhase\endlink 
00004  * object that is used to set up multiphase equilibrium problems (see \ref equilfunctions).
00005  */
00006 
00007 /*
00008  *  $Date: 2010-01-12 16:23:08 -0500 (Tue, 12 Jan 2010) $
00009  *  $Revision: 373 $
00010  */
00011 
00012 //  Copyright 2004  California Institute of Technology
00013 
00014 
00015 #ifndef CT_MULTIPHASE_H
00016 #define CT_MULTIPHASE_H
00017 
00018 #include "ct_defs.h"
00019 #include "DenseMatrix.h"
00020 #include "ThermoPhase.h"
00021 
00022 namespace Cantera {
00023 
00024   //! A class for multiphase mixtures. The mixture can contain any
00025   //! number of phases of any type. 
00026   /*!
00027    *   This object is the basic tool used by Cantera for use in
00028    *   Multiphase equilibrium calculations. 
00029    *
00030    *   It is a container for a set of phases. Each phase has a 
00031    *   given number of kmoles. Therefore, MultiPhase may be considered
00032    *   an "extrinsic" thermodynamic object, in contrast to the ThermoPhase
00033    *   object, which is an "intrinsic" thermodynamic object.
00034    *
00035    *   MultiPhase may be considered to be "upstream" of the ThermoPhase
00036    *   objects in the sense that setting a property within MultiPhase,
00037    *   such as temperature, pressure, or species mole number,
00038    *   affects the underlying ThermoPhase object, but not the 
00039    *   other way around. 
00040    *
00041    *    All phases have the same
00042    *    temperature and pressure, and a specified number of moles for
00043    *    each phase.
00044    *    The phases do not need to have the same elements. For example,
00045    *    a mixture might consist of a gaseous phase with elements (H,
00046    *    C, O, N), a solid carbon phase containing only element C,
00047    *    etc. A master element set will be constructed for the mixture
00048    *    that is the intersection of the elements of each phase.
00049    *
00050    *   Below, reference is made to global species and global elements.
00051    *   These refer to the collective species and elements encompassing
00052    *   all of the phases tracked by the object. 
00053    *
00054    *   The global element list kept by this object is an
00055    *   intersection of the element lists of all the phases that
00056    *   comprise the MultiPhase.
00057    *
00058    *   The global species list kept by this object is a
00059    *   concatenated list of all of the species in all the phases that
00060    *   comprise the MultiPhase. The ordering of species is contiguous
00061    *   with respect to the phase id.
00062    *
00063    *  @ingroup equilfunctions
00064    */
00065   class MultiPhase {
00066 
00067   public:
00068 
00069     //! Shorthand for an index variable that can't be negative
00070     typedef size_t       index_t;
00071 
00072     //! Shorthand for a ThermoPhase
00073     typedef ThermoPhase  phase_t;
00074 
00075     //! shorthand for a 2D matrix
00076     typedef DenseMatrix  array_t;
00077 
00078     //! Shorthand for a vector of pointers to ThermoPhase's
00079     typedef std::vector<phase_t*> phase_list;
00080 
00081     //! Constructor.
00082     /*!
00083      *   The constructor takes no arguments, since
00084      *   phases are added using method addPhase().
00085      */
00086     MultiPhase();
00087 
00088     //! Destructor.
00089     /*!
00090      *  Does nothing. Class MultiPhase does not take
00091      *  "ownership" (i.e. responsibility for destroying) the
00092      *   phase objects.
00093      */
00094     virtual ~MultiPhase() {}
00095 
00096     //! Add a vector of phases to the mixture
00097     /*!
00098      * See the single addPhases command. This just does a bunch of phases
00099      * at one time
00100      *   @param phases Vector of pointers to phases
00101      *   @param phaseMoles Vector of mole numbers in each phase (kmol)
00102      */
00103     void addPhases(phase_list& phases, const vector_fp& phaseMoles);
00104 
00105     //! Add all phases present in 'mix' to this mixture.
00106     /*!
00107      *  @param mix  Add all of the phases in another MultiPhase
00108      *              object to the current object.
00109      */
00110     void addPhases(MultiPhase& mix);
00111 
00112     //! Add a phase to the mixture.
00113     /*!
00114      *  This function must be called befure the init() function is called,
00115      *  which serves to freeze the MultiPhase.
00116      *
00117      *  @param p pointer to the phase object
00118      *  @param moles total number of moles of all species in this phase
00119      */
00120     void addPhase(phase_t* p, doublereal moles);
00121 
00122     /// Number of elements.
00123     int nElements() const { return int(m_nel); }
00124 
00125     //! Returns the string name of the global element \a m.
00126     /*!
00127      *  @param m index of the global element
00128      */
00129     std::string elementName(int m) const;
00130 
00131     //! Returns the index of the element with name \a name.
00132     /*!
00133      * @param name   String name of the global element
00134      */
00135     int elementIndex(std::string name) const;
00136 
00137     //! Number of species, summed over all phases.
00138     int nSpecies() const { return int(m_nsp); }
00139 
00140     //! Name of species with global index \a kGlob
00141     /*!
00142      * @param kGlob   global species index
00143      */
00144     std::string speciesName(const int kGlob) const;
00145 
00146     //! Returns the Number of atoms of global element \a mGlob in
00147     //! global species \a kGlob.
00148     /*!
00149      * @param kGlob   global species index
00150      * @param mGlob   global element index
00151      * @return        returns the number of atoms.
00152      */
00153     doublereal nAtoms(const int kGlob, const int mGlob) const;
00154 
00155     /// Returns the global Species mole fractions.
00156     /*!
00157      *   Write the array of species mole
00158      *   fractions into array \c x. The mole fractions are
00159      *   normalized to sum to one in each phase.
00160      *
00161      *    @param x  vector of mole fractions.
00162      *              Length = number of global species.
00163      */
00164     void getMoleFractions(doublereal* const x) const;
00165 
00166     //! Process phases and build atomic composition array. 
00167     /*!This method
00168      *  must be called after all phases are added, before doing
00169      *  anything else with the mixture. After init() has been called,
00170      *  no more phases may be added.
00171      */
00172     void init();
00173 
00174     //! Returns the name of the n'th phase
00175     /*!
00176      *   @param iph  phase Index
00177      */
00178     std::string phaseName(const index_t iph) const;
00179 
00180     //! Returns the index, given the phase name
00181     /*!
00182      * @param pName Name of the phase
00183      *
00184      * @return returns the index. A value of -1 means
00185      *         the phase isn't in the object.
00186      */
00187     int phaseIndex(const std::string &pName) const;
00188 
00189     //! Return the number of moles in phase n.
00190     /*!
00191      * @param n  Index of the phase.
00192      */
00193     doublereal phaseMoles(const index_t n) const;
00194 
00195     //! Set the number of moles of phase with index n.
00196     /*!
00197      * @param n     Index of the phase
00198      * @param moles Number of moles in the phase (kmol)
00199      */
00200     void setPhaseMoles(const index_t n, const doublereal moles);
00201 
00202     /// Return a %ThermoPhase reference to phase n.
00203     /*! The state of phase n is
00204      *  also updated to match the state stored locally in the
00205      *  mixture object.
00206      *
00207      * @param n  Phase Index
00208      *
00209      * @return   Reference to the %ThermoPhase object for the phase
00210      */
00211     phase_t& phase(index_t n);
00212 
00213     //! Returns the moles of global species \c k.
00214     /*!
00215      * Returns the moles of global species k.
00216      * units = kmol
00217      *
00218      * @param kGlob   Global species index k
00219      */
00220     doublereal speciesMoles(index_t kGlob) const;
00221 
00222     //! Return the global index of the species belonging to phase number \c p
00223     //! with local index \c k within the phase.
00224     /*!
00225      * Returns the index of the global species
00226      *
00227      * @param k local index of the species within the phase
00228      * @param p index of the phase
00229      */
00230     int speciesIndex(index_t k, index_t p) const {
00231       return m_spstart[p] + k;
00232     }
00233 
00234     //! Return the global index of the species belonging to phase name \c phaseName
00235     //! with species name \c speciesName
00236     /*!
00237      * Returns the index of the global species
00238      *
00239      * @param speciesName    Species Name
00240      * @param phaseName      Phase Name
00241      *
00242      * @return returns the global index
00243      *
00244      *  If the species or phase name is not recognized, this routine throws
00245      *  a CanteraError.
00246      */
00247     int speciesIndex(std::string speciesName, std::string phaseName);
00248 
00249     /// Minimum temperature for which all solution phases have
00250     /// valid thermo data. Stoichiometric phases are not
00251     /// considered, since they may have thermo data only valid for
00252     /// conditions for which they are stable.
00253     doublereal minTemp() const { return m_Tmin; }
00254 
00255     /// Maximum temperature for which all solution phases have
00256     /// valid thermo data. Stoichiometric phases are not
00257     /// considered, since they may have thermo data only valid for
00258     /// conditions for which they are stable.
00259     doublereal maxTemp() const { return m_Tmax; }
00260 
00261     /// Total charge (Coulombs).
00262     doublereal charge() const;
00263 
00264     /// Charge (Coulombs) of phase with index \a p.
00265     /*!
00266      * @param p     Phase Index
00267      */
00268     doublereal phaseCharge(index_t p) const;
00269 
00270     /// Total moles of global element \a m, summed over all phases.
00271     /*!
00272      * @param m   Index of the global element
00273      */
00274     doublereal elementMoles(index_t m) const;
00275 
00276     //!  Returns a vector of Chemical potentials.
00277     /*!
00278      *  Write into array \a mu the chemical
00279      *  potentials of all species [J/kmol]. The chemical
00280      *  potentials are related to the activities by
00281      *
00282      *  \f$
00283      *          \mu_k = \mu_k^0(T, P) + RT \ln a_k.
00284      *  \f$.
00285      *
00286      * @param mu Chemical potential vector.
00287      *           Length = num global species.
00288      *           Units = J/kmol.
00289      */
00290     void getChemPotentials(doublereal* mu) const;
00291 
00292     /// Returns a vector of Valid chemical potentials.
00293     /*!
00294      *   Write into array \a mu the
00295      *   chemical potentials of all species with thermo data valid
00296      *   for the current temperature [J/kmol]. For other species,
00297      *   set the chemical potential to the value \a not_mu. If \a
00298      *   standard is set to true, then the values returned are
00299      *   standard chemical potentials.
00300      *
00301      *   This method is designed for use in computing chemical
00302      *   equilibrium by Gibbs minimization. For solution phases (more
00303      *   than one species), this does the same thing as
00304      *   getChemPotentials. But for stoichiometric phases, this writes
00305      *   into array \a mu the user-specified value \a not_mu instead of
00306      *   the chemical potential if the temperature is outside the range
00307      *   for which the thermo data for the one species in the phase are
00308      *   valid. The need for this arises since many condensed phases
00309      *   have thermo data fit only for the temperature range for which
00310      *   they are stable. For example, in the NASA database, the fits
00311      *   for H2O(s) are only done up to 0 C, the fits for H2O(L) are
00312      *   only done from 0 C to 100 C, etc. Using the polynomial fits outside
00313      *   the range for which the fits were done can result in spurious
00314      *   chemical potentials, and can lead to condensed phases
00315      *   appearing when in fact they should be absent.
00316      *
00317      *   By setting \a not_mu to a large positive value, it is possible
00318      *   to force routines which seek to minimize the Gibbs free energy
00319      *   of the mixture to zero out any phases outside the temperature
00320      *   range for which their thermo data are valid.
00321      *
00322      * @param not_mu Value of the chemical potential to set
00323      *               species in phases, for which the thermo data 
00324      *               is not valid
00325      *
00326      * @param mu    Vector of chemical potentials
00327      *              length = Global species, units = J kmol-1
00328      *
00329      * @param standard  If this method is called with \a standard set to true, then
00330      *                  the composition-independent standard chemical potentials are
00331      *                  returned instead of the composition-dependent chemical
00332      *                  potentials.
00333      */
00334     void getValidChemPotentials(doublereal not_mu, doublereal* mu,
00335                                 bool standard = false) const;
00336 
00337     //! Temperature [K].
00338     doublereal temperature() const { return m_temp; }
00339 
00340     //! Set the mixture to a state of chemical equilibrium.
00341     /*!
00342      *    @param XY   Integer flag specifying properties to hold fixed.
00343      *    @param err  Error tolerance for \f$\Delta \mu/RT \f$ for
00344      *                all reactions. Also used as the relative error tolerance
00345      *                for the outer loop.
00346      *    @param maxsteps Maximum number of steps to take in solving
00347      *                    the fixed TP problem.
00348      *    @param maxiter Maximum number of "outer" iterations for
00349      *                   problems holding fixed something other than (T,P).
00350      *    @param loglevel Level of diagnostic output, written to a
00351      *                    file in HTML format.
00352      */
00353     doublereal equilibrate(int XY, doublereal err = 1.0e-9,
00354                            int maxsteps = 1000, int maxiter = 200, int loglevel = -99);
00355 
00356 
00357     /// Set the temperature [K].
00358     /*!
00359      * @param T   value of the temperature (Kelvin)
00360      */
00361     void setTemperature(const doublereal T);
00362 
00363     //! Set the state of the underlying ThermoPhase objects in one call
00364     /*!
00365      *   @param T      Temperature of the system (kelvin)
00366      *   @param Pres   pressure of the system (pascal)
00367      *                 (kmol)
00368      */
00369     void setState_TP(const doublereal T, const doublereal Pres);
00370 
00371     //! Set the state of the underlying ThermoPhase objects in one call
00372     /*!
00373      *   @param T      Temperature of the system (kelvin)
00374      *   @param Pres   pressure of the system (pascal)
00375      *   @param Moles  Vector of mole numbers of all the species in all the phases
00376      *                 (kmol)
00377      */
00378     void setState_TPMoles(const doublereal T, const doublereal Pres, const doublereal *Moles);
00379 
00380     /// Pressure [Pa].
00381     doublereal pressure() const {
00382       return m_press;
00383     }
00384 
00385     /// Volume [m^3].
00386     /*!
00387      * Returns the cummulative sum of the volumes of all the 
00388      * phases in the %MultiPhase.
00389      */
00390     doublereal volume() const;
00391 
00392     //! Set the pressure [Pa].
00393     /*!
00394      * @param P Set the pressure in the %MultiPhase object (Pa)
00395      */
00396     void setPressure(doublereal P) {
00397       m_press = P;
00398       updatePhases();
00399     }
00400 
00401     /// Enthalpy [J].
00402     doublereal enthalpy() const;
00403 
00404     /// Enthalpy [J].
00405     doublereal IntEnergy() const;
00406 
00407     /// Entropy [J/K].
00408     doublereal entropy() const;
00409 
00410     /// Gibbs function [J].
00411     doublereal gibbs() const;
00412 
00413     /// Heat capacity at constant pressure [J/K].
00414     doublereal cp() const;
00415 
00416     /// Number of phases.
00417     index_t nPhases() const {
00418       return m_np;
00419     }
00420 
00421     //! Return true is species \a kGlob is a species in a
00422     //! multicomponent solution phase.
00423     /*!
00424      * @param kGlob   index of the global species
00425      */
00426     bool solutionSpecies(index_t kGlob) const;
00427 
00428     //! Returns the phase index of the Kth "global" species
00429     /*!
00430      * @param kGlob Global species index.
00431      *
00432      * @return
00433      *     Returns the index of the owning phase.
00434      */
00435     int speciesPhaseIndex(const index_t kGlob) const;
00436 
00437     //! Returns the mole fraction of global species k
00438     /*!
00439      * @param kGlob Index of the global species.
00440      */
00441     doublereal moleFraction(const index_t kGlob) const;
00442 
00443     //! Set the Mole fractions of the nth phase
00444     /*!
00445      *  This function sets the mole fractions of the
00446      *  nth phase. Note, the mole number of the phase
00447      *  stays constant
00448      *
00449      * @param n    ID of the phase
00450      * @param x    Vector of input mole fractions.
00451      */
00452     void setPhaseMoleFractions(const index_t n, const doublereal* const x);
00453 
00454     //! Set the number numbers of species in the MultiPhase
00455     /*!
00456      *  @param xMap   CompositionMap of the species with 
00457      *                nonzero mole numbers
00458      *                units = kmol.
00459      */
00460     void setMolesByName(compositionMap& xMap);
00461 
00462     //! Set the Moles via a string containing their names.
00463     /*!
00464      * The string x is in the form of a composition map
00465      * Species which are not listed by name in the composition
00466      * map are set to zero.
00467      *
00468      * @param x string x in the form of a composition map
00469      *             where values are the moles of the species.
00470      */
00471     void setMolesByName(const std::string& x);
00472 
00473   
00474     //! Return a vector of global species mole numbers
00475     /*!
00476      *  Returns a vector of the number of moles of each species
00477      *  in the multiphase object.
00478      *
00479      * @param molNum Vector of doubles of length nSpecies
00480      *               containing the global mole numbers
00481      *               (kmol).
00482      */
00483     void getMoles(doublereal * molNum) const;
00484 
00485     //! Sets all of the global species mole numbers
00486     /*!
00487      *  Sets the number of moles of each species
00488      *  in the multiphase object.
00489      *
00490      * @param n    Vector of doubles of length nSpecies
00491      *             containing the global mole numbers
00492      *               (kmol).
00493      */
00494     void setMoles(const doublereal* n);
00495 
00496 
00497     //! Adds moles of a certain species to the mixture
00498     /*!
00499      *   @param indexS   Index of the species in the MultiPhase object
00500      *   @param addedMoles   Value of the moles that are added to the species.
00501      */
00502     void addSpeciesMoles(const int indexS, const doublereal addedMoles);
00503 
00504     //! Retrieves a vector of element abundances
00505     /*!
00506      * @param elemAbundances  Vector of element abundances
00507      * Length = number of elements in the MultiPhase object.
00508      * Index is the global element index
00509      * units is in kmol.
00510      */
00511     void getElemAbundances(doublereal * elemAbundances) const;
00512 
00513     //! Return true if the phase \a p has valid thermo data for
00514     //! the current temperature.
00515     /*!
00516      * @param p  Index of the phase.
00517      */
00518     bool tempOK(index_t p) const;
00519 
00520     // These methods are meant for internal use.
00521 
00522     //! Update the locally-stored composition within this object
00523     //! to match the current compositions of the phase objects.
00524     /*!
00525      *
00526      *  @deprecated  'update' is confusing within this context.
00527      *               Switching to the terminology 'uploadFrom'
00528      *               and 'downloadTo'. uploadFrom means to
00529      *               query the underlying ThermoPhase objects and 
00530      *               fill in the resulting information within
00531      *               this object. downloadTo means to take information
00532      *               from this object and put it into the underlying
00533      *               ThermoPhase objects.
00534      *               switch to uploadMoleFractionsFromPhases();
00535      */
00536     void updateMoleFractions();
00537 
00538     //! Update the locally-stored composition within this object
00539     //! to match the current compositions of the phase objects.
00540     /*!
00541      *    Query the underlying ThermoPhase objects for their moel
00542      *    fractions and fill in the mole fraction vector of this
00543      *    current object. Adjust element compositions within this
00544      *    object to match.
00545      *
00546      *    This is an upload operation in the sense that we are taking
00547      *    downstream information (ThermoPhase object info) and
00548      *    applying it to an upstream object (MultiPhase object).
00549      */
00550     void uploadMoleFractionsFromPhases();
00551 
00552   private:
00553 
00554     //! Set the states of the phase objects to the locally-stored
00555     //! state within this MultiPhase object.
00556     /*!
00557      *
00558      *  Note that if individual phases have T and P different
00559      *  than that stored locally, the phase T and P will be modified.
00560      *
00561      *    This is an download operation in the sense that we are taking
00562      *    upstream object information (MultiPhase object) and
00563      *    applying it to downstrean objects (ThermoPhase object information)
00564      *
00565      *    Therefore, the term, "update", is appropriate for a downstream
00566      *    operation.
00567      */
00568     void updatePhases() const;
00569 
00570     //! Calculate the element abundance vector
00571     void calcElemAbundances() const;
00572 
00573 
00574     //! Vector of the number of moles in each phase.
00575     /*!
00576      * Length = m_np, number of phases.
00577      */
00578     vector_fp m_moles;
00579 
00580     /**
00581      * Vector of the ThermoPhase Pointers.
00582      */
00583     std::vector<phase_t*> m_phase;
00584 
00585     //! Global Stoichiometric Coefficient array
00586     /*!
00587      *  This is a two dimensional array m_atoms(m, k). The first
00588      *  index is the global element index. The second index, k, is the
00589      *  global species index.
00590      *  The value is the number of atoms of type m in species k.
00591      */
00592     array_t m_atoms;
00593 
00594     /**
00595      * Locally storred vector of mole fractions of all species
00596      * comprising the MultiPhase object.
00597      */
00598     vector_fp m_moleFractions;
00599 
00600     //! Mapping between the global species number and the phase ID
00601     /*!
00602      *  m_spphase[kGlobal] = iPhase
00603      *  Length = number of global species
00604      */
00605     vector_int m_spphase;
00606 
00607     //! Vector of ints containing of first species index in the global list of species
00608     //! for each phase
00609     /*!
00610      *  kfirst = m_spstart[ip], kfirst is the index of the first species in the ip'th
00611      *                          phase.
00612      */
00613     vector_int m_spstart;
00614 
00615     //! String names of the global elements
00616     /*!
00617      *  This has a length equal to the number of global elements.
00618      */
00619     std::vector<std::string> m_enames;
00620 
00621     //! Atomic number of each element
00622     /*!
00623      *  This is the atomic number of each global element.
00624      */
00625     vector_int m_atomicNumber;
00626 
00627     //! Vector of species names in the problem
00628     /*!
00629      *   Vector is over all species defined in the object,
00630      *   the global species index.
00631      */
00632     std::vector<std::string> m_snames;
00633 
00634     //! Returns the global element index, given the element string name
00635     /*!
00636      * -> used in the construction. However, wonder if it needs to be global.
00637      */
00638     std::map<std::string, int> m_enamemap;
00639 
00640     /**
00641      *   Number of phases in the MultiPhase object
00642      */
00643     index_t  m_np;
00644 
00645     //! Current value of the temperature (kelvin)
00646     doublereal m_temp;
00647 
00648     //! Current value of the pressure (Pa)
00649     doublereal m_press;
00650 
00651     /**
00652      * Number of distinct elements in all of the phases
00653      */
00654     index_t m_nel;
00655     /**
00656      * Number of distinct species in all of the phases
00657      */
00658     index_t m_nsp;
00659 
00660     //! True if the init() routine has been called, and the MultiPhase frozen
00661     bool m_init;
00662 
00663     //! Global ID of the element corresponding to the electronic charge.
00664     /*!
00665      * If there is none, then this is equal to -1
00666      */
00667     int m_eloc;
00668 
00669     //! Vector of bools indicating whether temperatures are ok for phases.
00670     /*!
00671      * If the current temperature is outside the range of valid temperatures
00672      * for the phase thermodynamics, the phase flag is set to false.
00673      */
00674     mutable std::vector<bool> m_temp_OK;
00675 
00676     //! Minimum temperature for which thermo parameterizations are valid
00677     /*!
00678      *  Stoichiometric phases are ignored in this determination.
00679      *  units Kelvin
00680      */
00681     doublereal m_Tmin;
00682 
00683     //! Minimum temperature for which thermo parameterizations are valid
00684     /*!
00685      *  Stoichiometric phases are ignored in this determination.
00686      *  units Kelvin
00687      */
00688     doublereal m_Tmax;
00689 
00690     //! Vector of element abundances
00691     /*!
00692      *  m_elemAbundances[mGlobal] = kmol of element mGlobal summed over all
00693      *      species in all phases.
00694      */
00695     mutable vector_fp m_elemAbundances;
00696   };
00697 
00698   //! Function to output a MultiPhase description to a stream
00699   /*!
00700    *  Writes out a description of the contents of each phase of the 
00701    *  MultiPhase using the report function.
00702    *
00703    *  @param s ostream
00704    *  @param x  Reference to a MultiPhase
00705    *  @return returns a reference to the ostream
00706    */
00707   inline std::ostream& operator<<(std::ostream& s, Cantera::MultiPhase& x) {
00708     size_t ip;
00709     for (ip = 0; ip < x.nPhases(); ip++) {
00710       if (x.phase(ip).name() != "") {
00711         s << "*************** " << x.phase(ip).name() << " *****************" << std::endl;
00712       }
00713       else {
00714         s << "*************** Phase " << ip << " *****************" << std::endl;
00715       }
00716       s << "Moles: " << x.phaseMoles(ip) << std::endl;
00717 
00718       s << report(x.phase(ip)) << std::endl;
00719     }
00720     return s;
00721   }
00722 
00723   //!  Choose the optimum basis of species for the equilibrium calculations.
00724   /*!
00725    * This is done by 
00726    * choosing the species with the largest mole fraction 
00727    * not currently a linear combination of the previous components. 
00728    * Then, calculate the stoichiometric coefficient matrix for that 
00729    * basis. 
00730    *
00731    * Calculates the identity of the component species in the mechanism. 
00732    * Rearranges the solution data to put the component data at the 
00733    * front of the species list. 
00734    *
00735    * Then, calculates SC(J,I) the formation reactions for all noncomponent 
00736    * species in the mechanism. 
00737    *
00738    * Input 
00739    * --------- 
00740    * @param mphase   Pointer to the multiphase object. Contains the 
00741    *                 species mole fractions, which are used to pick the
00742    *                 current optimal species component basis.
00743    * @param orderVectorElements
00744    *                 Order vector for the elements. The element rows
00745    *                 in the formula matrix are
00746    *                 rearranged according to this vector.
00747    * @param orderVectorSpecies
00748    *                 Order vector for the species. The species are
00749    *                 rearranged according to this formula. The first
00750    *                 nCompoments of this vector contain the calculated
00751    *                 species components on exit.
00752    * @param doFormRxn  If true, the routine calculates the formation
00753    *                 reaction matrix based on the calculated 
00754    *                 component species. If false, this step is skipped.
00755    * 
00756    * Output 
00757    * --------- 
00758    *  @param usedZeroedSpecies = If true, then a species with a zero concentration
00759    *                     was used as a component. The problem may be
00760    *                     converged.
00761    * @param formRxnMatrix 
00762    *
00763    * @return      Returns the number of components.
00764    *
00765    *  @ingroup equilfunctions
00766    */
00767   int BasisOptimize( int *usedZeroedSpecies, bool doFormRxn,
00768                      MultiPhase *mphase, vector_int & orderVectorSpecies,
00769                      vector_int & orderVectorElements,
00770                      vector_fp & formRxnMatrix);
00771 
00772   //!   This subroutine handles the potential rearrangement of the constraint
00773   //!   equations represented by the Formula Matrix. 
00774   /*!
00775    *    Rearrangement is only
00776    *    necessary when the number of components is less than the number of
00777    *    elements. For this case, some constraints can never be satisfied 
00778    *    exactly, because the range space represented by the Formula
00779    *    Matrix of the components can't span the extra space. These 
00780    *    constraints, which are out of the range space of the component
00781    *    Formula matrix entries, are migrated to the back of the Formula
00782    *    matrix.
00783    *
00784    *    A prototypical example is an extra element column in 
00785    *    FormulaMatrix[], 
00786    *    which is identically zero. For example, let's say that argon is
00787    *    has an element column in FormulaMatrix[], but no species in the 
00788    *    mechanism
00789    *    actually contains argon. Then, nc < ne. Unless the entry for
00790    *    desired element abundance vector for Ar is zero, then this
00791    *    element abundance constraint can never be satisfied. The 
00792    *    constraint vector is not in the range space of the formula
00793    *    matrix.
00794    *    Also, without perturbation
00795    *    of FormulaMatrix[], BasisOptimize[] would produce a zero pivot 
00796    *    because the matrix
00797    *    would be singular (unless the argon element column was already the
00798    *    last column of  FormulaMatrix[]. 
00799    *       This routine borrows heavily from BasisOptimize algorithm. It 
00800    *    finds nc constraints which span the range space of the Component
00801    *    Formula matrix, and assigns them as the first nc components in the
00802    *    formular matrix. This guarrantees that BasisOptimize has a
00803    *    nonsingular matrix to invert.
00804    *  input 
00805    *    @param nComponents  Number of components calculated previously.
00806    *
00807    *    @param elementAbundances  Current value of the element abundances
00808    *
00809    *    @param mphase  Input pointer to a MultiPhase object
00810    *
00811    *    @param orderVectorSpecies input vector containing the ordering
00812    *                of the global species in mphase. This is used
00813    *                to extract the component basis of the mphase object.
00814    *
00815    *  output
00816    *     @param orderVectorElements Ouput vector containing the order
00817    *                      of the elements that is necessary for 
00818    *                      calculation of the formula matrix. 
00819    *
00820    *  @ingroup equilfunctions
00821    */
00822   int ElemRearrange(int nComponents,  const vector_fp & elementAbundances,
00823                     MultiPhase *mphase,
00824                     vector_int & orderVectorSpecies,
00825                     vector_int & orderVectorElements);
00826 
00827 #ifdef DEBUG_MODE
00828   //! External int that is used to turn on debug printing for the
00829   //! BasisOptimze program.
00830   /*!
00831    *   Set this to 1 if you want debug printing from BasisOptimize.
00832    */
00833   extern int BasisOptimize_print_lvl;
00834 #endif
00835 }
00836 
00837 #endif
Generated by  doxygen 1.6.3