Index
> Projects and Ideas > multiprecision calculator for win64Goto page 1, 2 Next |
| Author |
|
|
tthsqe 07 Mar 2013, 04:41
features:
+ all of the stuff from the previous 32 bit version + number of unlimited size via N[-,digits] (GMP) + memoziation of functions for machine integers + flexible and unified input language + able to return multiple values + uses asymptotically fastest known algorithms for constants and inverse functions + better interface with special characters + limited array support + SSA-based optimizer
Last edited by tthsqe on 10 Aug 2013, 15:21; edited 12 times in total |
|||||||||||
07 Mar 2013, 04:41 |
|
|
revolution 07 Mar 2013, 07:51
It looks like your precision of PI is too low to support the 48 digit output.
|
|||
|
07 Mar 2013, 07:51 |
|
|
hopcode 07 Mar 2013, 11:55
revolution wrote: It looks like your precision of PI is too low to support the 48 digit output. last method is lookup tables. choosing the method is not trivial in all cases. _________________ ⠓⠕⠏⠉⠕⠙⠑ |
|||
|
07 Mar 2013, 11:55 |
|
|
tthsqe 25 Mar 2013, 06:44
Quote: It looks like your precision of PI is too low to support the 48 digit output. it calculates all of those digits - the window on the right is not intended to display all of them Quote: in other cases, using Chebyshev is ok for at least 11 digits approximation on operations with functions not quite sure what this is referring to - the readme contains descriptions of the algorithms to compute stuff |
|||
|
25 Mar 2013, 06:44 |
|
|
bitRAKE 25 Mar 2013, 12:16
Awesome work.
N[Sqrt[2],69213] crashed in gmp.__gmp_doprnt_mpf2 N[Sqrt[2],69212] doesn't though. 
General limitation of N[] it seems, due to GMP.dll? Still reading the source. Is there a way to copy the output window? Because this already does about 90% of what I've been using other calculators for. I'd like to generate some number and partition it into dq/dd/etc.. pieces for inclusion in code. I'll try to break it some more... [edit]Every win64 OS comes with the Cambia Math Font - which is quite nice. My first example: Code: S[n]:=If[n==0,4,S[n-1]^2-2]; M[p]:=2^p-1; v=7;N[S[v-2]/M[v],20] _________________ ¯\(°_o)/¯ “languages are not safe - uses can be” Bjarne Stroustrup |
|||
|
25 Mar 2013, 12:16 |
|
|
hopcode 25 Mar 2013, 16:12
good work
i think GMP may be avoided in some way, because calculations are already done FPU on big numbers. the degree of a minimax polynomials is bound directly to the required accuracy of a function. Chebyshev may be used in all cases of bare functions, except sqrt and arcsin. a general reference here http://www.research.scea.com/gdc2003/fast-math-functions_p2.pdf and the same as textual, browse down to "A Note On Convergence" http://www.research.scea.com/research/pdfs/RGREENfastermath_GDC02.pdf _________________ ⠓⠕⠏⠉⠕⠙⠑ |
|||
|
25 Mar 2013, 16:12 |
|
|
tthsqe 25 Mar 2013, 17:10
Ok - I fixed the issue that bitrake found and some others.
What are some good formats for outputting the numbers in hex? |
|||
|
25 Mar 2013, 17:10 |
|
|
bitRAKE 29 Jun 2013, 06:14
FASM is very flexible, so if you output hex as 0xF... then it's also good for C/++. You probably don't want to include a whole string library to support general partitioning and concatenation.
p.s. I miss the delete key.  |
|||
|
29 Jun 2013, 06:14 |
|
|
tthsqe 10 Jul 2013, 07:09
So for hex string of integers, I think it is clear how it should be done:
Code: HexString[-3^200] dq 0xA40500E155074F5F, 0x7C1309091B5851DD, 0x028C2681BB89F949, 0x3DE56C85890CBCD0, 0xE02A79C3C14FB961 The result is in 2-complement and Little endian form. For hex strings of real numbers, I am thinking of just extending the 80-bit (http://en.wikipedia.org/wiki/Extended_precision) format by making a longer mantissa. So, the exponent and sign still has the same meaning, and the radix point is between the 62nd and 63rd bit of the most significant limb. The 63rd bit of the most significant limb should be 1. EX: Code: HexString[-p\^2] dq 0x08566a3fe0d0a228, 0x56e26cd9808c1ac7, 0x9de9e64df22ef2d2, 3, -1 would indicate that Code: -Pi^2 = - 2^4 * 0.9de9e64df22ef2d256e26cd9808c1ac708566a3fe0d0a228........ Again, the exponent on the 2 is 4 because the implied radix point in the representation is on the rhs of the 63rd bit of the most sig limb. This kind of representation is nice when rolling your own multiprecision arithmetic functions on numbers of a fixed size. |
|||
|
10 Jul 2013, 07:09 |
|
|
bitRAKE 15 Jul 2013, 01:21
Had some confusion and a crash:
Code: n=255;dn=90/n;a={0}; S[x]:=n*Sin[x]; dS[x]:=S[(x+1)*dn]-S[x*dn]; Do[a[[i]]=Floor[dS[i]],{i,0,n}] Not sure where the syntax error is? Sin probably takes radians rather than degrees. |
|||
|
15 Jul 2013, 01:21 |
|
|
tthsqe 15 Jul 2013, 02:17
All that is behavior that I expected.
When you do a={}, that is invalid syntax because the List function has an arg count of zero. When you do a={0}, you are only guaranteed that 32bytes(one entry) are allocated for a. No bounds checking are done for statements such aa a[[i]]=x, so if n is high enough, you will get seg fault. Code: n=256; dn=90/(n); a = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; S[x]:=n*Sin[x]; dS[x]:=S[(x+1)*dn]-S[x*dn]; Do[a[[i]]=Floor[dS[i]],{i,0,n-1}]; a {88,77,57,29,-2,-32,-59,-79,-89,-88,-77,-56,-28,3,33,60,79,88,87,75,53,26,-6,-36,-62,-81,-90,-87,-75,-53,-25,7,37,63,81,89,86,73,50,22,-10,-40,-65,-82,-90,-86,-72,-50,-21,11,40,65,82,89,85,70,47,18,-14,-43,-68,-84,-90,-85,-70,-46,-17,15,44,68,84,89,83,68,43,14,-17,-47,-70,-85,-90,-84,-67,-43,-13,18,47,71,85,89,82,65,40,10,-21,-50,-73,-87,-90,-82,-64,-39,-9,22,51,73,86,89,80,62,36,6,-25,-53,-75,-88,-89,-80,-62,-35,-5,26,54,75,87,88,78,59,32,2,-29,-57,-77,-88,-89,-79,-59,-32,-1,30,57,77,88,88,77,56,29,-2,-33,-60,-79,-89,-88,-76,-56,-28,3,34,60,79,88,87,74,53,25,-6,-37,-63,-81,-90,-87,-74,-52,-24,7,37,63,81,89,86,72,50,21,-10,-40,-65,-83,-90,-86,-72,-49,-20,11,41,66,82,89,85,70,46,17,-14,-44,-68,-84,-90,-85,-69,-46,-16,15,45,68,84,89,83,67,43,13,-18,-47,-71,-86,-90,-83,-67,-42,-12,19,48,71,85,89,82,64,39,9,-22,-51,-73,-87,-90,-82,-64,-38,-8,23,51,73,86,89,80,62,36,5,-26} With the source I just posted (slight type in previous), the above does what you want, but there are some obvious problems. There are two possibilities: 1. I could implement a CreateArray function which allocates space for an array of specified size 2. I could implement the missing Table operator In either case the proposed code would be 1. Code: n=255; dn=90/n; a = CreateArray[n]; S[x]:=n*Sin[x]; dS[x]:=S[(x+1)*dn]-S[x*dn]; Do[a[[i]]=Floor[dS[i]],{i,0,n-1}]; a 2. Code: n=255; dn=90/n; S[x]:=n*Sin[x]; dS[x]:=S[(x+1)*dn]-S[x*dn]; a = Table[Floor[dS[i]],{i,0,n-1}] Which is better? They are both quite doable... maybe both? And yes, I do not use degrees. |
|||
|
15 Jul 2013, 02:17 |
|
|
bitRAKE 15 Jul 2013, 03:44
That makes sense, no dynamic array. Although the routine you posted also crashes. Think I found the error (Evaluator.inc, line 132):
Code: .ARRAY: movzx edx,word[rbx+12] xor ecx,ecx shl edx,5 mov r8d,MEM_COMMIT mov r9d,PAGE_READWRITE call [VirtualAlloc] mov rdi,rax mov rsi,qword[rbx+0] my_pinsrq xmm6,rax,0,xmm0 mov ecx,[rbx+8] shl ecx,4 rep movsq jmp .ret *Looks like I might have had an older version of the code. _________________ ¯\(°_o)/¯ “languages are not safe - uses can be” Bjarne Stroustrup |
|||
|
15 Jul 2013, 03:44 |
|
|
bitRAKE 15 Jul 2013, 05:33
Okay, so I tried another method (weird character is Pi):
Code: n=255; dn=90/(n); S[x]:=Floor[n*Sin[x * §/180]]; dS[x]:=S[(x+1)*dn]-S[x*dn]; a=0;i=0; While[j=dS[i], If[j==2,a=a+Power[2,i]]; i=i+1 ]; b=i; While[i<256, If[dS[i]==1,a=a+Power[2,i]]; i=i+1 ]; HexString[a] The first loop does execute, but not the second. _________________ ¯\(°_o)/¯ “languages are not safe - uses can be” Bjarne Stroustrup |
|||
|
15 Jul 2013, 05:33 |
|
|
tthsqe 15 Jul 2013, 08:19
Man, i guess you should not be very confident in the bug-freenes of the program
I'll have a look at that and see what the issue is. I don't think that Code: While[j=dS[i], If[j==2,a=a+Power[2,i]]; i=i+1 ]; is doing what you want. Nor do I know what you are trying to do. It is just like Mathematica - the first argument of While should return a boole, and yours is calculating dS[i], saving dS[i] to j, then returning the value of dS[i]. It is true that I mapped booles onto machine integers (true -> -1, false->0), but While expects that its first argument, say x, is a machine integer, and says that x is true depending on whether the lowest byte of x is nonzero or not. PS: dynamic array allocation and the Table function are almost ready, if this is what you are trying to replicate... So I ran Code: n = 255; dn = 90/(n); S[x_] := Floor[n*Sin[x*Pi/180]]; dS[x_] := S[(x + 1)*dn] - S[x*dn]; a = 0; i = 0; While[BitAnd[j = dS[i], 255] =!= 0, If[j == 2, a = a + Power[2, i]]; i = i + 1]; b = i; While[i < 256, If[dS[i] == 1, a = a + Power[2, i]]; i = i + 1]; IntegerString[a, 16] in mathematica and the result was Code: 40010444494a95aadb6dddf7e00000020208444492494a552aaaaaab556ad5aa is this your intended number? |
|||
|
15 Jul 2013, 08:19 |
|
|
bitRAKE 15 Jul 2013, 08:54
DOS386, posted a link for a sine table generator. Since sine is an increasing function on [0-90] degrees, a sub-table of 256 bits are generated which indicate if 0,1, or 2 should be added to get the next table value. I wanted to reproduce that sub-table of bits. We need to know when the first repeat entry (zero delta) happens to have the code switch from adding 1,2 to adding 0,1.
So, the first loop stops when dS[i]==0 and does appear to do that, although the value of (i) is 156, and should be 209, iirc. I've changed the While function to j !=\ 0, and it seems to perform the same. dS[i] only returns 0,1,2. The second loop doesn't work, and I'm still debugging. My confidence only waned slightly, as I'm still learning what's what. The table equates to (in hex): 55ab56aad5555554aa5292492222104040000007efbbb6db55a9529222208003 ...but there is an extra set bit at the top. Actually, 156 is correct for the first loop as there are 100 bits set - for a total of 256. _________________ ¯\(°_o)/¯ “languages are not safe - uses can be” Bjarne Stroustrup Last edited by bitRAKE on 15 Jul 2013, 09:09; edited 2 times in total |
|||
|
15 Jul 2013, 08:54 |
|
|
tthsqe 15 Jul 2013, 09:04
Ok, so if the argument of the first while is (j = dS[i]) != 0, the calculator and mathematica both terminate with
a=683266112844956550622064541659919209898 i=156 so I think 156 is correct. It's really strange - the second While should not even be doing such a conversion. |
|||
|
15 Jul 2013, 09:04 |
|
|
bitRAKE 15 Jul 2013, 09:13
tthsqe wrote: in mathematica and the result was _________________ ¯\(°_o)/¯ “languages are not safe - uses can be” Bjarne Stroustrup |
|||
|
15 Jul 2013, 09:13 |
|
|
tthsqe 15 Jul 2013, 09:46
Ok, since the mistake was trivial I'll just have you correct it in your source, since mine is currently busy.
line 83 or so of printer.inc should read Code: Print_B_INT: test ecx,ecx jz .float push rbp mov rbp,rsp and rsp,-16 invoke gmp_sprintf,rdi,'%Zd',rsi add rdi,rax mov rsp,rbp pop rbp ret .float: mov rax,'Big Int ' stosq mov rcx,rsi call Convert_B_INT_to_M_FP ; not B_REAL_to_M_FP fstp st1 call PrintFloat ret Nevermind this: the newer posted version has this corrected |
|||
|
15 Jul 2013, 09:46 |
|
|
tthsqe 09 Aug 2013, 10:16
bitrake, I have seen you with an integer square root algorithm before, Could you post it here? It should work as follows:
input: positive integer x (3 qword limbs) return: s = Floor[Sqrt[x]] and r = x-s^2 Actually, if r is not zero, I do not need r nor do I need s. |
|||
|
09 Aug 2013, 10:16 |
|
| Goto page 1, 2 Next < Last Thread | Next Thread > |
Forum Rules:
|