[7449] | 1 | """Class Parallel_shallow_water_domain - |
---|
[3185] | 2 | 2D triangular domains for finite-volume computations of |
---|
| 3 | the shallow water equation, with extra structures to allow |
---|
[7449] | 4 | communication between other Parallel_domains and itself |
---|
[3185] | 5 | |
---|
| 6 | This module contains a specialisation of class Domain |
---|
| 7 | from module shallow_water.py |
---|
| 8 | |
---|
| 9 | Ole Nielsen, Stephen Roberts, Duncan Gray, Christopher Zoppou |
---|
| 10 | Geoscience Australia, 2004-2005 |
---|
[3580] | 11 | |
---|
[3185] | 12 | """ |
---|
| 13 | |
---|
[7877] | 14 | from anuga import Domain |
---|
[3429] | 15 | |
---|
[8107] | 16 | from anuga_parallel.parallel_generic_communications import * |
---|
[8224] | 17 | from anuga.abstract_2d_finite_volumes.neighbour_mesh import Mesh |
---|
[7400] | 18 | |
---|
| 19 | import numpy as num |
---|
| 20 | |
---|
[3185] | 21 | |
---|
[8260] | 22 | #Import matplotlib |
---|
[3427] | 23 | |
---|
[8260] | 24 | |
---|
[8262] | 25 | |
---|
[7449] | 26 | class 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 | |
---|