Changeset 7352


Ignore:
Timestamp:
Aug 11, 2009, 4:38:46 PM (15 years ago)
Author:
steve
Message:

Updating cross section

File:
1 edited

Legend:

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

    r7350 r7352  
    557557    # @note 'polyline' may contain multiple sections allowing complex shapes.
    558558    # @note Assume absolute UTM coordinates.
     559    def new_get_energy_through_cross_section(self, polyline,
     560                                         kind='total',
     561                                         verbose=False):
     562        """Obtain average energy head [m] across specified cross section.
     563
     564        Inputs:
     565            polyline: Representation of desired cross section - it may contain
     566                      multiple sections allowing for complex shapes. Assume
     567                      absolute UTM coordinates.
     568                      Format [[x0, y0], [x1, y1], ...]
     569            kind:     Select which energy to compute.
     570                      Options are 'specific' and 'total' (default)
     571
     572        Output:
     573            E: Average energy [m] across given segments for all stored times.
     574
     575        The average velocity is computed for each triangle intersected by
     576        the polyline and averaged weighted by segment lengths.
     577
     578        The typical usage of this function would be to get average energy of
     579        flow in a channel, and the polyline would then be a cross section
     580        perpendicular to the flow.
     581
     582        #FIXME (Ole) - need name for this energy reflecting that its dimension
     583        is [m].
     584        """
     585
     586
     587
     588        cross_section = Cross_section(self, polyline, verbose)
     589
     590        return cross_section.get_energy_through_cross_section()
     591
     592
     593    ##
     594    # @brief
     595    # @param polyline Representation of desired cross section.
     596    # @param kind Select energy type to compute ('specific' or 'total').
     597    # @param verbose True if this method is to be verbose.
     598    # @note 'polyline' may contain multiple sections allowing complex shapes.
     599    # @note Assume absolute UTM coordinates.
    559600    def get_energy_through_cross_section(self, polyline,
    560601                                         kind='total',
     
    646687        return average_energy
    647688
     689
    648690    ##
    649691    # @brief Run integrity checks on shallow water domain.
     
    679721        else:
    680722            distribute_using_vertex_limiter(self)
     723
     724
    681725
    682726    ##
     
    27742818 
    27752819
     2820    ##
     2821    # @brief calculate current energy flow through cross section
     2822    def get_energy_through_cross_section(self):
     2823        """Obtain average energy head [m] across specified cross section.
     2824
     2825        Output:
     2826            E: Average energy [m] across given segments for all stored times.
     2827
     2828        The average velocity is computed for each triangle intersected by
     2829        the polyline and averaged weighted by segment lengths.
     2830
     2831        The typical usage of this function would be to get average energy of
     2832        flow in a channel, and the polyline would then be a cross section
     2833        perpendicular to the flow.
     2834
     2835        #FIXME (Ole) - need name for this energy reflecting that its dimension
     2836        is [m].
     2837        """
     2838
     2839        # Get interpolated values
     2840        stage = self.domain.get_quantity('stage')
     2841        elevation = self.domain.get_quantity('elevation')
     2842        xmomentum = self.domain.get_quantity('xmomentum')
     2843        ymomentum = self.domain.get_quantity('ymomentum')
     2844
     2845        w = stage.get_values(interpolation_points=midpoints, use_cache=True)
     2846        z = elevation.get_values(interpolation_points=midpoints, use_cache=True)
     2847        uh = xmomentum.get_values(interpolation_points=midpoints,
     2848                                  use_cache=True)
     2849        vh = ymomentum.get_values(interpolation_points=midpoints,
     2850                                  use_cache=True)
     2851        h = w-z                # Depth
     2852
     2853        # Compute total length of polyline for use with weighted averages
     2854        total_line_length = 0.0
     2855        for segment in segments:
     2856            total_line_length += segment.length
     2857
     2858        # Compute and sum flows across each segment
     2859        average_energy = 0.0
     2860        for i in range(len(w)):
     2861            # Average velocity across this segment
     2862            if h[i] > epsilon:
     2863                # Use protection against degenerate velocities
     2864                u = uh[i]/(h[i] + h0/h[i])
     2865                v = vh[i]/(h[i] + h0/h[i])
     2866            else:
     2867                u = v = 0.0
     2868
     2869            speed_squared = u*u + v*v
     2870            kinetic_energy = 0.5*speed_squared/g
     2871
     2872            if kind == 'specific':
     2873                segment_energy = h[i] + kinetic_energy
     2874            elif kind == 'total':
     2875                segment_energy = w[i] + kinetic_energy
     2876            else:
     2877                msg = 'Energy kind must be either "specific" or "total".'
     2878                msg += ' I got %s' %kind
     2879
     2880            # Add to weighted average
     2881            weigth = segments[i].length/total_line_length
     2882            average_energy += segment_energy*weigth
     2883
     2884        return average_energy
     2885
    27762886
    27772887
Note: See TracChangeset for help on using the changeset viewer.