Ignore:
Timestamp:
Nov 27, 2009, 9:31:34 AM (13 years ago)
Author:
steve
Message:

Committing a version of shallow_water_balanced which passes it unit test
using a version of edge limiting which doesn't limit boundary edges. THis
is useful to allow linear functions to be reconstructed.

Had to play with the transmissive BC and use centroid values which is
set via domain set_centroid_transmissive_bc

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/abstract_2d_finite_volumes/quantity_ext.c

    r7519 r7573  
    339339  int i, k, k2, k3, k6;
    340340  long n;
    341   double qmin, qmax, qn, qc;
     341  double qmin, qmax, qn, qc, sign;
    342342  double dq, dqa[3], phi, r;
    343343 
     
    361361    }
    362362   
     363    sign = 0.0;
     364    if (qmin > 0) {
     365      sign = 1.0;
     366    } else if (qmax < 0) {
     367      sign = -1.0;
     368    }
     369
    363370    phi = 1.0;
    364     for (i=0; i<3; i++) {   
    365       r = 1.0;
    366      
    367       dq = edge_values[k3+i] - qc;     //Delta between edge and centroid values
    368       dqa[i] = dq;                      //Save dq for use in updating vertex values
    369      
    370       if (dq > 0.0) r = (qmax - qc)/dq;
    371       if (dq < 0.0) r = (qmin - qc)/dq;     
    372      
    373      
    374       phi = min( min(r*beta, 1.0), phi);   
     371    for (i=0; i<3; i++) { 
     372      dq = edge_values[k3+i] - qc;      //Delta between edge and centroid values
     373      dqa[i] = dq;                      //Save dq for use in updating vertex values 
     374     
     375      // FIXME SR 20091125 This caused problems in shallow_water_balanced
     376      // commenting out problem
     377      // Just limit non boundary edges so that we can reconstruct a linear function
     378      if (neighbours[k3+i] >= 0) {
     379        r = 1.0;
     380     
     381        if (dq > 0.0) r = (qmax - qc)/dq;
     382        if (dq < 0.0) r = (qmin - qc)/dq;     
     383           
     384        phi = min( min(r*beta, 1.0), phi);
     385        }
     386
     387      if (neighbours[k3+i] < 0) {
     388        r = 1.0;
     389     
     390        if (dq > 0.0 && sign == -1.0 ) r = (0.0 - qc)/dq;
     391        if (dq < 0.0 && sign ==  1.0 ) r = (0.0 - qc)/dq;     
     392           
     393        phi = min( min(r*beta, 1.0), phi);
     394        }
     395   
    375396    }
    376397   
Note: See TracChangeset for help on using the changeset viewer.