Changeset 7860 for anuga_work/development/2010-projects/anuga_1d/base/generic_domain.py

Timestamp:

Jun 18, 2010, 5:49:57 PM (13 years ago)

Author:

steve

Message:

Continuing to numpy the for loops

File:

• anuga_work/development/2010-projects/anuga_1d/base/generic_domain.py (modified) (4 diffs)

anuga_work/development/2010-projects/anuga_1d/base/generic_domain.py

@@ -40 +40 @@
         self.set_timestepping_method(timestepping_method)
 
-        self.wet_nodes = numpy.zeros((N,2), numpy.int) # should this be here
+
+        # work array for special vertices (like wet or dry cells
+        self.special_vertices = numpy.zeros((N,2), numpy.int) # should this be here
 
         #Allocate space for neighbour and boundary structures
@@ -58 +60 @@
 
         self.normals = numpy.zeros((N, 2), numpy.float)
-        
-        for i in range(N):
-            xl = self.coordinates[i]
-            xr = self.coordinates[i+1]
-            self.vertices[i,0] = xl
-            self.vertices[i,1] = xr
-
-            centroid = (xl+xr)/2.0
-            self.centroids[i] = centroid
-
-            msg = 'Coordinates should be ordered, smallest to largest'
-            assert xr>xl, msg
-            
-            #The normal vectors
-            # - point outward from each edge
-            # - are orthogonal to the edge
-            # - have unit length
-            # - Are enumerated by left vertex then right vertex normals
-
-            nl = -1.0
-            nr =  1.0
-            self.normals[i,:] = [nl, nr]
-
-            self.areas[i] = (xr-xl)
-            
-# #         print 'N', N
-# #         print 'Centroid', self.centroids
-# #         print 'Areas', self.areas
-# #         print 'Vertex_Coordinates', self.vertices
-            
-            #Initialise Neighbours (-1 means that it is a boundary neighbour)
-            self.neighbours[i, :] = [-1, -1]
-            #Initialise edge ids of neighbours
-            #Initialise vertex ids of neighbours
-            #In case of boundaries this slot is not used
-            #self.neighbour_edges[i, :] = [-1, -1]
-            self.neighbour_vertices[i, :] = [-1, -1]
+
+
+        xl = self.coordinates[0:-1]
+        xr = self.coordinates[1: ]
+
+        self.vertices[:,0] = xl
+        self.vertices[:,1] = xr
+
+        self.centroids[:] = (xl+xr)/2.0
+
+        msg = 'Coordinates should be ordered, smallest to largest'
+        assert (xr>xl).all(), msg
+
+        #The normal vectors
+        # - point outward from each edge
+        # - are orthogonal to the edge
+        # - have unit length
+        # - Are enumerated by left vertex then right vertex normals
+
+
+        self.normals[:,0] = -1.0
+        self.normals[:,1] =  1.0
+
+        self.areas[:] = (xr-xl)
+
+        #Initialise Neighbours (-1 means that it is a boundary neighbour)
+
+        self.neighbours[:, :] = -1
+
+        #Initialise edge ids of neighbours
+        self.neighbour_vertices[:, :] = -1
+
+
+#        for i in range(N):
+#            xl = self.coordinates[i]
+#            xr = self.coordinates[i+1]
+#            self.vertices[i,0] = xl
+#            self.vertices[i,1] = xr
+#
+#            centroid = (xl+xr)/2.0
+#            self.centroids[i] = centroid
+#
+#            msg = 'Coordinates should be ordered, smallest to largest'
+#            assert xr>xl, msg
+#
+#            #The normal vectors
+#            # - point outward from each edge
+#            # - are orthogonal to the edge
+#            # - have unit length
+#            # - Are enumerated by left vertex then right vertex normals
+#
+#            nl = -1.0
+#            nr =  1.0
+#            self.normals[i,:] = [nl, nr]
+#
+#            self.areas[i] = (xr-xl)
+#
+## #         print 'N', N
+## #         print 'Centroid', self.centroids
+## #         print 'Areas', self.areas
+## #         print 'Vertex_Coordinates', self.vertices
+#
+#            #Initialise Neighbours (-1 means that it is a boundary neighbour)
+#            self.neighbours[i, :] = [-1, -1]
+#            #Initialise edge ids of neighbours
+#            #Initialise vertex ids of neighbours
+#            #In case of boundaries this slot is not used
+#            #self.neighbour_edges[i, :] = [-1, -1]
+#            self.neighbour_vertices[i, :] = [-1, -1]
 
         self.build_vertexlist()
@@ -630 +664 @@
             return
 
+        if timestepping_method in [1, 2, 3]:
+            self.timestepping_method = ['euler', 'rk2', 'rk3'][timestepping_method-1]
+            return
+
         msg = '%s is an incorrect timestepping type'% timestepping_method
         raise Exception, msg
@@ -814 +852 @@
         return self.filename
 
+
+    def get_spatial_order(self):
+        return self.order
+
+    def set_spatial_order(self, order):
+
+        possible_orders = [1,2]
+
+        if order in possible_orders:
+            self.order = order
+            return
+
+        msg = '%s is an incorrect spatial order.\n'% limiter
+        msg += 'Possible orders are: '+ ", ".join(["%s" % el for el in possible_orders])
+        raise Exception, msg
+    
+
     def get_limiter(self):
         return self.limiter