Markus -
I think any bias-free formula is OK as long as the threshold for
pressure is reasonable. This is, of course, problem dependent.
I am generally against using floor-values in the code as they are (1)
problem-dependent, and (2) mask some real problems. One can almost
always find a reasonable (unbiased, right order of magnitued) estimate
for the limiting value based on physical or numerical
arguments. 1.e-10 or 1.e+10 does not really mean much when taken out
of context. Hardwiring -10/+10 and -2/+2, as you suggested, may not be
the best solution. Tying those constant to epsiln may work better.
Best wishes -
Tomek
-- On Wed, Jul 16, 2003 at 01:33:44PM +0100, Markus Gross wrote: > Hi again! > > Realizing that my fix was not much better than the original for any non > constructed case, i.e. when rounding error may be present, may I suggest the > following: > > !this is the original: > !!$ if (dp(i) .LT. 0.e0) then > !!$ scrch2(i) = scrch3(i+1) > !!$ else > !!$ scrch2(i) = scrch3(i-1) > !!$ endif > > !my "new" suggestion: > if (dp(i) .LT. -10.e0) THEN > scrch2(i) = scrch3(i+1) > ELSE > if (dp(i) .GT. 10.e0) THEN > scrch2(i) = scrch3(i-1) > ELSE > weight1 = 0.5d0*dabs(derf(dp(i)+2.d0)-1.d0) > weight2 = 0.5d0*dabs(derf(dp(i)-2.d0)+1.d0) > scrch2(i) = scrch3(i+1) * weight1 +& > scrch3(i-1) * weight2 +& > scrch3(i ) * (1.d0-(weight1+weight2)) > END IF > END IF > > This of course assumes that a -2.d0< dp(i) < 2.d0 is effectively a dp(i)=0 > which for my case is true. Modifications for other "small" values is of course > straight forward. > > Again, what do you think? A weibull may of course be faster than the erf... > > Regards, > > Markus. > >--
This archive was generated by hypermail 2b30 : Wed Jul 16 2003 - 09:52:04 CDT