Changeset 8281


Ignore:
Timestamp:
Dec 12, 2011, 5:11:47 PM (8 years ago)
Author:
steve
Message:

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

Location:
trunk/anuga_core/source
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/generic_domain.py

    r8235 r8281  
    165165        self.fractional_step_operators = []
    166166
     167
     168
     169        # by default domain is not parallel
     170        self.parallel = False
     171
     172        self.number_of_global_triangles = self.number_of_triangles
     173        self.number_of_global_nodes = self.number_of_nodes
     174
    167175        # Setup the ghost cell communication
    168176        if full_send_dict is None:
  • trunk/anuga_core/source/anuga/file/sww.py

    r8278 r8281  
    183183                                        points_georeference=\
    184184                                            domain.geo_reference)
     185
     186
     187        if domain.parallel:
     188            self.writer.store_parallel_data(fid,
     189                                        domain.number_of_global_triangles,
     190                                        domain.number_of_global_nodes,
     191                                        domain.tri_full_flag,
     192                                        domain.tri_l2g,
     193                                        domain.node_l2g)
    185194
    186195
     
    676685        outfile.variables['y'][:] = y #- geo_ref.get_yllcorner()
    677686        outfile.variables['volumes'][:] = volumes.astype(num.int32) #On Opteron 64
     687
     688
     689    def store_parallel_data(self,
     690                            outfile,
     691                            number_of_global_triangles,
     692                            number_of_global_nodes,
     693                            tri_full_flag = None,
     694                            tri_l2g = None,
     695                            node_l2g = None,
     696                            sww_precision=netcdf_float32,
     697                            verbose=False):
     698
     699
     700         # dimension definitions
     701        #outfile.createDimension('number_of_volumes', number_of_volumes)
     702        #outfile.createDimension('number_of_vertices', 3)
     703        #outfile.createDimension('numbers_in_range', 2)
     704
     705        print 'store parallel data'
     706        outfile.number_of_global_triangles = number_of_global_triangles
     707        outfile.number_of_global_nodes = number_of_global_nodes
     708
     709        # variable definitions
     710        outfile.createVariable('tri_l2g',  netcdf_int, ('number_of_volumes',))
     711        outfile.createVariable('node_l2g', netcdf_int, ('number_of_points',))
     712        outfile.createVariable('tri_full_flag', netcdf_int, ('number_of_volumes',))
     713
     714        print tri_l2g.shape
     715        print tri_l2g
     716        print outfile.variables['tri_l2g'].shape
     717
     718        outfile.variables['tri_l2g'][:] = tri_l2g.astype(num.int32)
     719
     720        print node_l2g.shape
     721        print node_l2g
     722        print outfile.variables['node_l2g'].shape
     723
     724        outfile.variables['node_l2g'][:] = node_l2g.astype(num.int32)
     725
     726        print tri_full_flag.shape
     727        print tri_full_flag
     728        print outfile.variables['tri_full_flag'].shape
     729
     730        outfile.variables['tri_full_flag'][:] = tri_full_flag.astype(num.int32)
    678731
    679732
  • trunk/anuga_core/source/anuga/utilities/sww_merge.py

    r8278 r8281  
    1414
    1515    output = domain_global_name+".sww"
    16     swwfiles = [ domain_global_name+"_P"+str(v)+"_"+str(np)+".sww" for v in range(np)]
     16    swwfiles = [ domain_global_name+"_P"+str(np)+"_"+str(v)+".sww" for v in range(np)]
    1717
    1818    _sww_merge(swwfiles, output, verbose)
     
    132132
    133133
     134
     135    from anuga.coordinate_transforms.geo_reference import Geo_reference
     136    geo_reference = Geo_reference()
    134137   
    135    
    136     sww.store_triangulation(fido, points, out_tris)
     138    sww.store_triangulation(fido, points, out_tris, points_georeference=geo_reference)
    137139
    138140    fido.order      = order
  • trunk/anuga_core/source/anuga/utilities/test_file_utils.py

    r8124 r8281  
    77from anuga.utilities.file_utils import copy_code_files, get_all_swwfiles
    88from anuga.utilities.file_utils import del_dir
    9 from sww_merge import sww_merge
     9from sww_merge import sww_merge, _sww_merge
    1010
    1111
     
    139139               
    140140        outfile = 'test_out.sww'
    141         sww_merge(['test1.sww', 'test2.sww'], outfile)
     141        _sww_merge(['test1.sww', 'test2.sww'], outfile)
    142142        self.failUnless(os.access(outfile, os.F_OK)) 
    143143       
  • trunk/anuga_core/source/anuga_parallel/distribute_mesh.py

    r8260 r8281  
    147147        del npart
    148148
    149         # Sometimes (usu. on x86_64), partMeshNodal returnes an array of zero
     149        # Sometimes (usu. on x86_64), partMeshNodal returns an array of zero
    150150        # dimensional arrays. Correct this.
    151151        if type(epart[0]) == num.ndarray:
     
    808808#########################################################
    809809
    810 def build_local_GA(nodes, triangles, boundaries, tri_index):
     810def build_local_GA(nodes, triangles, boundaries, tri_map):
    811811
    812812    Nnodes =len(nodes)
     
    823823        if nodes[i][0] > NGlobal:
    824824            NGlobal = nodes[i][0]
    825     index = num.zeros(int(NGlobal)+1, num.int)
    826     num.put(index, num.take(nodes, (0,), 1).astype(num.int), \
     825
     826    node_map = -1*num.ones(int(NGlobal)+1, num.int)
     827
     828    num.put(node_map, num.take(nodes, (0,), 1).astype(num.int), \
    827829        num.arange(Nnodes))
    828830       
     
    830832
    831833    GAtriangles = num.zeros((Ntriangles, 3), num.int)
    832     GAtriangles[:,0] = num.take(index, triangles[:,0])
    833     GAtriangles[:,1] = num.take(index, triangles[:,1])
    834     GAtriangles[:,2] = num.take(index, triangles[:,2])
     834    GAtriangles[:,0] = num.take(node_map, triangles[:,0])
     835    GAtriangles[:,1] = num.take(node_map, triangles[:,1])
     836    GAtriangles[:,2] = num.take(node_map, triangles[:,2])
    835837
    836838    # Change the triangle numbering in the boundaries
     
    838840    GAboundaries = {}
    839841    for b in boundaries:
    840         GAboundaries[tri_index[b[0]], b[1]] = boundaries[b]
     842        GAboundaries[tri_map[b[0]], b[1]] = boundaries[b]
    841843       
    842     del (index)
    843    
    844     return GAnodes, GAtriangles, GAboundaries
     844   
     845    return GAnodes, GAtriangles, GAboundaries, node_map
    845846
    846847
     
    868869#########################################################
    869870
    870 def build_local_commun(index, ghostc, fullc, nproc):
     871def build_local_commun(tri_map, ghostc, fullc, nproc):
    871872
    872873    # Initialise
     
    886887            ghost_recv[c] = [0, 0]
    887888            ghost_recv[c][1] = d
    888             ghost_recv[c][0] = num.take(index, d)
     889            ghost_recv[c][0] = num.take(tri_map, d)
    889890           
    890891    # Build a temporary copy of the full_send dictionary
     
    899900                tmp_send[neigh] = []
    900901            tmp_send[neigh].append([global_id, \
    901                                     index[global_id]])
     902                                    tri_map[global_id]])
    902903
    903904    # Extract the full send information and put it in the form
     
    958959        if id > NGlobal:
    959960            NGlobal = id
    960     index = num.zeros(int(NGlobal)+1, num.int)
    961     index[lower_t:upper_t]=num.arange(upper_t-lower_t)
     961    #index = num.zeros(int(NGlobal)+1, num.int)
     962    tri_map = -1*num.ones(int(NGlobal)+1, num.int)
     963    tri_map[lower_t:upper_t]=num.arange(upper_t-lower_t)
    962964    for i in range(len(submesh["ghost_triangles"])):
    963         index[submesh["ghost_triangles"][i][0]] = i+upper_t-lower_t
     965        tri_map[submesh["ghost_triangles"][i][0]] = i+upper_t-lower_t
    964966   
    965967    # Change the node numbering (and update the numbering in the
    966968    # triangles)
    967969
    968     [GAnodes, GAtriangles, GAboundary] = build_local_GA(nodes, triangles, boundaries, index)
     970    [GAnodes, GAtriangles, GAboundary, node_map] = \
     971    build_local_GA(nodes, triangles, boundaries, tri_map)
    969972
    970973    # Extract the local quantities
     
    984987    fcommun = submesh["full_commun"]
    985988    [ghost_rec, full_send] = \
    986                 build_local_commun(index, gcommun, fcommun, nproc)
    987 
    988     # Clean up before exiting
    989 
    990     del(index)
     989                build_local_commun(tri_map, gcommun, fcommun, nproc)
     990
    991991
    992992    return GAnodes, GAtriangles, GAboundary, quantities, ghost_rec, \
    993            full_send
     993           full_send, tri_map, node_map
    994994
    995995
     
    12841284    # datastructure
    12851285
    1286     [GAnodes, GAtriangles, boundary, quantities, ghost_rec, full_send] = \
     1286    [GAnodes, GAtriangles, boundary, quantities, \
     1287     ghost_rec, full_send, tri_map, node_map] = \
    12871288              build_local_mesh(submesh_cell, lower_t, upper_t, \
    12881289                               numproc)
     
    12901291    return GAnodes, GAtriangles, boundary, quantities,\
    12911292           ghost_rec, full_send,\
    1292            number_of_full_nodes, number_of_full_triangles
     1293           number_of_full_nodes, number_of_full_triangles, tri_map, node_map
     1294         
    12931295
    12941296
     
    13271329
    13281330    numprocs = len(triangles_per_proc)
    1329     points, vertices, boundary, quantities, ghost_recv_dict, full_send_dict = \
     1331    points, vertices, boundary, quantities, ghost_recv_dict, \
     1332            full_send_dict, tri_map, node_map = \
    13301333            build_local_mesh(submesh_cell, 0, triangles_per_proc[0], numprocs)
    13311334
    13321335
    13331336    return  points, vertices, boundary, quantities, ghost_recv_dict, \
    1334            full_send_dict
    1335 
    1336 
    1337 
    1338 
     1337           full_send_dict, tri_map, node_map
     1338           
     1339
     1340
     1341
     1342
  • trunk/anuga_core/source/anuga_parallel/parallel_api.py

    r8272 r8281  
    9999                ghost_recv_dict, full_send_dict,\
    100100                number_of_full_nodes, number_of_full_triangles,\
    101                 s2p_map, p2s_map =\
     101                s2p_map, p2s_map, tri_map, node_map =\
    102102                distribute_mesh(domain, verbose=verbose)
     103
     104        number_of_global_triangles = len(tri_map)
     105        number_of_global_nodes = len(node_map)
     106
     107        # Extract l2g maps
     108        tri_l2g  = extract_l2g_map(tri_map)
     109        node_l2g = extract_l2g_map(node_map)
    103110
    104111        # Send serial to parallel (s2p) and parallel to serial (p2s) triangle mapping to proc 1 .. numprocs
     
    115122        points, vertices, boundary, quantities,\
    116123                ghost_recv_dict, full_send_dict,\
    117                 number_of_full_nodes, number_of_full_triangles =\
     124                number_of_full_nodes, number_of_full_triangles, \
     125                tri_map, node_map =\
    118126                rec_submesh(0, verbose)
    119127
     128        number_of_global_triangles = len(tri_map)
     129        number_of_global_nodes = len(node_map)
     130
     131        # Extract l2g maps
     132        tri_l2g  = extract_l2g_map(tri_map)
     133        node_l2g = extract_l2g_map(node_map)
     134       
    120135        # Recieve serial to parallel (s2p) and parallel to serial (p2s) triangle mapping
    121136        s2p_map = receive(0)
     
    136151                             number_of_full_triangles=number_of_full_triangles,
    137152                             geo_reference=georef,
    138                              tri_map = s2p_map,
    139                              inv_tri_map = p2s_map) ## jj added this
     153                             number_of_global_triangles = number_of_global_triangles,
     154                             number_of_global_nodes = number_of_global_nodes,
     155                             s2p_map = s2p_map,
     156                             p2s_map = p2s_map, ## jj added this
     157                             tri_l2g = tri_l2g, ## SR added this
     158                             node_l2g = node_l2g)
    140159
    141160    #------------------------------------------------------------------------
     
    177196    # Subdivide the mesh
    178197    if verbose: print 'Subdivide mesh'
    179     nodes, triangles, boundary, triangles_per_proc, quantities, s2p_map, p2s_map = \
     198    nodes, triangles, boundary, triangles_per_proc, quantities, \
     199           s2p_map, p2s_map = \
    180200           pmesh_divide_metis_with_map(domain, numprocs)
    181201
    182202    #PETE: s2p_map (maps serial domain triangles to parallel domain triangles)
    183     #p2_map (maps parallel domain triangles to domain triangles)
     203    #      sp2_map (maps parallel domain triangles to domain triangles)
    184204
    185205
     
    204224
    205225    # Build the local mesh for processor 0
    206     points, vertices, boundary, quantities, ghost_recv_dict, full_send_dict =\
     226    points, vertices, boundary, quantities, \
     227            ghost_recv_dict, full_send_dict, tri_map, node_map =\
    207228              extract_hostmesh(submesh, triangles_per_proc)
    208229
     
    240261    return points, vertices, boundary, quantities,\
    241262           ghost_recv_dict, full_send_dict,\
    242            number_of_full_nodes, number_of_full_triangles, s2p_map, p2s_map
    243    
    244 
    245 
    246 
     263           number_of_full_nodes, number_of_full_triangles, \
     264           s2p_map, p2s_map, tri_map, node_map
     265   
     266
     267
     268def extract_l2g_map(map):
     269    # Extract l2g_map
     270
     271    import numpy as num
     272   
     273    b = num.arange(len(map))
     274
     275    l_ids = num.extract(map>-1,map)
     276    g_ids = num.extract(map>-1,b)
     277
     278#    print len(g_ids)
     279#    print len(l_ids)
     280#    print l_ids
     281
     282    l2g = num.zeros_like(g_ids)
     283    l2g[l_ids] = g_ids
     284
     285    return l2g
     286
  • trunk/anuga_core/source/anuga_parallel/parallel_shallow_water.py

    r8272 r8281  
    3333                 number_of_full_triangles=None,
    3434                 geo_reference=None,
    35                  tri_map=None,
    36                  inv_tri_map=None): #jj added this
     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
    3741
    3842        Domain.__init__(self,
     
    4852                        geo_reference=geo_reference) #jj added this
    4953       
     54
     55
     56        self.parallel = True
     57
    5058        # PETE: Find the number of full nodes and full triangles, this is a temporary fix
    5159        # until the bug with get_number_of_full_[nodes|triangles]() is fixed.
     
    6472
    6573        self.global_name = 'domain'
    66         self.tri_map = tri_map
    67         self.inv_tri_map = inv_tri_map
     74
     75        self.number_of_global_triangles=number_of_global_triangles
     76        self.number_of_global_nodes = number_of_global_nodes
     77
     78        self.s2p_map = s2p_map
     79        self.p2s_map = p2s_map
     80        self.tri_l2g = tri_l2g
     81        self.node_l2g = node_l2g
    6882
    6983
     
    7892
    7993        # Call parents method with processor number attached.
    80         Domain.set_name(self, name + '_P%d_%d' %(self.processor, self.numproc))
     94        Domain.set_name(self, name + '_P%d_%d' %(self.numproc, self.processor))
    8195
    8296
  • trunk/anuga_core/source/anuga_parallel/run_parallel_sw_merimbula.py

    r8277 r8281  
    3535
    3636
    37 from anuga_parallel import distribute, myid, numprocs, finalize
     37from anuga_parallel import distribute, myid, numprocs, finalize, barrier
    3838
    3939
     
    4242#--------------------------------------------------------------------------
    4343
    44 mesh_filename = "merimbula_10785_1.tsh" ; x0 = 756000.0 ; x1 = 756500.0
     44#mesh_filename = "merimbula_10785_1.tsh" ; x0 = 756000.0 ; x1 = 756500.0
    4545#mesh_filename = "merimbula_43200.tsh"   ; x0 = 756000.0 ; x1 = 756500.0
    46 #mesh_filename = "test-100.tsh" ; x0 = 0.25 ; x1 = 0.5
     46mesh_filename = "test-100.tsh" ; x0 = 0.25 ; x1 = 0.5
    4747yieldstep = 20
    4848finaltime = 500
     
    7070    domain = create_domain_from_file(mesh_filename)
    7171    domain.set_quantity('stage', Set_Stage(x0, x1, 2.0))
     72
     73    print domain.statistics()
     74    print domain.get_extent()
     75    print domain.get_extent(absolute=True)
     76    print domain.geo_reference
    7277else:
    7378    domain = None
     
    8691#domain.set_CFL(0.7)
    8792#domain.set_beta(1.5)
    88 domain.set_name('meribula')
     93domain.set_name('merimbula')
     94
     95for p in range(numprocs):
     96    if myid == p:
     97        print 'P%d'%p
     98        print domain.get_extent()
     99        print domain.get_extent(absolute=True)
     100        print domain.geo_reference
     101        #print domain.tri_map
     102        #print domain.inv_tri_map
     103        print domain.tri_l2g
     104        print domain.node_l2g
     105    else:
     106        pass
     107   
     108    barrier()
     109
    89110
    90111#------------------------------------------------------------------------------
  • trunk/anuga_core/source/anuga_parallel/test_parallel_distribute_mesh.py

    r8115 r8281  
    195195        # Test extract_hostmesh
    196196        #----------------------------------------------------------------------------------
    197         points, vertices, boundary, quantities, ghost_recv_dict, full_send_dict  =\
     197        points, vertices, boundary, quantities, \
     198                ghost_recv_dict, full_send_dict, tri_map, node_map  =\
    198199        extract_hostmesh(submesh, triangles_per_proc)
    199200
     
    222223        # Test rec_submesh
    223224        #----------------------------------------------------------------------------------
    224         points, vertices, boundary, quantities, ghost_recv_dict, full_send_dict, no_full_nodes, no_full_trigs = rec_submesh(0, verbose=False)   
     225        points, vertices, boundary, quantities, \
     226                ghost_recv_dict, full_send_dict, \
     227                no_full_nodes, no_full_trigs, tri_map, node_map  = \
     228                rec_submesh(0, verbose=False)   
    225229
    226230        if myid == 1:
     
    235239            true_full_send =  {0: [num.array([0, 1, 3, 4, 5]), num.array([ 5,  6,  8,  9, 10])], 2: [num.array([0, 1]), num.array([5, 6])]}
    236240
     241            true_tri_map  = num.array([ 6,  7,  8, -1,  9,  0, 1,  2,  3,  4,  5, 10, 11])
     242
     243            true_node_map = num.array([ 0,  1,  2,  7,  3,  4, -1,  8,  9,  5,  6, 10, 11])
     244
     245            assert_(num.allclose(tri_map,   true_tri_map))
     246            assert_(num.allclose(node_map,   true_node_map))
    237247            assert_(num.allclose(points,   true_points))
    238248            assert_(num.allclose(vertices, true_vertices))
     
    249259            true_vertices =  [[1, 5, 0], [1, 6, 2], [3, 6, 1], [4, 6, 3], [2, 6, 4], [2, 5, 1], [2, 10, 8], [4, 10, 2], [9, 10, 4], [0, 5, 7], [7, 5, 2]]
    250260
    251 
    252261            true_ghost_recv =   {0: [num.array([5, 6, 7, 8]), num.array([0, 1, 2, 3])], 1: [num.array([ 9, 10]), num.array([5, 6])]}
    253262
    254263            true_full_send =   {0: [num.array([0, 1, 2, 3, 4]), num.array([11, 12, 13, 14, 15])], 1: [num.array([0, 1]), num.array([11, 12])]}
    255264
    256 
     265            true_tri_map  = num.array([5, 6, 7, 8, -1, 9, 10, -1, -1, -1, -1, 0, 1, 2, 3, 4, -1])
     266
     267            true_node_map = num.array([ 0,  7, -1,  1,  2,  8 , 3,  4,  9,  5, -1,  6, 10])
     268
     269            assert_(num.allclose(tri_map,   true_tri_map))
     270            assert_(num.allclose(node_map,   true_node_map))
    257271            assert_(num.allclose(points,   true_points))
    258272            assert_(num.allclose(vertices, true_vertices))
Note: See TracChangeset for help on using the changeset viewer.