Changeset 6143


Ignore:
Timestamp:
Jan 13, 2009, 11:12:52 AM (10 years ago)
Author:
ole
Message:

Simplified flow protection at the inlet and wrote better diagnostics.

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

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/culvert_flows/Test_Culvert_Flat_Water_Lev.py

    r6059 r6143  
    2121     Transmissive_boundary, Time_boundary
    2222
    23 from anuga.culvert_flows.culvert_class import Culvert_flow_rating, Culvert_flow_energy
     23from anuga.culvert_flows.culvert_class import Culvert_flow
    2424from anuga.culvert_flows.culvert_routines import boyd_generalised_culvert_model
    2525     
     
    104104
    105105
    106 culvert_rating = Culvert_flow_rating(domain,
     106culvert_rating = Culvert_flow(domain,
    107107                       culvert_description_filename='example_rating_curve.csv',
    108108                       end_point0=[9.0, 2.5],
     
    111111
    112112
    113 culvert_energy = Culvert_flow_energy(domain,
     113culvert_energy = Culvert_flow(domain,
    114114                       label='Culvert No. 1',
    115115                       description='This culvert is a test unit 1.2m Wide by 0.75m High',   
  • anuga_core/source/anuga/culvert_flows/culvert_class.py

    r6142 r6143  
    2626class Above_interval(Exception): pass
    2727
     28# FIXME(Ole): Take a good hard look at logging here
    2829
    2930
     
    404405        """Adjust Q downwards depending on available water at inlet
    405406        """
     407   
     408        if delta_t < epsilon:
     409            # No need to adjust if time step is very small or zero
     410            # In this case the possible flow will be very large
     411            # anyway.
     412            return Q
    406413       
    407414        # Short hands
     
    410417        time = domain.get_time()
    411418        I = self.inlet
    412 
    413        
    414         # Find triangle with the smallest available flow
    415         max_Q = sys.maxint
    416         if delta_t > 0:
    417             for i in I.exchange_indices:
    418                 stage = dq['stage'].get_values(location='centroids',
    419                                                indices=[i])[0]
    420                 elevation = dq['elevation'].get_values(location='centroids',
    421                                                        indices=[i])[0]       
    422                 depth = stage-elevation
    423                 area = domain.areas[i]
    424 
    425                 # Possible rate based on this triangle
    426                 Q_possible = depth*I.exchange_area/delta_t
    427                                
    428                 # Use triangle with least depth to determine total flow
    429                 max_Q = min(max_Q, Q_possible)
    430 
    431            
     419        idx = I.exchange_indices   
     420
     421        # Find triangle with the smallest depth
     422        stage = dq['stage'].get_values(location='centroids',
     423                                               indices=[idx])
     424        elevation = dq['elevation'].get_values(location='centroids',
     425                                               indices=[idx])       
     426        depth = stage-elevation
     427        min_depth = min(depth.flat)
     428
     429        # Compute possible flow for exchange region based on
     430        # triangle with smallest depth
     431        max_Q = min_depth*I.exchange_area/delta_t       
     432       
     433        # Calculate the minimum in absolute terms of
     434        # the requsted flow and the possible flow
    432435        Q_reduced = sign(Q)*min(abs(Q), abs(max_Q))
    433436       
    434437        if abs(Q_reduced) < abs(Q):
    435             msg = '%.2fs: Computed extraction for this time interval (Q*dt) is ' % time
    436             msg += 'greater than current volume (V) at inlet.\n'
    437             msg += ' Q will be reduced from %.2f m^3/s to %.2f m^3/s.' % (Q, Q_reduced)
     438            msg = '%.2fs: Requested flow is ' % time
     439            msg += 'greater than what is supported by the smallest '
     440            msg += 'depth at inlet exchange area:\n        '
     441            msg += 'h_min*inlet_area/delta_t = %.2f*%.2f/%.2f '\
     442                % (min_depth, I.exchange_area, delta_t)
     443            msg += ' = %.2f m^3/s\n        ' % Q_reduced
     444            msg += 'Q will be reduced from %.2f m^3/s to %.2f m^3/s.' % (Q, Q_reduced)
    438445            if self.verbose is True:
    439446                print msg
     
    521528            # Flow will be purely controlled by uphill outlet face
    522529            if self.verbose is True:
    523                 print '%.2f - WARNING: Flow is running uphill.' % time
     530                print '%.2fs - WARNING: Flow is running uphill.' % time
    524531           
    525532        if hasattr(self, 'log_filename'):
Note: See TracChangeset for help on using the changeset viewer.