[cig-commits] r1298 - in trunk/aspect: include/aspect/material_model source source/compositional_initial_conditions source/material_model source/postprocess source/postprocess/visualization source/simulator

dannberg at dealii.org dannberg at dealii.org
Fri Oct 19 08:41:38 PDT 2012


Author: dannberg
Date: 2012-10-19 09:41:38 -0600 (Fri, 19 Oct 2012)
New Revision: 1298

Modified:
   trunk/aspect/include/aspect/material_model/interface.h
   trunk/aspect/source/adiabatic_conditions.cc
   trunk/aspect/source/compositional_initial_conditions/interface.cc
   trunk/aspect/source/material_model/interface.cc
   trunk/aspect/source/material_model/simple.cc
   trunk/aspect/source/postprocess/heat_flux_statistics.cc
   trunk/aspect/source/postprocess/table_heat_flux_statistics.cc
   trunk/aspect/source/postprocess/visualization/density.cc
   trunk/aspect/source/postprocess/visualization/friction_heating.cc
   trunk/aspect/source/postprocess/visualization/seismic_vp.cc
   trunk/aspect/source/postprocess/visualization/seismic_vs.cc
   trunk/aspect/source/postprocess/visualization/specific_heat.cc
   trunk/aspect/source/postprocess/visualization/thermal_expansivity.cc
   trunk/aspect/source/postprocess/visualization/thermodynamic_phase.cc
   trunk/aspect/source/postprocess/visualization/viscosity.cc
   trunk/aspect/source/postprocess/visualization/viscosity_ratio.cc
   trunk/aspect/source/simulator/assembly.cc
   trunk/aspect/source/simulator/core.cc
   trunk/aspect/source/simulator/helper_functions.cc
Log:
indent

Modified: trunk/aspect/include/aspect/material_model/interface.h
===================================================================
--- trunk/aspect/include/aspect/material_model/interface.h	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/include/aspect/material_model/interface.h	2012-10-19 15:41:38 UTC (rev 1298)
@@ -420,7 +420,7 @@
           std::vector<double> pressure;
           std::vector<std::vector<double> > composition;
           std::vector<SymmetricTensor<2,dim> > strain_rate;
-       };
+        };
 
         struct MaterialModelOutputs
         {
@@ -435,7 +435,7 @@
           bool is_compressible;
         };
 
-        virtual void compute_parameters(const MaterialModelInputs & in, MaterialModelOutputs & out);
+        virtual void compute_parameters(const MaterialModelInputs &in, MaterialModelOutputs &out);
 
         /**
          * @name Functions used in dealing with run-time parameters

Modified: trunk/aspect/source/adiabatic_conditions.cc
===================================================================
--- trunk/aspect/source/adiabatic_conditions.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/adiabatic_conditions.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -66,7 +66,7 @@
         //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 c=0;c<n_compositional_fields;++c)
+        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

Modified: trunk/aspect/source/compositional_initial_conditions/interface.cc
===================================================================
--- trunk/aspect/source/compositional_initial_conditions/interface.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/compositional_initial_conditions/interface.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -109,10 +109,10 @@
           }
           prm.leave_subsection ();
 
-	  Interface<dim> *plugin = std_cxx1x::get<dim>(registered_plugins).create_plugin (model_name, prm);
-	  plugin->initialize (geometry_model);
-	  return plugin;
-	}
+          Interface<dim> *plugin = std_cxx1x::get<dim>(registered_plugins).create_plugin (model_name, prm);
+          plugin->initialize (geometry_model);
+          return plugin;
+        }
     }
 
 

Modified: trunk/aspect/source/material_model/interface.cc
===================================================================
--- trunk/aspect/source/material_model/interface.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/material_model/interface.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -254,41 +254,41 @@
 
     template <int dim>
     Interface<dim>::MaterialModelInputs::MaterialModelInputs(unsigned int n_points, unsigned int n_comp)
-      {
-        position.resize(n_points);
-        temperature.resize(n_points);
-        pressure.resize(n_points);
-        composition.resize(n_points);
-        for(unsigned int i=0;i<n_points;++i)
-          composition[i].resize(n_comp);
-        strain_rate.resize(n_points);
-      }
+    {
+      position.resize(n_points);
+      temperature.resize(n_points);
+      pressure.resize(n_points);
+      composition.resize(n_points);
+      for (unsigned int i=0; i<n_points; ++i)
+        composition[i].resize(n_comp);
+      strain_rate.resize(n_points);
+    }
 
     template <int dim>
     Interface<dim>::MaterialModelOutputs::MaterialModelOutputs(unsigned int n_points)
-      {
-        viscosities.resize(n_points);
-        densities.resize(n_points);
-        thermal_expansion_coefficients.resize(n_points);
-        specific_heat.resize(n_points);
-        thermal_conductivities.resize(n_points);
-        compressibilities.resize(n_points);
-      }
+    {
+      viscosities.resize(n_points);
+      densities.resize(n_points);
+      thermal_expansion_coefficients.resize(n_points);
+      specific_heat.resize(n_points);
+      thermal_conductivities.resize(n_points);
+      compressibilities.resize(n_points);
+    }
 
     template <int dim>
     void
     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.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.is_compressible = is_compressible();
-      }
+        {
+          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.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.is_compressible = is_compressible();
+        }
     }
 
 

Modified: trunk/aspect/source/material_model/simple.cc
===================================================================
--- trunk/aspect/source/material_model/simple.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/material_model/simple.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -43,11 +43,11 @@
       const double temperature_dependence = std::max(std::min(std::exp(-thermal_viscosity_exponent*delta_temp/reference_T),1e2),1e-2);
 
       return temperature_dependence * eta;
-/*      return (this->n_compositional_fields()>0
-          ?
-          (6.5*composition[0]+1) * eta
-          :
-          eta);*/
+      /*      return (this->n_compositional_fields()>0
+                ?
+                (6.5*composition[0]+1) * eta
+                :
+                eta);*/
     }
 
 
@@ -122,11 +122,11 @@
              const Point<dim> &) const
     {
       return (this->n_compositional_fields()>0
-          ?
-          100.0 * compositional_fields[0] + reference_rho *
+              ?
+              100.0 * compositional_fields[0] + reference_rho *
               (1 - thermal_alpha * (temperature - reference_T))
-          :
-          reference_rho * (1 - thermal_alpha * (temperature - reference_T)));
+              :
+              reference_rho * (1 - thermal_alpha * (temperature - reference_T)));
     }
 
 

Modified: trunk/aspect/source/postprocess/heat_flux_statistics.cc
===================================================================
--- trunk/aspect/source/postprocess/heat_flux_statistics.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/heat_flux_statistics.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -55,10 +55,10 @@
       const FEValuesExtractors::Scalar temperature (dim+1);
       std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-      for (unsigned int c=0;c<this->n_compositional_fields();++c)
+      for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
         {
-        const FEValuesExtractors::Scalar temp(dim+2+c);
-        compositional_fields.push_back(temp);
+          const FEValuesExtractors::Scalar temp(dim+2+c);
+          compositional_fields.push_back(temp);
         }
 
       std::vector<Tensor<1,dim> > temperature_gradients (quadrature_formula.size());
@@ -94,7 +94,7 @@
                                                                  temperature_values);
                 fe_face_values[pressure].get_function_values (this->get_solution(),
                                                               pressure_values);
-                for(unsigned int c=0;c<this->n_compositional_fields();++c)
+                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]);
 
@@ -102,8 +102,8 @@
                 for (unsigned int q=0; q<fe_face_values.n_quadrature_points; ++q)
                   {
                     this->get_composition_values_at_q_point (composition_values,
-                                                           q,
-                                                           composition_values_at_q_point);
+                                                             q,
+                                                             composition_values_at_q_point);
 
                     const double thermal_conductivity
                       = this->get_material_model().thermal_conductivity(temperature_values[q],

Modified: trunk/aspect/source/postprocess/table_heat_flux_statistics.cc
===================================================================
--- trunk/aspect/source/postprocess/table_heat_flux_statistics.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/table_heat_flux_statistics.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -54,10 +54,10 @@
       const FEValuesExtractors::Scalar temperature (dim+1);
       std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-      for (unsigned int c=0;c<this->n_compositional_fields();++c)
+      for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
         {
-        const FEValuesExtractors::Scalar temp(dim+2+c);
-        compositional_fields.push_back(temp);
+          const FEValuesExtractors::Scalar temp(dim+2+c);
+          compositional_fields.push_back(temp);
         }
 
       std::vector<Tensor<1,dim> > temperature_gradients (quadrature_formula.size());
@@ -99,7 +99,7 @@
                                                                  temperature_values);
                 fe_face_values[pressure].get_function_values (this->get_solution(),
                                                               pressure_values);
-                for(unsigned int c=0;c<this->n_compositional_fields();++c)
+                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]);
 
@@ -107,8 +107,8 @@
                 for (unsigned int q=0; q<fe_face_values.n_quadrature_points; ++q)
                   {
                     this->get_composition_values_at_q_point (composition_values,
-                                                           q,
-                                                           composition_values_at_q_point);
+                                                             q,
+                                                             composition_values_at_q_point);
 
                     const double thermal_conductivity
                       = this->get_material_model().thermal_conductivity(temperature_values[q],

Modified: trunk/aspect/source/postprocess/visualization/density.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/density.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/density.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -63,7 +63,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/friction_heating.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/friction_heating.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/friction_heating.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -77,7 +77,7 @@
                  strain_rate);
 
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/seismic_vp.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/seismic_vp.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/seismic_vp.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -62,7 +62,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/seismic_vs.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/seismic_vs.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/seismic_vs.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -62,7 +62,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/specific_heat.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/specific_heat.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/specific_heat.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -62,7 +62,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/thermal_expansivity.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/thermal_expansivity.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/thermal_expansivity.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -63,7 +63,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/thermodynamic_phase.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/thermodynamic_phase.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/thermodynamic_phase.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -62,7 +62,7 @@
             const double pressure    = uh[q][dim];
             const double temperature = uh[q][dim+1];
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/postprocess/visualization/viscosity.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/viscosity.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/viscosity.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -77,7 +77,7 @@
                  strain_rate);
 
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int i=0;i<this->n_compositional_fields();++i)
+            for (unsigned int i=0; i<this->n_compositional_fields(); ++i)
               composition[i] = uh[q][dim+2+i];
 
             computed_quantities[q](0) = this->get_material_model().viscosity(temperature,

Modified: trunk/aspect/source/postprocess/visualization/viscosity_ratio.cc
===================================================================
--- trunk/aspect/source/postprocess/visualization/viscosity_ratio.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/postprocess/visualization/viscosity_ratio.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -77,7 +77,7 @@
                  strain_rate);
 
             std::vector<double> composition(this->n_compositional_fields());
-            for (unsigned int c=0;c<this->n_compositional_fields();++c)
+            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,

Modified: trunk/aspect/source/simulator/assembly.cc
===================================================================
--- trunk/aspect/source/simulator/assembly.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/simulator/assembly.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -265,7 +265,7 @@
           old_temperature_laplacians(quadrature.size()),
           old_old_temperature_laplacians(quadrature.size()),
           old_composition_values(n_compositional_fields,
-                                     std::vector<double>(quadrature.size())),
+                                 std::vector<double>(quadrature.size())),
           old_old_composition_values(n_compositional_fields,
                                      std::vector<double>(quadrature.size())),
           current_temperature_values(quadrature.size()),
@@ -673,21 +673,22 @@
                                       const std::vector<Point<dim> >     &evaluation_points,
                                       double                             &max_residual,
                                       double                             &max_velocity) const
-    {
+  {
 
     const unsigned int n_q_points = old_temperature.size();
 
     typename MaterialModel::Interface<dim>::MaterialModelInputs inputs (n_q_points, parameters.n_compositional_fields);
     typename MaterialModel::Interface<dim>::MaterialModelOutputs outputs (n_q_points);
 
-    for (unsigned int q=0; q<n_q_points; ++q) {
-      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 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;
-    }
+    for (unsigned int q=0; q<n_q_points; ++q)
+      {
+        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 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;
+      }
 
     material_model->compute_parameters(inputs,outputs);
 
@@ -770,7 +771,7 @@
         max_residual = std::max (residual,        max_residual);
         max_velocity = std::max (std::sqrt (u*u), max_velocity);
       }
-    }
+  }
 
   template <int dim>
   void
@@ -787,7 +788,7 @@
                                       const unsigned int                  composition_index,
                                       double                             &max_residual,
                                       double                             &max_velocity) const
-    {
+{
 
     const unsigned int n_q_points = old_composition.size();
 
@@ -807,8 +808,8 @@
         const double kappa = 0;
 
         const double kappa_Delta_C = kappa
-                                 * (old_composition_laplacians[q] +
-                                    old_old_composition_laplacians[q]) / 2;
+                                     * (old_composition_laplacians[q] +
+                                        old_old_composition_laplacians[q]) / 2;
 
         double residual
           = std::abs(dC_dt + u_grad_C - kappa_Delta_C);
@@ -818,7 +819,7 @@
         max_residual = std::max (residual,        max_residual);
         max_velocity = std::max (std::sqrt (u*u), max_velocity);
       }
-    }
+  }
 
 
 
@@ -855,32 +856,33 @@
 
     AssertIndexRange(index,parameters.n_compositional_fields+1);
 
-    if(index == 0) {
-      //make sure that all arguments we need for computing the residual are passed
-      Assert (old_strain_rates.size() > 0 && old_old_strain_rates.size() > 0
-              && old_pressure.size() > 0 && old_old_pressure.size() > 0,
-              ExcMessage ("Not enough parameters to calculate artificial viscosity "
-                          "for the temperature equation."));
+    if (index == 0)
+      {
+        //make sure that all arguments we need for computing the residual are passed
+        Assert (old_strain_rates.size() > 0 && old_old_strain_rates.size() > 0
+                && old_pressure.size() > 0 && old_old_pressure.size() > 0,
+                ExcMessage ("Not enough parameters to calculate artificial viscosity "
+                            "for the temperature equation."));
 
-      compute_temperature_system_residual(old_field,
-                                          old_old_field,
-                                          old_field_grads,
-                                          old_old_field_grads,
-                                          old_field_laplacians,
-                                          old_old_field_laplacians,
-                                          old_velocity_values,
-                                          old_old_velocity_values,
-                                          old_strain_rates,
-                                          old_old_strain_rates,
-                                          old_pressure,
-                                          old_old_pressure,
-                                          old_composition,
-                                          old_old_composition,
-                                          average_field,
-                                          evaluation_points,
-                                          max_residual,
-                                          max_velocity);
-    }
+        compute_temperature_system_residual(old_field,
+                                            old_old_field,
+                                            old_field_grads,
+                                            old_old_field_grads,
+                                            old_field_laplacians,
+                                            old_old_field_laplacians,
+                                            old_velocity_values,
+                                            old_old_velocity_values,
+                                            old_strain_rates,
+                                            old_old_strain_rates,
+                                            old_pressure,
+                                            old_old_pressure,
+                                            old_composition,
+                                            old_old_composition,
+                                            average_field,
+                                            evaluation_points,
+                                            max_residual,
+                                            max_velocity);
+      }
     else
       compute_composition_system_residual(old_field,
                                           old_old_field,
@@ -901,7 +903,7 @@
       // we don't have sensible timesteps during the first two iterations
       return max_viscosity;
     else
-    {
+      {
         Assert (old_time_step > 0, ExcInternalError());
 
         double entropy_viscosity;
@@ -937,33 +939,33 @@
 
     unsigned int n_q_points = material_model_inputs.temperature.size();
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      input_composition.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        input_composition.push_back(temp);
       }
 
     input_finite_element_values[input_temperature].get_function_values (current_linearization_point,
                                                                         material_model_inputs.temperature);
     input_finite_element_values[input_pressure].get_function_values(current_linearization_point,
                                                                     material_model_inputs.pressure);
-    if(compute_strainrate)
+    if (compute_strainrate)
       input_finite_element_values[input_velocities].get_function_symmetric_gradients(current_linearization_point,
                                                                                      material_model_inputs.strain_rate);
 
     // the values of the compositional fields are stored as blockvectors for each field
     // we have to extract them in this structure
     std::vector<std::vector<double> > composition_values (parameters.n_compositional_fields,
-        std::vector<double> (n_q_points));
+                                                          std::vector<double> (n_q_points));
 
-    for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    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]);
+                                                                            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 c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int q=0; q<n_q_points; ++q)
+      for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
         material_model_inputs.composition[q][c] = composition_values[c][q];
 
   }
@@ -1345,10 +1347,10 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     const unsigned int dofs_per_cell = scratch.finite_element_values.get_fe().dofs_per_cell;
@@ -1392,12 +1394,13 @@
     scratch.finite_element_values[velocities].get_function_values(current_linearization_point,
                                                                   scratch.current_velocity_values);
 
-    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]);
-    }
+    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]);
+      }
 
 
     // TODO: Compute artificial viscosity once per timestep instead of each time
@@ -1464,10 +1467,10 @@
         const double viscosity            = scratch.material_model_outputs.viscosities[q];
         const bool is_compressible        = scratch.material_model_outputs.is_compressible;
         const double compressibility      = (is_compressible
-                                            ?
-                                            scratch.material_model_outputs.compressibilities[q]
-                                            :
-                                            std::numeric_limits<double>::quiet_NaN() );
+                                             ?
+                                             scratch.material_model_outputs.compressibilities[q]
+                                             :
+                                             std::numeric_limits<double>::quiet_NaN() );
 
         const Tensor<1,dim>
         gravity = gravity_model->gravity_vector (scratch.finite_element_values.quadrature_point(q));

Modified: trunk/aspect/source/simulator/core.cc
===================================================================
--- trunk/aspect/source/simulator/core.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/simulator/core.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -794,17 +794,17 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> composition;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      composition.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        composition.push_back(temp);
       }
 
     //Velocity|Temperature|Normalized density and temperature|Weighted density and temperature|Density c_p temperature
 
     // compute density error
     if (parameters.refinement_strategy != "Temperature" && parameters.refinement_strategy != "Velocity"
-                                                        && parameters.refinement_strategy != "Composition")
+        && parameters.refinement_strategy != "Composition")
       {
         bool lookup_rho_c_p_T = (parameters.refinement_strategy == "Density c_p temperature");
 
@@ -831,9 +831,9 @@
         // the values of the compositional fields are stored as blockvectors for each field
         // we have to extract them in this structure
         std::vector<std::vector<double> > prelim_composition_values (parameters.n_compositional_fields,
-            std::vector<double> (quadrature.size()));
+                                                                     std::vector<double> (quadrature.size()));
         std::vector<std::vector<double> > composition_values (quadrature.size(),
-            std::vector<double> (parameters.n_compositional_fields));
+                                                              std::vector<double> (parameters.n_compositional_fields));
 
 
         typename DoFHandler<dim>::active_cell_iterator
@@ -847,13 +847,13 @@
                                                        pressure_values);
               fe_values[temperature].get_function_values (solution,
                                                           temperature_values);
-              for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+              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 c=0;c<parameters.n_compositional_fields;++c)
+              for (unsigned int q=0; q<quadrature.size(); ++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);
@@ -944,7 +944,7 @@
     // compute the errors for composition solution
     if (parameters.refinement_strategy == "Composition")
       {
-        for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+        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+c] = true;
@@ -961,7 +961,7 @@
       }
     else
       {
-        for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+        for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
           estimated_error_per_cell_C[c] = 0;
       }
 
@@ -1027,7 +1027,7 @@
 
           for (unsigned int i=0; i<estimated_error_per_cell.size(); ++i)
             estimated_error_per_cell(i)
-	      = estimated_error_per_cell_T(i)*(1.0+estimated_error_per_cell_rho(i));
+              = estimated_error_per_cell_T(i)*(1.0+estimated_error_per_cell_rho(i));
         }
       else if (parameters.refinement_strategy == "Density c_p temperature")
         {
@@ -1188,7 +1188,7 @@
                   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
+                    = 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);
                 }
 

Modified: trunk/aspect/source/simulator/helper_functions.cc
===================================================================
--- trunk/aspect/source/simulator/helper_functions.cc	2012-10-19 15:30:58 UTC (rev 1297)
+++ trunk/aspect/source/simulator/helper_functions.cc	2012-10-19 15:41:38 UTC (rev 1298)
@@ -259,7 +259,7 @@
                                                               const unsigned int                      q,
                                                               std::vector<double>                    &composition_values_at_q_point) const
   {
-    for(unsigned int k=0; k < composition_values_at_q_point.size(); ++k)
+    for (unsigned int k=0; k < composition_values_at_q_point.size(); ++k)
       composition_values_at_q_point[k] = composition_values[k][q];
   }
 
@@ -627,10 +627,10 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     std::vector<SymmetricTensor<2,dim> > strain_rates(n_q_points);
@@ -655,7 +655,7 @@
                                                       temperature_values);
           fe_values[velocities].get_function_symmetric_gradients (this->solution,
                                                                   strain_rates);
-          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             fe_values[compositional_fields[c]].get_function_values(this->solution,
                                                                    composition_values[c]);
 
@@ -829,10 +829,10 @@
     const FEValuesExtractors::Scalar pressure (dim);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     std::vector<double> pressure_values(n_q_points);
@@ -851,7 +851,7 @@
           fe_values.reinit (cell);
           fe_values[pressure].get_function_values (this->solution,
                                                    pressure_values);
-          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          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)
@@ -910,10 +910,10 @@
     const FEValuesExtractors::Scalar pressure (dim);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     std::vector<double> pressure_values(n_q_points);
@@ -932,7 +932,7 @@
           fe_values.reinit (cell);
           fe_values[pressure].get_function_values (this->solution,
                                                    pressure_values);
-          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          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)
@@ -991,10 +991,10 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     std::vector<double> pressure_values(n_q_points);
@@ -1015,7 +1015,7 @@
                                                    pressure_values);
           fe_values[temperature].get_function_values (this->solution,
                                                       temperature_values);
-          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             fe_values[compositional_fields[c]].get_function_values(this->solution,
                                                                    composition_values[c]);
 
@@ -1066,10 +1066,10 @@
     const FEValuesExtractors::Scalar temperature (dim+1);
     std::vector<FEValuesExtractors::Scalar> compositional_fields;
 
-    for (unsigned int c=0;c<parameters.n_compositional_fields;++c)
+    for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
       {
-      const FEValuesExtractors::Scalar temp(dim+2+c);
-      compositional_fields.push_back(temp);
+        const FEValuesExtractors::Scalar temp(dim+2+c);
+        compositional_fields.push_back(temp);
       }
 
     std::vector<double> pressure_values(n_q_points);
@@ -1091,7 +1091,7 @@
                                                    pressure_values);
           fe_values[temperature].get_function_values (this->solution,
                                                       temperature_values);
-          for(unsigned int c=0;c<parameters.n_compositional_fields;++c)
+          for (unsigned int c=0; c<parameters.n_compositional_fields; ++c)
             fe_values[compositional_fields[c]].get_function_values(this->solution,
                                                                    composition_values[c]);
 



More information about the CIG-COMMITS mailing list