Changeset 8260

Timestamp:

Nov 30, 2011, 7:57:30 PM (13 years ago)

Author:

janes

Message:

Adding anuga.error_api.py

Location:

trunk/anuga_core/source

Files:

• anuga/error_api.py (added)
• anuga_parallel/distribute_mesh.py (modified) (5 diffs)
• anuga_parallel/parallel_api.py (modified) (6 diffs)
• anuga_parallel/parallel_shallow_water.py (modified) (4 diffs)
• anuga_parallel/test_parallel_frac_op.py (modified) (5 diffs)

trunk/anuga_core/source/anuga_parallel/distribute_mesh.py

@@ -101 +101 @@
 
 def pmesh_divide_metis(domain, n_procs):
+    # Wrapper for old pmesh_divide_metis which does not return tri_index or r_tri_index
+    nodes, ttriangles, boundary, triangles_per_proc, quantities, tri_index, r_tri_index = pmesh_divide_metis_helper(domain, n_procs)
+    return nodes, ttriangles, boundary, triangles_per_proc, quantities
+
+def pmesh_divide_metis_with_map(domain, n_procs):
+    return pmesh_divide_metis_helper(domain, n_procs)
+
+def pmesh_divide_metis_helper(domain, n_procs):
     
     # Initialise the lists
@@ -110 +118 @@
     
     tri_list = []
+
+    # Serial to Parallel and Parallel to Serial Triangle index maps
+    tri_index = {}
+    r_tri_index = {} # reverse tri index, parallel to serial triangle index mapping
     
     # List indexed by processor of cumulative total of triangles allocated.
@@ -148 +160 @@
         # (local to the processor)
         
-        tri_index = {}
         triangles = []        
         for i in range(n_tri):
@@ -154 +165 @@
             tri_list[epart[i]].append(domain.triangles[i])
             tri_index[i] = ([epart[i], len(tri_list[epart[i]]) - 1])
+            r_tri_index[epart[i], len(tri_list[epart[i]]) - 1] = i
         
         # Order the triangle list so that all of the triangles belonging
@@ -202 +214 @@
     for i in range(len(triangles)):
         ttriangles[i] = triangles[i]
-    
-    return nodes, ttriangles, boundary, triangles_per_proc, quantities
+
+    #return nodes, ttriangles, boundary, triangles_per_proc, quantities
+    
+    return nodes, ttriangles, boundary, triangles_per_proc, quantities, tri_index, r_tri_index
 
 

trunk/anuga_core/source/anuga_parallel/parallel_api.py

@@ -20 +20 @@
     # Mesh partitioning using Metis
     from anuga_parallel.distribute_mesh import build_submesh
-    from anuga_parallel.distribute_mesh import pmesh_divide_metis
+    from anuga_parallel.distribute_mesh import pmesh_divide_metis_with_map
 
     from anuga_parallel.parallel_shallow_water import Parallel_domain
@@ -97 +97 @@
         points, vertices, boundary, quantities,\
                 ghost_recv_dict, full_send_dict,\
-                number_of_full_nodes, number_of_full_triangles =\
+                number_of_full_nodes, number_of_full_triangles,\
+                s2p_map, p2s_map =\
                 distribute_mesh(domain, verbose=verbose)
 
+        # Send serial to parallel (s2p) and parallel to serial (p2s) triangle mapping to proc 1 .. numprocs
+        for p in range(1, numprocs):
+            send(s2p_map, p)
+            send(p2s_map, p)
 
         if verbose: print 'Communication done'
@@ -112 +117 @@
                 rec_submesh(0, verbose)
 
+        # Recieve serial to parallel (s2p) and parallel to serial (p2s) triangle mapping
+        s2p_map = receive(0)
+        p2s_map = receive(0)
+
 
     #------------------------------------------------------------------------
@@ -125 +134 @@
                              number_of_full_nodes=number_of_full_nodes,
                              number_of_full_triangles=number_of_full_triangles,
-                             geo_reference=georef) ## jj added this
+                             geo_reference=georef,
+                             tri_map = s2p_map,
+                             inv_tri_map = p2s_map) ## jj added this
 
     #------------------------------------------------------------------------
@@ -165 +176 @@
     # Subdivide the mesh
     if verbose: print 'Subdivide mesh'
-    nodes, triangles, boundary, triangles_per_proc, quantities = \
-           pmesh_divide_metis(domain, numprocs)
-
-
+    nodes, triangles, boundary, triangles_per_proc, quantities, s2p_map, p2s_map = \
+           pmesh_divide_metis_with_map(domain, numprocs)
+
+    #PETE: s2p_map (maps serial domain triangles to parallel domain triangles)
+    #p2_map (maps parallel domain triangles to domain triangles)
 
 
@@ -227 +239 @@
     return points, vertices, boundary, quantities,\
            ghost_recv_dict, full_send_dict,\
-           number_of_full_nodes, number_of_full_triangles
-    
-
-
-
+           number_of_full_nodes, number_of_full_triangles, s2p_map, p2s_map
+    
+
+
+

trunk/anuga_core/source/anuga_parallel/parallel_shallow_water.py

@@ -20 +20 @@
 
 
+#Import matplotlib
+
+import matplotlib
+matplotlib.use('Agg')
+import matplotlib.pyplot as plt
+import matplotlib.tri as tri
 
 class Parallel_domain(Domain):
@@ -29 +35 @@
                  number_of_full_nodes=None,
                  number_of_full_triangles=None,
-                 geo_reference=None): #jj added this
+                 geo_reference=None,
+                 tri_map=None,
+                 inv_tri_map=None): #jj added this
 
         Domain.__init__(self,
@@ -57 +65 @@
 
         setup_buffers(self)
+
+        self.tri_map = tri_map
+        self.inv_tri_map = inv_tri_map
 
 
@@ -122 +133 @@
         return N[:self.number_of_full_nodes_tmp,:]
 
+    def get_tri_map(self):
+        return self.tri_map
+
+    def get_inv_tri_map(self):
+        return self.inv_tri_map
+
+    '''
+    Outputs domain triangulation, full triangles are shown in green while ghost triangles are shown in blue.
+    The default filename is "domain.png"
+    '''
+    def dump_triangulation(self, filename="domain.png"):
+        # Get vertex coordinates, partition full and ghost triangles based on self.tri_full_flag
+        vertices = self.get_vertex_coordinates()
+        full_mask = num.repeat(self.tri_full_flag == 1, 3)
+        ghost_mask = num.repeat(self.tri_full_flag == 0, 3)
+        
+        myid = pypar.rank()
+        numprocs = pypar.size()
+
+        if myid == 0:
+            fx = {}
+            fy = {}
+            gx = {}
+            gy = {}
+
+            # Proc 0 gathers full and ghost nodes from self and other processors
+            fx[0] = vertices[full_mask,0]
+            fy[0] = vertices[full_mask,1]
+            gx[0] = vertices[ghost_mask,0]
+            gy[0] = vertices[ghost_mask,1]
+            
+            for i in range(1,numprocs):
+                fx[i] = pypar.receive(i)
+                fy[i] = pypar.receive(i)
+                gx[i] = pypar.receive(i)
+                gy[i] = pypar.receive(i)
+
+            # Plot full triangles
+            for i in range(0, numprocs):
+                n = int(len(fx[i])/3)
+                            
+                triang = num.array(range(0,3*n))
+                triang.shape = (n, 3)
+                plt.triplot(fx[i], fy[i], triang, 'g-')
+
+            # Plot ghost triangles
+            for i in range(0, numprocs):
+                n = int(len(gx[i])/3)
+                            
+                triang = num.array(range(0,3*n))
+                triang.shape = (n, 3)
+                plt.triplot(gx[i], gy[i], triang, 'b--')
+
+            # Save triangulation to location pointed by filename
+            plt.savefig(filename)
+
+        else:
+            # Proc 1..numprocs send full and ghost triangles to Proc 0
+            pypar.send(vertices[full_mask,0], 0)
+            pypar.send(vertices[full_mask,1], 0)
+            pypar.send(vertices[ghost_mask,0], 0)
+            pypar.send(vertices[ghost_mask,1], 0)
+

trunk/anuga_core/source/anuga_parallel/test_parallel_frac_op.py

@@ -105 +105 @@
 ##-----------------------------------------------------------------------
 
-    if parallel: domain = distribute(domain)
+    if parallel:
+        domain = distribute(domain)
+        domain.dump_triangulation("frac_op_domain.png")
     
 
@@ -147 +149 @@
 
     ################ Define Fractional Operators ##########################
+
+    inlet0 = None
+    inlet1 = None
+    boyd_box0 = None
     
     inlet0 = Inlet_operator(domain, line0, Q0, debug = False)
@@ -152 +158 @@
     
     # Enquiry point [ 19.    2.5] is contained in two domains in 4 proc case
+    
     boyd_box0 = Boyd_box_operator(domain,
                                   end_points=[[9.0, 2.5],[19.0, 2.5]],
@@ -161 +168 @@
                                   manning=0.013,
                                   verbose=False, debug = False)
-
+        
     if inlet0 is not None: inlet0.print_statistics()
     if inlet1 is not None: inlet1.print_statistics()
@@ -202 +209 @@
     ##-----------------------------------------------------------------------
 
-    for t in domain.evolve(yieldstep = 1.0, finaltime = 38):
+    for t in domain.evolve(yieldstep = 0.1, finaltime = 38):
         domain.write_time()