Garth Wells avatar Garth Wells committed 48545ae

Minor fixes.

Comments (0)

Files changed (6)

dolfin/geometry/Point.h

   };
 
   /// Multiplication with scalar
-  inline Point operator*(double a, const Point& p) { return p*a; }
+  inline Point operator*(double a, const Point& p)
+  { return p*a; }
 
   inline std::ostream& operator<<(std::ostream& stream, const Point& point)
   { stream << point.str(false); return stream; }

dolfin/mesh/BoundaryComputation.cpp

                  "Unknown boundary type (%d)", type.c_str());
   }
 
-  const std::size_t this_process = MPI::process_number();
+  // Get my MPI process rank and number of MPI processes
+  const std::size_t my_rank = MPI::process_number();
+  const std::size_t num_processes = MPI::num_processes();
 
   // Open boundary mesh for editing
   const std::size_t D = mesh.topology().dim();
         {
           const std::size_t local_mesh_index = v->index();
 
-          if (boundary_vertices.find(local_mesh_index) == boundary_vertices.end())
+          if (boundary_vertices.find(local_mesh_index)
+              == boundary_vertices.end())
           {
             const std::size_t local_boundary_index = num_boundary_vertices;
             boundary_vertices[local_mesh_index] = local_boundary_index;
 
             // Determine "owner" of global_mesh_index
-            std::size_t owner = this_process;
+            std::size_t owner = my_rank;
 
             std::map<unsigned int, std::set<unsigned int> >::const_iterator
               other_processes_it = shared_vertices.find(local_mesh_index);
-
             if (other_processes_it != shared_vertices.end() && D > 1)
             {
-              const std::set<unsigned int> &other_processes = other_processes_it->second;
-              const std::size_t min_process = *std::min_element(other_processes.begin(), other_processes.end());
-              boundary.topology().shared_entities(0)[local_boundary_index] = other_processes;
+              const std::set<unsigned int>& other_processes
+                = other_processes_it->second;
+              const std::size_t min_process
+                = *std::min_element(other_processes.begin(),
+                                    other_processes.end());
+              boundary.topology().shared_entities(0)[local_boundary_index]
+                = other_processes;
 
               // FIXME: More sophisticated ownership determination
               if (min_process < owner)
                 owner = min_process;
             }
-            const std::size_t global_mesh_index = mesh.topology().global_indices(0)[local_mesh_index];
+            const std::size_t global_mesh_index
+              = mesh.topology().global_indices(0)[local_mesh_index];
             global_index_owner[global_mesh_index] = owner;
 
             // Update counts
-            if (owner == this_process)
+            if (owner == my_rank)
               num_owned_vertices++;
             num_boundary_vertices++;
           }
     }
   }
 
-  // Initiate boundary topology 
+  // Initiate boundary topology
   boundary.topology().init_global(0, MPI::sum(num_owned_vertices));
   boundary.topology().init_global(D-1, MPI::sum(num_boundary_cells));
 
   if (num_boundary_vertices > 0)
     vertex_map.init(boundary, 0, num_boundary_vertices);
   std::map<std::size_t, std::size_t>::const_iterator it;
-  for (it=boundary_vertices.begin(); it != boundary_vertices.end(); ++it)
+  for (it = boundary_vertices.begin(); it != boundary_vertices.end(); ++it)
     vertex_map[it->second] = it->first;
 
-  // Get vertex ownership distribution, and find index to start global numbering from
-  std::vector<std::size_t> ownership_distribution(MPI::num_processes());
+  // Get vertex ownership distribution, and find index to start global
+  // numbering from
+  std::vector<std::size_t> ownership_distribution(num_processes);
   MPI::all_gather(num_owned_vertices, ownership_distribution);
   std::size_t start_index = 0;
-  for (std::size_t j=0; j<MPI::process_number(); j++)
+  for (std::size_t j = 0; j < my_rank; j++)
     start_index += ownership_distribution[j];
 
-  // Set global indices of owned vertices, request global indices for vertices owned elsewhere
+  // Set global indices of owned vertices, request global indices for
+  // vertices owned elsewhere
   std::map<std::size_t, std::size_t> global_indices;
-  std::vector< std::vector<std::size_t> > request_global_indices(MPI::num_processes());
+  std::vector<std::vector<std::size_t> > request_global_indices(num_processes);
 
   std::size_t current_index = start_index;
-  for (std::size_t local_boundary_index=0; local_boundary_index<num_boundary_vertices; local_boundary_index++)
+  for (std::size_t local_boundary_index = 0;
+       local_boundary_index<num_boundary_vertices; local_boundary_index++)
   {
     const std::size_t local_mesh_index = vertex_map[local_boundary_index];
-    const std::size_t global_mesh_index = mesh.topology().global_indices(0)[local_mesh_index];
+    const std::size_t global_mesh_index
+      = mesh.topology().global_indices(0)[local_mesh_index];
 
     const std::size_t owner = global_index_owner[global_mesh_index];
-    if (owner != this_process)
+    if (owner != my_rank)
       request_global_indices[owner].push_back(global_mesh_index);
     else
       global_indices[global_mesh_index] = current_index++;
   }
 
   // Send and receive requests from other processes
-  std::vector< std::vector<std::size_t> > global_index_requests(MPI::num_processes());
+  std::vector<std::vector<std::size_t> > global_index_requests(num_processes);
   MPI::all_to_all(request_global_indices, global_index_requests);
 
   // Find response to requests of global indices
-  std::vector< std::vector<std::size_t> > respond_global_indices(MPI::num_processes());  
-  for (std::size_t i=0; i<MPI::num_processes(); i++)
+  std::vector<std::vector<std::size_t> > respond_global_indices(num_processes);
+  for (std::size_t i = 0; i < num_processes; i++)
   {
     const std::size_t N = global_index_requests[i].size();
     respond_global_indices[i].resize(N);
 
-    for (std::size_t j=0; j < N; j++)
-      respond_global_indices[i][j] = global_indices[global_index_requests[i][j]];  
+    for (std::size_t j = 0; j < N; j++)
+      respond_global_indices[i][j]
+        = global_indices[global_index_requests[i][j]];
   }
 
-  // Scatter responses back to requesting processes  
-  std::vector< std::vector<std::size_t> > global_index_responses(MPI::num_processes());
+  // Scatter responses back to requesting processes
+  std::vector<std::vector<std::size_t> > global_index_responses(num_processes);
   MPI::all_to_all(respond_global_indices, global_index_responses);
 
   // Update global_indices
-  for (std::size_t i=0; i<MPI::num_processes(); i++)
+  for (std::size_t i = 0; i < num_processes; i++)
   {
     const std::size_t N = global_index_responses[i].size();
     // Check that responses are the same size as the requests made
-    dolfin_assert(global_index_responses[i].size() == request_global_indices[i].size());
-    for (std::size_t j=0; j<N; j++)
+    dolfin_assert(global_index_responses[i].size()
+                  == request_global_indices[i].size());
+    for (std::size_t j = 0; j < N; j++)
     {
       const std::size_t global_mesh_index = request_global_indices[i][j];
       const std::size_t global_boundary_index = global_index_responses[i][j];
   }
 
   // Create vertices
-  for (std::size_t local_boundary_index=0; local_boundary_index < num_boundary_vertices; local_boundary_index++)
+  for (std::size_t local_boundary_index = 0;
+       local_boundary_index < num_boundary_vertices; local_boundary_index++)
   {
     const std::size_t local_mesh_index = vertex_map[local_boundary_index];
-    const std::size_t global_mesh_index = mesh.topology().global_indices(0)[local_mesh_index];
+    const std::size_t global_mesh_index
+      = mesh.topology().global_indices(0)[local_mesh_index];
     const std::size_t global_boundary_index = global_indices[global_mesh_index];
 
     Vertex v(mesh, local_mesh_index);
 
-    editor.add_vertex_global(local_boundary_index, global_boundary_index, v.point());
+    editor.add_vertex_global(local_boundary_index, global_boundary_index,
+                             v.point());
   }
 
   // Find global index to start cell numbering from for current process
-  std::vector<std::size_t> cell_distribution(MPI::num_processes());
+  std::vector<std::size_t> cell_distribution(num_processes);
   MPI::all_gather(num_boundary_cells, cell_distribution);
   std::size_t start_cell_index = 0;
-  for (std::size_t i=0; i<MPI::process_number(); i++)
+  for (std::size_t i = 0; i < my_rank; i++)
     start_cell_index += cell_distribution[i];
 
   // Create cells (facets) and map between boundary mesh cells and facets parent
   MeshFunction<std::size_t>& cell_map = boundary.entity_map(D - 1);
   if (num_boundary_cells > 0)
     cell_map.init(boundary, D - 1, num_boundary_cells);
-  std::vector<std::size_t> cell(boundary.type().num_vertices(boundary.topology().dim()));
+  std::vector<std::size_t>
+    cell(boundary.type().num_vertices(boundary.topology().dim()));
   std::size_t current_cell = 0;
   for (FacetIterator f(mesh); !f.end(); ++f)
   {
       for (std::size_t i = 0; i < cell.size(); i++)
         cell[i] = boundary_vertices[vertices[i]];
 
-      // Reorder vertices so facet is right-oriented w.r.t. facet normal
+      // Reorder vertices so facet is right-oriented w.r.t. facet
+      // normal
       reorder(cell, *f);
 
       // Create mapping from boundary cell to mesh facet if requested
   // ordering from destroying the orientation of facets accomplished
   // by calling reorder() below.
   editor.close(false);
-
 }
 //-----------------------------------------------------------------------------
 void BoundaryComputation::reorder(std::vector<std::size_t>& vertices,

dolfin/mesh/BoundaryMesh.cpp

 using namespace dolfin;
 
 //-----------------------------------------------------------------------------
-BoundaryMesh::BoundaryMesh(const Mesh& mesh, std::string type, bool order) : Mesh()
+BoundaryMesh::BoundaryMesh(const Mesh& mesh, std::string type, bool order)
+  : Mesh()
 {
   // Create boundary mesh
   BoundaryComputation::compute_boundary(mesh, type, *this);

dolfin/mesh/DynamicMeshEditor.h

     ~DynamicMeshEditor();
 
     /// Open mesh of given cell type, topological and geometrical dimension
-    void open(Mesh& mesh, CellType::Type type, std::size_t tdim, std::size_t gdim);
+    void open(Mesh& mesh, CellType::Type type, std::size_t tdim,
+              std::size_t gdim);
 
     /// Open mesh of given cell type, topological and geometrical dimension
     void open(Mesh& mesh, std::string type, std::size_t tdim, std::size_t gdim);
     void add_cell(std::size_t c, std::size_t v0, std::size_t v1);
 
     /// Add cell (triangle) with given vertices
-    void add_cell(std::size_t c, std::size_t v0, std::size_t v1, std::size_t v2);
+    void add_cell(std::size_t c,  std::size_t v0, std::size_t v1,
+                  std::size_t v2);
 
     /// Add cell (tetrahedron) with given vertices
-    void add_cell(std::size_t c, std::size_t v0, std::size_t v1, std::size_t v2, std::size_t v3);
+    void add_cell(std::size_t c, std::size_t v0, std::size_t v1,
+                  std::size_t v2, std::size_t v3);
 
     /// Close mesh, finish editing, and order entities locally
     void close(bool order=false);

dolfin/mesh/MeshTransformation.cpp

                                 const Point& p)
 {
   // Compute angle (radians)
-  const double theta = angle / 180.0 * DOLFIN_PI;
+  const double theta = angle/180.0*DOLFIN_PI;
 
   // Get coordinates of point
   const double* c = p.coordinates();
 
     // Set up rotation matrix
     const double S00 = cos(theta); const double S01 = -sin(theta);
-    const double S10 = sin(theta); const double S11 = cos(theta);
+    const double S10 = sin(theta); const double S11 =  cos(theta);
 
     // Rotate all points
     MeshGeometry& geometry = mesh.geometry();
   {
     // Set up 2D rotation matrix
     const double S00 = cos(theta); const double S01 = -sin(theta);
-    const double S10 = sin(theta); const double S11 = cos(theta);
+    const double S10 = sin(theta); const double S11 =  cos(theta);
 
     // Initialize 3D rotation matrix to identity matrix
     double R00 = 1.0; double R01 = 0.0; double R02 = 0.0;

dolfin/mesh/MeshTransformation.h

 {
 
    class Mesh;
+   class Point;
 
   /// This class implements various transformations of the coordinates
   /// of a mesh.
     ///         The coordinate axis around which to rotate the mesh.
     ///     point (_Point_)
     ///         The point around which to rotate the mesh.
-    static void rotate(Mesh& mesh, double angle, std::size_t axis, const Point& p);
+    static void rotate(Mesh& mesh, double angle, std::size_t axis,
+                       const Point& p);
 
   };
 
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for ProjectModifiedEvent.java.
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.