PureBasic Forums - English

Posted: **Sat Jan 05, 2008 8:00 pm**

Here is the procedure for division , seems to work fine, (so far) and there is a small test code where i use multiplication again to verify the division result

would like some comments on how to make it faster,

I guess one straight approach is to do away with strings, will be trying to port it into your code format in the coming week,
most of the code is self explantory, i return the remainder of the division via n2$

Code: Select all

Procedure.s Divide (myn1$,myn2$)
Protected sizemyn2=Len(myn2$)

Protected result
;this procedure divides myn1$ by myn2$ i.e myn1$/myn2$ and returns the quotient
result=Compare (myn1$,myn2$) 
myn1$=n1$;we need to reassign these to the padded values 
myn2$=n2$ 
Protected result$="" 
;Protected carry=0 
Protected sizemyn1=Len(myn1$) 

Protected temp_divident$=""
Protected temp_divisor_multiple$=""
Protected quotient$=""
Protected previous_remainder$=Left(myn1$,sizemyn2)
Protected temp_result
;PrintN ("myn2$:"+myn2$)

For i= sizemyn2 To sizemyn1
If i=sizemyn2
temp_divident$=previous_remainder$
Else
 ;we start with the part of myn1$ as long as myn2$
 temp_divident$=previous_remainder$+Mid(myn1$,i,1)
EndIf 
 ;PrintN ("divident:"+temp_divident$)
 temp_result=Compare(temp_divident$,myn2$)
 If (temp_result=0)Or (temp_result=-1) ; if temp_divident$>=myn2$
  For m=1 To 9
   temp_divisor_multiple$=Multiply (myn2$,Str(m))
   temp_result1=Compare(temp_divident$,temp_divisor_multiple$)
   ;PrintN ("m="+Str(m))
   If temp_result1=1
   m=m-1
   ;PrintN ("breaking..")
   Break 
   EndIf
  Next m
 Else
 m=0
 EndIf ; temp_result=0 or temp_result=-1 ends here
 
 If m>9
 m=9
 EndIf
 
 ;PrintN ("m="+Str(m))
 temp_divisor_multiple$=Multiply (myn2$,Str(m))
 
 quotient$=quotient$+Str(m)
 previous_remainder$=Subtract(temp_divident$,temp_divisor_multiple$)
 ;myn1$=LSet("",i,"0")+Right(myn1$,sizemyn1-i)
 ;PrintN ("temp_divisor"+temp_divisor_multiple$)
 ;PrintN("quotient$="+quotient$+";previous_remainder:"+previous_remainder$)
 

 Next i
n2$=previous_remainder$
result$=quotient$
;remove leading zeros
For i=1 To Len(result$)
 If Mid(result$,i,1)="0" And (flag=0)
  cnt=cnt+1
 Else
 Break
EndIf
Next i


result$=Right(result$,Len(result$)-cnt)
ProcedureReturn result$
EndProcedure
;compare function returns 0 if myn1>n2 and 1 if n2>n1 
;it pads the smaller number with leading 0s
;it removes any leading 0s in the original input
;the numbers obtained after a compare are suitable for any mathematical operation now




;num1$="1689087077234463463463474574575057957235732325357997563937563839566393745749395475748386339394756383937336383836363635715723233331243"
;num2$="500057383736353638340584638384739393030846363838373636363763636363636363638393930393735392847569475"
OpenConsole()

;PrintN ("num1$=" +num1$ + ";num2$=" + num2$)
;result=Compare (num1$,num2$)
;Result$= Add(n1$,n2$)
;PrintN ("Result::"+Str(result))
;result$=Multiply(num1$,num2$)
;PrintN ("Multiply::"+result$)
;Input()
;num2$=num1$
;num1$=result$

;num1=Len(num1$)
;num2=Len(num2$)
;PrintN (num1$+"::"+num2$)
;result$=Add(num1$,num2$)
;PrintN ("addition::"+result$)
;PrintN (Str(Val("00100")))

;fact$="1"
;For i=1 To 100
;PrintN(fact$)
;fact$=Multiply(fact$,Str(i))
;Next i

;PrintN("100!="+fact$)
num2$="53530583058305830581913467"
num1$="800970970593757057055000683553746232356"
result$=Divide(num1$,num2$)
PrintN ("Division:"+result$)
remainder$=n2$
PrintN ("Remainder:"+remainder$)
result1$=Multiply(result$,num2$)
result2$=Add(result1$,remainder$)
PrintN ("Multiply:"+result2$)
Input()
CloseConsole()

Posted: **Sat Jan 05, 2008 8:39 pm**

Hi,

fine, that integer math is back in the forum - I'm sure there is somewhere a clever person who is able to implement a fast division routine

Just one more code for doing fast integer calculations:
http://www.purebasic.fr/english/viewtop ... hlight=vnf

The routines above are very (!) fast - e.g. 2000! is done in less than 1/10sec

Code: Select all

fak.VNF
PutNumber(@"1",fak)
n=2000
While n
	VNFmuw(fak,n,fak)
	n-1
Wend
Debug Number(fak)

_____
PS changed font style of 'fast' to bold

Posted: **Sun Jan 06, 2008 6:16 am**

I'm getting some errors in this, but it might be bacause I've swapped out your subtract and multiply functions.

What happens when you use

num2$="80"
num1$="800"

I get answers like

Division:099
Remainder:90992720

perhaps I've broken it. I kept your compare function but I want to get rid of it as its a perf neck

I was working on something similar but it wasn't working out so I'm going to get some ideas from this, it works a lot better than mine

Posted: **Sun Jan 06, 2008 6:56 am**

80 and 800 work fine in the given code

the big thing to look for is when num2$ > longint

that is when the fun starts
example:
num2$="800000000000000000000000000000000000000000000000"
num1$="8000000000000000000000000000000000000000000000000033534543634645745757457000000000000000000000000000000000"

Posted: **Sun Jan 06, 2008 7:03 am**

adding negative numbers to your code is no big deal,
simply check at the time of input and make a truth table:
eg if it is add, the - sign makes it subtract , if both are - then it is add with a - in the answer
mul and div too will work the same way, no issues

the reason i came across compare() and all that is that this was "code as you think" way, not at all optimized

dont bother, my ideal version of the best search algo in a sorted tree is random jumps, not even binary or tert serach which most people use

Posted: **Sun Jan 06, 2008 11:29 am**

Finally!

There's not much info out there on large milti length long division algorithms but I came across this paper (PDF) with lots of details about the topic and some algorthms and ideas

http://citeseer.ist.psu.edu/cache/paper ... isited.pdf

Posted: **Sun Jan 06, 2008 3:23 pm**

Okay, I got my divide working. since it requires multiply, add and subtract too I put it all together to us as an include file.

There's still lots of room for optimisation of course. I decided to try a recusive algorithm to keep it simple but I guess it could risk overloading the stack

. I have had a couple of intermittant errors when the debugger is on and I'm suspecting that the overhead there is blowing the stack. It will need to recurse more if the length of the two numbers are very different (like dividing a 100 digit number with a 2 digit number) but I've seen it recurse up to 5 or 6000 deep and not crash.

I'm not happy with the speed yet but it works which is a start.

I've chucked some code on the end so you can copy-paste and run it

Code: Select all


Structure Bytes 
    Byte.b[0] 
EndStructure

Declare.s String(Char.s, Length.l)
Declare.l IsBiggerInt(BiggerNum.s, SmallerNum.s )
Declare.s BigDivide(DivThis.s, byThis.s)
Declare.s BigMultiply(Num1.s,Num2.s)
Declare.s BigSubtract(Num1.s, Num2.s) ; Answer = Num1 - Num2 (negative output not handled yet)
Declare.s BigAdd(Num1.s, Num2.s)
Declare.s BigSquared(Num.s)
Declare.s BigPow(Num.s, pow.l)


Procedure.s BigSquared(Num.s)

    carry.l = 0

    Num1Len = Len(num)
    Num2len = Num1Len
    
    AnswerLen = Num1len + Num2len ;+1
    
    Dim BCD1.l(Answerlen)
    Dim BCD2.l(Answerlen)
    Dim BCDA.l(Answerlen)
  
    ;Get Data out of string
    *BytePtr1.Bytes = @Num
    *BytePtr2.Bytes = @Num
    
    For i = 1 To Num1len 
        ;bcd1(i) = Val(Mid(Num1,Num1len - i +1 ,1))
        bcd1(i) = *BytePtr1\Byte[Num1len - i] -48
        bcd2(i) = bcd1(i)
    Next
    
    ;multiply
    For L1.l = 1 To Num1len    ;lower number
        carry = 0
        Digit.l = L1;1      
        For L2.l = 1 To Num2len  ;upper number
            BCDA(Digit) + (BCD1(L1) * BCD2(L2)) + carry
            carry = Int(BCDA(Digit) / 10)  ;int()
            BCDA(Digit) = BCDA(Digit) % 10
            Digit + 1    
        Next
                              
        If carry
            BCDA(Digit) = Carry
        EndIf                               
    Next
      
    ; set size
    Answer.s = ""

    For i = AnswerLen To 1 Step -1
        If i = answerlen And bcda(i) = 0
            ;skip trailing 0
        Else
            Answer = Answer + Str(BCDA(i))
        EndIf  
    Next

    ProcedureReturn Answer
    
EndProcedure

Procedure.s BigAdd(Num1.s, Num2.s)

    Carry.b = 0
    Num1Len = Len(num1)
    Num2Len = Len(num2)
    If Num1len > Num2Len
        AnswerLen = Num1len +1
    Else
        AnswerLen = Num2len + 1
    EndIf
    
    Dim BCD1.b(Answerlen)
    Dim BCD2.b(Answerlen)
    Dim BCDA.b(Answerlen)
    
    ;Get Data out of string
    
    For i = 1 To Num1len 
        bcd1(i) = Val(Mid(Num1,Num1len - i +1 ,1))
    Next

    For i = 1 To Num2len 
        bcd2(i) = Val(Mid(Num2,Num2len - i +1,1))
    Next    
       
    ;add
    For i = 1 To Answerlen   
        BCDA(i) = bcd1(i) + bcd2(i) + Carry
        If BCDA(i) > 9
            BCDA(i) - 10
            carry = 1
        Else
            carry = 0
        EndIf
    Next
    
    Answer.s = ""

    For i = AnswerLen To 1 Step -1
        If i = answerlen And bcda(i) = 0
            ;skip trailing 0
        Else
            Answer = Answer + Str(BCDA(i))
        EndIf
        
    Next
    
    ProcedureReturn Answer
    
EndProcedure

Procedure.s BigSubtract(Num1.s, Num2.s) ; Answer = Num1 - Num2 (negative output not handled yet)

    Carry.b = 0
    Num1Len = Len(num1)
    Num2Len = Len(num2)
    If Num1len > Num2Len
        AnswerLen = Num1len +1
    Else
        AnswerLen = Num2len + 1
    EndIf
    
    Dim BCD1.b(Answerlen)
    Dim BCD2.b(Answerlen)
    Dim BCDA.b(Answerlen)
    
    ;Get Data out of string
    
    For i = 1 To Num1len 
        bcd1(i) = Val(Mid(Num1,Num1len - i +1 ,1))
    Next

    For i = 1 To Num2len 
        bcd2(i) = Val(Mid(Num2,Num2len - i +1,1))
    Next   
       
    ;add
    For i = 1 To Answerlen   
        BCDA(i) = bcd1(i) - bcd2(i) - Carry
        If BCDA(i) < 0
            BCDA(i) + 10
            carry = 1
        Else
            carry = 0
        EndIf
    Next
    
    Answer.s = ""
    
    StringStart.l = #True

    For i = AnswerLen To 1 Step -1
        If StringStart = #True And bcda(i) = 0  ;If i = answerlen And bcda(i) = 0
            ;skip trailing 0
        Else
            Answer = Answer + Str(BCDA(i))
            StringStart = #False
        EndIf        
    Next    
   
    ProcedureReturn Answer

EndProcedure

Procedure.s BigMultiply(Num1.s,Num2.s)

    carry.l = 0

    Num1Len = Len(num1)
    Num2Len = Len(num2)
    
    AnswerLen = Num1len + Num2len ;+1
    
    Dim BCD1.b(Answerlen)
    Dim BCD2.b(Answerlen)
    Dim BCDA.c(Answerlen)
    
    ;Get Data out of string
    *BytePtr1.Bytes = @Num1
    *BytePtr2.Bytes = @Num2 
    
    For i = 1 To Num1len 
        bcd1(i) = *BytePtr1\Byte[Num1len - i] -48
    Next

    For i = 1 To Num2len 
        bcd2(i) = *BytePtr2\Byte[Num2len - i] -48
    Next

    ;multiply
    
    For L1.l = 1 To Num1len    ;lower number
        carry = 0
        Digit.l = L1;1      
        For L2.l = 1 To Num2len  ;upper number
            BCDA(Digit) + (BCD1(L1) * BCD2(L2)) + carry
            carry = BCDA(Digit) / 10  ;int()
            BCDA(Digit) = BCDA(Digit) % 10
            Digit + 1    
        Next
                              
        If carry
            BCDA(Digit) = Carry
        EndIf                               
    Next
      
    ; set size
   Answer.s = RSet(Chr(0), AnswerLen,Chr(0))

    *AnsPtr.Bytes = @Answer

    For i = AnswerLen To 1 Step -1
        *AnsPtr\byte[AnswerLen - i ] = BCDA(i) + 48
    Next
    
    If Left(Answer,1) = "0"
        Answer = Right(answer,AnswerLen-1)
    EndIf

    ProcedureReturn Answer
    
EndProcedure

Procedure.s BigPow(Num.s, pow.l)

    If pow = 1 ;
        ProcedureReturn num
    ElseIf pow = 2
        ProcedureReturn BigMultiply(Num,num)
    ElseIf pow = 3
        ProcedureReturn BigMultiply(Num, BigPow(Num,2))
    Else
        
        
        If pow % 2 = 0
            ProcedureReturn BigMultiply(BigPow(Num,pow / 2), BigPow(Num,pow / 2))  
        Else
            ProcedureReturn BigMultiply(BigMultiply(BigPow(Num,pow / 2), BigPow(Num,pow / 2)),num)  
        EndIf
    EndIf
    
EndProcedure

Procedure.l IsBiggerInt(BiggerNum.s, SmallerNum.s ) ;returns true if first int is bigger

    Biglen.l = Len(BiggerNum)
    SmallLen.l = Len(SmallerNum)
    
    If Biglen > SmallLen
        ProcedureReturn #True
    ElseIf Biglen < SmallLen
        ProcedureReturn #False
    Else ;same len, string compare
        If BiggerNum > SmallerNum
            ProcedureReturn #True
        Else
            ProcedureReturn #False
        EndIf       
    EndIf 

EndProcedure

Procedure.s BigDivide(DivThis.s, byThis.s)  ; DivThis / byThis (First param should be bigger or result is zero)
   
    Define Result.s , divisorInc.l
    
    DivThislen.l = Len(DivThis)
    byThisLen.l = Len(byThis)
    
    NewDivisor.s = byThis
    NewDivLen.l = byThisLen
    
    If IsBiggerInt(byThis, DivThis)
        ProcedureReturn "0"
    EndIf
    
    DifLen = DivThislen - NewDivLen

    If diflen > 0 
    
        divisorInc = Val(Left(DivThis,1)) / Val(Left(byThis,1)) -1
        If divisorInc = < 1
            divisorInc = 1
        EndIf
                
        NewDivisor = BigMultiply(NewDivisor, Str(divisorInc) + string("0",DifLen-1))       
        NewDivLen = Len(NewDivisor)
        Result =  BigAdd(Result, Str(divisorInc) + string("0",DifLen-1))
        ; subtract calced piece
        DivThis = BigSubtract(DivThis, NewDivisor)
        DivThislen = Len(DivThis) 
        
        result = bigadd(result,BigDivide(DivThis.s, byThis.s)) 
    Else
    
        For i = 2 To 10 
            If IsBiggerInt(BigMultiply(Str(i),byThis),DivThis)
                result = Str(i-1)
                Break
            EndIf    
        Next
    
    EndIf
   
    ProcedureReturn(result)


EndProcedure

Procedure.s String(Char.s, Length.l)

    If Len(Char) = 1 
        TempStr.s = Space(Length)
        ReplaceString(TempStr," ",Char,2) 
        ProcedureReturn TempStr
    Else
        ProcedureReturn ""
    EndIf

EndProcedure



; ============= Test code from here ===========

DisableDebugger

OpenConsole()

    Divide.s = "26098234092834075093854092383452340987340295703295702398573290457987043295783240957239032475320948799" 
    By.s     = "9309478523409785234085234095780329457987394752309870987"

    PrintN("Divide   " + divide)
    PrintN("By       " + by)
    PrintN("")
    PrintN(BigDivide(Divide, by )) 
    Input()


CloseConsole()

Posted: **Mon Jan 07, 2008 6:58 am**

well yep it seems very similar to mine, except i think the array instead of a string may make is faster
! cool thanks! I'll see what else I can add to this!

Posted: **Mon Jan 07, 2008 7:01 am**

how come this never broke for you:
If BiggerNum > SmallerNum
ProcedureReturn #True

i mean no reason why 0 pads etc cant break it

Posted: **Mon Jan 07, 2008 7:45 am**

Why would it break?

It only does this if the strings are of the same size, if they aren't then it takes the bigger one.

It works because:

1. My other functions always strip leading zeros (it's a check at the end of them)

2. I don't support negative numbers with a leading dash "-"

it may yet come back and bite me though.

I have an idea for an optimisation of my div function that I will try out later this week, for many sums the difference could be quite dramatic.

Posted: **Mon Jan 07, 2008 8:24 am**

what exactly are you trying to optimize? i was trying to do away with the compare..but no luck so far

Posted: **Mon Jan 07, 2008 11:54 am**

The "trick" for optimizing is to collect ciphers into groups so it can be handled faster...

Thats easy for add/sub/mul, for example:
123456 x 123456 can be done with doing 36 multiplications (1x1, 1x2, 1x3,... 6x6) and quite a lot of additions - which costs a lot of time.
But if grouped, let's say into pieces with three ciphers, you only have to do 4 multiplications: 123x123, 123x456, 456x123 and 456x456 - that should be already about 10 times faster!

I played around to find a large number of cipher which can be handled without problems, using quads at the beginning but returned to long variables which are able to work with numbers up to +/-2.000.000.000 - this allows to collect even 9 ciphers to a single group and have space for the "overflow" when adding two of these numbers.
So that gives an enormous boost to the calculating procedures and I'm quite happy about this

The most important problem (at least for me) is to find a routine to divide numbers using these groups!

All I was able to do for now is programming a (simple) division routine, which is used when the divisor is smaller than 1x10^10...

Michael

Posted: **Mon Jan 07, 2008 11:58 am**

well if i could go back to the old defense labs I used to be in back in 2000 i am sure I could get a x86 or a plain microcontroller to handle this code (nevermind the strings, they are just for my ease since i am coding at a much higher level)
it wont take much to throw a standard hardware at this and get a number

Posted: **Mon Jan 07, 2008 11:59 am**

Paul, if you could help figure out how to return the remainder via your code, it would be great!

Posted: **Mon Jan 07, 2008 1:44 pm**

I think this is what you mean.

Please note that I have updated the BigSubtract() function in the full code above as it was leaving trailing zeros in some cases that were causing problems.

Code: Select all


Procedure.s BigMod(DivThis.s, byThis.s)

    Divisor.s = BigDivide(DivThis, byThis)
    Multiplied.s = bigMultiply(byThis,Divisor)

    ProcedureReturn bigSubtract(DivThis,Multiplied)

EndProcedure

I don't know if there are shortcuts to calc mod or not (My math aith that good). This way, as I update the bigdivide function bigmod will get faster too!

I'm going to write some test code to work with quads to compare my functions with quad range random numbers in long loops to help make sure there are no situations where these are returning the wrong answer, if in the quad space they are correct they should scale upward okay

PureBasic Forums - English

generic module to add , subtract, mul 2 numbers of any size

division

negative numbers