source: trunk/anuga_core/source/anuga_parallel/parallel_shallow_water.py @ 8281

Last change on this file since 8281 was 8281, checked in by steve, 12 years ago

added in some parallel data to allow sww_merge to create a merged sww file

File size: 7.5 KB
RevLine 
[7449]1"""Class Parallel_shallow_water_domain -
[3185]22D triangular domains for finite-volume computations of
3the shallow water equation, with extra structures to allow
[7449]4communication between other Parallel_domains and itself
[3185]5
6This module contains a specialisation of class Domain
7from module shallow_water.py
8
9Ole Nielsen, Stephen Roberts, Duncan Gray, Christopher Zoppou
10Geoscience Australia, 2004-2005
[3580]11
[3185]12"""
13
[7877]14from anuga import Domain
[3429]15
[8107]16from anuga_parallel.parallel_generic_communications import *
[8224]17from anuga.abstract_2d_finite_volumes.neighbour_mesh import Mesh
[7400]18
19import numpy as num
20
[3185]21
[8260]22#Import matplotlib
[3427]23
[8260]24
[8262]25
[7449]26class Parallel_domain(Domain):
[3185]27
[3926]28    def __init__(self, coordinates, vertices,
29                 boundary=None,
30                 full_send_dict=None,
31                 ghost_recv_dict=None,
[8209]32                 number_of_full_nodes=None,
33                 number_of_full_triangles=None,
[8260]34                 geo_reference=None,
[8281]35                 number_of_global_triangles=None, ## SR added this
36                 number_of_global_nodes= None, ## SR added this
37                 s2p_map=None,
38                 p2s_map=None, #jj added this
39                 tri_l2g = None, ## SR added this
40                 node_l2g = None): ## SR added this
[3185]41
42        Domain.__init__(self,
43                        coordinates,
44                        vertices,
45                        boundary,
46                        full_send_dict=full_send_dict,
47                        ghost_recv_dict=ghost_recv_dict,
48                        processor=pypar.rank(),
[3926]49                        numproc=pypar.size(),
[8209]50                        number_of_full_nodes=number_of_full_nodes,
51                        number_of_full_triangles=number_of_full_triangles,
[8114]52                        geo_reference=geo_reference) #jj added this
[8224]53       
[8281]54
55
56        self.parallel = True
57
[8224]58        # PETE: Find the number of full nodes and full triangles, this is a temporary fix
59        # until the bug with get_number_of_full_[nodes|triangles]() is fixed.
[3185]60
[8224]61        if number_of_full_nodes is not None:
62            self.number_of_full_nodes_tmp = number_of_full_nodes
63        else:
64            self.number_of_full_nodes_tmp = get_number_of_nodes()
65
66        if number_of_full_triangles is not None:
67            self.number_of_full_triangles_tmp = number_of_full_triangles
68        else:
69            self.number_of_full_triangles_tmp = get_number_of_triangles()
70
[8107]71        setup_buffers(self)
[3185]72
[8272]73        self.global_name = 'domain'
[3185]74
[8281]75        self.number_of_global_triangles=number_of_global_triangles
76        self.number_of_global_nodes = number_of_global_nodes
[8260]77
[8281]78        self.s2p_map = s2p_map
79        self.p2s_map = p2s_map
80        self.tri_l2g = tri_l2g
81        self.node_l2g = node_l2g
82
83
[3185]84    def set_name(self, name):
85        """Assign name based on processor number
86        """
87
[3893]88        if name.endswith('.sww'):
89            name = name[:-4]
90
[8272]91        self.global_name = name
92
[3185]93        # Call parents method with processor number attached.
[8281]94        Domain.set_name(self, name + '_P%d_%d' %(self.numproc, self.processor))
[3185]95
96
[8272]97    def get_global_name(self):
98
99        return self.global_name
100
101
[3185]102    def update_timestep(self, yieldstep, finaltime):
103        """Calculate local timestep
104        """
105
[8107]106        communicate_flux_timestep(self, yieldstep, finaltime)
[3185]107
108        Domain.update_timestep(self, yieldstep, finaltime)
109
[3431]110
[8115]111
[3185]112    def update_ghosts(self):
113
114        # We must send the information from the full cells and
115        # receive the information for the ghost cells
116
117
[8107]118        communicate_ghosts(self)
[8224]119
120    def apply_fractional_steps(self):
121
122        for operator in self.fractional_step_operators:
123            operator()
124
125        # PETE: Make sure that there are no deadlocks here
126
127        self.update_ghosts()
128
[8272]129
130
131    def sww_merge(self, verbose=False):
132
133        if self.processor == 0 and self.numproc > 1:
134            import anuga.utilities.sww_merge as merge
135           
136            merge.sww_merge(self.get_global_name(),self.numproc,verbose)
137
138
[8224]139# =======================================================================
140# PETE: NEW METHODS FOR FOR PARALLEL STRUCTURES. Note that we assume the
141# first "number_of_full_[nodes|triangles]" are full [nodes|triangles]
142# For full triangles it is possible to enquire self.tri_full_flag == True
143# =======================================================================
144
145    def get_number_of_full_triangles(self, *args, **kwargs):
146        return self.number_of_full_triangles_tmp
147
148    def get_full_centroid_coordinates(self, *args, **kwargs):
149        C = self.mesh.get_centroid_coordinates(*args, **kwargs)
150        return C[:self.number_of_full_triangles_tmp, :]
151
152    def get_full_vertex_coordinates(self, *args, **kwargs):
153        V = self.mesh.get_vertex_coordinates(*args, **kwargs)
154        return V[:3*self.number_of_full_triangles_tmp,:]
155
156    def get_full_triangles(self, *args, **kwargs):
157        T = self.mesh.get_triangles(*args, **kwargs)
158        return T[:self.number_of_full_triangles_tmp,:]
159
160    def get_full_nodes(self, *args, **kwargs):
161        N = self.mesh.get_nodes(*args, **kwargs)
162        return N[:self.number_of_full_nodes_tmp,:]
163
[8260]164    def get_tri_map(self):
165        return self.tri_map
166
167    def get_inv_tri_map(self):
168        return self.inv_tri_map
169
170    '''
171    Outputs domain triangulation, full triangles are shown in green while ghost triangles are shown in blue.
172    The default filename is "domain.png"
173    '''
174    def dump_triangulation(self, filename="domain.png"):
175        # Get vertex coordinates, partition full and ghost triangles based on self.tri_full_flag
[8262]176
177        try:
178            import matplotlib
179            matplotlib.use('Agg')
180            import matplotlib.pyplot as plt
181            import matplotlib.tri as tri
182        except:
183            print "Couldn't import module from matplotlib, probably you need to update matplotlib"
184            raise
185
[8260]186        vertices = self.get_vertex_coordinates()
187        full_mask = num.repeat(self.tri_full_flag == 1, 3)
188        ghost_mask = num.repeat(self.tri_full_flag == 0, 3)
189       
190        myid = pypar.rank()
191        numprocs = pypar.size()
192
193        if myid == 0:
194            fx = {}
195            fy = {}
196            gx = {}
197            gy = {}
198
199            # Proc 0 gathers full and ghost nodes from self and other processors
200            fx[0] = vertices[full_mask,0]
201            fy[0] = vertices[full_mask,1]
202            gx[0] = vertices[ghost_mask,0]
203            gy[0] = vertices[ghost_mask,1]
204           
205            for i in range(1,numprocs):
206                fx[i] = pypar.receive(i)
207                fy[i] = pypar.receive(i)
208                gx[i] = pypar.receive(i)
209                gy[i] = pypar.receive(i)
210
211            # Plot full triangles
212            for i in range(0, numprocs):
213                n = int(len(fx[i])/3)
214                           
215                triang = num.array(range(0,3*n))
216                triang.shape = (n, 3)
217                plt.triplot(fx[i], fy[i], triang, 'g-')
218
219            # Plot ghost triangles
220            for i in range(0, numprocs):
221                n = int(len(gx[i])/3)
222                           
223                triang = num.array(range(0,3*n))
224                triang.shape = (n, 3)
225                plt.triplot(gx[i], gy[i], triang, 'b--')
226
227            # Save triangulation to location pointed by filename
228            plt.savefig(filename)
229
230        else:
231            # Proc 1..numprocs send full and ghost triangles to Proc 0
232            pypar.send(vertices[full_mask,0], 0)
233            pypar.send(vertices[full_mask,1], 0)
234            pypar.send(vertices[ghost_mask,0], 0)
235            pypar.send(vertices[ghost_mask,1], 0)
[8272]236           
[8260]237
Note: See TracBrowser for help on using the repository browser.