]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update library to use new cross_product signature
authorMatthias Maier <tamiko@43-1.org>
Tue, 15 Sep 2015 03:37:00 +0000 (22:37 -0500)
committerMatthias Maier <tamiko@43-1.org>
Tue, 15 Sep 2015 17:04:05 +0000 (12:04 -0500)
include/deal.II/numerics/vector_tools.templates.h
source/base/geometry_info.cc
source/fe/mapping_fe_field.cc
source/fe/mapping_q_generic.cc
source/grid/tria_accessor.cc
source/grid/tria_boundary.cc

index 7737b179bad5884d83da43dd2547c285d046eda9..99c45a06b7856ac2a84a6b417e856b206a51fdec 100644 (file)
@@ -2875,8 +2875,8 @@ namespace VectorTools
 
           // now compute the
           // two normals
-          cross_product (orthonormals[0], vector, tmp);
-          cross_product (orthonormals[1], vector, orthonormals[0]);
+          orthonormals[0] = cross_product (vector, tmp);
+          orthonormals[1] = cross_product (vector, orthonormals[0]);
 
           break;
         }
@@ -4343,9 +4343,9 @@ namespace VectorTools
           // associated with this face. We also must include the residuals from the
           // shape funcations associated with edges.
           Tensor<1, dim> tmp;
-          Tensor<1, dim> cross_product_i,
-                 cross_product_j,
-                 cross_product_rhs;
+          Tensor<1, dim> cross_product_i;
+          Tensor<1, dim> cross_product_j;
+          Tensor<1, dim> cross_product_rhs;
 
           // Loop to construct face linear system.
           for (unsigned int q_point = 0;
@@ -4394,33 +4394,25 @@ namespace VectorTools
                   const unsigned int j_face_idx = associated_face_dof_to_face_dof[j];
                   const unsigned int cell_j = fe.face_to_cell_index (j_face_idx, face);
 
-                  cross_product(cross_product_j,
-                                normal_vector,
-                                fe_face_values[vec].value(cell_j, q_point));
+                  cross_product_j =
+                    cross_product(normal_vector,
+                                  fe_face_values[vec].value(cell_j, q_point));
 
                   for (unsigned int i = 0; i < associated_face_dofs; ++i)
                     {
                       const unsigned int i_face_idx = associated_face_dof_to_face_dof[i];
                       const unsigned int cell_i = fe.face_to_cell_index (i_face_idx, face);
-                      cross_product(cross_product_i,
-                                    normal_vector,
-                                    fe_face_values[vec].value(cell_i, q_point));
-
-                      face_matrix(i,j)
-                      += fe_face_values.JxW (q_point)
-                         *cross_product_i
-                         *cross_product_j;
+                      cross_product_i =
+                        cross_product(normal_vector,
+                                      fe_face_values[vec].value(cell_i, q_point));
 
+                      face_matrix(i, j) += fe_face_values.JxW(q_point) *
+                                           cross_product_i * cross_product_j;
                     }
                   // compute rhs
-                  cross_product(cross_product_rhs,
-                                normal_vector,
-                                tmp);
-                  face_rhs(j)
-                  += fe_face_values.JxW (q_point)
-                     *cross_product_rhs
-                     *cross_product_j;
-
+                  cross_product_rhs = cross_product(normal_vector, tmp);
+                  face_rhs(j) += fe_face_values.JxW(q_point) *
+                                 cross_product_rhs * cross_product_j;
                 }
             }
 
@@ -5727,7 +5719,7 @@ namespace VectorTools
                     // we get here only for dim==3, but at least one isn't
                     // quite smart enough to notice this and warns when
                     // compiling the function in 2d
-                    cross_product (tangent, normals[0], normals[dim-2]);
+                    tangent = cross_product (normals[0], normals[dim-2]);
                     break;
                   default:
                     Assert (false, ExcNotImplemented());
index da75f2ba316693867903e02f92ce0605d777dcbe..db0b24e10dee16218f1135218e1fcbccc09f1425 100644 (file)
@@ -1820,10 +1820,7 @@ namespace internal
     Tensor<1,3>
     wedge_product (const Tensor<1,3> (&derivative)[2])
     {
-      Tensor<1,3> result;
-      cross_product (result, derivative[0], derivative[1]);
-
-      return result;
+      result cross_product (derivative[0], derivative[1]);
     }
 
 
index 12429b41547f65f1189dc3807a7eddd70baee95d..fcb2eb26af840ec308f9dc873931b8fe6bb78b11 100644 (file)
@@ -1043,10 +1043,11 @@ namespace internal
                                                       -1 : +1);
                   break;
                 case 2:
-                  cross_product (output_data.boundary_forms[i], data.aux[0][i]);
+                  output_data.boundary_forms[i] = cross_product(data.aux[0][i]);
                   break;
                 case 3:
-                  cross_product (output_data.boundary_forms[i], data.aux[0][i], data.aux[1][i]);
+                  output_data.boundary_forms[i] =
+                    cross_product(data.aux[0][i], data.aux[1][i]);
                   break;
                 default:
                   Assert(false, ExcNotImplemented());
@@ -1075,17 +1076,17 @@ namespace internal
 
                 if (dim==2)
                   {
-                    Tensor<1,spacedim> cell_normal;
                     const DerivativeForm<1,spacedim,dim> DX_t =
                       data.contravariant[point].transpose();
-                    cross_product(cell_normal,DX_t[0],DX_t[1]);
+
+                    Tensor<1, spacedim> cell_normal =
+                      cross_product(DX_t[0], DX_t[1]);
                     cell_normal /= cell_normal.norm();
 
                     // then compute the face normal from the face tangent
                     // and the cell normal:
-                    cross_product (output_data.boundary_forms[point],
-                                   data.aux[0][point], cell_normal);
-
+                    output_data.boundary_forms[point] =
+                      cross_product(data.aux[0][point], cell_normal);
                   }
 
               }
@@ -1291,9 +1292,11 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                                 ExcMessage("There is no cell normal in codim 2."));
 
                         if (dim==1)
-                          cross_product(output_data.normal_vectors[point], -DX_t[0]);
+                          output_data.normal_vectors[point] =
+                            cross_product(-DX_t[0]);
                         else //dim == 2
-                          cross_product(output_data.normal_vectors[point],DX_t[0],DX_t[1]);
+                          output_data.normal_vectors[point] =
+                            cross_product(DX_t[0], DX_t[1]);
 
                         output_data.normal_vectors[point] /= output_data.normal_vectors[point].norm();
 
index d37b6d71a53a6b60dbbc66d1e9e9feeafc82532a..155e6a23a938eae8d052cb833d2af9b8e22c83aa 100644 (file)
@@ -1501,10 +1501,11 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                         Assert( codim==1 , ExcMessage("There is no cell normal in codim 2."));
 
                         if (dim==1)
-                          cross_product(output_data.normal_vectors[point],
-                                        -DX_t[0]);
+                          output_data.normal_vectors[point] =
+                            cross_product(-DX_t[0]);
                         else //dim == 2
-                          cross_product(output_data.normal_vectors[point],DX_t[0],DX_t[1]);
+                          output_data.normal_vectors[point] =
+                            cross_product(DX_t[0], DX_t[1]);
 
                         output_data.normal_vectors[point] /= output_data.normal_vectors[point].norm();
 
@@ -1644,10 +1645,12 @@ namespace internal
                                                         -1 : +1);
                     break;
                   case 2:
-                    cross_product (output_data.boundary_forms[i], data.aux[0][i]);
+                    output_data.boundary_forms[i] =
+                      cross_product(data.aux[0][i]);
                     break;
                   case 3:
-                    cross_product (output_data.boundary_forms[i], data.aux[0][i], data.aux[1][i]);
+                    output_data.boundary_forms[i] =
+                      cross_product(data.aux[0][i], data.aux[1][i]);
                     break;
                   default:
                     Assert(false, ExcNotImplemented());
@@ -1675,16 +1678,17 @@ namespace internal
 
                   if (dim==2)
                     {
-                      Tensor<1,spacedim> cell_normal;
                       const DerivativeForm<1,spacedim,dim> DX_t =
                         data.contravariant[point].transpose();
-                      cross_product(cell_normal,DX_t[0],DX_t[1]);
+
+                      Tensor<1, spacedim> cell_normal =
+                        cross_product(DX_t[0], DX_t[1]);
                       cell_normal /= cell_normal.norm();
 
                       // then compute the face normal from the face tangent
                       // and the cell normal:
-                      cross_product (output_data.boundary_forms[point],
-                                     data.aux[0][point], cell_normal);
+                      output_data.boundary_forms[point] =
+                        cross_product(data.aux[0][point], cell_normal);
                     }
                 }
             }
index bb430edc849908488689fd7991f8377182c43f6c..361924f54821e22a814388e01f13bc37944575b1 100644 (file)
@@ -981,8 +981,7 @@ namespace
     const Tensor<1,3> v01 = accessor.vertex(1) - accessor.vertex(0);
     const Tensor<1,3> v02 = accessor.vertex(2) - accessor.vertex(0);
 
-    Tensor<1,3> normal;
-    cross_product(normal, v01, v02);
+    Tensor<1,3> normal = cross_product(v01, v02);
 
     const Tensor<1,3> v03 = accessor.vertex(3) - accessor.vertex(0);
 
@@ -1004,8 +1003,7 @@ namespace
     // the face is planar. then its area is 1/2 of the norm of the
     // cross product of the two diagonals
     const Tensor<1,3> v12 = accessor.vertex(2) - accessor.vertex(1);
-    Tensor<1,3> twice_area;
-    cross_product(twice_area, v03, v12);
+    Tensor<1,3> twice_area = cross_product(v03, v12);
     return 0.5 * twice_area.norm();
   }
 
index e68d058cd586370bbf70b6ccb37e2818ea977a3c..132c33381d0eff613bf96ac321de726589981d5f 100644 (file)
@@ -491,8 +491,7 @@ namespace internal
     Tensor<1,2>
     normalized_alternating_product (const Tensor<1,2> (&basis_vectors)[1])
     {
-      Tensor<1,2> tmp;
-      cross_product (tmp, basis_vectors[0]);
+      Tensor<1,2> tmp = cross_product (basis_vectors[0]);
       return tmp/tmp.norm();
     }
 
@@ -513,8 +512,7 @@ namespace internal
     Tensor<1,3>
     normalized_alternating_product (const Tensor<1,3> (&basis_vectors)[2])
     {
-      Tensor<1,3> tmp;
-      cross_product (tmp, basis_vectors[0], basis_vectors[1]);
+      Tensor<1,3> tmp = cross_product (basis_vectors[0], basis_vectors[1]);
       return tmp/tmp.norm();
     }
 
@@ -689,11 +687,8 @@ get_normals_at_vertices (const Triangulation<3>::face_iterator &face,
   { {1,2},{3,0},{0,3},{2,1}};
   for (unsigned int vertex=0; vertex<vertices_per_face; ++vertex)
     {
-      // first define the two tangent
-      // vectors at the vertex by
-      // using the two lines
-      // radiating away from this
-      // vertex
+      // first define the two tangent vectors at the vertex by using the
+      // two lines radiating away from this vertex
       const Tensor<1,3> tangents[2]
         = { face->vertex(neighboring_vertices[vertex][0])
             - face->vertex(vertex),
@@ -701,13 +696,9 @@ get_normals_at_vertices (const Triangulation<3>::face_iterator &face,
             - face->vertex(vertex)
           };
 
-      // then compute the normal by
-      // taking the cross
-      // product. since the normal is
-      // not required to be
-      // normalized, no problem here
-      cross_product (face_vertex_normals[vertex],
-                     tangents[0], tangents[1]);
+      // then compute the normal by taking the cross product. since the
+      // normal is not required to be normalized, no problem here
+      face_vertex_normals[vertex] = cross_product(tangents[0], tangents[1]);
     };
 }
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.