[cig-commits] commit 1891 by buerg to /var/svn/dealii/aspect

[dealii.demon at gmail.com]

Revision 1891

Guard composition contraint.

U   trunk/aspire/include/aspect/simulator.h
U   trunk/aspire/source/material_model/cantera.cc
U   trunk/aspire/source/simulator/core.cc
U   trunk/aspire/source/simulator/helper_functions.cc


http://www.dealii.org/websvn/revision.php?repname=Aspect+Repository&path=%2F&rev=1891&peg=1891

Diff:
Modified: trunk/aspire/include/aspect/simulator.h
===================================================================
--- trunk/aspire/include/aspect/simulator.h	2013-09-13 09:32:35 UTC (rev 1890)
+++ trunk/aspire/include/aspect/simulator.h	2013-09-13 19:45:22 UTC (rev 1891)
@@ -701,6 +701,8 @@
        */
       void denormalize_pressure(LinearAlgebra::BlockVector &vector);
       
+      void normalize_compositions ();
+      
       /**
        * Compute the maximal velocity throughout the domain. This is needed
        * to compute the size of the time step.

Modified: trunk/aspire/source/material_model/cantera.cc
===================================================================
--- trunk/aspire/source/material_model/cantera.cc	2013-09-13 09:32:35 UTC (rev 1890)
+++ trunk/aspire/source/material_model/cantera.cc	2013-09-13 19:45:22 UTC (rev 1891)
@@ -59,11 +59,11 @@
           
           for (unsigned int i = 0; i < n_species; ++i)
             {
-              out.compositional_sources[q][i] = molar_mass[i] * net[i];
+              out.compositional_sources[q][i] = molar_mass[i] * net[i] / out.densities[q];
               out.thermal_sources[q] += enthalpy_RT[i] * out.compositional_sources[q][i];
             }
           
-          out.thermal_sources[q] *= -1.0;
+          out.thermal_sources[q] *= -1.0 * ideal_gas_mix->enthalpy_mass ();
           out.thermal_expansion_coefficients[q] = ideal_gas_mix->thermalExpansionCoeff ();
           out.viscosities[q] = transport->viscosity ();
         }

Modified: trunk/aspire/source/simulator/core.cc
===================================================================
--- trunk/aspire/source/simulator/core.cc	2013-09-13 09:32:35 UTC (rev 1890)
+++ trunk/aspire/source/simulator/core.cc	2013-09-13 19:45:22 UTC (rev 1891)
@@ -1047,6 +1047,7 @@
               solve_composition (c);
             }
 
+          normalize_compositions ();
           break;
         }
         

Modified: trunk/aspire/source/simulator/helper_functions.cc
===================================================================
--- trunk/aspire/source/simulator/helper_functions.cc	2013-09-13 09:32:35 UTC (rev 1890)
+++ trunk/aspire/source/simulator/helper_functions.cc	2013-09-13 19:45:22 UTC (rev 1891)
@@ -80,18 +80,18 @@
               + (parameters.n_compositional_fields > 0 ?
                  parameters.n_compositional_fields * current_constraints_composition[0].memory_consumption () : 
                  0)) / mb << std::endl
-          << "* 3 Matrices " << (system_matrix_velocity.memory_consumption ()
-                                 + system_matrix_temperature.memory_consumption ()
-                                 + system_matrix_composition.memory_consumption ()) / mb << std::endl
+          << "* 3 Matrices " << (/*system_matrix_velocity.memory_consumption ()
+                                 +*/ system_matrix_temperature.memory_consumption ()
+                                /* + system_matrix_composition.memory_consumption ()*/) / mb << std::endl
           << "* 15 Vectors " << 5 * (solution_velocity.memory_consumption ()
                                      + solution_temperature.memory_consumption ()
                                      + solution_composition.memory_consumption ()) / mb << std::endl
-          << "* preconditioner " << (system_preconditioner_matrix.memory_consumption()
+          << "* preconditioner "// << (system_preconditioner_matrix.memory_consumption()
                                      /*+ Amg_preconditioner->memory_consumption()*/
                                      /*+Mp_preconditioner->memory_consumption()
-                                                                      +T_preconditioner->memory_consumption()*/)/mb
+                                                                      +T_preconditioner->memory_consumption()*///)/mb
           << std::endl
-          << "  - matrix " << system_preconditioner_matrix.memory_consumption()/mb << std::endl
+          << "  - matrix "/* << system_preconditioner_matrix.memory_consumption()/mb*/ << std::endl
           //<< "  - prec vel " << Amg_preconditioner->memory_consumption()/mb << std::endl
           << "  - prec mass " << 0/*Mp_preconditioner->memory_consumption()/mb*/ << std::endl
           << "  - prec T " << 0/*T_preconditioner->memory_consumption()/mb*/ << std::endl
@@ -549,6 +549,34 @@
 
     vector.block(1).add(correction, pressure_shape_function_integrals.block(1));
   }
+  
+  template <int dim>
+  void Simulator<dim>::normalize_compositions ()
+  {
+    for (unsigned int i = 0; i < solution_composition.block (0).size (); ++i)
+    {
+      double sum = 0.0;
+      
+      for (unsigned int c = 0; c < parameters.n_compositional_fields; ++c)
+      {
+        if (solution_composition.block (c) (i) < 0.0)
+          solution_composition.block (c) (i) = 0.0;
+        
+        else
+          if (solution_composition.block (c) (i) > 1.0)
+            solution_composition.block (c) (i) = 1.0;
+        
+        sum += solution_composition.block (c) (i);
+      }
+      
+      solution_composition.compress (VectorOperation::add);
+      
+      for (unsigned int c = 0; c < parameters.n_compositional_fields; ++c)
+        solution_composition.block (c) (i) /= sum;
+    }
+    
+    
+  }
 }
 
 // explicit instantiation of the functions we implement in this file
@@ -557,6 +585,7 @@
 #define INSTANTIATE(dim) \
   template void Simulator<dim>::normalize_pressure(LinearAlgebra::BlockVector &vector); \
   template void Simulator<dim>::denormalize_pressure(LinearAlgebra::BlockVector &vector); \
+  template void Simulator<dim>::normalize_compositions (); \
   template double Simulator<dim>::get_maximal_velocity (const LinearAlgebra::BlockVector &solution) const; \
   template std::pair<double,double> Simulator<dim>::get_extrapolated_composition_range (const unsigned int) const; \
   template std::pair<double,double> Simulator<dim>::get_extrapolated_temperature_range () const; \