Stephen Paulger / ttmath

Library for big numbers from http://www.ttmath.org/

Dependents:   PIDHeater82 Conceptcontroller_v_1_0 AlarmClockApp COG4050_adxl355_tilt ... more

TTMath is a small library which allows one to perform arithmetic operations with big unsigned integer, big signed integer and big floating point numbers. It provides standard mathematical operations like adding, subtracting, multiplying, dividing.

TTMath is BSD Licensed (new/modified BSD)

For more information about ttmath see http://www.ttmath.org/

ttmathuint_noasm.h@0:04a9f72bbca7, 2013-07-30 (annotated)

Committer:
stevep
Date:
Tue Jul 30 18:43:48 2013 +0000
Revision:
0:04a9f72bbca7
v0.9.3 of ttmath

Who changed what in which revision?

UserRevisionLine numberNew contents of line
stevep 0:04a9f72bbca7 1 /*
stevep 0:04a9f72bbca7 2 * This file is a part of TTMath Bignum Library
stevep 0:04a9f72bbca7 3 * and is distributed under the (new) BSD licence.
stevep 0:04a9f72bbca7 4 * Author: Tomasz Sowa <t.sowa@ttmath.org>
stevep 0:04a9f72bbca7 5 */
stevep 0:04a9f72bbca7 6
stevep 0:04a9f72bbca7 7 /*
stevep 0:04a9f72bbca7 8 * Copyright (c) 2006-2010, Tomasz Sowa
stevep 0:04a9f72bbca7 9 * All rights reserved.
stevep 0:04a9f72bbca7 10 *
stevep 0:04a9f72bbca7 11 * Redistribution and use in source and binary forms, with or without
stevep 0:04a9f72bbca7 12 * modification, are permitted provided that the following conditions are met:
stevep 0:04a9f72bbca7 13 *
stevep 0:04a9f72bbca7 14 * * Redistributions of source code must retain the above copyright notice,
stevep 0:04a9f72bbca7 15 * this list of conditions and the following disclaimer.
stevep 0:04a9f72bbca7 16 *
stevep 0:04a9f72bbca7 17 * * Redistributions in binary form must reproduce the above copyright
stevep 0:04a9f72bbca7 18 * notice, this list of conditions and the following disclaimer in the
stevep 0:04a9f72bbca7 19 * documentation and/or other materials provided with the distribution.
stevep 0:04a9f72bbca7 20 *
stevep 0:04a9f72bbca7 21 * * Neither the name Tomasz Sowa nor the names of contributors to this
stevep 0:04a9f72bbca7 22 * project may be used to endorse or promote products derived
stevep 0:04a9f72bbca7 23 * from this software without specific prior written permission.
stevep 0:04a9f72bbca7 24 *
stevep 0:04a9f72bbca7 25 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
stevep 0:04a9f72bbca7 26 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
stevep 0:04a9f72bbca7 27 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
stevep 0:04a9f72bbca7 28 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
stevep 0:04a9f72bbca7 29 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
stevep 0:04a9f72bbca7 30 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
stevep 0:04a9f72bbca7 31 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
stevep 0:04a9f72bbca7 32 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
stevep 0:04a9f72bbca7 33 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
stevep 0:04a9f72bbca7 34 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
stevep 0:04a9f72bbca7 35 * THE POSSIBILITY OF SUCH DAMAGE.
stevep 0:04a9f72bbca7 36 */
stevep 0:04a9f72bbca7 37
stevep 0:04a9f72bbca7 38 #ifndef headerfilettmathuint_noasm
stevep 0:04a9f72bbca7 39 #define headerfilettmathuint_noasm
stevep 0:04a9f72bbca7 40
stevep 0:04a9f72bbca7 41
stevep 0:04a9f72bbca7 42 #ifdef TTMATH_NOASM
stevep 0:04a9f72bbca7 43
stevep 0:04a9f72bbca7 44 /*!
stevep 0:04a9f72bbca7 45 \file ttmathuint_noasm.h
stevep 0:04a9f72bbca7 46 \brief template class UInt<uint> with methods without any assembler code
stevep 0:04a9f72bbca7 47
stevep 0:04a9f72bbca7 48 this file is included at the end of ttmathuint.h
stevep 0:04a9f72bbca7 49 */
stevep 0:04a9f72bbca7 50
stevep 0:04a9f72bbca7 51
stevep 0:04a9f72bbca7 52 namespace ttmath
stevep 0:04a9f72bbca7 53 {
stevep 0:04a9f72bbca7 54
stevep 0:04a9f72bbca7 55 /*!
stevep 0:04a9f72bbca7 56 returning the string represents the currect type of the library
stevep 0:04a9f72bbca7 57 we have following types:
stevep 0:04a9f72bbca7 58 asm_vc_32 - with asm code designed for Microsoft Visual C++ (32 bits)
stevep 0:04a9f72bbca7 59 asm_gcc_32 - with asm code designed for GCC (32 bits)
stevep 0:04a9f72bbca7 60 asm_vc_64 - with asm for VC (64 bit)
stevep 0:04a9f72bbca7 61 asm_gcc_64 - with asm for GCC (64 bit)
stevep 0:04a9f72bbca7 62 no_asm_32 - pure C++ version (32 bit) - without any asm code
stevep 0:04a9f72bbca7 63 no_asm_64 - pure C++ version (64 bit) - without any asm code
stevep 0:04a9f72bbca7 64 */
stevep 0:04a9f72bbca7 65 template<uint value_size>
stevep 0:04a9f72bbca7 66 const char * UInt<value_size>::LibTypeStr()
stevep 0:04a9f72bbca7 67 {
stevep 0:04a9f72bbca7 68 #ifdef TTMATH_PLATFORM32
stevep 0:04a9f72bbca7 69 static const char info[] = "no_asm_32";
stevep 0:04a9f72bbca7 70 #endif
stevep 0:04a9f72bbca7 71
stevep 0:04a9f72bbca7 72 #ifdef TTMATH_PLATFORM64
stevep 0:04a9f72bbca7 73 static const char info[] = "no_asm_64";
stevep 0:04a9f72bbca7 74 #endif
stevep 0:04a9f72bbca7 75
stevep 0:04a9f72bbca7 76 return info;
stevep 0:04a9f72bbca7 77 }
stevep 0:04a9f72bbca7 78
stevep 0:04a9f72bbca7 79
stevep 0:04a9f72bbca7 80 /*!
stevep 0:04a9f72bbca7 81 returning the currect type of the library
stevep 0:04a9f72bbca7 82 */
stevep 0:04a9f72bbca7 83 template<uint value_size>
stevep 0:04a9f72bbca7 84 LibTypeCode UInt<value_size>::LibType()
stevep 0:04a9f72bbca7 85 {
stevep 0:04a9f72bbca7 86 #ifdef TTMATH_PLATFORM32
stevep 0:04a9f72bbca7 87 LibTypeCode info = no_asm_32;
stevep 0:04a9f72bbca7 88 #endif
stevep 0:04a9f72bbca7 89
stevep 0:04a9f72bbca7 90 #ifdef TTMATH_PLATFORM64
stevep 0:04a9f72bbca7 91 LibTypeCode info = no_asm_64;
stevep 0:04a9f72bbca7 92 #endif
stevep 0:04a9f72bbca7 93
stevep 0:04a9f72bbca7 94 return info;
stevep 0:04a9f72bbca7 95 }
stevep 0:04a9f72bbca7 96
stevep 0:04a9f72bbca7 97
stevep 0:04a9f72bbca7 98 /*!
stevep 0:04a9f72bbca7 99 this method adds two words together
stevep 0:04a9f72bbca7 100 returns carry
stevep 0:04a9f72bbca7 101
stevep 0:04a9f72bbca7 102 this method is created only when TTMATH_NOASM macro is defined
stevep 0:04a9f72bbca7 103 */
stevep 0:04a9f72bbca7 104 template<uint value_size>
stevep 0:04a9f72bbca7 105 uint UInt<value_size>::AddTwoWords(uint a, uint b, uint carry, uint * result)
stevep 0:04a9f72bbca7 106 {
stevep 0:04a9f72bbca7 107 uint temp;
stevep 0:04a9f72bbca7 108
stevep 0:04a9f72bbca7 109 if( carry == 0 )
stevep 0:04a9f72bbca7 110 {
stevep 0:04a9f72bbca7 111 temp = a + b;
stevep 0:04a9f72bbca7 112
stevep 0:04a9f72bbca7 113 if( temp < a )
stevep 0:04a9f72bbca7 114 carry = 1;
stevep 0:04a9f72bbca7 115 }
stevep 0:04a9f72bbca7 116 else
stevep 0:04a9f72bbca7 117 {
stevep 0:04a9f72bbca7 118 carry = 1;
stevep 0:04a9f72bbca7 119 temp = a + b + carry;
stevep 0:04a9f72bbca7 120
stevep 0:04a9f72bbca7 121 if( temp > a ) // !(temp<=a)
stevep 0:04a9f72bbca7 122 carry = 0;
stevep 0:04a9f72bbca7 123 }
stevep 0:04a9f72bbca7 124
stevep 0:04a9f72bbca7 125 *result = temp;
stevep 0:04a9f72bbca7 126
stevep 0:04a9f72bbca7 127 return carry;
stevep 0:04a9f72bbca7 128 }
stevep 0:04a9f72bbca7 129
stevep 0:04a9f72bbca7 130
stevep 0:04a9f72bbca7 131
stevep 0:04a9f72bbca7 132 /*!
stevep 0:04a9f72bbca7 133 this method adding ss2 to the this and adding carry if it's defined
stevep 0:04a9f72bbca7 134 (this = this + ss2 + c)
stevep 0:04a9f72bbca7 135
stevep 0:04a9f72bbca7 136 c must be zero or one (might be a bigger value than 1)
stevep 0:04a9f72bbca7 137 function returns carry (1) (if it was)
stevep 0:04a9f72bbca7 138 */
stevep 0:04a9f72bbca7 139
stevep 0:04a9f72bbca7 140 template<uint value_size>
stevep 0:04a9f72bbca7 141 uint UInt<value_size>::Add(const UInt<value_size> & ss2, uint c)
stevep 0:04a9f72bbca7 142 {
stevep 0:04a9f72bbca7 143 uint i;
stevep 0:04a9f72bbca7 144
stevep 0:04a9f72bbca7 145 for(i=0 ; i<value_size ; ++i)
stevep 0:04a9f72bbca7 146 c = AddTwoWords(table[i], ss2.table[i], c, &table[i]);
stevep 0:04a9f72bbca7 147
stevep 0:04a9f72bbca7 148 TTMATH_LOGC("UInt::Add", c)
stevep 0:04a9f72bbca7 149
stevep 0:04a9f72bbca7 150 return c;
stevep 0:04a9f72bbca7 151 }
stevep 0:04a9f72bbca7 152
stevep 0:04a9f72bbca7 153
stevep 0:04a9f72bbca7 154 /*!
stevep 0:04a9f72bbca7 155 this method adds one word (at a specific position)
stevep 0:04a9f72bbca7 156 and returns a carry (if it was)
stevep 0:04a9f72bbca7 157
stevep 0:04a9f72bbca7 158 if we've got (value_size=3):
stevep 0:04a9f72bbca7 159 table[0] = 10;
stevep 0:04a9f72bbca7 160 table[1] = 30;
stevep 0:04a9f72bbca7 161 table[2] = 5;
stevep 0:04a9f72bbca7 162 and we call:
stevep 0:04a9f72bbca7 163 AddInt(2,1)
stevep 0:04a9f72bbca7 164 then it'll be:
stevep 0:04a9f72bbca7 165 table[0] = 10;
stevep 0:04a9f72bbca7 166 table[1] = 30 + 2;
stevep 0:04a9f72bbca7 167 table[2] = 5;
stevep 0:04a9f72bbca7 168
stevep 0:04a9f72bbca7 169 of course if there was a carry from table[2] it would be returned
stevep 0:04a9f72bbca7 170 */
stevep 0:04a9f72bbca7 171 template<uint value_size>
stevep 0:04a9f72bbca7 172 uint UInt<value_size>::AddInt(uint value, uint index)
stevep 0:04a9f72bbca7 173 {
stevep 0:04a9f72bbca7 174 uint i, c;
stevep 0:04a9f72bbca7 175
stevep 0:04a9f72bbca7 176 TTMATH_ASSERT( index < value_size )
stevep 0:04a9f72bbca7 177
stevep 0:04a9f72bbca7 178
stevep 0:04a9f72bbca7 179 c = AddTwoWords(table[index], value, 0, &table[index]);
stevep 0:04a9f72bbca7 180
stevep 0:04a9f72bbca7 181 for(i=index+1 ; i<value_size && c ; ++i)
stevep 0:04a9f72bbca7 182 c = AddTwoWords(table[i], 0, c, &table[i]);
stevep 0:04a9f72bbca7 183
stevep 0:04a9f72bbca7 184 TTMATH_LOGC("UInt::AddInt", c)
stevep 0:04a9f72bbca7 185
stevep 0:04a9f72bbca7 186 return c;
stevep 0:04a9f72bbca7 187 }
stevep 0:04a9f72bbca7 188
stevep 0:04a9f72bbca7 189
stevep 0:04a9f72bbca7 190
stevep 0:04a9f72bbca7 191
stevep 0:04a9f72bbca7 192
stevep 0:04a9f72bbca7 193 /*!
stevep 0:04a9f72bbca7 194 this method adds only two unsigned words to the existing value
stevep 0:04a9f72bbca7 195 and these words begin on the 'index' position
stevep 0:04a9f72bbca7 196 (it's used in the multiplication algorithm 2)
stevep 0:04a9f72bbca7 197
stevep 0:04a9f72bbca7 198 index should be equal or smaller than value_size-2 (index <= value_size-2)
stevep 0:04a9f72bbca7 199 x1 - lower word, x2 - higher word
stevep 0:04a9f72bbca7 200
stevep 0:04a9f72bbca7 201 for example if we've got value_size equal 4 and:
stevep 0:04a9f72bbca7 202 table[0] = 3
stevep 0:04a9f72bbca7 203 table[1] = 4
stevep 0:04a9f72bbca7 204 table[2] = 5
stevep 0:04a9f72bbca7 205 table[3] = 6
stevep 0:04a9f72bbca7 206 then let
stevep 0:04a9f72bbca7 207 x1 = 10
stevep 0:04a9f72bbca7 208 x2 = 20
stevep 0:04a9f72bbca7 209 and
stevep 0:04a9f72bbca7 210 index = 1
stevep 0:04a9f72bbca7 211
stevep 0:04a9f72bbca7 212 the result of this method will be:
stevep 0:04a9f72bbca7 213 table[0] = 3
stevep 0:04a9f72bbca7 214 table[1] = 4 + x1 = 14
stevep 0:04a9f72bbca7 215 table[2] = 5 + x2 = 25
stevep 0:04a9f72bbca7 216 table[3] = 6
stevep 0:04a9f72bbca7 217
stevep 0:04a9f72bbca7 218 and no carry at the end of table[3]
stevep 0:04a9f72bbca7 219
stevep 0:04a9f72bbca7 220 (of course if there was a carry in table[2](5+20) then
stevep 0:04a9f72bbca7 221 this carry would be passed to the table[3] etc.)
stevep 0:04a9f72bbca7 222 */
stevep 0:04a9f72bbca7 223 template<uint value_size>
stevep 0:04a9f72bbca7 224 uint UInt<value_size>::AddTwoInts(uint x2, uint x1, uint index)
stevep 0:04a9f72bbca7 225 {
stevep 0:04a9f72bbca7 226 uint i, c;
stevep 0:04a9f72bbca7 227
stevep 0:04a9f72bbca7 228 TTMATH_ASSERT( index < value_size - 1 )
stevep 0:04a9f72bbca7 229
stevep 0:04a9f72bbca7 230
stevep 0:04a9f72bbca7 231 c = AddTwoWords(table[index], x1, 0, &table[index]);
stevep 0:04a9f72bbca7 232 c = AddTwoWords(table[index+1], x2, c, &table[index+1]);
stevep 0:04a9f72bbca7 233
stevep 0:04a9f72bbca7 234 for(i=index+2 ; i<value_size && c ; ++i)
stevep 0:04a9f72bbca7 235 c = AddTwoWords(table[i], 0, c, &table[i]);
stevep 0:04a9f72bbca7 236
stevep 0:04a9f72bbca7 237 TTMATH_LOGC("UInt::AddTwoInts", c)
stevep 0:04a9f72bbca7 238
stevep 0:04a9f72bbca7 239 return c;
stevep 0:04a9f72bbca7 240 }
stevep 0:04a9f72bbca7 241
stevep 0:04a9f72bbca7 242
stevep 0:04a9f72bbca7 243
stevep 0:04a9f72bbca7 244 /*!
stevep 0:04a9f72bbca7 245 this static method addes one vector to the other
stevep 0:04a9f72bbca7 246 'ss1' is larger in size or equal to 'ss2'
stevep 0:04a9f72bbca7 247
stevep 0:04a9f72bbca7 248 ss1 points to the first (larger) vector
stevep 0:04a9f72bbca7 249 ss2 points to the second vector
stevep 0:04a9f72bbca7 250 ss1_size - size of the ss1 (and size of the result too)
stevep 0:04a9f72bbca7 251 ss2_size - size of the ss2
stevep 0:04a9f72bbca7 252 result - is the result vector (which has size the same as ss1: ss1_size)
stevep 0:04a9f72bbca7 253
stevep 0:04a9f72bbca7 254 Example: ss1_size is 5, ss2_size is 3
stevep 0:04a9f72bbca7 255 ss1: ss2: result (output):
stevep 0:04a9f72bbca7 256 5 1 5+1
stevep 0:04a9f72bbca7 257 4 3 4+3
stevep 0:04a9f72bbca7 258 2 7 2+7
stevep 0:04a9f72bbca7 259 6 6
stevep 0:04a9f72bbca7 260 9 9
stevep 0:04a9f72bbca7 261 of course the carry is propagated and will be returned from the last item
stevep 0:04a9f72bbca7 262 (this method is used by the Karatsuba multiplication algorithm)
stevep 0:04a9f72bbca7 263 */
stevep 0:04a9f72bbca7 264 template<uint value_size>
stevep 0:04a9f72bbca7 265 uint UInt<value_size>::AddVector(const uint * ss1, const uint * ss2, uint ss1_size, uint ss2_size, uint * result)
stevep 0:04a9f72bbca7 266 {
stevep 0:04a9f72bbca7 267 uint i, c = 0;
stevep 0:04a9f72bbca7 268
stevep 0:04a9f72bbca7 269 TTMATH_ASSERT( ss1_size >= ss2_size )
stevep 0:04a9f72bbca7 270
stevep 0:04a9f72bbca7 271 for(i=0 ; i<ss2_size ; ++i)
stevep 0:04a9f72bbca7 272 c = AddTwoWords(ss1[i], ss2[i], c, &result[i]);
stevep 0:04a9f72bbca7 273
stevep 0:04a9f72bbca7 274 for( ; i<ss1_size ; ++i)
stevep 0:04a9f72bbca7 275 c = AddTwoWords(ss1[i], 0, c, &result[i]);
stevep 0:04a9f72bbca7 276
stevep 0:04a9f72bbca7 277 TTMATH_VECTOR_LOGC("UInt::AddVector", c, result, ss1_size)
stevep 0:04a9f72bbca7 278
stevep 0:04a9f72bbca7 279 return c;
stevep 0:04a9f72bbca7 280 }
stevep 0:04a9f72bbca7 281
stevep 0:04a9f72bbca7 282
stevep 0:04a9f72bbca7 283
stevep 0:04a9f72bbca7 284
stevep 0:04a9f72bbca7 285 /*!
stevep 0:04a9f72bbca7 286 this method subtractes one word from the other
stevep 0:04a9f72bbca7 287 returns carry
stevep 0:04a9f72bbca7 288
stevep 0:04a9f72bbca7 289 this method is created only when TTMATH_NOASM macro is defined
stevep 0:04a9f72bbca7 290 */
stevep 0:04a9f72bbca7 291 template<uint value_size>
stevep 0:04a9f72bbca7 292 uint UInt<value_size>::SubTwoWords(uint a, uint b, uint carry, uint * result)
stevep 0:04a9f72bbca7 293 {
stevep 0:04a9f72bbca7 294 if( carry == 0 )
stevep 0:04a9f72bbca7 295 {
stevep 0:04a9f72bbca7 296 *result = a - b;
stevep 0:04a9f72bbca7 297
stevep 0:04a9f72bbca7 298 if( a < b )
stevep 0:04a9f72bbca7 299 carry = 1;
stevep 0:04a9f72bbca7 300 }
stevep 0:04a9f72bbca7 301 else
stevep 0:04a9f72bbca7 302 {
stevep 0:04a9f72bbca7 303 carry = 1;
stevep 0:04a9f72bbca7 304 *result = a - b - carry;
stevep 0:04a9f72bbca7 305
stevep 0:04a9f72bbca7 306 if( a > b ) // !(a <= b )
stevep 0:04a9f72bbca7 307 carry = 0;
stevep 0:04a9f72bbca7 308 }
stevep 0:04a9f72bbca7 309
stevep 0:04a9f72bbca7 310 return carry;
stevep 0:04a9f72bbca7 311 }
stevep 0:04a9f72bbca7 312
stevep 0:04a9f72bbca7 313
stevep 0:04a9f72bbca7 314
stevep 0:04a9f72bbca7 315
stevep 0:04a9f72bbca7 316 /*!
stevep 0:04a9f72bbca7 317 this method's subtracting ss2 from the 'this' and subtracting
stevep 0:04a9f72bbca7 318 carry if it has been defined
stevep 0:04a9f72bbca7 319 (this = this - ss2 - c)
stevep 0:04a9f72bbca7 320
stevep 0:04a9f72bbca7 321 c must be zero or one (might be a bigger value than 1)
stevep 0:04a9f72bbca7 322 function returns carry (1) (if it was)
stevep 0:04a9f72bbca7 323 */
stevep 0:04a9f72bbca7 324 template<uint value_size>
stevep 0:04a9f72bbca7 325 uint UInt<value_size>::Sub(const UInt<value_size> & ss2, uint c)
stevep 0:04a9f72bbca7 326 {
stevep 0:04a9f72bbca7 327 uint i;
stevep 0:04a9f72bbca7 328
stevep 0:04a9f72bbca7 329 for(i=0 ; i<value_size ; ++i)
stevep 0:04a9f72bbca7 330 c = SubTwoWords(table[i], ss2.table[i], c, &table[i]);
stevep 0:04a9f72bbca7 331
stevep 0:04a9f72bbca7 332 TTMATH_LOGC("UInt::Sub", c)
stevep 0:04a9f72bbca7 333
stevep 0:04a9f72bbca7 334 return c;
stevep 0:04a9f72bbca7 335 }
stevep 0:04a9f72bbca7 336
stevep 0:04a9f72bbca7 337
stevep 0:04a9f72bbca7 338
stevep 0:04a9f72bbca7 339
stevep 0:04a9f72bbca7 340 /*!
stevep 0:04a9f72bbca7 341 this method subtracts one word (at a specific position)
stevep 0:04a9f72bbca7 342 and returns a carry (if it was)
stevep 0:04a9f72bbca7 343
stevep 0:04a9f72bbca7 344 if we've got (value_size=3):
stevep 0:04a9f72bbca7 345 table[0] = 10;
stevep 0:04a9f72bbca7 346 table[1] = 30;
stevep 0:04a9f72bbca7 347 table[2] = 5;
stevep 0:04a9f72bbca7 348 and we call:
stevep 0:04a9f72bbca7 349 SubInt(2,1)
stevep 0:04a9f72bbca7 350 then it'll be:
stevep 0:04a9f72bbca7 351 table[0] = 10;
stevep 0:04a9f72bbca7 352 table[1] = 30 - 2;
stevep 0:04a9f72bbca7 353 table[2] = 5;
stevep 0:04a9f72bbca7 354
stevep 0:04a9f72bbca7 355 of course if there was a carry from table[2] it would be returned
stevep 0:04a9f72bbca7 356 */
stevep 0:04a9f72bbca7 357 template<uint value_size>
stevep 0:04a9f72bbca7 358 uint UInt<value_size>::SubInt(uint value, uint index)
stevep 0:04a9f72bbca7 359 {
stevep 0:04a9f72bbca7 360 uint i, c;
stevep 0:04a9f72bbca7 361
stevep 0:04a9f72bbca7 362 TTMATH_ASSERT( index < value_size )
stevep 0:04a9f72bbca7 363
stevep 0:04a9f72bbca7 364
stevep 0:04a9f72bbca7 365 c = SubTwoWords(table[index], value, 0, &table[index]);
stevep 0:04a9f72bbca7 366
stevep 0:04a9f72bbca7 367 for(i=index+1 ; i<value_size && c ; ++i)
stevep 0:04a9f72bbca7 368 c = SubTwoWords(table[i], 0, c, &table[i]);
stevep 0:04a9f72bbca7 369
stevep 0:04a9f72bbca7 370 TTMATH_LOGC("UInt::SubInt", c)
stevep 0:04a9f72bbca7 371
stevep 0:04a9f72bbca7 372 return c;
stevep 0:04a9f72bbca7 373 }
stevep 0:04a9f72bbca7 374
stevep 0:04a9f72bbca7 375
stevep 0:04a9f72bbca7 376 /*!
stevep 0:04a9f72bbca7 377 this static method subtractes one vector from the other
stevep 0:04a9f72bbca7 378 'ss1' is larger in size or equal to 'ss2'
stevep 0:04a9f72bbca7 379
stevep 0:04a9f72bbca7 380 ss1 points to the first (larger) vector
stevep 0:04a9f72bbca7 381 ss2 points to the second vector
stevep 0:04a9f72bbca7 382 ss1_size - size of the ss1 (and size of the result too)
stevep 0:04a9f72bbca7 383 ss2_size - size of the ss2
stevep 0:04a9f72bbca7 384 result - is the result vector (which has size the same as ss1: ss1_size)
stevep 0:04a9f72bbca7 385
stevep 0:04a9f72bbca7 386 Example: ss1_size is 5, ss2_size is 3
stevep 0:04a9f72bbca7 387 ss1: ss2: result (output):
stevep 0:04a9f72bbca7 388 5 1 5-1
stevep 0:04a9f72bbca7 389 4 3 4-3
stevep 0:04a9f72bbca7 390 2 7 2-7
stevep 0:04a9f72bbca7 391 6 6-1 (the borrow from previous item)
stevep 0:04a9f72bbca7 392 9 9
stevep 0:04a9f72bbca7 393 return (carry): 0
stevep 0:04a9f72bbca7 394 of course the carry (borrow) is propagated and will be returned from the last item
stevep 0:04a9f72bbca7 395 (this method is used by the Karatsuba multiplication algorithm)
stevep 0:04a9f72bbca7 396 */
stevep 0:04a9f72bbca7 397 template<uint value_size>
stevep 0:04a9f72bbca7 398 uint UInt<value_size>::SubVector(const uint * ss1, const uint * ss2, uint ss1_size, uint ss2_size, uint * result)
stevep 0:04a9f72bbca7 399 {
stevep 0:04a9f72bbca7 400 uint i, c = 0;
stevep 0:04a9f72bbca7 401
stevep 0:04a9f72bbca7 402 TTMATH_ASSERT( ss1_size >= ss2_size )
stevep 0:04a9f72bbca7 403
stevep 0:04a9f72bbca7 404 for(i=0 ; i<ss2_size ; ++i)
stevep 0:04a9f72bbca7 405 c = SubTwoWords(ss1[i], ss2[i], c, &result[i]);
stevep 0:04a9f72bbca7 406
stevep 0:04a9f72bbca7 407 for( ; i<ss1_size ; ++i)
stevep 0:04a9f72bbca7 408 c = SubTwoWords(ss1[i], 0, c, &result[i]);
stevep 0:04a9f72bbca7 409
stevep 0:04a9f72bbca7 410 TTMATH_VECTOR_LOGC("UInt::SubVector", c, result, ss1_size)
stevep 0:04a9f72bbca7 411
stevep 0:04a9f72bbca7 412 return c;
stevep 0:04a9f72bbca7 413 }
stevep 0:04a9f72bbca7 414
stevep 0:04a9f72bbca7 415
stevep 0:04a9f72bbca7 416
stevep 0:04a9f72bbca7 417 /*!
stevep 0:04a9f72bbca7 418 this method moves all bits into the left hand side
stevep 0:04a9f72bbca7 419 return value <- this <- c
stevep 0:04a9f72bbca7 420
stevep 0:04a9f72bbca7 421 the lowest *bit* will be held the 'c' and
stevep 0:04a9f72bbca7 422 the state of one additional bit (on the left hand side)
stevep 0:04a9f72bbca7 423 will be returned
stevep 0:04a9f72bbca7 424
stevep 0:04a9f72bbca7 425 for example:
stevep 0:04a9f72bbca7 426 let this is 001010000
stevep 0:04a9f72bbca7 427 after Rcl2_one(1) there'll be 010100001 and Rcl2_one returns 0
stevep 0:04a9f72bbca7 428 */
stevep 0:04a9f72bbca7 429 template<uint value_size>
stevep 0:04a9f72bbca7 430 uint UInt<value_size>::Rcl2_one(uint c)
stevep 0:04a9f72bbca7 431 {
stevep 0:04a9f72bbca7 432 uint i, new_c;
stevep 0:04a9f72bbca7 433
stevep 0:04a9f72bbca7 434 if( c != 0 )
stevep 0:04a9f72bbca7 435 c = 1;
stevep 0:04a9f72bbca7 436
stevep 0:04a9f72bbca7 437 for(i=0 ; i<value_size ; ++i)
stevep 0:04a9f72bbca7 438 {
stevep 0:04a9f72bbca7 439 new_c = (table[i] & TTMATH_UINT_HIGHEST_BIT) ? 1 : 0;
stevep 0:04a9f72bbca7 440 table[i] = (table[i] << 1) | c;
stevep 0:04a9f72bbca7 441 c = new_c;
stevep 0:04a9f72bbca7 442 }
stevep 0:04a9f72bbca7 443
stevep 0:04a9f72bbca7 444 TTMATH_LOGC("UInt::Rcl2_one", c)
stevep 0:04a9f72bbca7 445
stevep 0:04a9f72bbca7 446 return c;
stevep 0:04a9f72bbca7 447 }
stevep 0:04a9f72bbca7 448
stevep 0:04a9f72bbca7 449
stevep 0:04a9f72bbca7 450
stevep 0:04a9f72bbca7 451
stevep 0:04a9f72bbca7 452
stevep 0:04a9f72bbca7 453
stevep 0:04a9f72bbca7 454
stevep 0:04a9f72bbca7 455 /*!
stevep 0:04a9f72bbca7 456 this method moves all bits into the right hand side
stevep 0:04a9f72bbca7 457 c -> this -> return value
stevep 0:04a9f72bbca7 458
stevep 0:04a9f72bbca7 459 the highest *bit* will be held the 'c' and
stevep 0:04a9f72bbca7 460 the state of one additional bit (on the right hand side)
stevep 0:04a9f72bbca7 461 will be returned
stevep 0:04a9f72bbca7 462
stevep 0:04a9f72bbca7 463 for example:
stevep 0:04a9f72bbca7 464 let this is 000000010
stevep 0:04a9f72bbca7 465 after Rcr2_one(1) there'll be 100000001 and Rcr2_one returns 0
stevep 0:04a9f72bbca7 466 */
stevep 0:04a9f72bbca7 467 template<uint value_size>
stevep 0:04a9f72bbca7 468 uint UInt<value_size>::Rcr2_one(uint c)
stevep 0:04a9f72bbca7 469 {
stevep 0:04a9f72bbca7 470 sint i; // signed i
stevep 0:04a9f72bbca7 471 uint new_c;
stevep 0:04a9f72bbca7 472
stevep 0:04a9f72bbca7 473 if( c != 0 )
stevep 0:04a9f72bbca7 474 c = TTMATH_UINT_HIGHEST_BIT;
stevep 0:04a9f72bbca7 475
stevep 0:04a9f72bbca7 476 for(i=sint(value_size)-1 ; i>=0 ; --i)
stevep 0:04a9f72bbca7 477 {
stevep 0:04a9f72bbca7 478 new_c = (table[i] & 1) ? TTMATH_UINT_HIGHEST_BIT : 0;
stevep 0:04a9f72bbca7 479 table[i] = (table[i] >> 1) | c;
stevep 0:04a9f72bbca7 480 c = new_c;
stevep 0:04a9f72bbca7 481 }
stevep 0:04a9f72bbca7 482
stevep 0:04a9f72bbca7 483 c = (c != 0)? 1 : 0;
stevep 0:04a9f72bbca7 484
stevep 0:04a9f72bbca7 485 TTMATH_LOGC("UInt::Rcr2_one", c)
stevep 0:04a9f72bbca7 486
stevep 0:04a9f72bbca7 487 return c;
stevep 0:04a9f72bbca7 488 }
stevep 0:04a9f72bbca7 489
stevep 0:04a9f72bbca7 490
stevep 0:04a9f72bbca7 491
stevep 0:04a9f72bbca7 492
stevep 0:04a9f72bbca7 493 /*!
stevep 0:04a9f72bbca7 494 this method moves all bits into the left hand side
stevep 0:04a9f72bbca7 495 return value <- this <- c
stevep 0:04a9f72bbca7 496
stevep 0:04a9f72bbca7 497 the lowest *bits* will be held the 'c' and
stevep 0:04a9f72bbca7 498 the state of one additional bit (on the left hand side)
stevep 0:04a9f72bbca7 499 will be returned
stevep 0:04a9f72bbca7 500
stevep 0:04a9f72bbca7 501 for example:
stevep 0:04a9f72bbca7 502 let this is 001010000
stevep 0:04a9f72bbca7 503 after Rcl2(3, 1) there'll be 010000111 and Rcl2 returns 1
stevep 0:04a9f72bbca7 504 */
stevep 0:04a9f72bbca7 505 template<uint value_size>
stevep 0:04a9f72bbca7 506 uint UInt<value_size>::Rcl2(uint bits, uint c)
stevep 0:04a9f72bbca7 507 {
stevep 0:04a9f72bbca7 508 TTMATH_ASSERT( bits>0 && bits<TTMATH_BITS_PER_UINT )
stevep 0:04a9f72bbca7 509
stevep 0:04a9f72bbca7 510 uint move = TTMATH_BITS_PER_UINT - bits;
stevep 0:04a9f72bbca7 511 uint i, new_c;
stevep 0:04a9f72bbca7 512
stevep 0:04a9f72bbca7 513 if( c != 0 )
stevep 0:04a9f72bbca7 514 c = TTMATH_UINT_MAX_VALUE >> move;
stevep 0:04a9f72bbca7 515
stevep 0:04a9f72bbca7 516 for(i=0 ; i<value_size ; ++i)
stevep 0:04a9f72bbca7 517 {
stevep 0:04a9f72bbca7 518 new_c = table[i] >> move;
stevep 0:04a9f72bbca7 519 table[i] = (table[i] << bits) | c;
stevep 0:04a9f72bbca7 520 c = new_c;
stevep 0:04a9f72bbca7 521 }
stevep 0:04a9f72bbca7 522
stevep 0:04a9f72bbca7 523 TTMATH_LOGC("UInt::Rcl2", (c & 1))
stevep 0:04a9f72bbca7 524
stevep 0:04a9f72bbca7 525 return (c & 1);
stevep 0:04a9f72bbca7 526 }
stevep 0:04a9f72bbca7 527
stevep 0:04a9f72bbca7 528
stevep 0:04a9f72bbca7 529
stevep 0:04a9f72bbca7 530
stevep 0:04a9f72bbca7 531 /*!
stevep 0:04a9f72bbca7 532 this method moves all bits into the right hand side
stevep 0:04a9f72bbca7 533 C -> this -> return value
stevep 0:04a9f72bbca7 534
stevep 0:04a9f72bbca7 535 the highest *bits* will be held the 'c' and
stevep 0:04a9f72bbca7 536 the state of one additional bit (on the right hand side)
stevep 0:04a9f72bbca7 537 will be returned
stevep 0:04a9f72bbca7 538
stevep 0:04a9f72bbca7 539 for example:
stevep 0:04a9f72bbca7 540 let this is 000000010
stevep 0:04a9f72bbca7 541 after Rcr2(2, 1) there'll be 110000000 and Rcr2 returns 1
stevep 0:04a9f72bbca7 542 */
stevep 0:04a9f72bbca7 543 template<uint value_size>
stevep 0:04a9f72bbca7 544 uint UInt<value_size>::Rcr2(uint bits, uint c)
stevep 0:04a9f72bbca7 545 {
stevep 0:04a9f72bbca7 546 TTMATH_ASSERT( bits>0 && bits<TTMATH_BITS_PER_UINT )
stevep 0:04a9f72bbca7 547
stevep 0:04a9f72bbca7 548 uint move = TTMATH_BITS_PER_UINT - bits;
stevep 0:04a9f72bbca7 549 sint i; // signed
stevep 0:04a9f72bbca7 550 uint new_c;
stevep 0:04a9f72bbca7 551
stevep 0:04a9f72bbca7 552 if( c != 0 )
stevep 0:04a9f72bbca7 553 c = TTMATH_UINT_MAX_VALUE << move;
stevep 0:04a9f72bbca7 554
stevep 0:04a9f72bbca7 555 for(i=value_size-1 ; i>=0 ; --i)
stevep 0:04a9f72bbca7 556 {
stevep 0:04a9f72bbca7 557 new_c = table[i] << move;
stevep 0:04a9f72bbca7 558 table[i] = (table[i] >> bits) | c;
stevep 0:04a9f72bbca7 559 c = new_c;
stevep 0:04a9f72bbca7 560 }
stevep 0:04a9f72bbca7 561
stevep 0:04a9f72bbca7 562 c = (c & TTMATH_UINT_HIGHEST_BIT) ? 1 : 0;
stevep 0:04a9f72bbca7 563
stevep 0:04a9f72bbca7 564 TTMATH_LOGC("UInt::Rcr2", c)
stevep 0:04a9f72bbca7 565
stevep 0:04a9f72bbca7 566 return c;
stevep 0:04a9f72bbca7 567 }
stevep 0:04a9f72bbca7 568
stevep 0:04a9f72bbca7 569
stevep 0:04a9f72bbca7 570
stevep 0:04a9f72bbca7 571
stevep 0:04a9f72bbca7 572 /*!
stevep 0:04a9f72bbca7 573 this method returns the number of the highest set bit in x
stevep 0:04a9f72bbca7 574 if the 'x' is zero this method returns '-1'
stevep 0:04a9f72bbca7 575 */
stevep 0:04a9f72bbca7 576 template<uint value_size>
stevep 0:04a9f72bbca7 577 sint UInt<value_size>::FindLeadingBitInWord(uint x)
stevep 0:04a9f72bbca7 578 {
stevep 0:04a9f72bbca7 579 if( x == 0 )
stevep 0:04a9f72bbca7 580 return -1;
stevep 0:04a9f72bbca7 581
stevep 0:04a9f72bbca7 582 uint bit = TTMATH_BITS_PER_UINT - 1;
stevep 0:04a9f72bbca7 583
stevep 0:04a9f72bbca7 584 while( (x & TTMATH_UINT_HIGHEST_BIT) == 0 )
stevep 0:04a9f72bbca7 585 {
stevep 0:04a9f72bbca7 586 x = x << 1;
stevep 0:04a9f72bbca7 587 --bit;
stevep 0:04a9f72bbca7 588 }
stevep 0:04a9f72bbca7 589
stevep 0:04a9f72bbca7 590 return bit;
stevep 0:04a9f72bbca7 591 }
stevep 0:04a9f72bbca7 592
stevep 0:04a9f72bbca7 593
stevep 0:04a9f72bbca7 594
stevep 0:04a9f72bbca7 595 /*!
stevep 0:04a9f72bbca7 596 this method returns the number of the highest set bit in x
stevep 0:04a9f72bbca7 597 if the 'x' is zero this method returns '-1'
stevep 0:04a9f72bbca7 598 */
stevep 0:04a9f72bbca7 599 template<uint value_size>
stevep 0:04a9f72bbca7 600 sint UInt<value_size>::FindLowestBitInWord(uint x)
stevep 0:04a9f72bbca7 601 {
stevep 0:04a9f72bbca7 602 if( x == 0 )
stevep 0:04a9f72bbca7 603 return -1;
stevep 0:04a9f72bbca7 604
stevep 0:04a9f72bbca7 605 uint bit = 0;
stevep 0:04a9f72bbca7 606
stevep 0:04a9f72bbca7 607 while( (x & 1) == 0 )
stevep 0:04a9f72bbca7 608 {
stevep 0:04a9f72bbca7 609 x = x >> 1;
stevep 0:04a9f72bbca7 610 ++bit;
stevep 0:04a9f72bbca7 611 }
stevep 0:04a9f72bbca7 612
stevep 0:04a9f72bbca7 613 return bit;
stevep 0:04a9f72bbca7 614 }
stevep 0:04a9f72bbca7 615
stevep 0:04a9f72bbca7 616
stevep 0:04a9f72bbca7 617
stevep 0:04a9f72bbca7 618 /*!
stevep 0:04a9f72bbca7 619 this method sets a special bit in the 'value'
stevep 0:04a9f72bbca7 620 and returns the last state of the bit (zero or one)
stevep 0:04a9f72bbca7 621
stevep 0:04a9f72bbca7 622 bit is from <0,TTMATH_BITS_PER_UINT-1>
stevep 0:04a9f72bbca7 623
stevep 0:04a9f72bbca7 624 e.g.
stevep 0:04a9f72bbca7 625 uint x = 100;
stevep 0:04a9f72bbca7 626 uint bit = SetBitInWord(x, 3);
stevep 0:04a9f72bbca7 627 now: x = 108 and bit = 0
stevep 0:04a9f72bbca7 628 */
stevep 0:04a9f72bbca7 629 template<uint value_size>
stevep 0:04a9f72bbca7 630 uint UInt<value_size>::SetBitInWord(uint & value, uint bit)
stevep 0:04a9f72bbca7 631 {
stevep 0:04a9f72bbca7 632 TTMATH_ASSERT( bit < TTMATH_BITS_PER_UINT )
stevep 0:04a9f72bbca7 633
stevep 0:04a9f72bbca7 634 uint mask = 1;
stevep 0:04a9f72bbca7 635
stevep 0:04a9f72bbca7 636 if( bit > 0 )
stevep 0:04a9f72bbca7 637 mask = mask << bit;
stevep 0:04a9f72bbca7 638
stevep 0:04a9f72bbca7 639 uint last = value & mask;
stevep 0:04a9f72bbca7 640 value = value | mask;
stevep 0:04a9f72bbca7 641
stevep 0:04a9f72bbca7 642 return (last != 0) ? 1 : 0;
stevep 0:04a9f72bbca7 643 }
stevep 0:04a9f72bbca7 644
stevep 0:04a9f72bbca7 645
stevep 0:04a9f72bbca7 646
stevep 0:04a9f72bbca7 647
stevep 0:04a9f72bbca7 648
stevep 0:04a9f72bbca7 649
stevep 0:04a9f72bbca7 650 /*!
stevep 0:04a9f72bbca7 651 *
stevep 0:04a9f72bbca7 652 * Multiplication
stevep 0:04a9f72bbca7 653 *
stevep 0:04a9f72bbca7 654 *
stevep 0:04a9f72bbca7 655 */
stevep 0:04a9f72bbca7 656
stevep 0:04a9f72bbca7 657
stevep 0:04a9f72bbca7 658 /*!
stevep 0:04a9f72bbca7 659 multiplication: result_high:result_low = a * b
stevep 0:04a9f72bbca7 660 result_high - higher word of the result
stevep 0:04a9f72bbca7 661 result_low - lower word of the result
stevep 0:04a9f72bbca7 662
stevep 0:04a9f72bbca7 663 this methos never returns a carry
stevep 0:04a9f72bbca7 664 this method is used in the second version of the multiplication algorithms
stevep 0:04a9f72bbca7 665 */
stevep 0:04a9f72bbca7 666 template<uint value_size>
stevep 0:04a9f72bbca7 667 void UInt<value_size>::MulTwoWords(uint a, uint b, uint * result_high, uint * result_low)
stevep 0:04a9f72bbca7 668 {
stevep 0:04a9f72bbca7 669 #ifdef TTMATH_PLATFORM32
stevep 0:04a9f72bbca7 670
stevep 0:04a9f72bbca7 671 /*
stevep 0:04a9f72bbca7 672 on 32bit platforms we have defined 'unsigned long long int' type known as 'ulint' in ttmath namespace
stevep 0:04a9f72bbca7 673 this type has 64 bits, then we're using only one multiplication: 32bit * 32bit = 64bit
stevep 0:04a9f72bbca7 674 */
stevep 0:04a9f72bbca7 675
stevep 0:04a9f72bbca7 676 union uint_
stevep 0:04a9f72bbca7 677 {
stevep 0:04a9f72bbca7 678 struct
stevep 0:04a9f72bbca7 679 {
stevep 0:04a9f72bbca7 680 uint low; // 32 bits
stevep 0:04a9f72bbca7 681 uint high; // 32 bits
stevep 0:04a9f72bbca7 682 } u_;
stevep 0:04a9f72bbca7 683
stevep 0:04a9f72bbca7 684 ulint u; // 64 bits
stevep 0:04a9f72bbca7 685 } res;
stevep 0:04a9f72bbca7 686
stevep 0:04a9f72bbca7 687 res.u = ulint(a) * ulint(b); // multiply two 32bit words, the result has 64 bits
stevep 0:04a9f72bbca7 688
stevep 0:04a9f72bbca7 689 *result_high = res.u_.high;
stevep 0:04a9f72bbca7 690 *result_low = res.u_.low;
stevep 0:04a9f72bbca7 691
stevep 0:04a9f72bbca7 692 #else
stevep 0:04a9f72bbca7 693
stevep 0:04a9f72bbca7 694 /*
stevep 0:04a9f72bbca7 695 64 bits platforms
stevep 0:04a9f72bbca7 696
stevep 0:04a9f72bbca7 697 we don't have a native type which has 128 bits
stevep 0:04a9f72bbca7 698 then we're splitting 'a' and 'b' to 4 parts (high and low halves)
stevep 0:04a9f72bbca7 699 and using 4 multiplications (with additions and carry correctness)
stevep 0:04a9f72bbca7 700 */
stevep 0:04a9f72bbca7 701
stevep 0:04a9f72bbca7 702 uint_ a_;
stevep 0:04a9f72bbca7 703 uint_ b_;
stevep 0:04a9f72bbca7 704 uint_ res_high1, res_high2;
stevep 0:04a9f72bbca7 705 uint_ res_low1, res_low2;
stevep 0:04a9f72bbca7 706
stevep 0:04a9f72bbca7 707 a_.u = a;
stevep 0:04a9f72bbca7 708 b_.u = b;
stevep 0:04a9f72bbca7 709
stevep 0:04a9f72bbca7 710 /*
stevep 0:04a9f72bbca7 711 the multiplication is as follows (schoolbook algorithm with O(n^2) ):
stevep 0:04a9f72bbca7 712
stevep 0:04a9f72bbca7 713 32 bits 32 bits
stevep 0:04a9f72bbca7 714
stevep 0:04a9f72bbca7 715 +--------------------------------+
stevep 0:04a9f72bbca7 716 | a_.u_.high | a_.u_.low |
stevep 0:04a9f72bbca7 717 +--------------------------------+
stevep 0:04a9f72bbca7 718 | b_.u_.high | b_.u_.low |
stevep 0:04a9f72bbca7 719 +--------------------------------+--------------------------------+
stevep 0:04a9f72bbca7 720 | res_high1.u | res_low1.u |
stevep 0:04a9f72bbca7 721 +--------------------------------+--------------------------------+
stevep 0:04a9f72bbca7 722 | res_high2.u | res_low2.u |
stevep 0:04a9f72bbca7 723 +--------------------------------+--------------------------------+
stevep 0:04a9f72bbca7 724
stevep 0:04a9f72bbca7 725 64 bits 64 bits
stevep 0:04a9f72bbca7 726 */
stevep 0:04a9f72bbca7 727
stevep 0:04a9f72bbca7 728
stevep 0:04a9f72bbca7 729 uint_ temp;
stevep 0:04a9f72bbca7 730
stevep 0:04a9f72bbca7 731 res_low1.u = uint(b_.u_.low) * uint(a_.u_.low);
stevep 0:04a9f72bbca7 732
stevep 0:04a9f72bbca7 733 temp.u = uint(res_low1.u_.high) + uint(b_.u_.low) * uint(a_.u_.high);
stevep 0:04a9f72bbca7 734 res_low1.u_.high = temp.u_.low;
stevep 0:04a9f72bbca7 735 res_high1.u_.low = temp.u_.high;
stevep 0:04a9f72bbca7 736 res_high1.u_.high = 0;
stevep 0:04a9f72bbca7 737
stevep 0:04a9f72bbca7 738 res_low2.u_.low = 0;
stevep 0:04a9f72bbca7 739 temp.u = uint(b_.u_.high) * uint(a_.u_.low);
stevep 0:04a9f72bbca7 740 res_low2.u_.high = temp.u_.low;
stevep 0:04a9f72bbca7 741
stevep 0:04a9f72bbca7 742 res_high2.u = uint(b_.u_.high) * uint(a_.u_.high) + uint(temp.u_.high);
stevep 0:04a9f72bbca7 743
stevep 0:04a9f72bbca7 744 uint c = AddTwoWords(res_low1.u, res_low2.u, 0, &res_low2.u);
stevep 0:04a9f72bbca7 745 AddTwoWords(res_high1.u, res_high2.u, c, &res_high2.u); // there is no carry from here
stevep 0:04a9f72bbca7 746
stevep 0:04a9f72bbca7 747 *result_high = res_high2.u;
stevep 0:04a9f72bbca7 748 *result_low = res_low2.u;
stevep 0:04a9f72bbca7 749
stevep 0:04a9f72bbca7 750 #endif
stevep 0:04a9f72bbca7 751 }
stevep 0:04a9f72bbca7 752
stevep 0:04a9f72bbca7 753
stevep 0:04a9f72bbca7 754
stevep 0:04a9f72bbca7 755
stevep 0:04a9f72bbca7 756 /*!
stevep 0:04a9f72bbca7 757 *
stevep 0:04a9f72bbca7 758 * Division
stevep 0:04a9f72bbca7 759 *
stevep 0:04a9f72bbca7 760 *
stevep 0:04a9f72bbca7 761 */
stevep 0:04a9f72bbca7 762
stevep 0:04a9f72bbca7 763
stevep 0:04a9f72bbca7 764 /*!
stevep 0:04a9f72bbca7 765 this method calculates 64bits word a:b / 32bits c (a higher, b lower word)
stevep 0:04a9f72bbca7 766 r = a:b / c and rest - remainder
stevep 0:04a9f72bbca7 767
stevep 0:04a9f72bbca7 768 *
stevep 0:04a9f72bbca7 769 * WARNING:
stevep 0:04a9f72bbca7 770 * the c has to be suitably large for the result being keeped in one word,
stevep 0:04a9f72bbca7 771 * if c is equal zero there'll be a hardware interruption (0)
stevep 0:04a9f72bbca7 772 * and probably the end of your program
stevep 0:04a9f72bbca7 773 *
stevep 0:04a9f72bbca7 774 */
stevep 0:04a9f72bbca7 775 template<uint value_size>
stevep 0:04a9f72bbca7 776 void UInt<value_size>::DivTwoWords(uint a, uint b, uint c, uint * r, uint * rest)
stevep 0:04a9f72bbca7 777 {
stevep 0:04a9f72bbca7 778 // (a < c ) for the result to be one word
stevep 0:04a9f72bbca7 779 TTMATH_ASSERT( c != 0 && a < c )
stevep 0:04a9f72bbca7 780
stevep 0:04a9f72bbca7 781 #ifdef TTMATH_PLATFORM32
stevep 0:04a9f72bbca7 782
stevep 0:04a9f72bbca7 783 union
stevep 0:04a9f72bbca7 784 {
stevep 0:04a9f72bbca7 785 struct
stevep 0:04a9f72bbca7 786 {
stevep 0:04a9f72bbca7 787 uint low; // 32 bits
stevep 0:04a9f72bbca7 788 uint high; // 32 bits
stevep 0:04a9f72bbca7 789 } u_;
stevep 0:04a9f72bbca7 790
stevep 0:04a9f72bbca7 791 ulint u; // 64 bits
stevep 0:04a9f72bbca7 792 } ab;
stevep 0:04a9f72bbca7 793
stevep 0:04a9f72bbca7 794 ab.u_.high = a;
stevep 0:04a9f72bbca7 795 ab.u_.low = b;
stevep 0:04a9f72bbca7 796
stevep 0:04a9f72bbca7 797 *r = uint(ab.u / c);
stevep 0:04a9f72bbca7 798 *rest = uint(ab.u % c);
stevep 0:04a9f72bbca7 799
stevep 0:04a9f72bbca7 800 #else
stevep 0:04a9f72bbca7 801
stevep 0:04a9f72bbca7 802 uint_ c_;
stevep 0:04a9f72bbca7 803 c_.u = c;
stevep 0:04a9f72bbca7 804
stevep 0:04a9f72bbca7 805
stevep 0:04a9f72bbca7 806 if( a == 0 )
stevep 0:04a9f72bbca7 807 {
stevep 0:04a9f72bbca7 808 *r = b / c;
stevep 0:04a9f72bbca7 809 *rest = b % c;
stevep 0:04a9f72bbca7 810 }
stevep 0:04a9f72bbca7 811 else
stevep 0:04a9f72bbca7 812 if( c_.u_.high == 0 )
stevep 0:04a9f72bbca7 813 {
stevep 0:04a9f72bbca7 814 // higher half of 'c' is zero
stevep 0:04a9f72bbca7 815 // then higher half of 'a' is zero too (look at the asserts at the beginning - 'a' is smaller than 'c')
stevep 0:04a9f72bbca7 816 uint_ a_, b_, res_, temp1, temp2;
stevep 0:04a9f72bbca7 817
stevep 0:04a9f72bbca7 818 a_.u = a;
stevep 0:04a9f72bbca7 819 b_.u = b;
stevep 0:04a9f72bbca7 820
stevep 0:04a9f72bbca7 821 temp1.u_.high = a_.u_.low;
stevep 0:04a9f72bbca7 822 temp1.u_.low = b_.u_.high;
stevep 0:04a9f72bbca7 823
stevep 0:04a9f72bbca7 824 res_.u_.high = (unsigned int)(temp1.u / c);
stevep 0:04a9f72bbca7 825 temp2.u_.high = (unsigned int)(temp1.u % c);
stevep 0:04a9f72bbca7 826 temp2.u_.low = b_.u_.low;
stevep 0:04a9f72bbca7 827
stevep 0:04a9f72bbca7 828 res_.u_.low = (unsigned int)(temp2.u / c);
stevep 0:04a9f72bbca7 829 *rest = temp2.u % c;
stevep 0:04a9f72bbca7 830
stevep 0:04a9f72bbca7 831 *r = res_.u;
stevep 0:04a9f72bbca7 832 }
stevep 0:04a9f72bbca7 833 else
stevep 0:04a9f72bbca7 834 {
stevep 0:04a9f72bbca7 835 return DivTwoWords2(a, b, c, r, rest);
stevep 0:04a9f72bbca7 836 }
stevep 0:04a9f72bbca7 837
stevep 0:04a9f72bbca7 838 #endif
stevep 0:04a9f72bbca7 839 }
stevep 0:04a9f72bbca7 840
stevep 0:04a9f72bbca7 841
stevep 0:04a9f72bbca7 842 #ifdef TTMATH_PLATFORM64
stevep 0:04a9f72bbca7 843
stevep 0:04a9f72bbca7 844
stevep 0:04a9f72bbca7 845 /*!
stevep 0:04a9f72bbca7 846 this method is available only on 64bit platforms
stevep 0:04a9f72bbca7 847
stevep 0:04a9f72bbca7 848 the same algorithm like the third division algorithm in ttmathuint.h
stevep 0:04a9f72bbca7 849 but now with the radix=2^32
stevep 0:04a9f72bbca7 850 */
stevep 0:04a9f72bbca7 851 template<uint value_size>
stevep 0:04a9f72bbca7 852 void UInt<value_size>::DivTwoWords2(uint a, uint b, uint c, uint * r, uint * rest)
stevep 0:04a9f72bbca7 853 {
stevep 0:04a9f72bbca7 854 // a is not zero
stevep 0:04a9f72bbca7 855 // c_.u_.high is not zero
stevep 0:04a9f72bbca7 856
stevep 0:04a9f72bbca7 857 uint_ a_, b_, c_, u_, q_;
stevep 0:04a9f72bbca7 858 unsigned int u3; // 32 bit
stevep 0:04a9f72bbca7 859
stevep 0:04a9f72bbca7 860 a_.u = a;
stevep 0:04a9f72bbca7 861 b_.u = b;
stevep 0:04a9f72bbca7 862 c_.u = c;
stevep 0:04a9f72bbca7 863
stevep 0:04a9f72bbca7 864 // normalizing
stevep 0:04a9f72bbca7 865 uint d = DivTwoWordsNormalize(a_, b_, c_);
stevep 0:04a9f72bbca7 866
stevep 0:04a9f72bbca7 867 // loop from j=1 to j=0
stevep 0:04a9f72bbca7 868 // the first step (for j=2) is skipped because our result is only in one word,
stevep 0:04a9f72bbca7 869 // (first 'q' were 0 and nothing would be changed)
stevep 0:04a9f72bbca7 870 u_.u_.high = a_.u_.high;
stevep 0:04a9f72bbca7 871 u_.u_.low = a_.u_.low;
stevep 0:04a9f72bbca7 872 u3 = b_.u_.high;
stevep 0:04a9f72bbca7 873 q_.u_.high = DivTwoWordsCalculate(u_, u3, c_);
stevep 0:04a9f72bbca7 874 MultiplySubtract(u_, u3, q_.u_.high, c_);
stevep 0:04a9f72bbca7 875
stevep 0:04a9f72bbca7 876 u_.u_.high = u_.u_.low;
stevep 0:04a9f72bbca7 877 u_.u_.low = u3;
stevep 0:04a9f72bbca7 878 u3 = b_.u_.low;
stevep 0:04a9f72bbca7 879 q_.u_.low = DivTwoWordsCalculate(u_, u3, c_);
stevep 0:04a9f72bbca7 880 MultiplySubtract(u_, u3, q_.u_.low, c_);
stevep 0:04a9f72bbca7 881
stevep 0:04a9f72bbca7 882 *r = q_.u;
stevep 0:04a9f72bbca7 883
stevep 0:04a9f72bbca7 884 // unnormalizing for the remainder
stevep 0:04a9f72bbca7 885 u_.u_.high = u_.u_.low;
stevep 0:04a9f72bbca7 886 u_.u_.low = u3;
stevep 0:04a9f72bbca7 887 *rest = DivTwoWordsUnnormalize(u_.u, d);
stevep 0:04a9f72bbca7 888 }
stevep 0:04a9f72bbca7 889
stevep 0:04a9f72bbca7 890
stevep 0:04a9f72bbca7 891
stevep 0:04a9f72bbca7 892
stevep 0:04a9f72bbca7 893 template<uint value_size>
stevep 0:04a9f72bbca7 894 uint UInt<value_size>::DivTwoWordsNormalize(uint_ & a_, uint_ & b_, uint_ & c_)
stevep 0:04a9f72bbca7 895 {
stevep 0:04a9f72bbca7 896 uint d = 0;
stevep 0:04a9f72bbca7 897
stevep 0:04a9f72bbca7 898 for( ; (c_.u & TTMATH_UINT_HIGHEST_BIT) == 0 ; ++d )
stevep 0:04a9f72bbca7 899 {
stevep 0:04a9f72bbca7 900 c_.u = c_.u << 1;
stevep 0:04a9f72bbca7 901
stevep 0:04a9f72bbca7 902 uint bc = b_.u & TTMATH_UINT_HIGHEST_BIT; // carry from 'b'
stevep 0:04a9f72bbca7 903
stevep 0:04a9f72bbca7 904 b_.u = b_.u << 1;
stevep 0:04a9f72bbca7 905 a_.u = a_.u << 1; // carry bits from 'a' are simply skipped
stevep 0:04a9f72bbca7 906
stevep 0:04a9f72bbca7 907 if( bc )
stevep 0:04a9f72bbca7 908 a_.u = a_.u | 1;
stevep 0:04a9f72bbca7 909 }
stevep 0:04a9f72bbca7 910
stevep 0:04a9f72bbca7 911 return d;
stevep 0:04a9f72bbca7 912 }
stevep 0:04a9f72bbca7 913
stevep 0:04a9f72bbca7 914
stevep 0:04a9f72bbca7 915 template<uint value_size>
stevep 0:04a9f72bbca7 916 uint UInt<value_size>::DivTwoWordsUnnormalize(uint u, uint d)
stevep 0:04a9f72bbca7 917 {
stevep 0:04a9f72bbca7 918 if( d == 0 )
stevep 0:04a9f72bbca7 919 return u;
stevep 0:04a9f72bbca7 920
stevep 0:04a9f72bbca7 921 u = u >> d;
stevep 0:04a9f72bbca7 922
stevep 0:04a9f72bbca7 923 return u;
stevep 0:04a9f72bbca7 924 }
stevep 0:04a9f72bbca7 925
stevep 0:04a9f72bbca7 926
stevep 0:04a9f72bbca7 927 template<uint value_size>
stevep 0:04a9f72bbca7 928 unsigned int UInt<value_size>::DivTwoWordsCalculate(uint_ u_, unsigned int u3, uint_ v_)
stevep 0:04a9f72bbca7 929 {
stevep 0:04a9f72bbca7 930 bool next_test;
stevep 0:04a9f72bbca7 931 uint_ qp_, rp_, temp_;
stevep 0:04a9f72bbca7 932
stevep 0:04a9f72bbca7 933 qp_.u = u_.u / uint(v_.u_.high);
stevep 0:04a9f72bbca7 934 rp_.u = u_.u % uint(v_.u_.high);
stevep 0:04a9f72bbca7 935
stevep 0:04a9f72bbca7 936 TTMATH_ASSERT( qp_.u_.high==0 || qp_.u_.high==1 )
stevep 0:04a9f72bbca7 937
stevep 0:04a9f72bbca7 938 do
stevep 0:04a9f72bbca7 939 {
stevep 0:04a9f72bbca7 940 bool decrease = false;
stevep 0:04a9f72bbca7 941
stevep 0:04a9f72bbca7 942 if( qp_.u_.high == 1 )
stevep 0:04a9f72bbca7 943 decrease = true;
stevep 0:04a9f72bbca7 944 else
stevep 0:04a9f72bbca7 945 {
stevep 0:04a9f72bbca7 946 temp_.u_.high = rp_.u_.low;
stevep 0:04a9f72bbca7 947 temp_.u_.low = u3;
stevep 0:04a9f72bbca7 948
stevep 0:04a9f72bbca7 949 if( qp_.u * uint(v_.u_.low) > temp_.u )
stevep 0:04a9f72bbca7 950 decrease = true;
stevep 0:04a9f72bbca7 951 }
stevep 0:04a9f72bbca7 952
stevep 0:04a9f72bbca7 953 next_test = false;
stevep 0:04a9f72bbca7 954
stevep 0:04a9f72bbca7 955 if( decrease )
stevep 0:04a9f72bbca7 956 {
stevep 0:04a9f72bbca7 957 --qp_.u;
stevep 0:04a9f72bbca7 958 rp_.u += v_.u_.high;
stevep 0:04a9f72bbca7 959
stevep 0:04a9f72bbca7 960 if( rp_.u_.high == 0 )
stevep 0:04a9f72bbca7 961 next_test = true;
stevep 0:04a9f72bbca7 962 }
stevep 0:04a9f72bbca7 963 }
stevep 0:04a9f72bbca7 964 while( next_test );
stevep 0:04a9f72bbca7 965
stevep 0:04a9f72bbca7 966 return qp_.u_.low;
stevep 0:04a9f72bbca7 967 }
stevep 0:04a9f72bbca7 968
stevep 0:04a9f72bbca7 969
stevep 0:04a9f72bbca7 970 template<uint value_size>
stevep 0:04a9f72bbca7 971 void UInt<value_size>::MultiplySubtract(uint_ & u_, unsigned int & u3, unsigned int & q, uint_ v_)
stevep 0:04a9f72bbca7 972 {
stevep 0:04a9f72bbca7 973 uint_ temp_;
stevep 0:04a9f72bbca7 974
stevep 0:04a9f72bbca7 975 uint res_high;
stevep 0:04a9f72bbca7 976 uint res_low;
stevep 0:04a9f72bbca7 977
stevep 0:04a9f72bbca7 978 MulTwoWords(v_.u, q, &res_high, &res_low);
stevep 0:04a9f72bbca7 979
stevep 0:04a9f72bbca7 980 uint_ sub_res_high_;
stevep 0:04a9f72bbca7 981 uint_ sub_res_low_;
stevep 0:04a9f72bbca7 982
stevep 0:04a9f72bbca7 983 temp_.u_.high = u_.u_.low;
stevep 0:04a9f72bbca7 984 temp_.u_.low = u3;
stevep 0:04a9f72bbca7 985
stevep 0:04a9f72bbca7 986 uint c = SubTwoWords(temp_.u, res_low, 0, &sub_res_low_.u);
stevep 0:04a9f72bbca7 987
stevep 0:04a9f72bbca7 988 temp_.u_.high = 0;
stevep 0:04a9f72bbca7 989 temp_.u_.low = u_.u_.high;
stevep 0:04a9f72bbca7 990 c = SubTwoWords(temp_.u, res_high, c, &sub_res_high_.u);
stevep 0:04a9f72bbca7 991
stevep 0:04a9f72bbca7 992 if( c )
stevep 0:04a9f72bbca7 993 {
stevep 0:04a9f72bbca7 994 --q;
stevep 0:04a9f72bbca7 995
stevep 0:04a9f72bbca7 996 c = AddTwoWords(sub_res_low_.u, v_.u, 0, &sub_res_low_.u);
stevep 0:04a9f72bbca7 997 AddTwoWords(sub_res_high_.u, 0, c, &sub_res_high_.u);
stevep 0:04a9f72bbca7 998 }
stevep 0:04a9f72bbca7 999
stevep 0:04a9f72bbca7 1000 u_.u_.high = sub_res_high_.u_.low;
stevep 0:04a9f72bbca7 1001 u_.u_.low = sub_res_low_.u_.high;
stevep 0:04a9f72bbca7 1002 u3 = sub_res_low_.u_.low;
stevep 0:04a9f72bbca7 1003 }
stevep 0:04a9f72bbca7 1004
stevep 0:04a9f72bbca7 1005 #endif // #ifdef TTMATH_PLATFORM64
stevep 0:04a9f72bbca7 1006
stevep 0:04a9f72bbca7 1007
stevep 0:04a9f72bbca7 1008
stevep 0:04a9f72bbca7 1009 } //namespace
stevep 0:04a9f72bbca7 1010
stevep 0:04a9f72bbca7 1011
stevep 0:04a9f72bbca7 1012 #endif //ifdef TTMATH_NOASM
stevep 0:04a9f72bbca7 1013 #endif
stevep 0:04a9f72bbca7 1014
stevep 0:04a9f72bbca7 1015
stevep 0:04a9f72bbca7 1016
stevep 0:04a9f72bbca7 1017
stevep 0:04a9f72bbca7 1018