Hi everyone,<br> I believe the routines in GRHydro_SlopeLimiter.F90 are all either mislabelled or simply the wrong formulae. As a quick reminder, for the upwind slope ("dupw") and downwind slope ("dloc"), we expect to get the following results:<br>
<br>If dupw * dloc <= 0, delta=0 [slope limiters go flat near extrema, by design]<br><br>else<br><br>if (minmod)<br><br>delta = sign(dupw) * min( abs(dupw), abs(dloc))<br><br>else if (MC)<br><br>delta = sign(dupw) * min( 2*abs(dupw), 2*abs(dloc), abs(dupw+dloc)/2 )<br>
<br>
else if (superbee)<br><br>delta = sign(dupw) * max[ min(2*abs(dupw), dloc) , min(2*abs(dloc),dupw) ]<br><br>On the rather nice wikipedia page, <a href="http://en.wikipedia.org/wiki/Flux_limiter">http://en.wikipedia.org/wiki/Flux_limiter</a><br>
<br>they take variables r = dupw/dloc and phi = delta / dloc, but the results are the same.<br><br>Attached are a spreadsheet where I tried to reproduce what the current routine actually calculates for various inputs, along with the recommended code fix for GRHydro_TVDReconstruct.F90.<br>
<br>In summmary: the current "MC2" and "minmod" options work like they should. "MC1" is actually minmod with a floor for very small values. minmod2, minmod3, and superbee are all incorrect expressions.<br>
<br>Cheers,<br>Josh<br>