CCS C Software and Maintenance Offers
FAQFAQ   FAQForum Help   FAQOfficial CCS Support   SearchSearch  RegisterRegister 

ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

CCS does not monitor this forum on a regular basis.

Please do not post bug reports on this forum. Send them to CCS Technical Support

32 bit multiplication
Goto page 1, 2  Next
 
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion
View previous topic :: View next topic  
Author Message
stma



Joined: 16 Feb 2004
Posts: 26

View user's profile Send private message

32 bit multiplication
PostPosted: Wed Nov 08, 2006 5:00 am     Reply with quote

Hi,
Can anyone point me in the right direction for 32 bit multiplication help.
Obviously any 32bit multiplication will result in a 64 bit result.
Does some form of scaling need to be performed to produce a sensible result or is it split into two 32bit variables.
Thanks
Steve
Ttelmah
Guest







PostPosted: Wed Nov 08, 2006 5:51 am     Reply with quote

CCS, does not include any support for 64bit.
The existing 32bit multiplication routines, will simply overflow if the result is larger than 32bits.
You need to look at some of the standard libraries from MicroChip for example, and translate them yourself into a form that will work on CCS. What you would probably do is declare a structure containing 8 bytes, and 'typedef' this to be a int64. I did this in the past for an int24 type.

Best Wishes
stma



Joined: 16 Feb 2004
Posts: 26

View user's profile Send private message

PostPosted: Wed Nov 08, 2006 7:59 am     Reply with quote

Thanks very much for the reply. Will have a look at doing it that way.
Just done some scribblings and need only to multiply two 18 bit numbers initially.
Has anyone any ideas on scaling.
Im not a mathematician so keep it simple!
Thanks
Steve
Ken Johnson



Joined: 23 Mar 2006
Posts: 197
Location: Lewisburg, WV

View user's profile Send private message

PostPosted: Wed Nov 08, 2006 10:01 am     Reply with quote

To multiply 2 18-bit numbers, I'd just do c = a * b, where a, b, c are 32-bit. Unsigned would be better than signed, if they're always positive numbers.

Ken
Ttelmah
Guest







PostPosted: Wed Nov 08, 2006 4:44 pm     Reply with quote

Just for anyone who wants to start thinking about multiplying 32bit by 32bit, to give a 64bit result, I have generated the following:
Code:

union bits64 {
   int8 b[8];
   int32 w[2];
};
typedef union bits64 int64;

int64 add64x64(int64 a,int64 b) {
   //simple 64 bit addition
   #asm
   movf b.b[0],0
   addwf a.b[0],1
   movf b.b[1],0
   addwfc a.b[1],1
   movf b.b[2],0
   addwfc a.b[2],1
   movf b.b[3],0
   addwfc a.b[3],1
   movf b.b[4],0
   addwfc a.b[4],1
   movf b.b[5],0
   addwfc a.b[5],1
   movf b.b[5],0
   addwfc a.b[6],1
   movf b.b[6],0
   addwfc a.b[6],1
   movf b.b[7],0
   addwfc a.b[7],1
   #endasm
   return a;
}

int64 mult32x32(int32 a,int32 b) {
   int64 temp1,temp2;
   int8 i;
   //zero output
   for (i=0;i<8;i++)
      temp2.b[i]=0;
   //Routine to add two 64bit values stored LSB first
   //Not trying to break any records for efficiency, so all in C, and not using
   //hardware multiply
   //start by moving the int32, into the low half of the int64
   temp1.w[0]=b;
   //and clear the top half
   for (i=4;i<8;i++)
      temp1.b[i]=0;
   //Now need to work through all 32 bits in the 'a' variable
   for (i=0;i<32;i++) {
      //if source bit in 'a' is one, perform addition
      if (bit_test(a,i))
         temp2=add64x64(temp1,temp2);
      //rotate the second value here   
      shift_left(&temp1,8,0);
   }
}

Called using:

test=mult32x32(123456789,456789);

Where test is another int64 variable, correctly generates the byte pattern:
F9 6F A5 2E 4A 33

Which is the right value for 56393703190521.

Best Wishes
Ttelmah
Guest







PostPosted: Thu Nov 09, 2006 3:08 am     Reply with quote

Ken Johnson wrote:
To multiply 2 18-bit numbers, I'd just do c = a * b, where a, b, c are 32-bit. Unsigned would be better than signed, if they're always positive numbers.

Ken

Unfortunately, that doesn't solve the problem.
Multiplying 2 18bit numbers, in the 'worst case', potentially gives a 36bit result. Won't fit into a int32.

Best Wishes
Ken Johnson



Joined: 23 Mar 2006
Posts: 197
Location: Lewisburg, WV

View user's profile Send private message

PostPosted: Thu Nov 09, 2006 8:07 am     Reply with quote

18 plus 18 is more than 32 ! right, of course.

I had a funny feeling when I posted that, but it didn't sink in until sometime last night, just before I went to sleep . . .

Sorry 'bout that!
Ken
Ttelmah
Guest







PostPosted: Thu Nov 09, 2006 8:24 am     Reply with quote

It is one of those things that is not obvious, till you think about it, and is made 'worse' by having as standard a '32bit' arithmetic library, which just throws away the bits about 32...
Generally, for most applications, 32bits is 'big enough', but it is a problem sometimes waiting to sneak out and catch you!.
I wrote 32bit libraries years ago, for the Z80, and for those, had to deal with generating a potentially 64bit result for the multiplication. In fact the 'flowchart', is what I implemented in the code I posted. This could be made fairly efficient, if the 'scratch' was made global, so the value didn't have to keep being passed backwards/forwards to the addition routine, but even as posted, it'll work quite well. Not using the hardware multiply function, makes it portable to 16x, and 18x chips.
The 'extra' needed, is the division routine, since this will be required to print the result if required. Not terribly hard, but requires the generation of a 64bit 'compare' as well.

Best Wishes
SherpaDoug



Joined: 07 Sep 2003
Posts: 1640
Location: Cape Cod Mass USA

View user's profile Send private message

PostPosted: Thu Nov 09, 2006 1:39 pm     Reply with quote

I am curious, what is the application that rerquires more than 32 bits of resolution? I can maybe imagine some things in astronomy, nuclear physics, and maybe international finances, but none of them would be handled by a PIC! Are you sure you couldn't make do with an 18X14 multiply?
_________________
The search for better is endless. Instead simply find very good and get the job done.
Ttelmah
Guest







PostPosted: Thu Nov 09, 2006 5:00 pm     Reply with quote

Ongoing to my earlier post, though not 'elegant', I have assembled a set of basic routines for 64bit unsigned integer arithmetic.
Code:

union bits64 {
   int8 b[8];
   int32 w[2];
};
typedef union bits64 int64;
int64 scratch64;
int1 BORROW=FALSE;
#bit CARRY=0xFD8.0

//Macro to add two 64bit values. Result into the first
#define M_add64x64(M_a,M_b) \
   #asm\
   movf M_b.b[0],0\
   addwf M_a.b[0],1\
   movf M_b.b[1],0\
   addwfc M_a.b[1],1\
   movf M_b.b[2],0\
   addwfc M_a.b[2],1\
   movf M_b.b[3],0\
   addwfc M_a.b[3],1\
   movf M_b.b[4],0\
   addwfc M_a.b[4],1\
   movf M_b.b[5],0\
   addwfc M_a.b[5],1\
   movf M_b.b[6],0\
   addwfc M_a.b[6],1\
   movf M_b.b[7],0\
   addwfc M_a.b[7],1\
   #endasm
//Macro as above, to subtract two 64bit values.
#define M_sub64x64(M_a,M_b) \
   #asm\
   movf M_b.b[0],0\
   subwf M_a.b[0],1\
   movf M_b.b[1],0\
   subwfb M_a.b[1],1\
   movf M_b.b[2],0\
   subwfb M_a.b[2],1\
   movf M_b.b[3],0\
   subwfb M_a.b[3],1\
   movf M_b.b[4],0\
   subwfb M_a.b[4],1\
   movf M_b.b[5],0\
   subwfb M_a.b[5],1\
   movf M_b.b[6],0\
   subwfb M_a.b[6],1\
   movf M_b.b[7],0\
   subwfb M_a.b[7],1\
   #endasm

//Macros to shift a 64bit value
#define M_shiftleft64(x) shift_left(&x,8,0)
#define M_shiftright64(x) shift_right(&x,8,0)
//Macro to zero a 64bit value
#define M_zero64(x) x.w[0]=0L;x.w[1]=0L
#define M_iszero64(x) ((x.w[0]==0L)&&(x.w[1]==0L))
#define bit64_set(x,i) x.b[i>>3]|=(1<<(i&7))

int64 cast32x64(int32 a) {
   //routine to convert a 32bit int to 64bit
   int64 temp;
   int8 i;
   temp.w[0]=a;
   //and clear the top half
   temp.w[1]=0L;
   return temp;
}

int64 I_A_LT_B(int64 a,int64 b) {
   //Internal routine to _compare_ two 64 bit values
   //Actually performs subtraction, and returns this, with the 'borrow'
   //flag in the global variable 'BORROW'. Performs a-b, hence flag is
   //set if A less than B
   M_sub64x64(a,b);
   if (CARRY==0) BORROW=TRUE;
   else BORROW=FALSE;
   return a;
}

int64 add64x64(int64 a,int64 b) {
   //simple 64 bit addition call - returns (a+b)
   M_add64x64(a,b);
   return a;
}

int64 sub64x64(int64 a,int64 b) {
   //64 bit subtraction call as above - returns (a-b)
   M_sub64x64(a,b);
   return a;
}

int64 mult32x32(int32 a,int32 b) {
   //Routine to multiply two 32bit values with a 64bit result
   int64 temp2;
   int8 i;
   //zero output
   M_zero64(temp2);
   //Not trying to break any records for efficiency, so all in C, and not using
   //hardware multiply
   //start by moving the int32, into the low half of the int64
   scratch64=cast32x64(b);
   //Now need to work through all 32 bits in the 'a' variable
   for (i=0;i<32;i++) {
      //if source bit in 'a' is one, perform addition
      if (bit_test(a,i))
         M_add64x64(temp2,scratch64);
      //rotate the second value here
      M_shiftleft64(scratch64);
   }
   return temp2;
}

int64 div64x64(int64 a,int64 b) {
   //Routine to divide 'a' by 'b' in 64bit arithmetic
   //returns with result, leaving remainder in the scratch64
   int8 bitno=0,ctr;
   int64 temp;
   if (M_iszero64(b)) {
      if (!M_iszero64(a)) {
         //need maximum result
         temp.w[0]=temp.w[1]=0xFFFFFFFF;
         M_zero64(scratch64);
         return temp;
      }
      //Else 0/0=1
      M_zero64(temp);
      M_zero64(scratch64);
      temp.b[0]=1;
      return temp;
      //return 1
   }
   M_zero64(temp);
   bitno=0;
   //Now position divisor to suit find top bit in A, and rotate b till it's
   //top bit is in the same position
   if (a.w[1]!=0) {
      for(ctr=31;!bit_test(a.w[1], ctr);ctr--) ;
      //Now rotate b till it's top bit matches.
      while (!bit_test(b.w[1], ctr)) {
         M_shiftleft64(b);
         ++bitno;
      }
   }
   else {
      for(ctr=31;!bit_test(a.w[0], ctr);ctr--) ;
      bitno=ctr;
      //Now rotate b till it's top bit matches.
      while (!bit_test(b.w[1], ctr)) {
         M_shiftleft64(b);
         ++bitno;
      }
   }     
   //bitno now stores how far b had to rotate
   while (bitno) {
      //Loop for bitno bits, performing subtract, shift & test
      scratch64=I_A_LT_B(a,b);
      if (!BORROW) {
         a=scratch64;
         //Set output bit if no borrow
         bit64_set(temp,bitno);
      }
      if (M_iszero64(a)) break;
      M_shiftright64(b);
      bitno--;
   }
   scratch64=a;
   return temp;
}

These are only for the PIC18 (for the PIC16, the location of the borrow bit would need to be recoded, and a set of macros would need to be added to simulate the add with carry, and subtract with borrow instruction which doesn't exist on these chips). There is also a 'caveat' (which is one that exists fairly frequently in the CCS code), that the cast32x64 routine, should always put it's result into a temporary variable, and not be directly used inside a function call. If you code:

test=div64x64(test,cast32x64(45678L));

The compiler will overwrite the top bit of the value handed to the subroutine (Ihave no idea why!), while if you code:

temp=cast32x64(45678L);
test=div64x64(test,temp);

It works fine.
I have slightly tested these on half a dozen 'long' arithmetic operations, so am fairly sure they are close. Smile
Because the remainder is left in scratch64, from any division, doing decimal output is fairly easy. You can just repeatedly divide by ten, output the resulting digit, and then carry on using the remainder of the last sum.
The multiplication, is now a lot more efficient.

Best Wishes
stma



Joined: 16 Feb 2004
Posts: 26

View user's profile Send private message

PostPosted: Fri Nov 10, 2006 6:13 am     Reply with quote

Thanks to all for the replies.
Might get away with 18X14 (or less).
However the 32x32 to 64 routines will come in extremely useful.

On a similar subject, where is a good resource to look into floating point maths in ccs c. Starting from scratch!
Thanks again
Steve
kkelis



Joined: 13 Dec 2006
Posts: 5

View user's profile Send private message

PostPosted: Mon Mar 10, 2014 4:45 pm     Reply with quote

I know this is a very old thread, but hopefully i will get some help. Using the above routines provided by Ttelmah i am trying to make a multiplication and then divide the product. So far i was able to do addition,subtractiom, multiplication or a division but one at a time. When i try to divide the multiplication product i get 0 result or the LCD just goes blank. I dont know if this is too much for the PIC to handle...Any help appreciated

Here is my code:

Code:


#include <18f452.h>
#use delay(clock=40000000)
#fuses HS,NOWDT,NOLVP,PROTECT,PUT,BROWNOUT
#include "LCD.c"


union bits64 {
   int8 b[8];
   int32 w[2];
};
typedef union bits64 int64;
int64 scratch64;
int1 BORROW=FALSE;
#bit CARRY=0xFD8.0

//Macro to add two 64bit values. Result into the first
#define M_add64x64(M_a,M_b) \
   #asm\
   movf M_b.b[0],0\
   addwf M_a.b[0],1\
   movf M_b.b[1],0\
   addwfc M_a.b[1],1\
   movf M_b.b[2],0\
   addwfc M_a.b[2],1\
   movf M_b.b[3],0\
   addwfc M_a.b[3],1\
   movf M_b.b[4],0\
   addwfc M_a.b[4],1\
   movf M_b.b[5],0\
   addwfc M_a.b[5],1\
   movf M_b.b[6],0\
   addwfc M_a.b[6],1\
   movf M_b.b[7],0\
   addwfc M_a.b[7],1\
   #endasm
//Macro as above, to subtract two 64bit values.
#define M_sub64x64(M_a,M_b) \
   #asm\
   movf M_b.b[0],0\
   subwf M_a.b[0],1\
   movf M_b.b[1],0\
   subwfb M_a.b[1],1\
   movf M_b.b[2],0\
   subwfb M_a.b[2],1\
   movf M_b.b[3],0\
   subwfb M_a.b[3],1\
   movf M_b.b[4],0\
   subwfb M_a.b[4],1\
   movf M_b.b[5],0\
   subwfb M_a.b[5],1\
   movf M_b.b[6],0\
   subwfb M_a.b[6],1\
   movf M_b.b[7],0\
   subwfb M_a.b[7],1\
   #endasm

//Macros to shift a 64bit value
#define M_shiftleft64(x) shift_left(&x,8,0)
#define M_shiftright64(x) shift_right(&x,8,0)
//Macro to zero a 64bit value
#define M_zero64(x) x.w[0]=0L;x.w[1]=0L
#define M_iszero64(x) ((x.w[0]==0L)&&(x.w[1]==0L))
#define bit64_set(x,i) x.b[i>>3]|=(1<<(i&7))

int64 cast32x64(int32 a) {
   //routine to convert a 32bit int to 64bit
   int64 temp;
   int8 i;
   temp.w[0]=a;
   //and clear the top half
   temp.w[1]=0L;
   return temp;
}

int64 I_A_LT_B(int64 a,int64 b) {
   //Internal routine to _compare_ two 64 bit values
   //Actually performs subtraction, and returns this, with the 'borrow'
   //flag in the global variable 'BORROW'. Performs a-b, hence flag is
   //set if A less than B
   M_sub64x64(a,b);
   if (CARRY==0) BORROW=TRUE;
   else BORROW=FALSE;
   return a;
}

int64 add64x64(int64 a,int64 b) {
   //simple 64 bit addition call - returns (a+b)
   M_add64x64(a,b);
   return a;
}

int64 sub64x64(int64 a,int64 b) {
   //64 bit subtraction call as above - returns (a-b)
   M_sub64x64(a,b);
   return a;
}

int64 mult32x32(int32 a,int32 b) {
   //Routine to multiply two 32bit values with a 64bit result
   int64 temp2;
   int8 i;
   //zero output
   M_zero64(temp2);
   //Not trying to break any records for efficiency, so all in C, and not using
   //hardware multiply
   //start by moving the int32, into the low half of the int64
   scratch64=cast32x64(b);
   //Now need to work through all 32 bits in the 'a' variable
   for (i=0;i<32;i++) {
      //if source bit in 'a' is one, perform addition
      if (bit_test(a,i))
         M_add64x64(temp2,scratch64);
      //rotate the second value here
      M_shiftleft64(scratch64);
   }
   return temp2;
}

int64 div64x64(int64 a,int64 b) {
   //Routine to divide 'a' by 'b' in 64bit arithmetic
   //returns with result, leaving remainder in the scratch64
   int8 bitno=0,ctr;
   int64 temp;
   if (M_iszero64(b)) {
      if (!M_iszero64(a)) {
         //need maximum result
         temp.w[0]=temp.w[1]=0xFFFFFFFF;
         M_zero64(scratch64);
         return temp;
      }
      //Else 0/0=1
      M_zero64(temp);
      M_zero64(scratch64);
      temp.b[0]=1;
      return temp;
      //return 1
   }
   M_zero64(temp);
   bitno=0;
   //Now position divisor to suit find top bit in A, and rotate b till it's
   //top bit is in the same position
   if (a.w[1]!=0) {
      for(ctr=31;!bit_test(a.w[1], ctr);ctr--) ;
      //Now rotate b till it's top bit matches.
      while (!bit_test(b.w[1], ctr)) {
         M_shiftleft64(b);
         ++bitno;
      }
   }
   else {
      for(ctr=31;!bit_test(a.w[0], ctr);ctr--) ;
      bitno=ctr;
      //Now rotate b till it's top bit matches.
      while (!bit_test(b.w[1], ctr)) {
         M_shiftleft64(b);
         ++bitno;
      }
   }     
   //bitno now stores how far b had to rotate
   while (bitno) {
      //Loop for bitno bits, performing subtract, shift & test
      scratch64=I_A_LT_B(a,b);
      if (!BORROW) {
         a=scratch64;
         //Set output bit if no borrow
         bit64_set(temp,bitno);
      }
      if (M_iszero64(a)) break;
      M_shiftright64(b);
      bitno--;
   }
   scratch64=a;
   return temp;
}
//////////////////////////////////////////////////////////////////////


void main()
{

lcd_init();

  for(;;)
   {
int64 test1,test2,b;
int32 fout,a;
unsigned int32 d;

a=34359738368;
fout=20000000;
b=1000000000;


test1=mult32x32(fout,a);

//test1=150000000;

test2=div64x64(test1,b);

d=test2;

printf(lcd_putc,"\f%lu",d);
delay_ms(100);

   }   
}
Ttelmah



Joined: 11 Mar 2010
Posts: 19549

View user's profile Send private message

PostPosted: Tue Mar 11, 2014 3:44 am     Reply with quote

The reason it is hanging, is an 'a' typed where a 'b' should be.
Code:

   else {
      for(ctr=31;!bit_test(b.w[0], ctr);ctr--) ;
      bitno=ctr;
      //Now rotate b till it's top bit matches.
      while (!bit_test(b.w[1], ctr)) {
         M_shiftleft64(b);
         ++bitno;
      }
   }     

for the 'else' section in the div64 routine.
As I said, 'close'. Thought I had checked this before posting...

However there are problems (it'll always give zero with your values).
Remember that the biggest number CCS can code directly, is an int32 at 4294967295. Your value into 'a', is above this. The compiler does not check for overflows like this....

Remember also that the compiler does not understand how to write 'into' an int64, so will not clear the high bytes for you. So to load a value, you have to write to the int32 'parts'. So (if you coded 'a' as an int64:
Code:

a.w[0]=0x0L;
a.w[1]=0x8;

Will write your 64bit value (0x800000000) into the int64.
While you can put in the smaller value into b, with:
Code:

b=cast32x64(1000000000);


The point is that the compiler knows nothing about 64 bit values, so does not know how to cast up/down to these. You need to use the supplied cast if you want to do this. Similarly for the returned values, you should access the 32bit values from the union to read these.

Best Wishes
Douglas Kennedy



Joined: 07 Sep 2003
Posts: 755
Location: Florida

View user's profile Send private message AIM Address

PostPosted: Tue Mar 11, 2014 7:15 am     Reply with quote

Well ,
I'll echo a point Sherpadoug made. Why do you need it?
The most accurate physical calculation from Quantum Mechanics can be notated with just 12 significant decimal places. Quantum Mechanics underpins every physical mechanism and measurement.
Now Mathematics in abstraction can use notation of unlimited precision.
If you count all the grains of sand on all the world's beaches you may need integer notation with more than 12 decimal places. With 58 decimal digits you may be able to count all the atoms in the visible universe. 64 bit binary notation has equivalency to 19 decimal significant digit notation.
Often a quest for precision is triggered by the misconception that a result in binary notation that is translated into decimal notation is inaccurate.
The translational difference doesn't exist for integers but is a natural consequence for fractions and irrational numbers.
1/10 is accurate in decimal but in binary notation it can only be approximated just as 1/3 can only be approximated in the decimal system.
asmboy



Joined: 20 Nov 2007
Posts: 2128
Location: albany ny

View user's profile Send private message AIM Address

PostPosted: Tue Mar 11, 2014 8:00 am     Reply with quote

Quote:

Why do you need it?

cuz it would be handy and faster than the way i do it now.

l example: equation for the ohmic value of
a synchronous(phase sensitive) detector"X,Y,Z"
instrument i'm working on right now:

ohms=(((fullscale_span_value*ADC_count)/2^16*pregain*postgain)*16bit_correction_ref)/frequency_based_correction_value


in an actual instrument the first term can easily be
fullscale*adc_count
100000*65535

and correction ref( in my particular case is 75500 and a typical corr_value is in the range of 50000 to 82000

depending on pre and post gain - fractional ohm results are possible,and all this handily overflows 32bits

thats why we have floats
Very Happy Very Happy Very Happy Very Happy Very Happy
Display posts from previous:   
Post new topic   Reply to topic    CCS Forum Index -> General CCS C Discussion All times are GMT - 6 Hours
Goto page 1, 2  Next
Page 1 of 2

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group