From 5cf9a9faf99c5b71854df288e1afb8c9e3952241 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 30 Jul 2010 19:04:37 +0000 Subject: [PATCH] Incremental update --- Cantera/src/equil/vcs_VolPhase.cpp | 39 ++++++++++++++---------- Cantera/src/equil/vcs_VolPhase.h | 28 ++++++++++++++--- Cantera/src/equil/vcs_phaseStability.cpp | 12 +++----- 3 files changed, 51 insertions(+), 28 deletions(-) diff --git a/Cantera/src/equil/vcs_VolPhase.cpp b/Cantera/src/equil/vcs_VolPhase.cpp index afd547563..19c93a2c4 100644 --- a/Cantera/src/equil/vcs_VolPhase.cpp +++ b/Cantera/src/equil/vcs_VolPhase.cpp @@ -51,6 +51,8 @@ namespace VCSnonideal { m_useCanteraCalls(false), TP_ptr(0), v_totalMoles(0.0), + creationMoleNumbers_(0), + creationGlobalRxnNumbers_(0), m_phiVarIndex(-1), m_totalVol(0.0), m_vcsStateStatus(VCS_STATECALC_OLD), @@ -110,7 +112,9 @@ namespace VCSnonideal { m_MFStartIndex(b.m_MFStartIndex), m_useCanteraCalls(b.m_useCanteraCalls), TP_ptr(b.TP_ptr), - v_totalMoles(b.v_totalMoles), + v_totalMoles(b.v_totalMoles), + creationMoleNumbers_(0), + creationGlobalRxnNumbers_(0), m_phiVarIndex(-1), m_totalVol(b.m_totalVol), m_vcsStateStatus(VCS_STATECALC_OLD), @@ -214,7 +218,9 @@ namespace VCSnonideal { v_totalMoles = b.v_totalMoles; Xmol = b.Xmol; - fractionCreationDelta_ = b.fractionCreationDelta_; + + creationMoleNumbers_ = b.creationMoleNumbers_; + creationGlobalRxnNumbers_ = b.creationGlobalRxnNumbers_; m_phi = b.m_phi; m_phiVarIndex = b.m_phiVarIndex; @@ -309,10 +315,12 @@ namespace VCSnonideal { } Xmol.resize(nspecies, 0.0); - fractionCreationDelta_.resize(nspecies, 0.0); + creationMoleNumbers_.resize(nspecies, 0.0); + creationGlobalRxnNumbers_.resize(nspecies, -1); for (int i = 0; i < nspecies; i++) { Xmol[i] = 1.0/nspecies; - fractionCreationDelta_[i] = 1.0/nspecies; + creationMoleNumbers_[i] = 1.0/nspecies; + creationGlobalRxnNumbers_[i] = IndSpecies[i] - m_numElemConstraints; } SS0ChemicalPotential.resize(nspecies, -1.0); @@ -675,7 +683,7 @@ namespace VCSnonideal { */ if (stateCalc == VCS_STATECALC_OLD) { if (v_totalMoles > 0.0) { - fractionCreationDelta_ = Xmol; + creationMoleNumbers_ = Xmol; } } @@ -1107,7 +1115,7 @@ namespace VCSnonideal { resize(VP_ID_, nsp, nelem, PhaseName.c_str()); } TP_ptr->getMoleFractions(VCS_DATA_PTR(Xmol)); - fractionCreationDelta_ = Xmol; + creationMoleNumbers_ = Xmol; _updateMoleFractionDependencies(); /* @@ -1161,20 +1169,16 @@ namespace VCSnonideal { void vcs_VolPhase::setCreationMoleNumbers(const double * const n_k, const std::vector &creationGlobalRxnNumbers) { - vcs_dcopy(VCS_DATA_PTR(fractionCreationDelta_), n_k, m_numSpecies); - // vcs_icopy(VCS_DATA_PTR(creationGlobalRxnNumbers_), VCS_DATA_PTR(creationGlobalRxnNumbers), m_numSpecies); + vcs_dcopy(VCS_DATA_PTR(creationMoleNumbers_), n_k, m_numSpecies); + vcs_icopy(VCS_DATA_PTR(creationGlobalRxnNumbers_), VCS_DATA_PTR(creationGlobalRxnNumbers), m_numSpecies); } - - // void vcs_VolPhase::setFractionCreationDeltas(const double * const F_k) { - // for (int k = 0; k < m_numSpecies; k++) { - // fractionCreationDelta_[k] = F_k[k]; - // } - //} /***************************************************************************/ - const std::vector & vcs_VolPhase::fractionCreationDeltas() const { - return fractionCreationDelta_; + const std::vector & vcs_VolPhase::creationMoleNumbers(std::vector &creationGlobalRxnNumbers) const { + creationGlobalRxnNumbers = creationGlobalRxnNumbers_; + return creationMoleNumbers_; } + /***************************************************************************/ // Sets the total moles in the phase @@ -1387,6 +1391,9 @@ namespace VCSnonideal { void vcs_VolPhase::setSpGlobalIndexVCS(const int spIndex, const int spGlobalIndex) { IndSpecies[spIndex] = spGlobalIndex; + if (spGlobalIndex >= m_numElemConstraints) { + creationGlobalRxnNumbers_[spIndex] = spGlobalIndex - m_numElemConstraints; + } } /**********************************************************************/ diff --git a/Cantera/src/equil/vcs_VolPhase.h b/Cantera/src/equil/vcs_VolPhase.h index ce763402f..9363a2068 100644 --- a/Cantera/src/equil/vcs_VolPhase.h +++ b/Cantera/src/equil/vcs_VolPhase.h @@ -425,11 +425,11 @@ namespace VCSnonideal { */ void setCreationMoleNumbers(const double * const n_k, const std::vector &creationGlobalRxnNumbers); - //! Sets the fractionCreationDelta's within the phase object + //! Return a const reference to the creationMoleNumbers storred in the object. /*! - * @param F_k Pointer to a vector of F_k's + * @return Returns a const reference to the vector of creationMoleNumbers */ - // void setFractionCreationDeltas( const double * const F_k); + const std::vector & creationMoleNumbers(std::vector &creationGlobalRxnNumbers) const; //! Return a const reference to the fractionCreationDeltas storred in the //! object. @@ -879,11 +879,31 @@ namespace VCSnonideal { //! in the phase std::vector Xmol; + //! Vector of current creationMoleNumbers_ + /*! + * These are the actual unknowns in the phase stability problem + */ + std::vector creationMoleNumbers_; + + //! Vector of creation global reaction numbers for the phase stability problem + /*! + * The phase stability problem requires a global reaction number for each + * species in the phase. Usually this is the krxn = kglob - M for species + * in the phase that are not components. For component species, the + * choice of the reaction is one which maximimes the chance that the phase + * pops into (or remains in) existence. + * The index here is the local phase species index. + * the value of the variable is the global vcs reaction number. Note, + * that the global reaction number will go out of order when the species positions + * are swapped. So, this number has to be recalculated. + */ + std::vector creationGlobalRxnNumbers_; + //! Vector of current fractionalCreationDeltas /*! * These are the actual unknowns in the problem */ - std::vector fractionCreationDelta_; + // std::vector fractionCreationDelta_; //! If the potential is a solution variable in VCS, it acts as a species. diff --git a/Cantera/src/equil/vcs_phaseStability.cpp b/Cantera/src/equil/vcs_phaseStability.cpp index c055c76e9..17cbf9a2c 100644 --- a/Cantera/src/equil/vcs_phaseStability.cpp +++ b/Cantera/src/equil/vcs_phaseStability.cpp @@ -405,10 +405,8 @@ namespace VCSnonideal { } else { vector fracDelta(Vphase->nSpecies()); vector X_est(Vphase->nSpecies()); - - // fracDelta = Vphase->creationMoleNumbers(creationGlobalRxnNumbers); - fracDelta = Vphase->fractionCreationDeltas(); - + fracDelta = Vphase->creationMoleNumbers(creationGlobalRxnNumbers); + double sumFrac = 0.0; for (k = 0; k < Vphase->nSpecies(); k++) { sumFrac += fracDelta[k]; @@ -536,8 +534,8 @@ namespace VCSnonideal { // Get the storred estimate for the composition of the phase if // it gets created - // fracDelta_new = Vphase->creationMoleNumbers(creationGlobalRxnNumbers); - fracDelta_new = Vphase->fractionCreationDeltas(); + fracDelta_new = Vphase->creationMoleNumbers(creationGlobalRxnNumbers); + bool oneIsComponent = false; std::vector componentList; @@ -793,8 +791,6 @@ namespace VCSnonideal { Vphase->setMoleFractionsState(0.0, VCS_DATA_PTR(X_est), VCS_STATECALC_PHASESTABILITY); Vphase->setCreationMoleNumbers(VCS_DATA_PTR(fracDelta_new), creationGlobalRxnNumbers); - //Vphase->setFractionCreationDeltas(VCS_DATA_PTR(fracDelta_new)); - }