[cig-commits] r1294 - in branches/active_compositions: include/aspect/material_model source source/material_model source/postprocess source/postprocess/visualization source/simulator

dannberg at dealii.org dannberg at dealii.org
Thu Oct 18 15:46:52 PDT 2012


Author: dannberg
Date: 2012-10-18 16:46:52 -0600 (Thu, 18 Oct 2012)
New Revision: 1294

Modified:
   branches/active_compositions/include/aspect/material_model/interface.h
   branches/active_compositions/source/adiabatic_conditions.cc
   branches/active_compositions/source/material_model/interface.cc
   branches/active_compositions/source/postprocess/depth_average.cc
   branches/active_compositions/source/postprocess/heat_flux_statistics.cc
   branches/active_compositions/source/postprocess/table_heat_flux_statistics.cc
   branches/active_compositions/source/postprocess/visualization.cc
   branches/active_compositions/source/postprocess/visualization/density.cc
   branches/active_compositions/source/postprocess/visualization/friction_heating.cc
   branches/active_compositions/source/postprocess/visualization/seismic_vp.cc
   branches/active_compositions/source/postprocess/visualization/seismic_vs.cc
   branches/active_compositions/source/postprocess/visualization/specific_heat.cc
   branches/active_compositions/source/postprocess/visualization/thermal_expansivity.cc
   branches/active_compositions/source/postprocess/visualization/thermodynamic_phase.cc
   branches/active_compositions/source/postprocess/visualization/viscosity_ratio.cc
   branches/active_compositions/source/simulator/assembly.cc
   branches/active_compositions/source/simulator/core.cc
   branches/active_compositions/source/simulator/helper_functions.cc
Log:
Remove Vp,Vs and thermodynamic phases from compute_parameters in the material model interface and use the index c in all loops counting the compositional fields

Modified: branches/active_compositions/include/aspect/material_model/interface.h
===================================================================
--- branches/active_compositions/include/aspect/material_model/interface.h	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/include/aspect/material_model/interface.h	2012-10-18 22:46:52 UTC (rev 1294)
@@ -429,16 +429,13 @@
           std::vector<double> viscosities;
           std::vector<double> densities;
           std::vector<double> thermal_expansion_coefficients;
-          std::vector<double> seismic_Vp;
-          std::vector<double> seismic_Vs;
           std::vector<double> specific_heat;
           std::vector<double> thermal_conductivities;
           std::vector<double> compressibilities;
-          std::vector<int> thermodynamic_phases;
           bool is_compressible;
         };
 
-        virtual void compute_parameters(MaterialModelInputs & in, MaterialModelOutputs & out);
+        virtual void compute_parameters(const MaterialModelInputs & in, MaterialModelOutputs & out);
 
         /**
          * @name Functions used in dealing with run-time parameters

Modified: branches/active_compositions/source/adiabatic_conditions.cc
===================================================================
--- branches/active_compositions/source/adiabatic_conditions.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/adiabatic_conditions.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -66,8 +66,8 @@
         //TODO: we look up the composition at the representative point, but we should
         // use averaged compositional values here. Right?
         std::vector<double> initial_composition(n_compositional_fields);
-        for (unsigned int k=0;k<n_compositional_fields;++k)
-          initial_composition[k] = compositional_initial_conditions.initial_composition(representative_point, k);
+        for (unsigned int c=0;c<n_compositional_fields;++c)
+          initial_composition[c] = compositional_initial_conditions.initial_composition(representative_point, c);
 
         // get material parameters and the magnitude of gravity. we assume
         // that gravity always points along the depth direction. this

Modified: branches/active_compositions/source/material_model/interface.cc
===================================================================
--- branches/active_compositions/source/material_model/interface.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/material_model/interface.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -270,29 +270,23 @@
         viscosities.resize(n_points);
         densities.resize(n_points);
         thermal_expansion_coefficients.resize(n_points);
-        seismic_Vp.resize(n_points);
-        seismic_Vs.resize(n_points);
         specific_heat.resize(n_points);
         thermal_conductivities.resize(n_points);
         compressibilities.resize(n_points);
-        thermodynamic_phases.resize(n_points);
       }
 
     template <int dim>
     void
-    Interface<dim>::compute_parameters(struct MaterialModelInputs &in, struct MaterialModelOutputs &out)
+    Interface<dim>::compute_parameters(const struct MaterialModelInputs &in, struct MaterialModelOutputs &out)
     {
       for (unsigned int i=0; i < in.temperature.size(); ++i)
       {
         out.viscosities[i]                    = viscosity                     (in.temperature[i], in.pressure[i], in.composition[i], in.strain_rate[i], in.position[i]);
         out.densities[i]                      = density                       (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
         out.thermal_expansion_coefficients[i] = thermal_expansion_coefficient (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
-        out.seismic_Vp[i]                     = seismic_Vp                    (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
-        out.seismic_Vs[i]                     = seismic_Vs                    (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
         out.specific_heat[i]                  = specific_heat                 (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
         out.thermal_conductivities[i]         = thermal_conductivity          (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
         out.compressibilities[i]              = compressibility               (in.temperature[i], in.pressure[i], in.composition[i], in.position[i]);
-        out.thermodynamic_phases[i]           = thermodynamic_phase           (in.temperature[i], in.pressure[i], in.composition[i]);
         out.is_compressible = is_compressible();
       }
     }

Modified: branches/active_compositions/source/postprocess/depth_average.cc
===================================================================
--- branches/active_compositions/source/postprocess/depth_average.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/depth_average.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -75,8 +75,8 @@
         // add temperature and the compositional fields that follow
         // it immediately
         this->get_depth_average_temperature(temp[i++]);
-        for (unsigned int k=0; k<this->n_compositional_fields(); ++k)
-          this->get_depth_average_composition(k, temp[i++]);
+        for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+          this->get_depth_average_composition(c, temp[i++]);
         this->get_adiabatic_conditions().get_adiabatic_temperature_profile(temp[i++],100);
         this->get_depth_average_velocity_magnitude(temp[i++]);
         this->get_depth_average_sinking_velocity(temp[i++]);
@@ -112,8 +112,8 @@
           const std::string filename = this->get_output_directory() + "depthaverage.plt";
           std::ofstream f (filename.c_str());
           f << "# time, depth, avg T";
-          for (unsigned int k=0; k<this->n_compositional_fields(); ++k)
-            f << ", avg C_" << k+1;
+          for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+            f << ", avg C_" << c+1;
           f << ", adiabatic T, velocity magnitude, avg sinking velocity, avg Vs, avg Vp, avg viscosity" << std::endl;
           f << std::scientific;
 

Modified: branches/active_compositions/source/postprocess/heat_flux_statistics.cc
===================================================================
--- branches/active_compositions/source/postprocess/heat_flux_statistics.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/heat_flux_statistics.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -55,9 +55,9 @@
       const FEValuesExtractors::Scalar temperature (dim+1);
       std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-      for (unsigned int q=0;q<this->n_compositional_fields();++q)
+      for (unsigned int c=0;c<this->n_compositional_fields();++c)
         {
-        const FEValuesExtractors::Scalar temp(dim+2+q);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
         compositional_fields.push_back(temp);
         }
 
@@ -94,9 +94,9 @@
                                                                  temperature_values);
                 fe_face_values[pressure].get_function_values (this->get_solution(),
                                                               pressure_values);
-                for(unsigned int i=0;i<this->n_compositional_fields();++i)
-                  fe_face_values[compositional_fields[i]].get_function_values(this->get_solution(),
-                                                                              composition_values[i]);
+                for(unsigned int c=0;c<this->n_compositional_fields();++c)
+                  fe_face_values[compositional_fields[c]].get_function_values(this->get_solution(),
+                                                                              composition_values[c]);
 
                 double local_normal_flux = 0;
                 for (unsigned int q=0; q<fe_face_values.n_quadrature_points; ++q)

Modified: branches/active_compositions/source/postprocess/table_heat_flux_statistics.cc
===================================================================
--- branches/active_compositions/source/postprocess/table_heat_flux_statistics.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/table_heat_flux_statistics.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -54,9 +54,9 @@
       const FEValuesExtractors::Scalar temperature (dim+1);
       std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-      for (unsigned int q=0;q<this->n_compositional_fields();++q)
+      for (unsigned int c=0;c<this->n_compositional_fields();++c)
         {
-        const FEValuesExtractors::Scalar temp(dim+2+q);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
         compositional_fields.push_back(temp);
         }
 
@@ -99,9 +99,9 @@
                                                                  temperature_values);
                 fe_face_values[pressure].get_function_values (this->get_solution(),
                                                               pressure_values);
-                for(unsigned int i=0;i<this->n_compositional_fields();++i)
-                  fe_face_values[compositional_fields[i]].get_function_values(this->get_solution(),
-                                                                              composition_values[i]);
+                for(unsigned int c=0;c<this->n_compositional_fields();++c)
+                  fe_face_values[compositional_fields[c]].get_function_values(this->get_solution(),
+                                                                              composition_values[c]);
 
                 double local_normal_flux = 0;
                 for (unsigned int q=0; q<fe_face_values.n_quadrature_points; ++q)

Modified: branches/active_compositions/source/postprocess/visualization/density.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/density.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/density.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -63,8 +63,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().density(temperature,
                                                                            pressure,

Modified: branches/active_compositions/source/postprocess/visualization/friction_heating.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/friction_heating.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/friction_heating.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -77,8 +77,8 @@
                  strain_rate);
 
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = 2 * this->get_material_model().viscosity(temperature,
                                                                                  pressure,

Modified: branches/active_compositions/source/postprocess/visualization/seismic_vp.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/seismic_vp.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/seismic_vp.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -62,8 +62,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().seismic_Vp(temperature,
                                                                               pressure,

Modified: branches/active_compositions/source/postprocess/visualization/seismic_vs.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/seismic_vs.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/seismic_vs.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -62,8 +62,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().seismic_Vs(temperature,
                                                                               pressure,

Modified: branches/active_compositions/source/postprocess/visualization/specific_heat.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/specific_heat.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/specific_heat.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -62,8 +62,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().specific_heat(temperature,
                                                                                  pressure,

Modified: branches/active_compositions/source/postprocess/visualization/thermal_expansivity.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/thermal_expansivity.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/thermal_expansivity.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -63,8 +63,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().thermal_expansion_coefficient(temperature,
                                         pressure,

Modified: branches/active_compositions/source/postprocess/visualization/thermodynamic_phase.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/thermodynamic_phase.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/thermodynamic_phase.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -62,8 +62,8 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().thermodynamic_phase(temperature,
                                                                                        pressure,

Modified: branches/active_compositions/source/postprocess/visualization/viscosity_ratio.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization/viscosity_ratio.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization/viscosity_ratio.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -77,8 +77,8 @@
                  strain_rate);
 
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
-              composition[i] = uh[q][dim+2+i];
+            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+              composition[c] = uh[q][dim+2+c];
 
             computed_quantities[q](0) = this->get_material_model().viscosity_ratio(temperature,
                                                                                    pressure,

Modified: branches/active_compositions/source/postprocess/visualization.cc
===================================================================
--- branches/active_compositions/source/postprocess/visualization.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/postprocess/visualization.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -129,8 +129,8 @@
       std::vector<std::string> solution_names (dim, "velocity");
       solution_names.push_back ("p");
       solution_names.push_back ("T");
-      for (unsigned int i=0; i<this->n_compositional_fields(); ++i)
-        solution_names.push_back ("C_" + boost::lexical_cast<std::string>(i+1));
+      for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+        solution_names.push_back ("C_" + boost::lexical_cast<std::string>(c+1));
 
 
       std::vector<DataComponentInterpretation::DataComponentInterpretation>
@@ -138,7 +138,7 @@
                       DataComponentInterpretation::component_is_part_of_vector);
       interpretation.push_back (DataComponentInterpretation::component_is_scalar);
       interpretation.push_back (DataComponentInterpretation::component_is_scalar);
-      for (unsigned int i=0; i<this->n_compositional_fields(); ++i)
+      for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
         interpretation.push_back (DataComponentInterpretation::component_is_scalar);
 
       data_out.add_data_vector (this->get_solution(),

Modified: branches/active_compositions/source/simulator/assembly.cc
===================================================================
--- branches/active_compositions/source/simulator/assembly.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/simulator/assembly.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -684,8 +684,8 @@
       inputs.temperature[q] = (old_temperature[q] + old_old_temperature[q]) / 2;
       inputs.position[q] = evaluation_points[q];
       inputs.pressure[q] = (old_pressure[q] + old_old_pressure[q]) / 2;
-      for (unsigned int i=0; i<parameters.n_compositional_fields; ++i)
-        inputs.composition[q][i] = (old_composition[i][q] + old_old_composition[i][q]) / 2;
+      for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
+        inputs.composition[q][c] = (old_composition[c][q] + old_old_composition[c][q]) / 2;
       inputs.strain_rate[q] = (old_strain_rates[q] + old_old_strain_rates[q]) / 2;
     }
 
@@ -937,9 +937,9 @@
 
     unsigned int n_q_points = material_model_inputs.temperature.size();
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       input_composition.push_back(temp);
       }
 
@@ -956,15 +956,15 @@
     std::vector<std::vector<double> > composition_values (parameters.n_compositional_fields,
         std::vector<double> (n_q_points));
 
-    for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-      input_finite_element_values[input_composition[i]].get_function_values(current_linearization_point,
-                                                                      composition_values[i]);
+    for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+      input_finite_element_values[input_composition[c]].get_function_values(current_linearization_point,
+                                                                      composition_values[c]);
 
     // then we copy these values to exchange the inner and outer vector, because for the material
     // model we need a vector with values of all the compositional fields for every quadrature point
     for(unsigned int q=0;q<n_q_points;++q)
-      for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-        material_model_inputs.composition[q][i] = composition_values[i][q];
+      for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+        material_model_inputs.composition[q][c] = composition_values[c][q];
 
   }
 
@@ -1345,9 +1345,9 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -1392,11 +1392,11 @@
     scratch.finite_element_values[velocities].get_function_values(current_linearization_point,
                                                                   scratch.current_velocity_values);
 
-    for(unsigned int q=0;q<parameters.n_compositional_fields;++q) {
-      scratch.finite_element_values[compositional_fields[q]].get_function_values(old_solution,
-                                                                                  scratch.old_composition_values[q]);
-      scratch.finite_element_values[compositional_fields[q]].get_function_values(old_old_solution,
-                                                                                  scratch.old_old_composition_values[q]);
+    for(unsigned int c=0;c<parameters.n_compositional_fields;++c) {
+      scratch.finite_element_values[compositional_fields[c]].get_function_values(old_solution,
+                                                                                  scratch.old_composition_values[c]);
+      scratch.finite_element_values[compositional_fields[c]].get_function_values(old_old_solution,
+                                                                                  scratch.old_old_composition_values[c]);
     }
 
 

Modified: branches/active_compositions/source/simulator/core.cc
===================================================================
--- branches/active_compositions/source/simulator/core.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/simulator/core.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -572,10 +572,10 @@
                        n_T = system_dofs_per_block[2];
     unsigned int       n_C_sum = 0;
     std::vector<unsigned int> n_C (parameters.n_compositional_fields+1);
-    for (unsigned int i=0; i<parameters.n_compositional_fields; ++i)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-        n_C[i] = system_dofs_per_block[i+3];
-        n_C_sum += n_C[i];
+        n_C[c] = system_dofs_per_block[c+3];
+        n_C_sum += n_C[c];
       }
 
 
@@ -603,8 +603,8 @@
           << n_u + n_p + n_T + n_C_sum
           << " (" << n_u << '+' << n_p << '+'<< n_T;
 
-    for (unsigned int i=0; i<parameters.n_compositional_fields; ++i)
-      pcout << '+' << n_C[i];
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
+      pcout << '+' << n_C[c];
 
     pcout <<')'
           << std::endl
@@ -632,13 +632,13 @@
       {
         unsigned int n_C_so_far = 0;
 
-        for (unsigned int i=0; i<parameters.n_compositional_fields; ++i)
+        for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
           {
             system_partitioning.push_back(system_index_set.get_view(n_u+n_p+n_T+n_C_so_far,
-                                                                    n_u+n_p+n_T+n_C_so_far+n_C[i]));
+                                                                    n_u+n_p+n_T+n_C_so_far+n_C[c]));
             system_relevant_partitioning.push_back(system_relevant_set.get_view(n_u+n_p+n_T+n_C_so_far,
-                                                                                n_u+n_p+n_T+n_C_so_far+n_C[i]));
-            n_C_so_far += n_C[i];
+                                                                                n_u+n_p+n_T+n_C_so_far+n_C[c]));
+            n_C_so_far += n_C[c];
           }
       }
     }
@@ -795,9 +795,9 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> composition;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       composition.push_back(temp);
       }
 
@@ -848,14 +848,14 @@
                                                        pressure_values);
               fe_values[temperature].get_function_values (solution,
                                                           temperature_values);
-              for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-                fe_values[composition[i]].get_function_values (solution,
-                                                               prelim_composition_values[i]);
+              for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+                fe_values[composition[c]].get_function_values (solution,
+                                                               prelim_composition_values[c]);
               // then we copy these values to exchange the inner and outer vector, because for the material
               // model we need a vector with values of all the compositional fields for every quadrature point
               for(unsigned int q=0;q<quadrature.size();++q)
-                for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-                  composition_values[q][i] = prelim_composition_values[i][q];
+                for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+                  composition_values[q][c] = prelim_composition_values[c][q];
 
               cell->get_dof_indices (local_dof_indices);
 
@@ -945,15 +945,15 @@
     // compute the errors for composition solution
     if (parameters.refinement_strategy == "Composition")
       {
-        for (unsigned int i=0;i<parameters.n_compositional_fields;++i)
+        for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
           {
             std::vector<bool> composition_component (dim+2+parameters.n_compositional_fields, false);
-            composition_component[dim+2+i] = true;
+            composition_component[dim+2+c] = true;
             KellyErrorEstimator<dim>::estimate (dof_handler,
                                                 QGauss<dim-1>(parameters.composition_degree+1),
                                                 typename FunctionMap<dim>::type(),
                                                 solution,
-                                                estimated_error_per_cell_C[i],
+                                                estimated_error_per_cell_C[c],
                                                 composition_component,
                                                 0,
                                                 0,
@@ -962,8 +962,8 @@
       }
     else
       {
-        for (unsigned int i=0;i<parameters.n_compositional_fields;++i)
-          estimated_error_per_cell_C[i] = 0;
+        for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          estimated_error_per_cell_C[c] = 0;
       }
 
     // compute the errors for the stokes solution
@@ -1039,9 +1039,9 @@
         {
           for (unsigned int i=0; i<estimated_error_per_cell.size(); ++i)
             estimated_error_per_cell(i) = 0;
-          for (unsigned int k=0; k<parameters.n_compositional_fields; ++k)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             for (unsigned int i=0; i<estimated_error_per_cell.size(); ++i)
-              estimated_error_per_cell(i) += estimated_error_per_cell_C[k](i);
+              estimated_error_per_cell(i) += estimated_error_per_cell_C[c](i);
         }
       else
         AssertThrow(false, ExcNotImplemented());
@@ -1148,12 +1148,12 @@
 
           current_linearization_point.block(2) = solution.block(2);
 
-          for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             {
-              assemble_composition_system (n);
-              build_composition_preconditioner(n);
-              solve_temperature_or_composition(1+n); // this is correct, 0 would be temperature
-              current_linearization_point.block(3+n) = solution.block(3+n);
+              assemble_composition_system (c);
+              build_composition_preconditioner(c);
+              solve_temperature_or_composition(1+c); // this is correct, 0 would be temperature
+              current_linearization_point.block(3+c) = solution.block(3+c);
             }
 
           assemble_stokes_system();
@@ -1184,13 +1184,13 @@
               rebuild_stokes_matrix = true;
               std::vector<double> composition_residual (parameters.n_compositional_fields,0);
 
-              for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
+              for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
                 {
-                  assemble_composition_system (n);
-                  build_composition_preconditioner(n);
-                  composition_residual[n]
-		    = solve_temperature_or_composition(1+n); // 1+n is correct, because 0 is for temperature
-                  current_linearization_point.block(3+n) = solution.block(3+n);
+                  assemble_composition_system (c);
+                  build_composition_preconditioner(c);
+                  composition_residual[c]
+		    = solve_temperature_or_composition(1+c); // 1+n is correct, because 0 is for temperature
+                  current_linearization_point.block(3+c) = solution.block(3+c);
                 }
 
               assemble_stokes_system();
@@ -1204,8 +1204,8 @@
               pcout << "   Nonlinear residuals: " << temperature_residual
                     << ", " << stokes_residual;
 
-              for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
-                pcout << ", " << composition_residual[n];
+              for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
+                pcout << ", " << composition_residual[c];
 
               pcout << std::endl
                     << std::endl;
@@ -1213,16 +1213,16 @@
               if (iteration == 0)
                 {
                   initial_temperature_residual = temperature_residual;
-                  for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
-                    initial_composition_residual[n] = composition_residual[n];
+                  for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
+                    initial_composition_residual[c] = composition_residual[c];
                   initial_stokes_residual      = stokes_residual;
                 }
               else
                 {
 //TODO: make this a parameter in the input file
                   double max = 0.0;
-                  for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
-                    max = std::max(composition_residual[n]/initial_composition_residual[n],max);
+                  for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
+                    max = std::max(composition_residual[c]/initial_composition_residual[c],max);
                   if (std::max(std::max(stokes_residual/initial_stokes_residual,
                                         temperature_residual/initial_temperature_residual),
                                max) < 1e-3)
@@ -1243,10 +1243,10 @@
           assemble_temperature_system ();
           solve_temperature_or_composition(0);
 
-          for (unsigned int n=0; n<parameters.n_compositional_fields; ++n)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             {
-              assemble_composition_system (n);
-              solve_temperature_or_composition(1+n); // 1+n is correct, because 0 is the temperature
+              assemble_composition_system (c);
+              solve_temperature_or_composition(1+c); // 1+n is correct, because 0 is the temperature
             }
 
           // ...and then iterate the solution

Modified: branches/active_compositions/source/simulator/helper_functions.cc
===================================================================
--- branches/active_compositions/source/simulator/helper_functions.cc	2012-10-18 22:12:21 UTC (rev 1293)
+++ branches/active_compositions/source/simulator/helper_functions.cc	2012-10-18 22:46:52 UTC (rev 1294)
@@ -627,9 +627,9 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -655,9 +655,9 @@
                                                       temperature_values);
           fe_values[velocities].get_function_symmetric_gradients (this->solution,
                                                                   strain_rates);
-          for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-            fe_values[compositional_fields[i]].get_function_values(this->solution,
-                                                                   composition_values[i]);
+          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+            fe_values[compositional_fields[c]].get_function_values(this->solution,
+                                                                   composition_values[c]);
 
           for (unsigned int q = 0; q < n_q_points; ++q)
             {
@@ -829,9 +829,9 @@
     const FEValuesExtractors::Scalar pressure (dim);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -851,9 +851,9 @@
           fe_values.reinit (cell);
           fe_values[pressure].get_function_values (this->solution,
                                                    pressure_values);
-          for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-            fe_values[compositional_fields[i]].get_function_values(this->solution,
-                                                                   composition_values[i]);
+          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+            fe_values[compositional_fields[c]].get_function_values(this->solution,
+                                                                   composition_values[c]);
           for (unsigned int q = 0; q < n_q_points; ++q)
             {
               const double depth = geometry_model->depth(fe_values.quadrature_point(q));
@@ -910,9 +910,9 @@
     const FEValuesExtractors::Scalar pressure (dim);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -932,9 +932,9 @@
           fe_values.reinit (cell);
           fe_values[pressure].get_function_values (this->solution,
                                                    pressure_values);
-          for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-            fe_values[compositional_fields[i]].get_function_values(this->solution,
-                                                                   composition_values[i]);
+          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+            fe_values[compositional_fields[c]].get_function_values(this->solution,
+                                                                   composition_values[c]);
           for (unsigned int q = 0; q < n_q_points; ++q)
             {
               const double depth = geometry_model->depth(fe_values.quadrature_point(q));
@@ -991,9 +991,9 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -1015,9 +1015,9 @@
                                                    pressure_values);
           fe_values[temperature].get_function_values (this->solution,
                                                       temperature_values);
-          for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-            fe_values[compositional_fields[i]].get_function_values(this->solution,
-                                                                   composition_values[i]);
+          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+            fe_values[compositional_fields[c]].get_function_values(this->solution,
+                                                                   composition_values[c]);
 
 
           extract_composition_values_at_q_point (composition_values,
@@ -1066,9 +1066,9 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int q=0;q<parameters.n_compositional_fields;++q)
+    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+q);
+      const FEValuesExtractors::Scalar temp(dim+2+c);
       compositional_fields.push_back(temp);
       }
 
@@ -1091,9 +1091,9 @@
                                                    pressure_values);
           fe_values[temperature].get_function_values (this->solution,
                                                       temperature_values);
-          for(unsigned int i=0;i<parameters.n_compositional_fields;++i)
-            fe_values[compositional_fields[i]].get_function_values(this->solution,
-                                                                   composition_values[i]);
+          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+            fe_values[compositional_fields[c]].get_function_values(this->solution,
+                                                                   composition_values[c]);
 
           extract_composition_values_at_q_point (composition_values,
                                                  0,



More information about the CIG-COMMITS mailing list