Changeset 6584


Ignore:
Timestamp:
Mar 23, 2009, 3:55:37 PM (16 years ago)
Author:
ole
Message:

Fixed MUX direction bug (again). Mux stores south as positive whereas ANUGA
regards north as positive. That means that the y component of the momentum
vector was swapped.
This bug was fixed in 2007 in changeset:4348 but had crept in again
when the mux reader was changed.
The new tests were also based on the wrong assumption of the north
direction being positive..

Location:
anuga_core/source/anuga/shallow_water
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/shallow_water/data_manager.py

    r6556 r6584  
    55425542    output:
    55435543      basename_out: name of sts file in which mux2 data is stored.
     5544     
     5545     
     5546     
     5547    NOTE: South is positive in mux files so sign of y-component of velocity is reverted
    55445548    """
    55455549
     
    57725776
    57735777            xmomentum[j,i] = ua * h
    5774             ymomentum[j,i] = va * h
     5778            ymomentum[j,i] = -va * h # South is positive in mux files
     5779
    57755780
    57765781    outfile.close()
  • anuga_core/source/anuga/shallow_water/test_data_manager.py

    r6556 r6584  
    55205520        na and va quantities will be the Easting values.
    55215521        Depth and ua will be the Northing value.
     5522       
     5523        The mux file format has south as positive so
     5524        this function will swap the sign for va. 
    55225525        """
    55235526
     
    55595562            quantities_init[0].append(this_ha) # HA
    55605563            quantities_init[1].append(this_ua) # UA
    5561             quantities_init[2].append(this_va) # VA
     5564            quantities_init[2].append(this_va) # VA 
    55625565               
    55635566        file_handle, base_name = tempfile.mkstemp("")
     
    60376040                quantities_init[2].append(num.ones(time_step_count,num.Float)*this_va) #
    60386041            else:
    6039                 quantities_init[2].append(va[i])           
     6042                quantities_init[2].append(-va[i]) # South is negative in MUX
    60406043
    60416044        file_handle, base_name = tempfile.mkstemp("")
     
    61596162        assert num.allclose(xvelocity,ua)
    61606163        msg='incorrect gauge va time series returned'
    6161         assert num.allclose(yvelocity,va)
     6164        assert num.allclose(yvelocity, -va)
    61626165
    61636166    def test_urs2sts_read_mux2_pyII(self):
     
    62116214
    62126215        msg='Incorrect gauge depths returned'
    6213         assert num.allclose(elevation,-depth),msg
     6216        assert num.allclose(elevation, -depth),msg
    62146217        msg='incorrect gauge height time series returned'
    6215         assert num.allclose(stage,ha)
     6218        assert num.allclose(stage, ha)
    62166219        msg='incorrect gauge ua time series returned'
    6217         assert num.allclose(xvelocity,ua)
     6220        assert num.allclose(xvelocity, ua)
    62186221        msg='incorrect gauge va time series returned'
    6219         assert num.allclose(yvelocity,va)
     6222        assert num.allclose(yvelocity, -va) # South is positive in MUX
    62206223
    62216224    def test_urs2sts_read_mux2_pyIII(self):
     
    62896292        assert num.allclose(xvelocity,ua)
    62906293        msg='incorrect gauge va time series returned'
    6291         assert num.allclose(yvelocity,va)
     6294        assert num.allclose(yvelocity, -va) # South is positive in mux
    62926295       
    62936296
     
    63916394                if j == 0: assert num.allclose(data[i][:parameters_index], ha0[permutation[i], :])
    63926395                if j == 1: assert num.allclose(data[i][:parameters_index], ua0[permutation[i], :])
    6393                 if j == 2: assert num.allclose(data[i][:parameters_index], va0[permutation[i], :])
     6396                if j == 2: assert num.allclose(data[i][:parameters_index], -va0[permutation[i], :])
    63946397       
    63956398
     
    65976600                        x = unpack('f',f.read(4))[0]
    65986601                        #print time, x, q_time[time, point_i]
     6602                        if q == 'VA': x = -x # South is positive in MUX
    65996603                        assert abs(q_time[time, point_i]-x)<epsilon
    66006604
     
    66706674                    #print 'v ', data[i][:parameters_index][8]                   
    66716675               
    6672                     assert num.allclose(data[i][:parameters_index], va1[permutation[i], :])
     6676                    # South is positive in MUX
     6677                    assert num.allclose(data[i][:parameters_index], -va1[permutation[i], :])
    66736678                   
    66746679       
     
    71927197                msg = 'urs stage is not equal to sts stage for for source %i and station %i' %(source_number,j)
    71937198                assert num.allclose(urs_stage[index_start_urs_z:index_end_urs_z],
    7194                                 sts_stage[index_start_stage:index_end_stage],
    7195                                 rtol=1.0e-6, atol=1.0e-5 ), msg                               
     7199                                    sts_stage[index_start_stage:index_end_stage],
     7200                                    rtol=1.0e-6, atol=1.0e-5 ), msg                               
    71967201                               
    71977202                # check that urs e velocity and sts xmomentum are the same
     
    72067211                msg = 'urs n velocity is not equivalent to sts y momentum for source %i and station %i' %(source_number,j)
    72077212                assert num.allclose(urs_n[index_start_urs_n:index_end_urs_n]*(urs_stage[index_start_urs_n:index_end_urs_n]-elevation[j]),
    7208                                 sts_ymom[index_start_y:index_end_y],
     7213                                -sts_ymom[index_start_y:index_end_y],
    72097214                                rtol=1.0e-5, atol=1.0e-4 ), msg
    72107215                                               
Note: See TracChangeset for help on using the changeset viewer.