Compute (inaccurate) square root (2012)

getroot.c

#include <stdio.h>
#include <stdlib.h>

typedef struct MATH_BIGNUMBER {
	unsigned int size_word;	//value size in unsigned short(word).
	unsigned short *value;	//pointer to value array in unsigned short.
} MATH_BigNumber;

typedef union MATH_BIGNUMBER_WORK32 {
	unsigned int dword;
	unsigned short word[2];
} MATH_BigNumber_Work32;

unsigned int GetSquareRoot_PutData(unsigned char c);
MATH_BigNumber *MATH_BigNumber_Initialise(unsigned int size_byte);
unsigned int MATH_BigNumber_Free(MATH_BigNumber *bignum);
unsigned int MATH_BigNumber_DumpHex(MATH_BigNumber *bignum);
unsigned int MATH_BigNumber_Copy(MATH_BigNumber *destination, MATH_BigNumber *source);
unsigned int MATH_BigNumber_PutString_Base10(MATH_BigNumber *bignum);
unsigned int MATH_BigNumber_Add_16Bit(MATH_BigNumber *bignum, unsigned short add);
unsigned short MATH_BigNumber_Divide_16Bit(MATH_BigNumber *bignum, unsigned short divisor);
unsigned int MATH_BigNumber_Multiple_16Bit(MATH_BigNumber *bignum, unsigned short times);
int MATH_BigNumber_Compare(MATH_BigNumber *bignum, MATH_BigNumber *target);
unsigned int MATH_BigNumber_Subtract(MATH_BigNumber *bignum, MATH_BigNumber *sub);

unsigned int put_counter;

int main(int argc, char *argv[])
{
	MATH_BigNumber *mainnum, *subnum, *temp;
	unsigned int i, j, option;
	unsigned char input[8], *answer;
	unsigned short root, digits, memory, root2;

	printf("<getroot.exe> Compute Square Root(Digit_by_digit)\n");
	option = 0;

	if(argc == 1){
		printf("\nInput number you want to compute root(less than 65536) >");
		fgets(input, sizeof(input), stdin);
		root = strtol(input, NULL, 10);
		printf("\nInput digit length you want(less than 65536) >");
		fgets(input, sizeof(input), stdin);
		digits = strtol(input, NULL, 10);
		printf("\nInput use memory size(less than 65536)(0:Auto) >");
		fgets(input, sizeof(input), stdin);
		memory = strtol(input, NULL, 10);
	} else if(argc < 4){
		printf("Usage:getroot.exe root digits memory_size\n");
		return 0;
	} else{
		root = strtol(argv[1], NULL, 10);
		digits = strtol(argv[2], NULL, 10);
		memory = strtol(argv[3], NULL, 10);
		if(argc >= 5){
			option = strtol(argv[4], NULL, 10);
		}
	}

	if(memory == 0){
		//memory = (((digits + 9) / 10) << 2) + 8;
		memory = (digits + 11) / 12 * 10;
	}

printf("Square root %d digits:%d memory:%dBytes\n", root, digits, memory);

	mainnum = MATH_BigNumber_Initialise(memory);
	if(option >= 1){
		for(i = 0; i < mainnum->size_word; i++){
			mainnum->value[i] = 0xffff;
		}
		printf("Max number:");
		MATH_BigNumber_PutString_Base10(mainnum);
		printf("\n");
		for(i = 0; i < mainnum->size_word; i++){
			mainnum->value[i] = 0x0000;
		}
	}
	subnum = MATH_BigNumber_Initialise(memory);
	temp = MATH_BigNumber_Initialise(memory);

//compute integer
	put_counter = 1;
	if(root >= 100){
		if(root >= 10000){
			root2 = root % 10000;
			root /= 10000;
			MATH_BigNumber_Add_16Bit(mainnum, root);
			for(j = 0; j < 11; j++){
				MATH_BigNumber_Copy(temp, subnum);
				MATH_BigNumber_Add_16Bit(temp, j);
				MATH_BigNumber_Multiple_16Bit(temp, j);
				if(option >= 2){
					printf("\t%d:", j);
					MATH_BigNumber_PutString_Base10(mainnum);
					printf(", ");
					MATH_BigNumber_PutString_Base10(temp);
					printf("\n");
				}
				if(MATH_BigNumber_Compare(mainnum, temp) == -1){
					j--;
					break;
				}
			}
			GetSquareRoot_PutData(j + 0x30);
			MATH_BigNumber_Copy(temp, subnum);
			MATH_BigNumber_Add_16Bit(temp, j);
			MATH_BigNumber_Multiple_16Bit(temp, j);
			MATH_BigNumber_Subtract(mainnum, temp);
			MATH_BigNumber_Multiple_16Bit(mainnum, 100);
			MATH_BigNumber_Copy(temp, subnum);
			MATH_BigNumber_Add_16Bit(temp, j * 2);
			MATH_BigNumber_Multiple_16Bit(temp, 10);
			MATH_BigNumber_Copy(subnum, temp);
			if(option >= 1){
				printf(" next:");
				MATH_BigNumber_PutString_Base10(mainnum);
				printf(", ");
				MATH_BigNumber_PutString_Base10(subnum);
				printf("\n");
			}
			root = root2;
		}
		root2 = root % 100;
		root /= 100;
		MATH_BigNumber_Add_16Bit(mainnum, root);
		for(j = 0; j < 11; j++){
			MATH_BigNumber_Copy(temp, subnum);
			MATH_BigNumber_Add_16Bit(temp, j);
			MATH_BigNumber_Multiple_16Bit(temp, j);
			if(option >= 2){
				printf("\t%d:", j);
				MATH_BigNumber_PutString_Base10(mainnum);
				printf(", ");
				MATH_BigNumber_PutString_Base10(temp);
				printf("\n");
			}
			if(MATH_BigNumber_Compare(mainnum, temp) == -1){
				j--;
				break;
			}
		}
		GetSquareRoot_PutData(j + 0x30);
		MATH_BigNumber_Copy(temp, subnum);
		MATH_BigNumber_Add_16Bit(temp, j);
		MATH_BigNumber_Multiple_16Bit(temp, j);
		MATH_BigNumber_Subtract(mainnum, temp);
		MATH_BigNumber_Multiple_16Bit(mainnum, 100);
		MATH_BigNumber_Copy(temp, subnum);
		MATH_BigNumber_Add_16Bit(temp, j * 2);
		MATH_BigNumber_Multiple_16Bit(temp, 10);
		MATH_BigNumber_Copy(subnum, temp);
		if(option >= 1){
			printf(" next:");
			MATH_BigNumber_PutString_Base10(mainnum);
			printf(", ");
			MATH_BigNumber_PutString_Base10(subnum);
			printf("\n");
		}
		root = root2;
	}
	MATH_BigNumber_Add_16Bit(mainnum, root);
	for(j = 0; j < 11; j++){
		MATH_BigNumber_Copy(temp, subnum);
		MATH_BigNumber_Add_16Bit(temp, j);
		MATH_BigNumber_Multiple_16Bit(temp, j);
		if(option >= 2){
			printf("\t%d:", j);
			MATH_BigNumber_PutString_Base10(mainnum);
			printf(", ");
			MATH_BigNumber_PutString_Base10(temp);
			printf("\n");
		}
		if(MATH_BigNumber_Compare(mainnum, temp) == -1){
			j--;
			break;
		}
	}
	GetSquareRoot_PutData(j + 0x30);
	MATH_BigNumber_Copy(temp, subnum);
	MATH_BigNumber_Add_16Bit(temp, j);
	MATH_BigNumber_Multiple_16Bit(temp, j);
	MATH_BigNumber_Subtract(mainnum, temp);
	MATH_BigNumber_Multiple_16Bit(mainnum, 100);
	MATH_BigNumber_Copy(temp, subnum);
	MATH_BigNumber_Add_16Bit(temp, j * 2);
	MATH_BigNumber_Multiple_16Bit(temp, 10);
	MATH_BigNumber_Copy(subnum, temp);
	if(option >= 1){
		printf(" next:");
		MATH_BigNumber_PutString_Base10(mainnum);
		printf(", ");
		MATH_BigNumber_PutString_Base10(subnum);
		printf("\n");
	}

	printf(".");
//compute decimal
	put_counter = 0;
	for(i = 0; i < (unsigned int)digits; i++){
		for(j = 0; j < 11; j++){
			MATH_BigNumber_Copy(temp, subnum);
			MATH_BigNumber_Add_16Bit(temp, j);
			MATH_BigNumber_Multiple_16Bit(temp, j);
			if(option >= 2){
				printf("\t%d:", j);
				MATH_BigNumber_PutString_Base10(mainnum);
				printf(", ");
				MATH_BigNumber_PutString_Base10(temp);
				printf("\n");
			}
			if(MATH_BigNumber_Compare(mainnum, temp) == -1){
				j--;
				break;
			}
		}
		GetSquareRoot_PutData(j + 0x30);
		MATH_BigNumber_Copy(temp, subnum);
		MATH_BigNumber_Add_16Bit(temp, j);
		MATH_BigNumber_Multiple_16Bit(temp, j);
		MATH_BigNumber_Subtract(mainnum, temp);
		MATH_BigNumber_Multiple_16Bit(mainnum, 100);
		MATH_BigNumber_Copy(temp, subnum);
		MATH_BigNumber_Add_16Bit(temp, j * 2);
		MATH_BigNumber_Multiple_16Bit(temp, 10);
		MATH_BigNumber_Copy(subnum, temp);
		if(option >= 1){
			printf(" next:");
			MATH_BigNumber_PutString_Base10(mainnum);
			printf(", ");
			MATH_BigNumber_PutString_Base10(subnum);
			printf("\n");
			if(option >= 3){
				MATH_BigNumber_DumpHex(mainnum);
				printf("\n");
			}
		}
	}
    
    printf("\n");

	return 0;
}

unsigned int GetSquareRoot_PutData(unsigned char c)
{
	if(put_counter == 0){
		fputc('\n', stdout);
		fputc('>', stdout);
		fputc('\t', stdout);
	}

	fputc(c, stdout);
	put_counter++;

	if(put_counter == 10 || put_counter == 20 || put_counter == 30 || put_counter == 40 || put_counter == 50){
		fputc(' ', stdout);
		if(put_counter == 50){
			put_counter = 0;
		}
	}
	return 0;
}

MATH_BigNumber *MATH_BigNumber_Initialise(unsigned int size_byte)
{
	MATH_BigNumber *bignum;
	unsigned int i;

	bignum = malloc(sizeof(MATH_BigNumber));
	if(bignum == 0){
		printf("\nMATH_BigNumber_Initialise:Can't allocate memory(bignum). Abort.\n");
		exit(EXIT_FAILURE);
	}

	bignum->size_word = (size_byte + 1) >> 1;
	bignum->value = malloc(bignum->size_word << 1);
	if(bignum->value == 0){
		printf("\nMATH_BigNumber_Initialise:Can't allocate memory(bignum->value). Abort.\n");
		exit(EXIT_FAILURE);
	}
	for(i = 0; i < bignum->size_word; i++){
		bignum->value[i] = 0;
	}

	return bignum;
}

unsigned int MATH_BigNumber_Free(MATH_BigNumber *bignum)
{
	free(bignum->value);
	free(bignum);

	return 0;
}

unsigned int MATH_BigNumber_DumpHex(MATH_BigNumber *bignum)
{
	unsigned int i;

	printf("0x");
	for(i = bignum->size_word ; i > 0; i--){
		printf("%04X", bignum->value[i - 1]);
	}

	return 0;
}

unsigned int MATH_BigNumber_Copy(MATH_BigNumber *destination, MATH_BigNumber *source)
{
	unsigned int i; 

	if(destination->size_word != source->size_word){
		printf("\nMATH_BigNumber_Copy:Not same size. Abort.\n");
		exit(EXIT_FAILURE);
	}

	for(i = 0; i < destination->size_word; i++){
		destination->value[i] = source->value[i];
	}

	return 0;
}

unsigned int MATH_BigNumber_PutString_Base10(MATH_BigNumber *bignum)
{
	unsigned char *s;
	unsigned int i, str_size;
	MATH_BigNumber *temp;

	str_size = (((bignum->size_word << 1) + 3) / 4 * 10) + 1;
	s = malloc(str_size);
	if(s == 0){
		printf("\nMATH_BigNumber_PutString_Base10:Can't allocate memory. Abort.\n");
		exit(EXIT_FAILURE);
	}


	s[str_size - 1] = 0x00;

	temp = MATH_BigNumber_Initialise(bignum->size_word << 1);
	MATH_BigNumber_Copy(temp, bignum);

	for(i = 0; i < str_size - 1; i++){
		s[str_size - 2 - i] = MATH_BigNumber_Divide_16Bit(temp, 10) + 0x30;
	}

	for(i = 0; i < str_size - 2; i++){
		if(s[i] != '0'){
			break;
		}
	}

	printf("%s", &s[i]);

	MATH_BigNumber_Free(temp);
	free(s);

	return 0;
}

unsigned int MATH_BigNumber_Add_16Bit(MATH_BigNumber *bignum, unsigned short add)
{
	MATH_BigNumber_Work32 temp;
	unsigned int temp2;
	unsigned int i;

	temp.word[1] = add;

	for(i = 0; i < bignum->size_word; i++){
		temp.word[0] = bignum->value[i];
		temp2 = temp.word[1];
		temp.word[1] = 0;
		temp.dword += temp2;
		bignum->value[i] = temp.word[0];
	}
	if(temp.word[1] != 0){
		printf("\nMATH_BigNumber_Add_16Bit:Overflow. Abort.\n");
		exit(EXIT_FAILURE);
	}
	return 0;
}

//divide [bignum] by [divisor]. answer is [bignum]. return value is residue.
unsigned short MATH_BigNumber_Divide_16Bit(MATH_BigNumber *bignum, unsigned short divisor)
{
	MATH_BigNumber_Work32 temp;
	unsigned int temp2;
	unsigned int i, max;

	if(divisor == 0){
		printf("\nMATH_BigNumber_Divide_16Bit:Divide by zero. Abort.\n");
		exit(EXIT_FAILURE);
	}

	max = 0;
	for(i = 0; i < bignum->size_word - 1; i++){
		if(bignum->value[i] != 0x0000){
			max = i;
		}
	}

	temp.word[1] = 0;
	for(i = max + 1; i > 0; i--){
		temp.word[0] = bignum->value[i - 1];
		temp2 = temp.dword % divisor;
		temp.dword /= divisor;
		bignum->value[i - 1] = temp.word[0];
		temp.word[1] = temp2;
	}

	return temp2;
}

unsigned int MATH_BigNumber_Multiple_16Bit(MATH_BigNumber *bignum, unsigned short times)
{
	MATH_BigNumber_Work32 temp;
	unsigned int temp2;
	unsigned int i;

	temp2 = 0;
	for(i = 0; i < bignum->size_word; i++){
		temp.word[1] = 0;
		temp.word[0] = bignum->value[i];
		temp.dword *= times;
		temp.dword += temp2;
		bignum->value[i] = temp.word[0];
		temp2 = temp.word[1];
	}
	if(temp2 != 0){
		printf("\nMATH_BigNumber_Multiple_16Bit:Overflow. Abort.\n");
		exit(EXIT_FAILURE);
	}
	return 0;
}

//bignum > target : 1
//bignum == target : 0
//bignum < target : -1
int MATH_BigNumber_Compare(MATH_BigNumber *bignum, MATH_BigNumber *target)
{
	unsigned int i, bignum_max, target_max;

	bignum_max = 0;
	for(i = 0; i < bignum->size_word - 1; i++){
		if(bignum->value[i] != 0x0000){
			bignum_max = i;
		}
	}

	target_max = 0;
	for(i = 0; i < target->size_word - 1; i++){
		if(target->value[i] != 0x0000){
			target_max = i;
		}
	}

	if(bignum_max > target_max){
		return 1;
	}
	if(bignum_max < target_max){
		return -1;
	}

	for(i = bignum_max + 1; i > 0; i--){
		if(bignum->value[i - 1] != target->value[i - 1]){
			if(bignum->value[i - 1] > target->value[i - 1]){
				return 1;
			}
			return -1;
		}
	}
	return 0;
}

//bignum -= sub
/*
unsigned int MATH_BigNumber_Subtract(MATH_BigNumber *bignum, MATH_BigNumber *sub)
{
	MATH_BigNumber_Work32 temp;
	unsigned int i, sub_max, borrow;
	int n;

	n = MATH_BigNumber_Compare(bignum, sub);
	if(n == -1){
		printf("\nMATH_BigNumber_Subtract:Underflow. Abort.\n");
		exit(EXIT_FAILURE);
	}
	if(n == 0){
		for(i = 0; i < bignum->size_word; i++){
			bignum->value[i] = 0;
		}
		return 0;
	}

printf("\n");
MATH_BigNumber_DumpHex(bignum);
printf("\n");
MATH_BigNumber_DumpHex(sub);
printf("\n");

	sub_max = 0;
	for(i = 0; i < sub->size_word - 1; i++){
		if(sub->value[i] != 0x0000){
			sub_max = i;
		}
	}

	borrow = 0;
	temp.word[1] = 0;
	for(i = 0; i <= sub_max; i++){
		temp.word[0] = bignum->value[i];
		if(temp.dword < sub->value[i] + borrow){
			temp.word[1] = 1;
			if(borrow == 1){
				temp.dword -= 1;
			}
			borrow = 1;
		} else{
			if(borrow == 1){
				temp.dword -= 1;
			}
			borrow = 0;
		}
		temp.dword -= sub->value[i];
		bignum->value[i] = temp.word[0];
	}
	for(; i < bignum->size_word; i++){
		if(borrow == 0){
			break;
		}
		temp.word[0] = bignum->value[i];
		if(temp.word[0] == 0){
			temp.word[1] = 1;
			temp.dword -= 1;
		}
		bignum->value[i] = temp.word[0];
	}
	return 0;
}
*/

unsigned int MATH_BigNumber_Subtract(MATH_BigNumber *bignum, MATH_BigNumber *sub)
{
	MATH_BigNumber_Work32 temp;
	unsigned int i, sub_max;
	int n;

	n = MATH_BigNumber_Compare(bignum, sub);
	if(n == -1){
		printf("\nMATH_BigNumber_Subtract:Underflow. Abort.\n");
		exit(EXIT_FAILURE);
	}
	if(n == 0){
		for(i = 0; i < bignum->size_word; i++){
			bignum->value[i] = 0;
		}
		return 0;
	}

	sub_max = 0;
	for(i = 0; i < sub->size_word - 1; i++){
		if(sub->value[i] != 0x0000){
			sub_max = i;
		}
	}

	temp.word[1] = 0;
	for(i = 0; i <= sub_max; i++){
		temp.word[0] = bignum->value[i];
		if(temp.dword == 0xffff){	//-1
			bignum->value[i + 1]--;
		}
		if(temp.dword < (unsigned int)sub->value[i]){
			bignum->value[i + 1]--;
			temp.word[1] = 1;
		}
		temp.dword -= sub->value[i];
		bignum->value[i] = temp.word[0];
	}
	for(; i < bignum->size_word; i++){
		if(bignum->value[i] == 0xffff){	//-1
			bignum->value[i + 1]--;
		} else{
			break;
		}
	}
	return 0;
}