Square root of a negative rral value

Square root of a negative rral value

Post by Terenc » Fri, 07 Nov 2008 08:33:03


I mentioned in another post that I was converting a mainframe Fortran
IV program to F77 (and F90).

I found that running the converted program with the original check
data set that went with this, crashes on trying to take the square
root of a negative variable (to get the standard deviation as a
measure for selecting later optimum adjustments).

I was utterly lost to discover why all the many computational results
checked exactly to some 6 places, EXCEPT for one standard error and
its significance (both unused later) and the program crashing on one
square root operation before terminating correctly with a final table.

Obviously the mainframe did not crash and the PC did on the same data!
After some days (!!) I inserted before the SQRT line equivalent to
this
std=SQRT(var)

with the line
IF (var.LT.0.d0) var = -var

Problem solved. (And I hope SQRT(0.0) gives a zero result!)

I deduce that mainframe runtime replaces the subject variable of a
square root operation with an absolute value and continues, wheras the
PC runtime informs the user of the error of trying to take the root of
a negative value and stops.

Aaagh!
 
 
 

Square root of a negative rral value

Post by Glen Herrm » Fri, 07 Nov 2008 08:59:38


I have seen this many times. Note that six places isn't enough,
if the result is negative in the seventh place.

When the standard deviation is zero, or very close to zero
it is possible for rounding in the variance to give a
negative result. Probably better to set it to zero
instead of the sqrt() of the absolute value, but it will be
very small in any case.

See http://www.yqcomputer.com/ #Definition_and_calculation

The popular standard deviation algorithm uses an identity that is
true mathematically, but not always numerically. That is, the
mean of the squares minus the square of the mean should be positive
but with some types of rounding, it isn't.

-- glen

 
 
 

Square root of a negative rral value

Post by Gordon San » Fri, 07 Nov 2008 09:54:21

On 2008-11-05 19:33:03 -0400, Terence < XXXX@XXXXX.COM > said:


Did you first try printing out the offending value? In floating point
so see small magnitudes. I seriously doubt that any mainframe or PC
would silently take the absolute value of a square root. Much more
likely that the rounding effects have changed.

Calculating a standard deviation that is small is an absolutely classical
demonstration of the computational instability of the moment based formulas
that filled statistics text books from the days of hand calculators. Much
better to do a two pass or provisional means algorithm. Trivial test examples
that are exact fits tend to demonstrate these problems.
 
 
 

Square root of a negative rral value

Post by Glen Herrm » Fri, 07 Nov 2008 10:30:08


(snip)


x87 might just leave the original value on the stack.

The design is that you would do the checking in software before
using the x87 floating point instructions, but when generated inline
the checking is often not done. I have seen it for sin(1e100),
I wouldn't be surprised for sqrt(-1.23).

But yes, different rounding is more likely. Note, for example
that S/360 (and successors) truncate the quotient on floating
point divide, where many others round. That could go either
way in the standard deviation case, though.

-- glen
 
 
 

Square root of a negative rral value

Post by Terenc » Fri, 07 Nov 2008 13:37:50

Well, yes, there was a 100 by100 matrix inversion, and a triple matrix
multiplication, and values wandered from about 0.4 up to 10**11, so it
might not be considered to be a stable calculation.
But the input is real-world data.

I was limited to REAL*8.
I accept the good point that 0.0 is a better default than the SQRT of
the absolute value, since the sample standard deviations SHOULD be
positive.

Anyway it works now, but I wonder just how many more pre-IBM360
algorithms have real-time problmes?? My best guide to "practical
methods" is the UK 1961 National Physica Laboratory White Paper on
Modern Computing Methods - often meaning hand-cranked calculators!
 
 
 

Square root of a negative rral value

Post by e p chandl » Fri, 07 Nov 2008 17:15:35


> std=SQRT(var) >> >> with the line >> IF (var.LT.0.d0) var = -var> >> > Problem solved. (And I hope SQRT(0.0) gives a zero result!)> >> > I deduce that mainframe runtime replaces the subject variable of > > square root operation with an absolute value and continues, wheras the> > PC runtime informs the user of the error of trying to take the root of> > a negative value and stops.> >> > Aaagh!

IIRC Forran I and II did this rather than crash. Likewise, early IBM-
PC Fortran was designed so that dividing by zero was not a fatal
error, instead the computation would march on.

-- E
 
 
 

Square root of a negative rral value

Post by e p chandl » Fri, 07 Nov 2008 17:27:35


A recent web search took me to a page from Cambridge about their EDSAC
1 and 2 computers. Amazing what they accomplished with vacuum tubes
and mercury delay lines. I bet that most people in the USA know little
about UK computing history (or jet aviation, etc.).

Interesting to see that programming was done in machine or assembly
code, with paper tape or mag tape.

- E
 
 
 

Square root of a negative rral value

Post by Terenc » Fri, 07 Nov 2008 18:59:49


Oh Yes!
I worked for a time on guidance for a large UK missile, at a London
office.
We used the Mercury computer for trajectory calculations.
First thing, was to walk into the computer with a box of double
triodes valves, unplug and replace any of the thousands with no
filaments visibly lit (while in standby mode).
Then exit, go to run mode, feed the next target paper tape and cross
fingers for a 20 minute no-fail run, collect the output tape, feed it
into the teleprinter and go off with the printout.
The "screen" was just an oscilloscope for the engineers.

A few years later later I was installing IBM 1620 (addition and
multiply was table look-up, division was a subtraction routine) paper
tape machines for those who couldn't afford a 704 (and what went with
it), and had come to the limit of plug-board capacity relay
calculators (though I did a lot of plug-board computations on water
control strategy for the Thames Water Board, using 1930-1950 data. Now
we do it with the Out-of-Kilter algorithm.
Remember Remingtons?.
 
 
 

Square root of a negative rral value

Post by user » Fri, 07 Nov 2008 21:13:04


I'd like to ask the same as Gordon, what happens when you try to insert a print
statement and look at the offending value. Is it a negative number that is
almost zero within roundoff error, or is it some totally ridiculous negative
number ?

You might even find that inserting a write statement masks the problem, in which
case something more serious could going wrong.
 
 
 

Square root of a negative rral value

Post by Gordon San » Fri, 07 Nov 2008 22:15:12

On 2008-11-06 00:37:50 -0400, Terence < XXXX@XXXXX.COM > said:


Some folks would take the view that the negative square root is
a hint that there are problems to be solved. Not shoved under the rug
or papered over or forgotten or otherwise ignored.


Why not take a look at Ake Bjorck's SIAM book on Least Squares. It
is much much much more modern than NPL's book on Modern Computing
Methods published by HMSO in the 1960s. A lot of its stuff is based
on fixed point with double length accumulation, the fi_2 mode.
 
 
 

Square root of a negative rral value

Post by Arjen Mark » Fri, 07 Nov 2008 22:40:48


> std=SQRT(var) >> >> with the line >> IF (var.LT.0.d0) var = -var> >> > Problem solved. (And I hope SQRT(0.0) gives a zero result!)> >> > I deduce that mainframe runtime replaces the subject variable of > > square root operation with an absolute value and continues, wheras the> > PC runtime informs the user of the error of trying to take the root of> > a negative value and stops.> >> > Aaagh!

As the others have said, most probably a rounding problem.
I came across the very same problem when a user tried to
determine the standard deviation of a set of equal numbers.

The "smart" one-pass algorithm I had used came up with a tiny
negative variance. And so the program crashed in the next
step.

As for this not happening in earlier computers: they used
a wide variety of implementations to do floating-point
computations. And so did early PC compilers (remember
the ones that did not have a mathematical coprocessor?).

I think I prefer the current, much more homogeneous
situation. But I can imagine your frustration.

Regards,

Arjen
 
 
 

Square root of a negative rral value

Post by nospa » Sat, 08 Nov 2008 02:29:23


Square roots of negative numbers can also happen for other reasons. The
usual formula for computing airspeed from instrumentation data involves
subtracting the static pressure from the total pressure and then taking
a square root. I've seen plenty of naive versions of such code that
regularly blow up when the aircraft is sitting on the ramp doing
preflight tests. The coders don't check for negative values because the
static pressure "can't" be greater than the total pressure (that and
because lots of coders are just plain careless about checking much of
anything). While the "true" static pressure can't be greater than the
total pressure, the measured data are another matter. Measured data are
imperfect. If the airplane is sitting on the ramp, it is quite easy for
the measured static pressure to be greater than the measured total
pressure, be it from sensor noise, imperfect calibration, or what.

I might add that the best answer in this case is zero; not the square
root of the absolute value. But one can't conclude that for all cases of
square roots of negatives. It depends on the particular application. For
computing airspeed, zero turns out to be the best answer.

--
Richard Maine | Good judgment comes from experience;
email: last name at domain . net | experience comes from bad judgment.
domain: summertriangle | -- Mark Twain
 
 
 

Square root of a negative rral value

Post by Gib Bogl » Sat, 08 Nov 2008 04:59:03


What happens if the plane is going backwards?
 
 
 

Square root of a negative rral value

Post by Paul van D » Sat, 08 Nov 2008 05:10:17


The pilot ejects?
 
 
 

Square root of a negative rral value

Post by nospa » Sat, 08 Nov 2008 05:13:15


Well, one might comment on that from two perspectives.

1. The definition. Airspeed is a magnitude. It doesn't have direction.
If you are going backwards, that's still a positive airspeed (with an
angle of attack of 180 degrees - or you could call it an angle of
sidelsliip of 180 degrees, depending on how you normalize those things).

2. The practical. It can, in fact, quite easily happen on the ground
that the airplane is going backwards (wind relative). All it takes is a
little tailwind with the airplane sitting still. But the usual
arrangement of sensors won't pick that up. The total pressure sensor is
only good for modest angles of attack and sideslip. There exist sensors
that can handle more extreme cases, but that's not the normal
arrangement. (Well, maybe on helicopters, but I don't "do" them.) With
normal sensors, you pretty much won't pick up anything useful when the
plane is going backwards. It will show as total pressure being equal to
static (i.e. zero airspeed) plus or minus some noise... and the noise
will get you into that negative square root territory.

Anyway, back to the somewhat vaguely Fortran-related part, the square
root of the absolute value would not be the right answer. Zero would be
about as well as you could do with normal instrumentats. If you happened
to have instruments suitable to the case, then you wouldn't have a
negative value.

--
Richard Maine | Good judgment comes from experience;
email: last name at domain . net | experience comes from bad judgment.
domain: summertriangle | -- Mark Twain