Ignore:
Timestamp:
Jun 18, 2010, 5:49:57 PM (12 years ago)
Author:
steve
Message:

Continuing to numpy the for loops

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/development/2010-projects/anuga_1d/base/generic_domain.py

    r7842 r7860  
    4040        self.set_timestepping_method(timestepping_method)
    4141
    42         self.wet_nodes = numpy.zeros((N,2), numpy.int) # should this be here
     42
     43        # work array for special vertices (like wet or dry cells
     44        self.special_vertices = numpy.zeros((N,2), numpy.int) # should this be here
    4345
    4446        #Allocate space for neighbour and boundary structures
     
    5860
    5961        self.normals = numpy.zeros((N, 2), numpy.float)
    60        
    61         for i in range(N):
    62             xl = self.coordinates[i]
    63             xr = self.coordinates[i+1]
    64             self.vertices[i,0] = xl
    65             self.vertices[i,1] = xr
    66 
    67             centroid = (xl+xr)/2.0
    68             self.centroids[i] = centroid
    69 
    70             msg = 'Coordinates should be ordered, smallest to largest'
    71             assert xr>xl, msg
    72            
    73             #The normal vectors
    74             # - point outward from each edge
    75             # - are orthogonal to the edge
    76             # - have unit length
    77             # - Are enumerated by left vertex then right vertex normals
    78 
    79             nl = -1.0
    80             nr =  1.0
    81             self.normals[i,:] = [nl, nr]
    82 
    83             self.areas[i] = (xr-xl)
    84            
    85 # #         print 'N', N
    86 # #         print 'Centroid', self.centroids
    87 # #         print 'Areas', self.areas
    88 # #         print 'Vertex_Coordinates', self.vertices
    89            
    90             #Initialise Neighbours (-1 means that it is a boundary neighbour)
    91             self.neighbours[i, :] = [-1, -1]
    92             #Initialise edge ids of neighbours
    93             #Initialise vertex ids of neighbours
    94             #In case of boundaries this slot is not used
    95             #self.neighbour_edges[i, :] = [-1, -1]
    96             self.neighbour_vertices[i, :] = [-1, -1]
     62
     63
     64        xl = self.coordinates[0:-1]
     65        xr = self.coordinates[1: ]
     66
     67        self.vertices[:,0] = xl
     68        self.vertices[:,1] = xr
     69
     70        self.centroids[:] = (xl+xr)/2.0
     71
     72        msg = 'Coordinates should be ordered, smallest to largest'
     73        assert (xr>xl).all(), msg
     74
     75        #The normal vectors
     76        # - point outward from each edge
     77        # - are orthogonal to the edge
     78        # - have unit length
     79        # - Are enumerated by left vertex then right vertex normals
     80
     81
     82        self.normals[:,0] = -1.0
     83        self.normals[:,1] =  1.0
     84
     85        self.areas[:] = (xr-xl)
     86
     87        #Initialise Neighbours (-1 means that it is a boundary neighbour)
     88
     89        self.neighbours[:, :] = -1
     90
     91        #Initialise edge ids of neighbours
     92        self.neighbour_vertices[:, :] = -1
     93
     94
     95#        for i in range(N):
     96#            xl = self.coordinates[i]
     97#            xr = self.coordinates[i+1]
     98#            self.vertices[i,0] = xl
     99#            self.vertices[i,1] = xr
     100#
     101#            centroid = (xl+xr)/2.0
     102#            self.centroids[i] = centroid
     103#
     104#            msg = 'Coordinates should be ordered, smallest to largest'
     105#            assert xr>xl, msg
     106#
     107#            #The normal vectors
     108#            # - point outward from each edge
     109#            # - are orthogonal to the edge
     110#            # - have unit length
     111#            # - Are enumerated by left vertex then right vertex normals
     112#
     113#            nl = -1.0
     114#            nr =  1.0
     115#            self.normals[i,:] = [nl, nr]
     116#
     117#            self.areas[i] = (xr-xl)
     118#
     119## #         print 'N', N
     120## #         print 'Centroid', self.centroids
     121## #         print 'Areas', self.areas
     122## #         print 'Vertex_Coordinates', self.vertices
     123#
     124#            #Initialise Neighbours (-1 means that it is a boundary neighbour)
     125#            self.neighbours[i, :] = [-1, -1]
     126#            #Initialise edge ids of neighbours
     127#            #Initialise vertex ids of neighbours
     128#            #In case of boundaries this slot is not used
     129#            #self.neighbour_edges[i, :] = [-1, -1]
     130#            self.neighbour_vertices[i, :] = [-1, -1]
    97131
    98132        self.build_vertexlist()
     
    630664            return
    631665
     666        if timestepping_method in [1, 2, 3]:
     667            self.timestepping_method = ['euler', 'rk2', 'rk3'][timestepping_method-1]
     668            return
     669
    632670        msg = '%s is an incorrect timestepping type'% timestepping_method
    633671        raise Exception, msg
     
    814852        return self.filename
    815853
     854
     855    def get_spatial_order(self):
     856        return self.order
     857
     858    def set_spatial_order(self, order):
     859
     860        possible_orders = [1,2]
     861
     862        if order in possible_orders:
     863            self.order = order
     864            return
     865
     866        msg = '%s is an incorrect spatial order.\n'% limiter
     867        msg += 'Possible orders are: '+ ", ".join(["%s" % el for el in possible_orders])
     868        raise Exception, msg
     869   
     870
    816871    def get_limiter(self):
    817872        return self.limiter
Note: See TracChangeset for help on using the changeset viewer.