[cig-commits] commit 2266 by buerg to /var/svn/dealii/aspect
dealii.demon at gmail.com
dealii.demon at gmail.com
Wed Jan 8 10:50:23 PST 2014
Revision 2266
Add simple integration for reaction.
U trunk/aspire/include/aspect/material_model/cantera.h
U trunk/aspire/source/material_model/cantera.cc
http://www.dealii.org/websvn/revision.php?repname=Aspect+Repository&path=%2F&rev=2266&peg=2266
Diff:
Modified: trunk/aspire/include/aspect/material_model/cantera.h
===================================================================
--- trunk/aspire/include/aspect/material_model/cantera.h 2013-12-28 21:07:07 UTC (rev 2265)
+++ trunk/aspire/include/aspect/material_model/cantera.h 2014-01-08 18:50:21 UTC (rev 2266)
@@ -28,9 +28,8 @@
#include <aspect/material_model/interface.h>
#include <aspect/simulator_access.h>
#include <cantera/IdealGasMix.h>
+#include <cantera/integrators.h>
#include <cantera/transport.h>
-#include <cantera/zeroD/IdealGasConstPressureReactor.h>
-#include <cantera/zeroD/ReactorNet.h>
namespace aspect
{
@@ -54,11 +53,6 @@
*/
/**
- * Constructor
- */
- Cantera ();
-
- /**
* Destructor
*/
~Cantera ();
@@ -121,7 +115,7 @@
MaterialModelOutputs &out) const;
virtual void react (const MaterialModelInputs &in,
const double time_step,
- MaterialModelOutputs &out) const;
+ MaterialModelOutputs &out);
/**
@@ -167,12 +161,12 @@
*/
private:
- ::Cantera::IdealGasConstPressureReactor* reactor;
::Cantera::IdealGasMix* ideal_gas_mix;
- ::Cantera::ReactorNet reactor_net;
::Cantera::Transport* transport;
+ double* composition_values;
double* enthalpy_RT;
- double* net;
+ double* mass_fractions;
+ double* net_production_rates;
double reference_eta;
double reference_pressure;
double reference_rho;
@@ -181,6 +175,8 @@
double reference_thermal_conductivity;
std::vector<double> molecular_weights;
std::vector<std::string> composition_names;
+ unsigned int n_reactions;
+ unsigned int n_species;
};
}
Modified: trunk/aspire/source/material_model/cantera.cc
===================================================================
--- trunk/aspire/source/material_model/cantera.cc 2013-12-28 21:07:07 UTC (rev 2265)
+++ trunk/aspire/source/material_model/cantera.cc 2014-01-08 18:50:21 UTC (rev 2266)
@@ -32,19 +32,14 @@
namespace MaterialModel
{
template <int dim>
- Cantera<dim>::Cantera () {
- reactor = new ::Cantera::IdealGasConstPressureReactor ();
- reactor_net.addReactor (reactor);
- }
-
- template <int dim>
Cantera<dim>::~Cantera () {
+ delete composition_values;
delete enthalpy_RT;
delete ideal_gas_mix;
- delete reactor;
+ delete net_production_rates;
delete transport;
}
-
+
template <int dim>
void
Cantera<dim>::evaluate (const MaterialModelInputs &in, MaterialModelOutputs &out) const
@@ -65,22 +60,46 @@
const double factor = -1.0 * ::Cantera::GasConstant * out.densities[q] * out.densities[q] * in.temperature[q] / reference_pressure;
- for (unsigned int i = 0; i < in.composition[0].size (); ++i)
+ for (unsigned int i = 0; i < n_species; ++i)
out.density_derivatives[q][i + 1] = factor / molecular_weights[i];
}
}
+
template <int dim>
void
Cantera<dim>::react (const MaterialModelInputs &in,
- const double time_step,
- MaterialModelOutputs &out) const
+ const double time_step,
+ MaterialModelOutputs &out)
{
Assert (in.composition[0].size () >= 1, ExcMessage ("Need compositional field!"));
for (unsigned int q = 0; q < out.viscosities.size (); ++q)
{
+ for (unsigned int i = 0; i < n_species; ++i)
+ composition_values[i] = in.composition[q][i];
+ for (unsigned int t = 0; t < 10; ++t)
+ {
+ ideal_gas_mix->setState_TPY (in.temperature[q], reference_pressure, composition_values);
+ ideal_gas_mix->getNetProductionRates (net_production_rates);
+
+ const double density = ideal_gas_mix->density ();
+
+ for (unsigned int i = 0;i < n_species; ++i)
+ composition_values[i] += 0.1 * time_step * net_production_rates[i] * molecular_weights[i] / density;
+ }
+
+ double sum = 0.0;
+
+ for (unsigned int i = 0; i < n_species; ++i)
+ {
+ out.compositional_increments[q][i] = composition_values[i] - in.composition[q][i];
+ sum += ::Cantera::GasConstant * in.temperature[q] * enthalpy_RT[i] * out.compositional_increments[q][i] / molecular_weights[i];
+ }
+
+ ideal_gas_mix->setState_TPY (in.temperature[q], reference_pressure, composition_values);
+ out.thermal_increments[q] = -1.0 * sum / ideal_gas_mix->cp_mass ();
}
}
@@ -229,18 +248,19 @@
reference_thermal_alpha = prm.get_double ("Reference thermal expansion coefficient");
reference_thermal_conductivity = prm.get_double ("Reference thermal conductivity");
ideal_gas_mix = new ::Cantera::IdealGasMix (prm.get ("XML file"));
-
- const unsigned int n_species = ideal_gas_mix->nSpecies ();
-
+ n_reactions = ideal_gas_mix->nReactions ();
+ n_species = ideal_gas_mix->nSpecies ();
composition_names.resize (n_species);
for (unsigned int i = 0; i < n_species; ++i)
composition_names[i] = ideal_gas_mix->speciesName (i);
+ composition_values = new double[n_species];
enthalpy_RT = new double[n_species];
+ net_production_rates = new double[n_species];
+ ideal_gas_mix->getEnthalpy_RT (enthalpy_RT);
molecular_weights.resize (n_species);
ideal_gas_mix->getMolecularWeights (molecular_weights);
- reactor->setThermoMgr (*ideal_gas_mix);
transport = newTransportMgr ("Mix", ideal_gas_mix);
}
prm.leave_subsection ();
@@ -252,9 +272,9 @@
void
Cantera<dim>::get_composition_names (std::vector<std::string>& names) const
{
- AssertDimension (names.size (), composition_names.size ());
+ AssertDimension (names.size (), n_species);
- for (unsigned int i = 0; i < composition_names.size (); ++i)
+ for (unsigned int i = 0; i < n_species; ++i)
names[i] = composition_names[i];
}
}
More information about the CIG-COMMITS
mailing list