R version 2.9.0 (2009-04-17)
Copyright (C) 2009 The R Foundation for Statistical Computing
ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

>  # weird.r             September 2009
>  # 
>  # some weirdness
>  #
>  # why not test equality with floating point numbers
>  # 
>    seven <- .7
>    two <- .2
>    one <- .1
>  # 
>  # surprised?
>    1 == (seven+two)+one
[1] FALSE
>    1 == seven+(two+one)
[1] TRUE
>    1 == (seven+one)+two
[1] TRUE
>  # 
>  # now look at representation
>    sprintf("%a",c(seven,two,one))       # three numbers
[1] "0x1.6666666666666p-1" "0x1.999999999999ap-3" "0x1.999999999999ap-4"
>    sprintf("%a",c(seven+two,one))       # partial sum
[1] "0x1.cccccccccccccp-1" "0x1.999999999999ap-4"
>    sprintf("%a",(seven+two)+one)        # why it's not one
[1] "0x1.fffffffffffffp-1"
>  #
>  # another from Kyle White
>    allbut970 <- c((959:969),(971:1023))
>    sum(2^allbut970)
[1] Inf
>    length(allbut970)
[1] 64
>  # reference -- this should be biggest number
>    biggest <- sum(2^(971:1023))
>    sprintf("%a",biggest)
[1] "0x1.fffffffffffffp+1023"
>  #
>    sum(2^allbut970[-1])
[1] 1.797693e+308
>    sum(2^allbut970[-1])+2^allbut970[1]
[1] 1.797693e+308
>  #
>  # ok, that's weird enough, what are they internally?
>    sprintf("%a",sum(2^allbut970[-1]))
[1] "0x1.fffffffffffffp+1023"
>    sprintf("%a",sum(2^allbut970[-1])+2^allbut970[1])
[1] "0x1.fffffffffffffp+1023"
>  #
>  # now look at internal view of cumulative sums
>    c1 <- sprintf("%a",cumsum(2^allbut970[1:63]))
>    c2 <- sprintf("%a",cumsum(2^allbut970[2:64]))
>    cbind(c1,c2)
      c1                        c2                       
 [1,] "0x1p+959"                "0x1p+960"               
 [2,] "0x1.8p+960"              "0x1.8p+961"             
 [3,] "0x1.cp+961"              "0x1.cp+962"             
 [4,] "0x1.ep+962"              "0x1.ep+963"             
 [5,] "0x1.fp+963"              "0x1.fp+964"             
 [6,] "0x1.f8p+964"             "0x1.f8p+965"            
 [7,] "0x1.fcp+965"             "0x1.fcp+966"            
 [8,] "0x1.fep+966"             "0x1.fep+967"            
 [9,] "0x1.ffp+967"             "0x1.ffp+968"            
[10,] "0x1.ff8p+968"            "0x1.ff8p+969"           
[11,] "0x1.ffcp+969"            "0x1.7fep+971"           
[12,] "0x1.7ffp+971"            "0x1.bffp+972"           
[13,] "0x1.bff8p+972"           "0x1.dff8p+973"          
[14,] "0x1.dffcp+973"           "0x1.effcp+974"          
[15,] "0x1.effep+974"           "0x1.f7fep+975"          
[16,] "0x1.f7ffp+975"           "0x1.fbffp+976"          
[17,] "0x1.fbff8p+976"          "0x1.fdff8p+977"         
[18,] "0x1.fdffcp+977"          "0x1.feffcp+978"         
[19,] "0x1.feffep+978"          "0x1.ff7fep+979"         
[20,] "0x1.ff7ffp+979"          "0x1.ffbffp+980"         
[21,] "0x1.ffbff8p+980"         "0x1.ffdff8p+981"        
[22,] "0x1.ffdffcp+981"         "0x1.ffeffcp+982"        
[23,] "0x1.ffeffep+982"         "0x1.fff7fep+983"        
[24,] "0x1.fff7ffp+983"         "0x1.fffbffp+984"        
[25,] "0x1.fffbff8p+984"        "0x1.fffdff8p+985"       
[26,] "0x1.fffdffcp+985"        "0x1.fffeffcp+986"       
[27,] "0x1.fffeffep+986"        "0x1.ffff7fep+987"       
[28,] "0x1.ffff7ffp+987"        "0x1.ffffbffp+988"       
[29,] "0x1.ffffbff8p+988"       "0x1.ffffdff8p+989"      
[30,] "0x1.ffffdffcp+989"       "0x1.ffffeffcp+990"      
[31,] "0x1.ffffeffep+990"       "0x1.fffff7fep+991"      
[32,] "0x1.fffff7ffp+991"       "0x1.fffffbffp+992"      
[33,] "0x1.fffffbff8p+992"      "0x1.fffffdff8p+993"     
[34,] "0x1.fffffdffcp+993"      "0x1.fffffeffcp+994"     
[35,] "0x1.fffffeffep+994"      "0x1.ffffff7fep+995"     
[36,] "0x1.ffffff7ffp+995"      "0x1.ffffffbffp+996"     
[37,] "0x1.ffffffbff8p+996"     "0x1.ffffffdff8p+997"    
[38,] "0x1.ffffffdffcp+997"     "0x1.ffffffeffcp+998"    
[39,] "0x1.ffffffeffep+998"     "0x1.fffffff7fep+999"    
[40,] "0x1.fffffff7ffp+999"     "0x1.fffffffbffp+1000"   
[41,] "0x1.fffffffbff8p+1000"   "0x1.fffffffdff8p+1001"  
[42,] "0x1.fffffffdffcp+1001"   "0x1.fffffffeffcp+1002"  
[43,] "0x1.fffffffeffep+1002"   "0x1.ffffffff7fep+1003"  
[44,] "0x1.ffffffff7ffp+1003"   "0x1.ffffffffbffp+1004"  
[45,] "0x1.ffffffffbff8p+1004"  "0x1.ffffffffdff8p+1005" 
[46,] "0x1.ffffffffdffcp+1005"  "0x1.ffffffffeffcp+1006" 
[47,] "0x1.ffffffffeffep+1006"  "0x1.fffffffff7fep+1007" 
[48,] "0x1.fffffffff7ffp+1007"  "0x1.fffffffffbffp+1008" 
[49,] "0x1.fffffffffbff8p+1008" "0x1.fffffffffdff8p+1009"
[50,] "0x1.fffffffffdffcp+1009" "0x1.fffffffffeffcp+1010"
[51,] "0x1.fffffffffeffep+1010" "0x1.ffffffffff7fep+1011"
[52,] "0x1.ffffffffff7ffp+1011" "0x1.ffffffffffbffp+1012"
[53,] "0x1.ffffffffffcp+1012"   "0x1.ffffffffffep+1013"  
[54,] "0x1.ffffffffffep+1013"   "0x1.fffffffffffp+1014"  
[55,] "0x1.fffffffffffp+1014"   "0x1.fffffffffff8p+1015" 
[56,] "0x1.fffffffffff8p+1015"  "0x1.fffffffffffcp+1016" 
[57,] "0x1.fffffffffffcp+1016"  "0x1.fffffffffffep+1017" 
[58,] "0x1.fffffffffffep+1017"  "0x1.ffffffffffffp+1018" 
[59,] "0x1.ffffffffffffp+1018"  "0x1.ffffffffffff8p+1019"
[60,] "0x1.ffffffffffff8p+1019" "0x1.ffffffffffffcp+1020"
[61,] "0x1.ffffffffffffcp+1020" "0x1.ffffffffffffep+1021"
[62,] "0x1.ffffffffffffep+1021" "0x1.fffffffffffffp+1022"
[63,] "0x1.fffffffffffffp+1022" "0x1.fffffffffffffp+1023"
>  #
>   rm(list=ls())
>   q()