Ignore:
Timestamp:
Sep 12, 2007, 1:05:19 PM (17 years ago)
Author:
ole
Message:

More small optimisations

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/shallow_water/shallow_water_ext.c

    r4730 r4731  
    280280  // Flux computation
    281281  denom = s_max-s_min;
    282   if (denom == 0.0) {  // FIXME (Ole): Try using epsilon here
     282  if (denom < epsilon) { // FIXME (Ole): Try using H0 here
    283283    for (i=0; i<3; i++) edgeflux[i] = 0.0;
    284284    *max_speed = 0.0;
     
    623623  double u, v, reduced_speed;
    624624
    625   // Protect against initesimal and negative heights
    626   for (k=0; k<N; k++) {
    627     hc = wc[k] - zc[k];
    628 
    629     if (hc < minimum_allowed_height) {
     625
     626  // Protect against initesimal and negative heights 
     627  if (maximum_allowed_speed < epsilon) {
     628    for (k=0; k<N; k++) {
     629      hc = wc[k] - zc[k];
     630
     631      if (hc < minimum_allowed_height) {
    630632       
    631       //Old code: Set momentum to zero and ensure h is non negative
    632       //xmomc[k] = 0.0;
    633       //ymomc[k] = 0.0;
    634       //if (hc <= 0.0) wc[k] = zc[k];
    635 
    636 
    637       //New code: Adjust momentum to guarantee speeds are physical
    638       //          ensure h is non negative
    639       //FIXME (Ole): This is only implemented in this C extension and
    640       //             has no Python equivalent
    641            
    642       if (hc <= 0.0) {
    643         wc[k] = zc[k];
     633        // Set momentum to zero and ensure h is non negative
    644634        xmomc[k] = 0.0;
    645635        ymomc[k] = 0.0;
    646       } else {
    647         //Reduce excessive speeds derived from division by small hc
    648         //FIXME (Ole): This may be unnecessary with new slope limiters
    649         //in effect.
    650        
    651         u = xmomc[k]/hc;
    652         if (fabs(u) > maximum_allowed_speed) {
    653           reduced_speed = maximum_allowed_speed * u/fabs(u);
    654           //printf("Speed (u) has been reduced from %.3f to %.3f\n",
    655           //     u, reduced_speed);
    656           xmomc[k] = reduced_speed * hc;
    657         }
    658 
    659         v = ymomc[k]/hc;
    660         if (fabs(v) > maximum_allowed_speed) {
    661           reduced_speed = maximum_allowed_speed * v/fabs(v);
    662           //printf("Speed (v) has been reduced from %.3f to %.3f\n",
    663           //     v, reduced_speed);
    664           ymomc[k] = reduced_speed * hc;
    665         }
     636        if (hc <= 0.0) wc[k] = zc[k];
     637      }
     638    }
     639  } else {
     640   
     641    // Protect against initesimal and negative heights
     642    for (k=0; k<N; k++) {
     643      hc = wc[k] - zc[k];
     644
     645      if (hc < minimum_allowed_height) {
     646
     647        //New code: Adjust momentum to guarantee speeds are physical
     648        //          ensure h is non negative
     649        //FIXME (Ole): This is only implemented in this C extension and
     650        //             has no Python equivalent
     651             
     652        if (hc <= 0.0) {
     653                wc[k] = zc[k];
     654        xmomc[k] = 0.0;
     655        ymomc[k] = 0.0;
     656        } else {
     657          //Reduce excessive speeds derived from division by small hc
     658        //FIXME (Ole): This may be unnecessary with new slope limiters
     659        //in effect.
     660         
     661          u = xmomc[k]/hc;
     662          if (fabs(u) > maximum_allowed_speed) {
     663            reduced_speed = maximum_allowed_speed * u/fabs(u);
     664            //printf("Speed (u) has been reduced from %.3f to %.3f\n",
     665            //   u, reduced_speed);
     666            xmomc[k] = reduced_speed * hc;
     667          }
     668
     669          v = ymomc[k]/hc;
     670          if (fabs(v) > maximum_allowed_speed) {
     671            reduced_speed = maximum_allowed_speed * v/fabs(v);
     672            //printf("Speed (v) has been reduced from %.3f to %.3f\n",
     673            //   v, reduced_speed);
     674            ymomc[k] = reduced_speed * hc;
     675          }
     676        }
    666677      }
    667678    }
     
    669680  return 0;
    670681}
    671 
    672682
    673683
     
    680690                              double cw) {
    681691
    682   //Assign windfield values to momentum updates
     692  // Assign windfield values to momentum updates
    683693
    684694  int k;
Note: See TracChangeset for help on using the changeset viewer.