Changeset 4827


Ignore:
Timestamp:
Nov 19, 2007, 11:39:08 AM (16 years ago)
Author:
ole
Message:

Better timestepping statistics including derived quantities such as velocity and Froude numbers.

Location:
anuga_core/source/anuga
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/abstract_2d_finite_volumes/domain.py

    r4805 r4827  
    648648       
    649649
     650        # qwidth determines the text field used for quantities
     651        qwidth = self.qwidth = 12
     652
     653
    650654        msg = ''
    651655        if self.min_timestep == self.max_timestep:
     
    712716           
    713717            # Find index of largest computed flux speed
    714             k = argmax(self.max_speed)
     718            k = self.k = argmax(self.max_speed)
    715719
    716720            x, y = self.get_centroid_coordinates()[k]
     
    727731                msg += '(timestep=%.6f)\n' %(0)     
    728732           
    729             # Report all quantity values at vertices
    730             msg += '    Quantity \t vertex values\t\t\t\t\t centroid values\n'
     733            # Report all quantity values at vertices, edges and centroid
     734            msg += '    Quantity'
     735            msg += '------------'
    731736            for name in self.quantities:
    732737                q = self.quantities[name]
    733738
    734739                V = q.get_values(location='vertices', indices=[k])[0]
     740                E = q.get_values(location='edges', indices=[k])[0]
    735741                C = q.get_values(location='centroids', indices=[k])                 
    736742               
    737                 s = '    %s:\t %.4f,\t %.4f,\t %.4f,\t %.4f\n'\
    738                     %(name, V[0], V[1], V[2], C[0])                 
     743                s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     744                     %(name.ljust(qwidth), V[0], V[1], V[2])
     745
     746                s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     747                     %(name.ljust(qwidth), E[0], E[1], E[2])
     748
     749                s += '    %s: centroid_value = %.4f\n'\
     750                     %(name.ljust(qwidth), C[0])                               
    739751
    740752                msg += s
  • anuga_core/source/anuga/shallow_water/shallow_water_domain.py

    r4769 r4827  
    461461        self.writer.store_timestep(name)
    462462
     463       
     464    def timestepping_statistics(self, track_speeds=False):
     465        """Return string with time stepping statistics for printing or logging
     466
     467        Optional boolean keyword track_speeds decides whether to report
     468        location of smallest timestep as well as a histogram and percentile
     469        report.
     470        """
     471
     472        from Numeric import sqrt
     473        from anuga.config import epsilon, g               
     474
     475
     476        # Call basic machinery from parent class
     477        msg = Generic_Domain.timestepping_statistics(self, track_speeds)
     478
     479        if track_speeds is True:
     480
     481            # qwidth determines the text field used for quantities
     482            qwidth = self.qwidth
     483       
     484            # Triangle with maximum speed
     485            k = self.k
     486
     487            # Report some derived quantities at vertices, edges and centroid
     488            # specific to the shallow water wave equation
     489
     490            z = self.quantities['elevation']
     491            w = self.quantities['stage']           
     492
     493            Vw = w.get_values(location='vertices', indices=[k])[0]
     494            Ew = w.get_values(location='edges', indices=[k])[0]
     495            Cw = w.get_values(location='centroids', indices=[k])
     496
     497            Vz = z.get_values(location='vertices', indices=[k])[0]
     498            Ez = z.get_values(location='edges', indices=[k])[0]
     499            Cz = z.get_values(location='centroids', indices=[k])                             
     500               
     501
     502            name = 'depth'
     503            Vh = Vw-Vz
     504            Eh = Ew-Ez
     505            Ch = Cw-Cz
     506           
     507            s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     508                 %(name.ljust(qwidth), Vh[0], Vh[1], Vh[2])
     509           
     510            s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     511                 %(name.ljust(qwidth), Eh[0], Eh[1], Eh[2])
     512           
     513            s += '    %s: centroid_value = %.4f\n'\
     514                 %(name.ljust(qwidth), Ch[0])                               
     515           
     516            msg += s
     517
     518            uh = self.quantities['xmomentum']
     519            vh = self.quantities['ymomentum']
     520
     521            Vuh = uh.get_values(location='vertices', indices=[k])[0]
     522            Euh = uh.get_values(location='edges', indices=[k])[0]
     523            Cuh = uh.get_values(location='centroids', indices=[k])
     524           
     525            Vvh = vh.get_values(location='vertices', indices=[k])[0]
     526            Evh = vh.get_values(location='edges', indices=[k])[0]
     527            Cvh = vh.get_values(location='centroids', indices=[k])
     528
     529            # Speeds in each direction
     530            Vu = Vuh/(Vh + epsilon)
     531            Eu = Euh/(Eh + epsilon)
     532            Cu = Cuh/(Ch + epsilon)           
     533            name = 'U'
     534            s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     535                 %(name.ljust(qwidth), Vu[0], Vu[1], Vu[2])
     536           
     537            s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     538                 %(name.ljust(qwidth), Eu[0], Eu[1], Eu[2])
     539           
     540            s += '    %s: centroid_value = %.4f\n'\
     541                 %(name.ljust(qwidth), Cu[0])                               
     542           
     543            msg += s
     544
     545            Vv = Vvh/(Vh + epsilon)
     546            Ev = Evh/(Eh + epsilon)
     547            Cv = Cvh/(Ch + epsilon)           
     548            name = 'V'
     549            s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     550                 %(name.ljust(qwidth), Vv[0], Vv[1], Vv[2])
     551           
     552            s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     553                 %(name.ljust(qwidth), Ev[0], Ev[1], Ev[2])
     554           
     555            s += '    %s: centroid_value = %.4f\n'\
     556                 %(name.ljust(qwidth), Cv[0])                               
     557           
     558            msg += s
     559
     560
     561            # Froude number in each direction
     562            name = 'Froude (x)'
     563            Vfx = Vu/(sqrt(g*Vh) + epsilon)
     564            Efx = Eu/(sqrt(g*Eh) + epsilon)
     565            Cfx = Cu/(sqrt(g*Ch) + epsilon)
     566           
     567            s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     568                 %(name.ljust(qwidth), Vfx[0], Vfx[1], Vfx[2])
     569           
     570            s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     571                 %(name.ljust(qwidth), Efx[0], Efx[1], Efx[2])
     572           
     573            s += '    %s: centroid_value = %.4f\n'\
     574                 %(name.ljust(qwidth), Cfx[0])                               
     575           
     576            msg += s
     577
     578
     579            name = 'Froude (y)'
     580            Vfy = Vv/(sqrt(g*Vh) + epsilon)
     581            Efy = Ev/(sqrt(g*Eh) + epsilon)
     582            Cfy = Cv/(sqrt(g*Ch) + epsilon)
     583           
     584            s  = '    %s: vertex_values =  %.4f,\t %.4f,\t %.4f\n'\
     585                 %(name.ljust(qwidth), Vfy[0], Vfy[1], Vfy[2])
     586           
     587            s += '    %s: edge_values =    %.4f,\t %.4f,\t %.4f\n'\
     588                 %(name.ljust(qwidth), Efy[0], Efy[1], Efy[2])
     589           
     590            s += '    %s: centroid_value = %.4f\n'\
     591                 %(name.ljust(qwidth), Cfy[0])                               
     592           
     593            msg += s           
     594
     595               
     596
     597        return msg
     598       
     599       
    463600
    464601#=============== End of class Shallow Water Domain ===============================
  • anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

    r4823 r4827  
    11721172        #--------------------------------------------------------------
    11731173        for t in domain.evolve(yieldstep = 0.1, finaltime = 3.0):
     1174            #print domain.timestepping_statistics(track_speeds=True)
    11741175            #domain.write_time()
    11751176            pass
Note: See TracChangeset for help on using the changeset viewer.