Cantera  3.0.0
Loading...
Searching...
No Matches
IonsFromNeutralVPSSTP.cpp
Go to the documentation of this file.
1/**
2 * @file IonsFromNeutralVPSSTP.cpp
3 * Definitions for the object which treats ionic liquids as made of ions as species
4 * even though the thermodynamics is obtained from the neutral molecule representation.
5 * (see @ref thermoprops
6 * and class @link Cantera::IonsFromNeutralVPSSTP IonsFromNeutralVPSSTP@endlink).
7 *
8 * Header file for a derived class of ThermoPhase that handles variable pressure
9 * standard state methods for calculating thermodynamic properties that are
10 * further based upon expressions for the excess Gibbs free energy expressed as
11 * a function of the mole fractions.
12 */
13
14// This file is part of Cantera. See License.txt in the top-level directory or
15// at https://cantera.org/license.txt for license and copyright information.
16
22
23#include <boost/algorithm/string.hpp>
24#include <fstream>
25
26namespace Cantera
27{
28
29IonsFromNeutralVPSSTP::IonsFromNeutralVPSSTP(const string& inputFile, const string& id_)
30{
31 warn_deprecated("class IonsFromNeutralVPSSTP", "To be removed after Cantera 3.0");
32 initThermoFile(inputFile, id_);
33}
34
35// ------------ Molar Thermodynamic Properties ----------------------
36
38{
39 getPartialMolarEnthalpies(m_work.data());
40 return mean_X(m_work);
41}
42
44{
45 getPartialMolarEntropies(m_work.data());
46 return mean_X(m_work);
47}
48
50{
51 getChemPotentials(m_work.data());
52 return mean_X(m_work);
53}
54
56{
57 getPartialMolarCp(m_work.data());
58 return mean_X(m_work);
59}
60
62{
63 // Need to revisit this, as it is wrong
64 getPartialMolarCp(m_work.data());
65 return mean_X(m_work);
66}
67
68// -- Activities, Standard States, Activity Concentrations -----------
69
71 vector<double>& charges, vector<size_t>& neutMolIndex) const
72{
73 coeffs = fm_neutralMolec_ions_;
74 charges = m_speciesCharge;
75 neutMolIndex = fm_invert_ionForNeutral;
76}
77
79{
80 // Update the activity coefficients
82
83 // take the exp of the internally stored coefficients.
84 for (size_t k = 0; k < m_kk; k++) {
85 ac[k] = exp(lnActCoeff_Scaled_[k]);
86 }
87}
88
89// --------- Partial Molar Properties of the Solution -------------
90
92{
93 size_t icat, jNeut;
94 double xx, fact2;
95
96 // Get the standard chemical potentials of neutral molecules
97 neutralMoleculePhase_->getStandardChemPotentials(muNeutralMolecule_.data());
98
99 switch (ionSolnType_) {
100 case cIonSolnType_PASSTHROUGH:
101 neutralMoleculePhase_->getChemPotentials(mu);
102 break;
103 case cIonSolnType_SINGLEANION:
104 neutralMoleculePhase_->getLnActivityCoefficients(lnActCoeff_NeutralMolecule_.data());
105 fact2 = 2.0 * RT() * log(2.0);
106
107 // Do the cation list
108 for (size_t k = 0; k < cationList_.size(); k++) {
109 // Get the id for the next cation
110 icat = cationList_[k];
111 jNeut = fm_invert_ionForNeutral[icat];
112 xx = std::max(SmallNumber, moleFractions_[icat]);
113 mu[icat] = muNeutralMolecule_[jNeut] + fact2 + RT() * (lnActCoeff_NeutralMolecule_[jNeut] + log(xx));
114 }
115
116 // Do the anion list
117 icat = anionList_[0];
118 jNeut = fm_invert_ionForNeutral[icat];
119 xx = std::max(SmallNumber, moleFractions_[icat]);
120 mu[icat] = RT() * log(xx);
121
122 // Do the list of neutral molecules
123 for (size_t k = 0; k < passThroughList_.size(); k++) {
124 icat = passThroughList_[k];
125 jNeut = fm_invert_ionForNeutral[icat];
126 xx = std::max(SmallNumber, moleFractions_[icat]);
127 mu[icat] = muNeutralMolecule_[jNeut] + RT() * (lnActCoeff_NeutralMolecule_[jNeut] + log(xx));
128 }
129 break;
130
131 case cIonSolnType_SINGLECATION:
132 throw CanteraError("IonsFromNeutralVPSSTP::getChemPotentials", "Unknown type");
133 case cIonSolnType_MULTICATIONANION:
134 throw CanteraError("IonsFromNeutralVPSSTP::getChemPotentials", "Unknown type");
135 default:
136 throw CanteraError("IonsFromNeutralVPSSTP::getChemPotentials", "Unknown type");
137 }
138}
139
141{
142 // Get the nondimensional standard state enthalpies
143 getEnthalpy_RT(hbar);
144
145 // dimensionalize it.
146 for (size_t k = 0; k < m_kk; k++) {
147 hbar[k] *= RT();
148 }
149
150 // Update the activity coefficients, This also update the internally stored
151 // molalities.
154 for (size_t k = 0; k < m_kk; k++) {
155 hbar[k] -= RT() * temperature() * dlnActCoeffdT_Scaled_[k];
156 }
157}
158
160{
161 // Get the nondimensional standard state entropies
162 getEntropy_R(sbar);
163
164 // Update the activity coefficients, This also update the internally stored
165 // molalities.
168
169 for (size_t k = 0; k < m_kk; k++) {
170 double xx = std::max(moleFractions_[k], SmallNumber);
171 sbar[k] += - lnActCoeff_Scaled_[k] -log(xx) - temperature() * dlnActCoeffdT_Scaled_[k];
172 }
173
174 // dimensionalize it.
175 for (size_t k = 0; k < m_kk; k++) {
176 sbar[k] *= GasConstant;
177 }
178}
179
180void IonsFromNeutralVPSSTP::getdlnActCoeffdlnX_diag(double* dlnActCoeffdlnX_diag) const
181{
184
185 for (size_t k = 0; k < m_kk; k++) {
186 dlnActCoeffdlnX_diag[k] = dlnActCoeffdlnX_diag_[k];
187 }
188}
189
190void IonsFromNeutralVPSSTP::getdlnActCoeffdlnN_diag(double* dlnActCoeffdlnN_diag) const
191{
194
195 for (size_t k = 0; k < m_kk; k++) {
196 dlnActCoeffdlnN_diag[k] = dlnActCoeffdlnN_diag_[k];
197 }
198}
199
200void IonsFromNeutralVPSSTP::getdlnActCoeffdlnN(const size_t ld, double* dlnActCoeffdlnN)
201{
204 double* data = & dlnActCoeffdlnN_(0,0);
205 for (size_t k = 0; k < m_kk; k++) {
206 for (size_t m = 0; m < m_kk; m++) {
207 dlnActCoeffdlnN[ld * k + m] = data[m_kk * k + m];
208 }
209 }
210}
211
213{
214 // This is a two phase process. First, we calculate the standard states
215 // within the neutral molecule phase.
216 neutralMoleculePhase_->setState_TP(temperature(), pressure());
217
218 // Calculate the partial molar volumes, and then the density of the fluid
220}
221
223{
224 // Download the neutral mole fraction vector into the vector,
225 // NeutralMolecMoleFractions_[]
226 neutralMoleculePhase_->getMoleFractions(NeutralMolecMoleFractions_.data());
227
228 // Zero the mole fractions
229 for (size_t k = 0; k < m_kk; k++) {
230 mf[k] = 0.0;
231 }
232
233 // Use the formula matrix to calculate the relative mole numbers.
234 for (size_t jNeut = 0; jNeut < numNeutralMoleculeSpecies_; jNeut++) {
235 for (size_t k = 0; k < m_kk; k++) {
236 double fmij = fm_neutralMolec_ions_[k + jNeut * m_kk];
237 mf[k] += fmij * NeutralMolecMoleFractions_[jNeut];
238 }
239 }
240
241 // Normalize the new mole fractions
242 double sum = 0.0;
243 for (size_t k = 0; k < m_kk; k++) {
244 sum += mf[k];
245 }
246 for (size_t k = 0; k < m_kk; k++) {
247 mf[k] /= sum;
248 }
249}
250
252{
253 size_t icat, jNeut;
254 double fmij;
255 double sum = 0.0;
256
257 // Zero the vector we are trying to find.
258 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
260 }
261 sum = -1.0;
262 for (size_t k = 0; k < m_kk; k++) {
263 sum += moleFractions_[k];
264 }
265 if (fabs(sum) > 1.0E-11) {
266 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions",
267 "molefracts don't sum to one: {}", sum);
268 }
269
270 switch (ionSolnType_) {
271 case cIonSolnType_PASSTHROUGH:
272 for (size_t k = 0; k < m_kk; k++) {
274 }
275 break;
276
277 case cIonSolnType_SINGLEANION:
278 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
280 }
281
282 for (size_t k = 0; k < cationList_.size(); k++) {
283 // Get the id for the next cation
284 icat = cationList_[k];
285 jNeut = fm_invert_ionForNeutral[icat];
286 if (jNeut != npos) {
287 fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
288 AssertTrace(fmij != 0.0);
289 NeutralMolecMoleFractions_[jNeut] += moleFractions_[icat] / fmij;
290 }
291 }
292
293 for (size_t k = 0; k < passThroughList_.size(); k++) {
294 icat = passThroughList_[k];
295 jNeut = fm_invert_ionForNeutral[icat];
296 fmij = fm_neutralMolec_ions_[ icat + jNeut * m_kk];
297 NeutralMolecMoleFractions_[jNeut] += moleFractions_[icat] / fmij;
298 }
299
300 for (size_t k = 0; k < m_kk; k++) {
302 }
303 for (jNeut = 0; jNeut < numNeutralMoleculeSpecies_; jNeut++) {
304 for (size_t k = 0; k < m_kk; k++) {
305 fmij = fm_neutralMolec_ions_[k + jNeut * m_kk];
307 }
308 }
309 for (size_t k = 0; k < m_kk; k++) {
310 if (fabs(moleFractionsTmp_[k]) > 1.0E-13) {
311 // Check to see if we have in fact found the inverse.
312 if (anionList_[0] != k) {
313 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions",
314 "neutral molecule calc error");
315 } else {
316 // For the single anion case, we will allow some slippage
317 if (fabs(moleFractionsTmp_[k]) > 1.0E-5) {
318 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions",
319 "neutral molecule calc error - anion");
320 }
321 }
322 }
323 }
324
325 // Normalize the Neutral Molecule mole fractions
326 sum = 0.0;
327 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
329 }
330 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
332 }
333 break;
334
335 case cIonSolnType_SINGLECATION:
336 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions", "Unknown type");
337 break;
338 case cIonSolnType_MULTICATIONANION:
339 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions", "Unknown type");
340 break;
341 default:
342 throw CanteraError("IonsFromNeutralVPSSTP::calcNeutralMoleculeMoleFractions", "Unknown type");
343 break;
344 }
345}
346
347void IonsFromNeutralVPSSTP::getNeutralMoleculeMoleGrads(const double* const dx, double* const dy) const
348{
349 double sumy, sumdy;
350
351 // check sum dx = 0
352 // Zero the vector we are trying to find.
353 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
354 y_[k] = 0.0;
355 dy[k] = 0.0;
356 }
357
358 switch (ionSolnType_) {
359
360 case cIonSolnType_PASSTHROUGH:
361 for (size_t k = 0; k < m_kk; k++) {
362 dy[k] = dx[k];
363 }
364 break;
365
366 case cIonSolnType_SINGLEANION:
367 for (size_t k = 0; k < cationList_.size(); k++) {
368 // Get the id for the next cation
369 size_t icat = cationList_[k];
370 size_t jNeut = fm_invert_ionForNeutral[icat];
371 if (jNeut != npos) {
372 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
373 AssertTrace(fmij != 0.0);
374 const double temp = 1.0/fmij;
375 dy[jNeut] += dx[icat] * temp;
376 y_[jNeut] += moleFractions_[icat] * temp;
377 }
378 }
379
380 for (size_t k = 0; k < passThroughList_.size(); k++) {
381 size_t icat = passThroughList_[k];
382 size_t jNeut = fm_invert_ionForNeutral[icat];
383 double fmij = fm_neutralMolec_ions_[ icat + jNeut * m_kk];
384 const double temp = 1.0/fmij;
385 dy[jNeut] += dx[icat] * temp;
386 y_[jNeut] += moleFractions_[icat] * temp;
387 }
388 // Normalize the Neutral Molecule mole fractions
389 sumy = 0.0;
390 sumdy = 0.0;
391 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
392 sumy += y_[k];
393 sumdy += dy[k];
394 }
395 sumy = 1.0 / sumy;
396 for (size_t k = 0; k < numNeutralMoleculeSpecies_; k++) {
397 dy[k] = dy[k] * sumy - y_[k]*sumdy*sumy*sumy;
398 }
399
400 break;
401
402 case cIonSolnType_SINGLECATION:
403 throw CanteraError("IonsFromNeutralVPSSTP::getNeutralMoleculeMoleGrads",
404 "Unknown type");
405 break;
406 case cIonSolnType_MULTICATIONANION:
407 throw CanteraError("IonsFromNeutralVPSSTP::getNeutralMoleculeMoleGrads",
408 "Unknown type");
409 break;
410 default:
411 throw CanteraError("IonsFromNeutralVPSSTP::getNeutralMoleculeMoleGrads",
412 "Unknown type");
413 break;
414 }
415}
416
418{
421 neutralMoleculePhase_->setMoleFractions(NeutralMolecMoleFractions_.data());
422}
423
424// ------------ Partial Molar Properties of the Solution ------------
425
426//! Return the factor overlap
427/*!
428 * @param elnamesVN
429 * @param elemVectorN
430 * @param nElementsN
431 * @param elnamesVI
432 * @param elemVectorI
433 * @param nElementsI
434 */
435static double factorOverlap(const vector<string>& elnamesVN ,
436 const vector<double>& elemVectorN,
437 const size_t nElementsN,
438 const vector<string>& elnamesVI ,
439 const vector<double>& elemVectorI,
440 const size_t nElementsI)
441{
442 double fMax = 1.0E100;
443 for (size_t mi = 0; mi < nElementsI; mi++) {
444 if (elnamesVI[mi] != "E" && elemVectorI[mi] > 1.0E-13) {
445 double eiNum = elemVectorI[mi];
446 for (size_t mn = 0; mn < nElementsN; mn++) {
447 if (elnamesVI[mi] == elnamesVN[mn]) {
448 if (elemVectorN[mn] <= 1.0E-13) {
449 return 0.0;
450 }
451 fMax = std::min(fMax, elemVectorN[mn]/eiNum);
452 }
453 }
454 }
455 }
456 return fMax;
457}
458
460 const AnyMap& rootNode)
461{
462 ThermoPhase::setParameters(phaseNode, rootNode);
463 m_rootNode = rootNode;
464}
465
467{
468 if (m_input.hasKey("neutral-phase")) {
469 string neutralName = m_input["neutral-phase"].asString();
470 const auto& slash = boost::ifind_last(neutralName, "/");
471 if (slash) {
472 string fileName(neutralName.begin(), slash.begin());
473 neutralName = string(slash.end(), neutralName.end());
474 AnyMap infile = AnyMap::fromYamlFile(fileName,
475 m_input.getString("__file__", ""));
476 AnyMap& phaseNode = infile["phases"].getMapWhere("name", neutralName);
477 setNeutralMoleculePhase(newThermo(phaseNode, infile));
478 } else {
479 AnyMap& phaseNode = m_rootNode["phases"].getMapWhere("name", neutralName);
480 setNeutralMoleculePhase(newThermo(phaseNode, m_rootNode));
481 }
482 }
483
485 throw CanteraError(
486 "IonsFromNeutralVPSSTP::initThermo",
487 "The neutral phase has not been initialized. Are you missing the "
488 "'neutral-phase' key?"
489 );
490 }
491
492 size_t nElementsN = neutralMoleculePhase_->nElements();
493 const vector<string>& elnamesVN = neutralMoleculePhase_->elementNames();
494 vector<double> elemVectorN(nElementsN);
495
496 size_t nElementsI = nElements();
497 const vector<string>& elnamesVI = elementNames();
498 vector<double> elemVectorI(nElementsI);
499
500 if (indexSpecialSpecies_ == npos) {
501 throw CanteraError(
502 "IonsFromNeutralVPSSTP::initThermo",
503 "No special-species were specified in the phase."
504 );
505 }
506 for (size_t m = 0; m < nElementsI; m++) {
507 elemVectorI[m] = nAtoms(indexSpecialSpecies_, m);
508 }
509
510 for (size_t jNeut = 0; jNeut < numNeutralMoleculeSpecies_; jNeut++) {
511 for (size_t m = 0; m < nElementsN; m++) {
512 elemVectorN[m] = neutralMoleculePhase_->nAtoms(jNeut, m);
513 }
514
515 double fac = factorOverlap(elnamesVN, elemVectorN, nElementsN,
516 elnamesVI ,elemVectorI, nElementsI);
517 if (fac > 0.0) {
518 for (size_t m = 0; m < nElementsN; m++) {
519 for (size_t mi = 0; mi < nElementsI; mi++) {
520 if (elnamesVN[m] == elnamesVI[mi]) {
521 elemVectorN[m] -= fac * elemVectorI[mi];
522 }
523
524 }
525 }
526 }
528
529 for (size_t k = 0; k < m_kk; k++) {
530 for (size_t m = 0; m < nElementsI; m++) {
531 elemVectorI[m] = nAtoms(k, m);
532 }
533 fac = factorOverlap(elnamesVN, elemVectorN, nElementsN,
534 elnamesVI ,elemVectorI, nElementsI);
535 if (fac > 0.0) {
536 for (size_t m = 0; m < nElementsN; m++) {
537 for (size_t mi = 0; mi < nElementsI; mi++) {
538 if (elnamesVN[m] == elnamesVI[mi]) {
539 elemVectorN[m] -= fac * elemVectorI[mi];
540 }
541 }
542 }
543 bool notTaken = true;
544 for (size_t iNeut = 0; iNeut < jNeut; iNeut++) {
545 if (fm_invert_ionForNeutral[k] == iNeut) {
546 notTaken = false;
547 }
548 }
549 if (notTaken) {
550 fm_invert_ionForNeutral[k] = jNeut;
551 } else {
552 throw CanteraError("IonsFromNeutralVPSSTP::initThermo",
553 "Simple formula matrix generation failed, one cation is shared between two salts");
554 }
555 }
556 fm_neutralMolec_ions_[k + jNeut * m_kk] += fac;
557 }
558
559 // Ok check the work
560 for (size_t m = 0; m < nElementsN; m++) {
561 if (fabs(elemVectorN[m]) > 1.0E-13) {
562 throw CanteraError("IonsFromNeutralVPSSTP::initThermo",
563 "Simple formula matrix generation failed");
564 }
565 }
566 }
567
569}
570
572{
575 phaseNode["neutral-phase"] = neutralMoleculePhase_->name();
576 }
577}
578
579void IonsFromNeutralVPSSTP::setNeutralMoleculePhase(shared_ptr<ThermoPhase> neutral)
580{
581 neutralMoleculePhase_ = neutral;
582 geThermo = dynamic_cast<GibbsExcessVPSSTP*>(neutralMoleculePhase_.get());
592 y_.resize(numNeutralMoleculeSpecies_, 0.0);
593 dlnActCoeff_NeutralMolecule_.resize(numNeutralMoleculeSpecies_, 0.0);
594 dX_NeutralMolecule_.resize(numNeutralMoleculeSpecies_, 0.0);
595 for (size_t k = 0; k < nSpecies(); k++) {
596 providePDSS(k)->setParent(this, k);
597 }
598}
599
600shared_ptr<ThermoPhase> IonsFromNeutralVPSSTP::getNeutralMoleculePhase()
601{
603}
604
605bool IonsFromNeutralVPSSTP::addSpecies(shared_ptr<Species> spec)
606{
607 bool added = GibbsExcessVPSSTP::addSpecies(spec);
608 if (added) {
609 moleFractions_.push_back(0.0);
610 moleFractionsTmp_.push_back(0.0);
611 m_work.push_back(0.0);
612 fm_invert_ionForNeutral.push_back(npos);
614
615 if (spec->charge > 0) {
616 cationList_.push_back(m_kk-1);
617 } else if (spec->charge < 0) {
618 anionList_.push_back(m_kk-1);
619 } else {
620 passThroughList_.push_back(m_kk-1);
621 }
622
623 if (spec->input.hasKey("equation-of-state")) {
624 auto& ss = spec->input["equation-of-state"].getMapWhere(
625 "model", "ions-from-neutral-molecule");
626 if (ss.getBool("special-species", false)) {
628 }
629 }
630 }
631 return added;
632}
633
635{
636 size_t icat, jNeut;
637 // Get the activity coefficients of the neutral molecules
638 neutralMoleculePhase_->getLnActivityCoefficients(lnActCoeff_NeutralMolecule_.data());
639
640 switch (ionSolnType_) {
641 case cIonSolnType_PASSTHROUGH:
642 break;
643 case cIonSolnType_SINGLEANION:
644 // Do the cation list
645 for (size_t k = 0; k < cationList_.size(); k++) {
646 // Get the id for the next cation
647 icat = cationList_[k];
648 jNeut = fm_invert_ionForNeutral[icat];
649 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
651 }
652
653 // Do the anion list
654 icat = anionList_[0];
655 jNeut = fm_invert_ionForNeutral[icat];
656 lnActCoeff_Scaled_[icat]= 0.0;
657
658 // Do the list of neutral molecules
659 for (size_t k = 0; k < passThroughList_.size(); k++) {
660 icat = passThroughList_[k];
661 jNeut = fm_invert_ionForNeutral[icat];
663 }
664 break;
665
666 case cIonSolnType_SINGLECATION:
667 throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type");
668 break;
669 case cIonSolnType_MULTICATIONANION:
670 throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type");
671 break;
672 default:
673 throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type");
674 break;
675 }
676}
677
678void IonsFromNeutralVPSSTP::getdlnActCoeffds(const double dTds, const double* const dXds,
679 double* dlnActCoeffds) const
680{
681 size_t icat, jNeut;
682 // Get the activity coefficients of the neutral molecules
683 if (!geThermo) {
684 for (size_t k = 0; k < m_kk; k++) {
685 dlnActCoeffds[k] = dXds[k] / moleFractions_[k];
686 }
687 return;
688 }
689
690 getNeutralMoleculeMoleGrads(dXds, dX_NeutralMolecule_.data());
691
692 // All mole fractions returned to normal
693 geThermo->getdlnActCoeffds(dTds, dX_NeutralMolecule_.data(), dlnActCoeff_NeutralMolecule_.data());
694
695 switch (ionSolnType_) {
696 case cIonSolnType_PASSTHROUGH:
697 break;
698 case cIonSolnType_SINGLEANION:
699 // Do the cation list
700 for (size_t k = 0; k < cationList_.size(); k++) {
701 // Get the id for the next cation
702 icat = cationList_[k];
703 jNeut = fm_invert_ionForNeutral[icat];
704 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
705 dlnActCoeffds[icat] = dlnActCoeff_NeutralMolecule_[jNeut]/fmij;
706 }
707
708 // Do the anion list
709 icat = anionList_[0];
710 jNeut = fm_invert_ionForNeutral[icat];
711 dlnActCoeffds[icat]= 0.0;
712
713 // Do the list of neutral molecules
714 for (size_t k = 0; k < passThroughList_.size(); k++) {
715 icat = passThroughList_[k];
716 jNeut = fm_invert_ionForNeutral[icat];
717 dlnActCoeffds[icat] = dlnActCoeff_NeutralMolecule_[jNeut];
718 }
719 break;
720
721 case cIonSolnType_SINGLECATION:
722 throw CanteraError("IonsFromNeutralVPSSTP::getdlnActCoeffds", "Unimplemented type");
723 break;
724 case cIonSolnType_MULTICATIONANION:
725 throw CanteraError("IonsFromNeutralVPSSTP::getdlnActCoeffds", "Unimplemented type");
726 break;
727 default:
728 throw CanteraError("IonsFromNeutralVPSSTP::getdlnActCoeffds", "Unimplemented type");
729 break;
730 }
731}
732
734{
735 size_t icat, jNeut;
736
737 // Get the activity coefficients of the neutral molecules
738 if (!geThermo) {
739 dlnActCoeffdT_Scaled_.assign(m_kk, 0.0);
740 return;
741 }
742
744
745 switch (ionSolnType_) {
746 case cIonSolnType_PASSTHROUGH:
747 break;
748 case cIonSolnType_SINGLEANION:
749 // Do the cation list
750 for (size_t k = 0; k < cationList_.size(); k++) {
751 //! Get the id for the next cation
752 icat = cationList_[k];
753 jNeut = fm_invert_ionForNeutral[icat];
754 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
756 }
757
758 // Do the anion list
759 icat = anionList_[0];
760 jNeut = fm_invert_ionForNeutral[icat];
761 dlnActCoeffdT_Scaled_[icat]= 0.0;
762
763 // Do the list of neutral molecules
764 for (size_t k = 0; k < passThroughList_.size(); k++) {
765 icat = passThroughList_[k];
766 jNeut = fm_invert_ionForNeutral[icat];
768 }
769 break;
770
771 case cIonSolnType_SINGLECATION:
772 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeffdT", "Unimplemented type");
773 break;
774 case cIonSolnType_MULTICATIONANION:
775 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeffdT", "Unimplemented type");
776 break;
777 default:
778 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeffdT", "Unimplemented type");
779 break;
780 }
781}
782
784{
785 size_t icat, jNeut;
786
787 // Get the activity coefficients of the neutral molecules
788 if (!geThermo) {
789 dlnActCoeffdlnX_diag_.assign(m_kk, 0.0);
790 return;
791 }
792
794
795 switch (ionSolnType_) {
796 case cIonSolnType_PASSTHROUGH:
797 break;
798 case cIonSolnType_SINGLEANION:
799 // Do the cation list
800 for (size_t k = 0; k < cationList_.size(); k++) {
801 // Get the id for the next cation
802 icat = cationList_[k];
803 jNeut = fm_invert_ionForNeutral[icat];
804 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
806 }
807
808 // Do the anion list
809 icat = anionList_[0];
810 jNeut = fm_invert_ionForNeutral[icat];
811 dlnActCoeffdlnX_diag_[icat]= 0.0;
812
813 // Do the list of neutral molecules
814 for (size_t k = 0; k < passThroughList_.size(); k++) {
815 icat = passThroughList_[k];
816 jNeut = fm_invert_ionForNeutral[icat];
818 }
819 break;
820
821 case cIonSolnType_SINGLECATION:
822 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnX_diag", "Unimplemented type");
823 break;
824 case cIonSolnType_MULTICATIONANION:
825 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnX_diag", "Unimplemented type");
826 break;
827 default:
828 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnX_diag", "Unimplemented type");
829 break;
830 }
831}
832
834{
835 size_t icat, jNeut;
836
837 // Get the activity coefficients of the neutral molecules
838 if (!geThermo) {
839 dlnActCoeffdlnN_diag_.assign(m_kk, 0.0);
840 return;
841 }
842
844
845 switch (ionSolnType_) {
846 case cIonSolnType_PASSTHROUGH:
847 break;
848 case cIonSolnType_SINGLEANION:
849 // Do the cation list
850 for (size_t k = 0; k < cationList_.size(); k++) {
851 // Get the id for the next cation
852 icat = cationList_[k];
853 jNeut = fm_invert_ionForNeutral[icat];
854 double fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk];
856 }
857
858 // Do the anion list
859 icat = anionList_[0];
860 jNeut = fm_invert_ionForNeutral[icat];
861 dlnActCoeffdlnN_diag_[icat]= 0.0;
862
863 // Do the list of neutral molecules
864 for (size_t k = 0; k < passThroughList_.size(); k++) {
865 icat = passThroughList_[k];
866 jNeut = fm_invert_ionForNeutral[icat];
868 }
869 break;
870
871 case cIonSolnType_SINGLECATION:
872 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN_diag", "Unimplemented type");
873 break;
874 case cIonSolnType_MULTICATIONANION:
875 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN_diag", "Unimplemented type");
876 break;
877 default:
878 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN_diag", "Unimplemented type");
879 break;
880 }
881}
882
884{
885 size_t kcat = 0, kNeut = 0, mcat = 0, mNeut = 0;
886 double fmij = 0.0;
888 // Get the activity coefficients of the neutral molecules
889 if (!geThermo) {
890 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN", "dynamic cast failed");
891 }
892 size_t nsp_ge = geThermo->nSpecies();
894
895 switch (ionSolnType_) {
896 case cIonSolnType_PASSTHROUGH:
897 break;
898 case cIonSolnType_SINGLEANION:
899 // Do the cation list
900 for (size_t k = 0; k < cationList_.size(); k++) {
901 for (size_t m = 0; m < cationList_.size(); m++) {
902 kcat = cationList_[k];
903
904 kNeut = fm_invert_ionForNeutral[kcat];
905 fmij = fm_neutralMolec_ions_[kcat + kNeut * m_kk];
907
908 mcat = cationList_[m];
909 mNeut = fm_invert_ionForNeutral[mcat];
910 double mfmij = fm_neutralMolec_ions_[mcat + mNeut * m_kk];
911
912 dlnActCoeffdlnN_(kcat,mcat) = dlnActCoeffdlnN_NeutralMolecule_(kNeut,mNeut) * mfmij / fmij;
913
914 }
915 for (size_t m = 0; m < passThroughList_.size(); m++) {
916 mcat = passThroughList_[m];
917 mNeut = fm_invert_ionForNeutral[mcat];
918 dlnActCoeffdlnN_(kcat, mcat) = dlnActCoeffdlnN_NeutralMolecule_(kNeut, mNeut) / fmij;
919 }
920 }
921
922 // Do the anion list -> anion activity coefficient is one
923 kcat = anionList_[0];
924 kNeut = fm_invert_ionForNeutral[kcat];
925 for (size_t k = 0; k < m_kk; k++) {
926 dlnActCoeffdlnN_(kcat, k) = 0.0;
927 dlnActCoeffdlnN_(k, kcat) = 0.0;
928 }
929
930 // Do the list of neutral molecules
931 for (size_t k = 0; k < passThroughList_.size(); k++) {
932 kcat = passThroughList_[k];
933 kNeut = fm_invert_ionForNeutral[kcat];
935
936 for (size_t m = 0; m < m_kk; m++) {
937 mcat = passThroughList_[m];
938 mNeut = fm_invert_ionForNeutral[mcat];
939 dlnActCoeffdlnN_(kcat, mcat) = dlnActCoeffdlnN_NeutralMolecule_(kNeut, mNeut);
940 }
941
942 for (size_t m = 0; m < cationList_.size(); m++) {
943 mcat = cationList_[m];
944 mNeut = fm_invert_ionForNeutral[mcat];
945 dlnActCoeffdlnN_(kcat, mcat) = dlnActCoeffdlnN_NeutralMolecule_(kNeut,mNeut);
946 }
947 }
948 break;
949
950 case cIonSolnType_SINGLECATION:
951 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN", "Unimplemented type");
952 break;
953 case cIonSolnType_MULTICATIONANION:
954 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN", "Unimplemented type");
955 break;
956 default:
957 throw CanteraError("IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN", "Unimplemented type");
958 break;
959 }
960}
961
962}
Header for intermediate ThermoPhase object for phases which consist of ions whose thermodynamics is c...
Declarations for the class PDSS_IonsFromNeutral ( which handles calculations for a single ion in a fl...
Declaration for class Cantera::Species.
Headers for the factory class that can create known ThermoPhase objects (see Thermodynamic Properties...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
bool hasKey(const string &key) const
Returns true if the map contains an item named key.
Definition AnyMap.cpp:1423
const string & getString(const string &key, const string &default_) const
If key exists, return it as a string, otherwise return default_.
Definition AnyMap.cpp:1530
static AnyMap fromYamlFile(const string &name, const string &parent_name="")
Create an AnyMap from a YAML file.
Definition AnyMap.cpp:1771
void zero()
Set all of the entries to zero.
Definition Array.h:143
virtual void resize(size_t n, size_t m, double v=0.0)
Resize the array, and fill the new entries with 'v'.
Definition Array.cpp:47
Base class for exceptions thrown by Cantera classes.
GibbsExcessVPSSTP is a derived class of ThermoPhase that handles variable pressure standard state met...
Array2D dlnActCoeffdlnN_
Storage for the current derivative values of the gradients with respect to logarithm of the species m...
vector< double > lnActCoeff_Scaled_
Storage for the current values of the activity coefficients of the species.
vector< double > dlnActCoeffdlnX_diag_
Storage for the current derivative values of the gradients with respect to logarithm of the mole frac...
vector< double > moleFractions_
Storage for the current values of the mole fractions of the species.
vector< double > dlnActCoeffdT_Scaled_
Storage for the current derivative values of the gradients with respect to temperature of the log of ...
virtual void getdlnActCoeffdT(double *dlnActCoeffdT) const
Get the array of temperature derivatives of the log activity coefficients.
void compositionChanged() override
Apply changes to the state which are needed after the composition changes.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
vector< double > dlnActCoeffdlnN_diag_
Storage for the current derivative values of the gradients with respect to logarithm of the mole frac...
void getdlnActCoeffdlnN(const size_t ld, double *const dlnActCoeffdlnN) override
Get the array of derivatives of the log activity coefficients with respect to the log of the species ...
vector< double > lnActCoeff_NeutralMolecule_
Storage vector for the neutral molecule ln activity coefficients.
vector< size_t > anionList_
List of the species in this ThermoPhase which are anion species.
void getdlnActCoeffds(const double dTds, const double *const dXds, double *dlnActCoeffds) const override
Get the change in activity coefficients wrt changes in state (temp, mole fraction,...
double enthalpy_mole() const override
Return the Molar enthalpy. Units: J/kmol.
AnyMap m_rootNode
Root node of the AnyMap which contains this phase definition.
vector< double > dlnActCoeffdlnX_diag_NeutralMolecule_
Storage vector for the neutral molecule d ln activity coefficients dX - diagonal component.
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
void getDissociationCoeffs(vector< double > &fm_neutralMolec_ions, vector< double > &charges, vector< size_t > &neutMolIndex) const
Get the Salt Dissociation Coefficients.
virtual void calcNeutralMoleculeMoleFractions() const
Calculate neutral molecule mole fractions.
virtual void calcIonMoleFractions(double *const mf) const
Calculate ion mole fractions from neutral molecule mole fractions.
vector< size_t > fm_invert_ionForNeutral
Mapping between ion species and neutral molecule for quick invert.
vector< double > muNeutralMolecule_
Storage vector for the neutral molecule chemical potentials.
void s_update_dlnActCoeff_dlnN_diag() const
Update the derivative of the log of the activity coefficients wrt log(number of moles) - diagonal com...
void getParameters(AnyMap &phaseNode) const override
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
vector< double > fm_neutralMolec_ions_
Formula Matrix for composition of neutral molecules in terms of the molecules in this ThermoPhase.
size_t indexSpecialSpecies_
Index of special species.
IonsFromNeutralVPSSTP(const string &inputFile="", const string &id="")
Construct an IonsFromNeutralVPSSTP object from an input file.
double cv_mole() const override
Molar heat capacity at constant volume. Units: J/kmol/K.
vector< double > dlnActCoeffdlnN_diag_NeutralMolecule_
Storage vector for the neutral molecule d ln activity coefficients dlnN.
vector< size_t > cationList_
List of the species in this ThermoPhase which are cation species.
Array2D dlnActCoeffdlnN_NeutralMolecule_
Storage vector for the neutral molecule d ln activity coefficients dlnN.
void s_update_dlnActCoeff_dlnN() const
Update the derivative of the log of the activity coefficients wrt log(number of moles) - diagonal com...
double entropy_mole() const override
Molar entropy. Units: J/kmol/K.
void calcDensity() override
Calculate the density of the mixture using the partial molar volumes and mole fractions as input.
shared_ptr< ThermoPhase > neutralMoleculePhase_
This is a pointer to the neutral Molecule Phase.
vector< double > dlnActCoeffdT_NeutralMolecule_
Storage vector for the neutral molecule d ln activity coefficients dT.
size_t numNeutralMoleculeSpecies_
Number of neutral molecule species.
double cp_mole() const override
Molar heat capacity at constant pressure. Units: J/kmol/K.
void compositionChanged() override
Apply changes to the state which are needed after the composition changes.
vector< double > moleFractionsTmp_
Temporary mole fraction vector.
void s_update_dlnActCoeff_dlnX_diag() const
Update the derivative of the log of the activity coefficients wrt log(mole fraction)
void s_update_lnActCoeff() const
Update the activity coefficients.
double gibbs_mole() const override
Molar Gibbs function. Units: J/kmol.
vector< size_t > passThroughList_
List of the species in this ThermoPhase which are passed through to the neutralMoleculePhase ThermoPh...
void getdlnActCoeffdlnX_diag(double *dlnActCoeffdlnX_diag) const override
Get the array of ln mole fraction derivatives of the log activity coefficients - diagonal component o...
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void getdlnActCoeffdlnN_diag(double *dlnActCoeffdlnN_diag) const override
Get the array of log species mole number derivatives of the log activity coefficients.
void getActivityCoefficients(double *ac) const override
Get the array of non-dimensional molar-based activity coefficients at the current solution temperatur...
void s_update_dlnActCoeffdT() const
Update the temperature derivative of the ln activity coefficients.
void setParameters(const AnyMap &phaseNode, const AnyMap &rootNode=AnyMap()) override
Set equation of state parameters from an AnyMap phase description.
vector< double > NeutralMolecMoleFractions_
Mole fractions using the Neutral Molecule Mole fraction basis.
void getNeutralMoleculeMoleGrads(const double *const dx, double *const dy) const
Calculate neutral molecule mole fractions.
IonSolnType_enumType ionSolnType_
Ion solution type.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies for the species in the mixture.
void getdlnActCoeffdlnN(const size_t ld, double *const dlnActCoeffdlnN) override
Get the array of derivatives of the log activity coefficients with respect to the log of the species ...
virtual void setParent(VPStandardStateTP *phase, size_t k)
Set the parent VPStandardStateTP object of this PDSS object.
Definition PDSS.h:411
void assignDensity(const double density_)
Set the internally stored constant density (kg/m^3) of the phase.
Definition Phase.cpp:718
size_t nSpecies() const
Returns the number of species in the phase.
Definition Phase.h:245
size_t m_kk
Number of species in the phase.
Definition Phase.h:947
double temperature() const
Temperature (K).
Definition Phase.h:662
double nAtoms(size_t k, size_t m) const
Number of atoms of element m in species k.
Definition Phase.cpp:103
size_t nElements() const
Number of elements.
Definition Phase.cpp:30
double mean_X(const double *const Q) const
Evaluate the mole-fraction-weighted mean of an array Q.
Definition Phase.cpp:737
vector< double > m_speciesCharge
Vector of species charges. length m_kk.
Definition Phase.h:958
virtual void setParameters(const AnyMap &phaseNode, const AnyMap &rootNode=AnyMap())
Set equation of state parameters from an AnyMap phase description.
virtual void getParameters(AnyMap &phaseNode) const
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
virtual void getdlnActCoeffdlnN_diag(double *dlnActCoeffdlnN_diag) const
Get the array of log species mole number derivatives of the log activity coefficients.
double RT() const
Return the Gas Constant multiplied by the current temperature.
virtual void getPartialMolarCp(double *cpbar) const
Return an array of partial molar heat capacities for the species in the mixture.
virtual void getdlnActCoeffds(const double dTds, const double *const dXds, double *dlnActCoeffds) const
Get the change in activity coefficients wrt changes in state (temp, mole fraction,...
virtual void initThermo()
Initialize the ThermoPhase object after all species have been set up.
void initThermoFile(const string &inputFile, const string &id)
Initialize a ThermoPhase object using an input file.
virtual void getdlnActCoeffdlnX_diag(double *dlnActCoeffdlnX_diag) const
Get the array of ln mole fraction derivatives of the log activity coefficients - diagonal component o...
AnyMap m_input
Data supplied via setParameters.
double pressure() const override
Returns the current pressure of the phase.
void getEntropy_R(double *sr) const override
Get the array of nondimensional Entropy functions for the standard state species at the current T and...
void getEnthalpy_RT(double *hrt) const override
Get the nondimensional Enthalpy functions for the species at their standard states at the current T a...
#define AssertTrace(expr)
Assertion must be true or an error is thrown.
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition ct_defs.h:120
shared_ptr< ThermoPhase > newThermo(const AnyMap &phaseNode, const AnyMap &rootNode)
Create a new ThermoPhase object and initialize it.
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const size_t npos
index returned by functions to indicate "no position"
Definition ct_defs.h:195
const vector< string > & elementNames()
Get a vector of the names of the elements defined in Cantera.
Definition Elements.cpp:227
const double SmallNumber
smallest number to compare to zero.
Definition ct_defs.h:158
static double factorOverlap(const vector< string > &elnamesVN, const vector< double > &elemVectorN, const size_t nElementsN, const vector< string > &elnamesVI, const vector< double > &elemVectorI, const size_t nElementsI)
Return the factor overlap.
void warn_deprecated(const string &source, const AnyBase &node, const string &message)
A deprecation warning for syntax in an input file.
Definition AnyMap.cpp:1926
Contains declarations for string manipulation functions within Cantera.