Changeset 5290


Ignore:
Timestamp:
May 7, 2008, 5:40:11 PM (16 years ago)
Author:
ole
Message:

Refactored options, so that the use of centroid_velocities is a separate
option to tight_slope_limiters.

Files:
7 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/config.py

    r5181 r5290  
    8888tight_slope_limiters = True
    8989
     90# Use centroid velocities to reconstruct momentum at vertices
     91# This option has a firts order flavour to it, but we still have second order
     92# reconstruction of stage and this option speeds ANUGA up
     93#
     94# This option must be used with tight_slope_limiters.
     95use_centroid_velocities = True
    9096
    9197
  • anuga_core/source/anuga/shallow_water/shallow_water_domain.py

    r5186 r5290  
    105105from anuga.config import optimised_gradient_limiter
    106106from anuga.config import use_edge_limiter
     107from anuga.config import use_centroid_velocities
     108
    107109
    108110#---------------------
     
    181183        # Limiters
    182184        self.use_edge_limiter = use_edge_limiter
    183 
    184185        self.optimised_gradient_limiter = optimised_gradient_limiter
    185                
     186        self.use_centroid_velocities = use_centroid_velocities
     187
    186188
    187189    def set_all_limiters(self, beta):
     
    741743              Ymom.vertex_values,
    742744              Elevation.vertex_values,
    743               int(domain.optimise_dry_cells))
     745              int(domain.optimise_dry_cells),
     746              int(domain.use_centroid_velocities))
    744747
    745748
     
    800803
    801804
    802     #Take bed elevation into account when water heights are small
     805    # Take bed elevation into account when water heights are small
    803806    balance_deep_and_shallow(domain)
    804807
    805     #Compute edge values by interpolation
     808    # Compute edge values by interpolation
    806809    for name in domain.conserved_quantities:
    807810        Q = domain.quantities[name]
  • anuga_core/source/anuga/shallow_water/shallow_water_ext.c

    r5238 r5290  
    144144                         
    145145  // Compute Froude number; v/sqrt(gh)
    146  
     146  // FIXME (Ole): Not currently in use
     147   
    147148  double froude_number;
    148149 
     
    549550                              double H0,
    550551                              int tight_slope_limiters,
     552                              int use_centroid_velocities,                           
    551553                              double alpha_balance) {
    552554
    553   int k, k3, i, excessive_froude_number=0;
     555  int k, k3, i; //, excessive_froude_number=0;
    554556
    555557  double dz, hmin, alpha, h_diff, hc_k;
    556558  double epsilon = 1.0e-6; // FIXME: Temporary measure
    557   double g = 9.81; // FIXME: Temporary measure
    558   double hv[3], h; // Depths at vertices
    559   double Fx, Fy; // Froude numbers
     559  double hv[3]; // Depths at vertices
     560  //double Fx, Fy; // Froude numbers
    560561  double uc, vc; // Centroid speeds
    561562
     
    701702    if (alpha < 1) {     
    702703      for (i=0; i<3; i++) {
     704     
    703705        // FIXME (Ole): Simplify when (if) hvbar gets retired       
    704706        if (beta_h > epsilon) {   
     
    709711
    710712        // Update momentum at vertices
    711 
    712 
    713         // Update momentum at vertices
    714         // Update momentum as a linear combination of
    715         // xmomc and ymomc (shallow) and momentum
    716         // from extrapolator xmomv and ymomv (deep).
    717 
    718 
    719         //xmomv[k3+i] = (1-alpha)*xmomc[k] + alpha*xmomv[k3+i];
    720         //ymomv[k3+i] = (1-alpha)*ymomc[k] + alpha*ymomv[k3+i];
    721 
    722         if (tight_slope_limiters == 1) {
    723           // FIXME(Ole): Here's what I think (as of 17 Nov 2007)
    724           // we need to do. Simple and efficient:
     713        if (use_centroid_velocities == 1) {
     714          // This is a simple, efficient and robust option
     715          // It uses first order approximation of velocities, but retains
     716          // the order used by stage.
    725717       
    726718          // Speeds at centroids
     
    738730          xmomv[k3+i] = uc*hv[i];
    739731          ymomv[k3+i] = vc*hv[i];
     732         
    740733        } else {
    741734          // Update momentum as a linear combination of
    742735          // xmomc and ymomc (shallow) and momentum
    743736          // from extrapolator xmomv and ymomv (deep).
    744           // FIXME (Ole): Is this really needed? Could we use the above
    745           // instead?
     737          // This assumes that values from xmomv and ymomv have
     738          // been established e.g. by the gradient limiter.
     739
    746740         
    747741          xmomv[k3+i] = (1-alpha)*xmomc[k] + alpha*xmomv[k3+i];
     
    749743       
    750744        }
    751       }
    752     }
    753        
    754 
    755        
    756     if (0) {  // FIXME(Ole): Disabled while testing balancing of velocities above
    757       //if (tight_slope_limiters == 1) {               
    758    
    759       // Ensure that the Froude number is kept realistic at vertices
    760       // FIXME (Ole): I think it could be used to adjust alpha down
    761       // whenever Fr is too large. Possible make sure it doesn't deviate
    762       // too much from the value at the centroid. I like this idea!
    763      
    764       // FIXME (Ole): currently only used with tights_SL
    765      
    766       // FIXME (Ole): may not be necessary now
    767 
    768       excessive_froude_number=0;   
    769       for (i=0; i<3; i++) {   
    770      
    771         // Recalculate depth at vertex i
    772         h = wv[k3+i] - zv[k3+i];
    773        
    774         Fx = compute_froude_number(xmomv[k3+i], h, g, epsilon);
    775         Fy = compute_froude_number(ymomv[k3+i], h, g, epsilon);
    776        
    777         if ( (fabs(Fx) > 100.0) || (fabs(Fy) > 100.0)) {
    778           // FIXME: use max_froude - or base it on centroid value of F
    779           // printf("Excessive Froude number detected: %f or %f\n", Fx, Fy);
    780           excessive_froude_number=1;         
    781         }
    782       }
    783      
    784       if (excessive_froude_number) {
    785    
    786         // printf("Adjusting momentum to first order.\n");
    787         // Go back to first order (alpha = 0) for this triangle
    788         for (i=0; i<3; i++) {         
    789           xmomv[k3+i] = xmomc[k];
    790           ymomv[k3+i] = ymomc[k];
    791          
    792           if (beta_h > epsilon) {         
    793             wv[k3+i] = zv[k3+i] + hvbar[k3+i];
    794           } else {
    795             wv[k3+i] = zv[k3+i] + hc_k;
    796           }
    797         }       
    798745      }
    799746    }
     
    10941041                                  double* ymom_vertex_values,
    10951042                                  double* elevation_vertex_values,
    1096                                   int optimise_dry_cells) {
     1043                                  int optimise_dry_cells,
     1044                                  int use_centroid_velocities) {
    10971045                                 
    10981046                                 
     
    13591307      // One internal neighbour and gradient is in direction of the neighbour's centroid
    13601308     
    1361       // Find the only internal neighbour
     1309      // Find the only internal neighbour (k1?)
    13621310      for (k2=k3;k2<k3+3;k2++){
    13631311        // Find internal neighbour of triangle k     
     
    13671315          break;
    13681316      }
     1317     
    13691318      if ((k2==k3+3)) {
    13701319        // If we didn't find an internal neighbour
     
    13961345      dy2=dx2*dy1;
    13971346      dx2*=dx1;
    1398      
    13991347     
    14001348     
     
    15611509  double beta_w, beta_w_dry, beta_uh, beta_uh_dry, beta_vh, beta_vh_dry;   
    15621510  double minimum_allowed_height, epsilon;
    1563   int optimise_dry_cells, number_of_elements, e;
     1511  int optimise_dry_cells, number_of_elements, e, use_centroid_velocities;
    15641512 
    15651513  // Provisional jumps from centroids to v'tices and safety factor re limiting
    15661514  // by which these jumps are limited
    15671515  // Convert Python arguments to C
    1568   if (!PyArg_ParseTuple(args, "OOOOOOOOOOOOOi",
     1516  if (!PyArg_ParseTuple(args, "OOOOOOOOOOOOOii",
    15691517                        &domain,
    15701518                        &surrogate_neighbours,
     
    15801528                        &ymom_vertex_values,
    15811529                        &elevation_vertex_values,
    1582                         &optimise_dry_cells)) {                 
     1530                        &optimise_dry_cells,
     1531                        &use_centroid_velocities)) {                   
    15831532                       
    15841533    PyErr_SetString(PyExc_RuntimeError,
     
    16251574                                   (double*) ymom_vertex_values -> data,
    16261575                                   (double*) elevation_vertex_values -> data,
    1627                                    optimise_dry_cells);
     1576                                   optimise_dry_cells,
     1577                                   use_centroid_velocities);
    16281578  if (e == -1) {
    16291579    // Use error string set inside computational routine
     
    22712221  double H0, beta_h;
    22722222
    2273   int N, tight_slope_limiters; //, err;
     2223  int N, tight_slope_limiters, use_centroid_velocities; //, err;
    22742224
    22752225  // Convert Python arguments to C
     
    23182268 
    23192269 
     2270  Tmp = PyObject_GetAttrString(domain, "use_centroid_velocities");
     2271  if (!Tmp) {
     2272    PyErr_SetString(PyExc_RuntimeError, "shallow_water_ext.c: balance_deep_and_shallow could not obtain object use_centroid_velocities from domain");
     2273    return NULL;
     2274  } 
     2275  use_centroid_velocities = PyInt_AsLong(Tmp);
     2276  Py_DECREF(Tmp);
     2277 
     2278 
     2279 
    23202280  N = wc -> dimensions[0];
    23212281  _balance_deep_and_shallow(N,
     
    23322292                            H0,
    23332293                            (int) tight_slope_limiters,
     2294                            (int) use_centroid_velocities,                         
    23342295                            alpha_balance);
    23352296
  • anuga_core/source/anuga/shallow_water/test_data_manager.py

    r5288 r5290  
    235235
    236236        self.domain.tight_slope_limiters = 0 # Backwards compatibility
     237        self.domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
     238       
    237239       
    238240        sww = get_dataobject(self.domain)       
     
    309311        #domain.tight_slope_limiters = 1
    310312        domain.tight_slope_limiters = 0 # Backwards compatibility
     313        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
     314       
    311315       
    312316        sww = get_dataobject(domain)
     
    73187322        # Look at sww file and see what happens when
    73197323        # domain.tight_slope_limiters = 1
    7320         domain.tight_slope_limiters = 0
     7324        domain.tight_slope_limiters = 0
     7325        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)       
    73217326       
    73227327        Br = Reflective_boundary(domain)
     
    74347439        # Create basic mesh (100m x 5m)
    74357440        width = 5
    7436         length = 100
     7441        length = 50
    74377442        t_end = 10
    7438         points, vertices, boundary = rectangular(length, width, 100, 5)
     7443        points, vertices, boundary = rectangular(length, width, 50, 5)
    74397444
    74407445        # Create shallow water domain
     
    75227527        # Create basic mesh (100m x 5m)
    75237528        width = 5
    7524         length = 100
     7529        length = 50
    75257530        t_end = 3
    75267531        points, vertices, boundary = rectangular(length, width,
     
    76517656        # Create basic mesh (100m x 5m)
    76527657        width = 5
    7653         length = 100
     7658        length = 50
    76547659        t_end = 1
    76557660        points, vertices, boundary = rectangular(length, width,
  • anuga_core/source/anuga/shallow_water/test_shallow_water_domain.py

    r5178 r5290  
    14371437        domain.H0 = 0 # Backwards compatibility (6/2/7)
    14381438        domain.beta_h = 0.2 # Backwards compatibility (14/2/7)
     1439        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
     1440       
    14391441
    14401442        #-----------------------------------------------------------------
     
    26572659        # FIXME (Ole): Need tests where this is commented out
    26582660        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                 
     2661        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
    26592662       
    26602663               
     
    38693872        # FIXME (Ole): Need tests where these two are commented out
    38703873        domain.H0 = 0        # Backwards compatibility (6/2/7)       
    3871         domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)         
     3874        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)
     3875        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)                 
    38723876
    38733877        #Bed-slope and friction
     
    39413945       
    39423946        # FIXME (Ole): Need tests where this is commented out
    3943         domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)         
     3947        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)
     3948        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)     
    39443949       
    39453950        #Bed-slope and friction at vertices (and interpolated elsewhere)
     
    40384043        # FIXME (Ole): Need tests where this is commented out
    40394044        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                 
    4040         domain.H0 = 0 # Backwards compatibility (6/2/7)       
     4045        domain.H0 = 0 # Backwards compatibility (6/2/7)
     4046        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
     4047       
    40414048
    40424049        #Bed-slope and friction at vertices (and interpolated elsewhere)
     
    41414148        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                 
    41424149        domain.H0 = 0 # Backwards compatibility (6/2/7)
     4150        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
     4151       
    41434152
    41444153        #Bed-slope and friction at vertices (and interpolated elsewhere)
     
    42404249        domain.H0 = 0        # Backwards compatibility (6/2/7)       
    42414250        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                 
     4251        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)
    42424252       
    42434253               
     
    45044514        # FIXME (Ole): Need tests where these two are commented out
    45054515        domain.H0 = 0        # Backwards compatibility (6/2/7)       
    4506         domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)                         
     4516        domain.tight_slope_limiters = 0 # Backwards compatibility (14/4/7)
     4517        domain.use_centroid_velocities = 0 # Backwards compatibility (7/5/8)       
    45074518       
    45084519
     
    56745685    suite = unittest.makeSuite(Test_Shallow_Water,'test')
    56755686
    5676     #suite = unittest.makeSuite(Test_Shallow_Water,'test_balance_deep_and_shallow_Froude')
    5677    
     5687    #suite = unittest.makeSuite(Test_Shallow_Water,'test_bedslope_problem_first_order_moresteps')   
    56785688    #suite = unittest.makeSuite(Test_Shallow_Water,'test_fitting_using_shallow_water_domain')   
    56795689    #suite = unittest.makeSuite(Test_Shallow_Water,'test_tight_slope_limiters')
  • anuga_validation/automated_validation_tests/fitting/validate_benchmark_fit.py

    r4919 r5290  
    77memory is measured here. What is actually happening is the mesh
    88generation is using less memory, in version 4917.  In version 4897
    9 mesh gen uses and frees up more memory, so phython asks for less extra
     9mesh gen uses and frees up more memory, so python asks for less extra
    1010memory in the fitting stage.  So overall version 4917 is an
    1111improvement.
     
    6565            self.assert_(time<time_standard*1.2)
    6666            mem_standard = 302404
    67             # 
     67            #
    6868            self.assert_(mem<mem_standard*1.2)
    6969           
     
    9090            time_standard = 14.2
    9191            self.assert_(time<time_standard*1.2)
    92             mem_standard = 21008 # this will need to be updated, see top notes
    93             self.assert_(mem<mem_standard*1.2)           
     92            mem_standard = 30000 # Updated by Ole 20080507
     93            msg = 'Memory used was %d, mem_standard is %d' %(mem, mem_standard)           
     94            assert mem < mem_standard*1.2, msg           
    9495
    9596        elif host.find('pc-31569') == 0: # DSG's PC
  • anuga_validation/automated_validation_tests/validate_all.py

    r4842 r5290  
    22"""
    33
    4 import os
     4import os, time
    55
    66
     
    2828
    2929
     30print
     31print '----------------------------------------------------------'
     32print 'Running all validation tests - this may take several hours'
     33print '----------------------------------------------------------'
    3034
    3135for path, filename in validation_dirs_and_files:
     
    3337
    3438    os.chdir(path)
    35     s = 'python %s' %( filename)
     39    s = 'python %s' %(filename)
     40    print s
     41    t0 = time.time()
    3642    os.system(s)
    37     # print 'current dir', os.getcwd()
     43    print 'That took %.2f seconds' %(time.time()-t0)
    3844   
    3945    # Back to parent directory
Note: See TracChangeset for help on using the changeset viewer.