floating point precision: is true double possible?

Forum 2005-2010 (read only) Software Development
#1

I'm writing a short program to follow the position of the sun, and the astronomical calculations require more significant digits than a single precision floating point can represent. On the arduino, as I understand, there is no real "double" type -- it exists, but has the same precision as a standard "float" (maybe 8 digits of precision).

Does an arduino library exist that allows higher floating point precision -- either arbitrary precision or a true double? I'm aware of the severe ram restriction and only need a few of these variables -- no huge matrices -- but the inability to do any floating-point math with precision higher than 8 digits is a very onerous limitation.

Thanks, -Jon

#2

How about fixed-point maths on "long long" types? If it's for astronomical calculations, it doesn't need to be fast.

#3

What do you plan to do with the results?

-j

#4

Yes, using long long integers might work in some parts of the code, but there are operations that need to be carried out as floating point (taking exponents, trig functions, etc.) or would involve significant code-contortions. If it's possible to just load a separate library and get bigger floats, I would prefer that approach to writing obfuscated mathematical code (i.e., trying to take the cosine of an integer and square it, etc.).

It's true that this doesn't need to be fast -- anything in the sky is going to make only one full revolution in 86400 seconds (1 day). I'm also reasonably confident that this is well within the capabilities of the Arduino, since -- although the required precision is high -- there are relatively few variables at any given point and the speed requirement is minimal.

Regarding what will be done with the results: that depends on whether I end up even getting good results -- this is my first adventure with Arduino. Eventually I'd like to drive a heliostat (servo- or stepper-controlled tracking mount), but that hasn't been built yet nor have parts even been purchased.

-Jon

#5

Eventually I'd like to drive a heliostat (servo- or stepper-controlled tracking mount)

This suggests to me that the precision you actually need will be determined by the resolution of your tracking hardware.

My knee-jerk reaction is that you don't really need double precision. But then, I don't know the formulas/values involved, either.

-j

#6

This suggests to me that the precision you actually need will be determined by the resolution of your tracking hardware.

My knee-jerk reaction is that you don't really need double precision. But then, I don't know the formulas/values involved, either.

Yes, the need for precision lies in the formulas. The equations are designed to be accurate for many years (otherwise the program would need to be continually re-written & updated), so the math contains terms that vary over many different timescales, from seconds to centuries.

As I mentioned originally, I can simplify some of this, combine terms, make approximations, etc., to make it work -- but this just seems like such an inelegant solution, compared to simply using an additional 4 bytes of ram for a proper double-float. I'm quite surprised that this, of all things, would be a sticking-point.

#7

I'm still learning the avr gcc C language that the Arduino uses, but can't one just include and link the double math library that is available?

http://www.nongnu.org/avr-libc/user-manual/group__avr__math.html

Lefty

#8

but this just seems like such an inelegant solution, compared to simply using an additional 4 bytes of ram for a proper double-float. I'm quite surprised that this, of all things, would be a sticking-point.

Welcome to the world of 8 bit microcontrollers. They're not designed for applications like yours.

And it isn't the 4 bytes of RAM so much as it is the silicon real estate to make the 32 bit paths and structures (adders, etc), and keep the cost of the processor at US$3.50.

I think I'd use integer math on a long long data type (64 bit int) and use seconds as the unit. A long long would give you about 2 billion centuries at 1s resolution. (This is assuming there isn't significant fractional data at 1s resolution).

-j

#9

can't one just include and link the double math library that is available?

Isn't it more complicated than that? How do you persuade the Arduino that sizeof (double) != sizeof (float)?

#10

I know this is possible because the TAPR orginzation did a controller that use a PIC chip in the past. It is not any more powerful than the Arduino.

One suggestion is the floating point external processor that can be hooked up

http://www.sparkfun.com/commerce/product_info.php?products_id=8129

You can talk to it via serial, I2C or SPI

This may or may not have the precision, but it does have the speed.

#11

the math contains terms that vary over many different timescales, from seconds to centuries.

I think I'm going to go back to my original thought: I really don't think you need double precision to do this. Floating point can handle very large numbers and very small numbers, and perform relatively accurate calculations even when very large and very small numbers are involved.

-j

#12

I really don't think you need double precision to do this.

The following is a very common and widely used formula for calculating distance between two sets of GPS coordinates. 

 dist = acos(sin(lat1)*sin(lat2)+cos(lat1)*cos(lat2)*cos(lon2-lon1));

Unfortunately, this formula can not be relied upon with only 32-bit floats. The challenge is that when you combine small and large scale values in calculations, you loose precision. Triangulation as above relates narrow distances to global coordinates several magnitudes larger and this is the challenge. Whether your required output is a single bit (true/false, left/right) or a high precision number really doesn't matter. You still need the intermediate calculations to be precise in order to determine correct output.

The following is a rewrite that work well with 32-bit floats.

  // alternative formulae (haversine)
  const float two=2.0;
  dist = two*asin(sqrt(square(sin((lat1-lat2)/two)) + 
         cos(lat1)*cos(lat2)*square(sin((lon2-lon1)/two))));

Debugging precision issues is difficult because an algorithm may work well for one set of coordinates, but fail miserably for another set. Realizing this and figuring out how to rewrite the formulas (if at all possible) to avoid rounding/precision issues is anything but trivial.

I'm no expert whatsoever in astronomical calculations, but I do respect that our planet earth is pretty small in comparison to our solar system at large. Even with double precision it is a formidable challenge to come up with formulas that avoid computational pitfalls due to rounding and loss of precision. GPS calculations in this context is like first grade math in comparison so I have no reason to doubt that double precision may indeed be required.

#13

Sun position code here: http://www.mowcius.co.uk/suntrackercode.html

Accuracy to 0.5deg at least. If you re-implemented the seconds code (which I removed as it's really not necessary) then you could probably get it slightly faster.

Any questions, just ask.

Mowcius

#14

Accuracy to 0.5deg at least

#15

I'm interetsed in the revised formula BenF posted - I'm using this & multiplying by Earth's diameter in m (6371000), though am not getting a valid result (getting 50~300m deviation when my GPS module is stationary, though I know the deviation is actually 5-10m). Is there any other conversion required?

#16

Earth's diameter in m (6371000)

That might be your issue

Earth's Diameter at the Equator: 7,926.28 miles (12,756.1 km)

Earth's Diameter at the Poles: 7,899.80 miles (12,713.5 km)

Mowcius

#17

Is there any other conversion required?

Arguments for latitude and longitude need to be passed as radians (not degrees). Could this be what you have missed?

The return value (also in radians) must then be multiplied with the earth radius. For this I use the the IUGG mean radius defined as 6371009 meters.

http://en.wikipedia.org/wiki/Earth_radius

#18

though am not getting a valid result (getting 50~300m

Have you checked your datum? For instance, the UK OS datum produces errors of about this magnitude, compared to the WGS datum used by GPS.