《北京邮电大学数值及符号计算实验报告》由会员分享，可在线阅读，更多相关《北京邮电大学数值及符号计算实验报告（17页珍藏版）》请在金锄头文库上搜索。
1、1 方程求根的牛顿迭代法实验报告方程求根的牛顿迭代法实验报告 目录目录 方程求根的牛顿迭代法实验报告 1 1 实验环境及要求：. 1 2 实验内容：. 1 3 实验步骤：. 2 3.1 浮点数的平方根 2 3.2 整数 c 的平方根 5 3.3 求解浮点数倒数 6 3.4 奇数 C 的逆 1 mod 223 . 8 4 实验结果及分析说明. 10 4.1 浮点数开方 10 4.2 整数开根号 10 4.3 浮点数求倒数 10 4.4 整数倒数 10 5 附录一（源代码）. 11 1 实验环境及要求：实验环境及要求： Intel i5 core，Ubuntu 12.04 X64 操作系统，使用
2、Linux c+ 嵌入 AT int s; double result; double p=0.4144142135623730951,q=0.5946699141100893; asm( “movq %2,%rcx;“ “movq $0x7ff0000000000000,%rax;“ “andq %rcx, %rax;“ “shrq $52,%rax;“ “subq $1023,%rax;“ “movq %rax,%rdx;“ “movq $0x800fffffffffffff,%rax;“ “andq %rax,%rcx;“ “movq $0x3ff0000000000000,%rax;“
3、 “orq %rax,%rcx;“ :“+c“(c),“=d“(s) :“m“(c); /step 2 /the process of Newton iteration. 3 times iterations reach to 64 bits precisiong asm( “movsd %3,%xmm0;“ “mulsd %1,%xmm0;“ “addsd %2,%xmm0;“ “movl $5,%ecx;“ “movq %3,%xmm1;“ 4 “movabsq $4611686018427387904,%rax;“ /the first operan is 2. thats the fo
4、rmat of double float “movq %rax,%xmm2;“ “start:;“ “cmpl $0,%ecx;“ “jz endloop;“ “decl %ecx;“ “movq %xmm1,%rax;“ “divsd %xmm0,%xmm1;“ “addsd %xmm1,%xmm0;“ “movq %rax,%xmm1;“ “divsd %xmm2,%xmm0;“ “jmp start;“ “endloop:“ “movq %xmm0,%rdx;“ :“+d“(result) :“m“(p),“m“(q),“m“(c); /step3 /add the multipicit
5、y s/2 asm( “movq %rbx,%rdx;“ “shrq $1,%rbx;“ “movq $0x800fffffffffffff,%rcx;“ “andq %rcx,%rax;“ “movq $0x3ff,%rcx;“ “addq %rbx,%rcx;“ “salq $52,%rcx;“ “orq %rcx,%rax;“ “andq $1,%rdx;“ “cmpl $0,%edx;“ “jz end3;“ “movabsq $4609047870845172685,%rcx;“ “movq %rax,%xmm1;“ “movq %rcx,%xmm0;“ “mulsd %xmm1,%
6、xmm0;“ “movq %xmm0,%rax;“ “end3:;“ “movq %rax,%rdx;“ :“+d“(result) :“a“(result),“b“(s); return result; 代码说明，代码使用的是 linux c+嵌套汇编，浮点运算采用的是 xmm 指令。代 码分为三个部分，一是提取 11 位的阶数，然后开根号，再将阶数返回。 5 3.2 整数整数 c 的平方根的平方根 其牛顿迭代式： 0 +1= 1 2 (+ ) 其中取0 迭代时若+1文件中的函数 sqrt() 结果为：3.87298334620742。 结果一样，实现正确。 同时使用了 linux(ubun
7、tu12.04)系统的提供的中的函数 colock_gettime()获 取运行时间，该函数可获取纳秒级别的时间。 如图，实验实现的浮点平方根求解 15 的运行时间为 146357ns。而 c+提供的函数为 4322ns。相差 33.86 倍。 4.2 整数开根号整数开根号 整数 15 的结果为 3，结果正确。 4.3 浮点数求倒数浮点数求倒数 浮点数 1.2 的倒数为0.833333333333333，结果正确。 4.4 整数倒数整数倒数 987 在 32 位机器中的倒数为 2945990739，由图只 2945990739*987=1 mod 2 32 11 5 附录一（源代码）附录一（源
8、代码） 源代码运行说明：本代码在 ubuntu12.04 x64 linux 系统上开发，Intel i5 cpu。 本代码编译命令为：g+ sourcefile.cpp Wall lrt o object 运行命令：./object 其中，librt.so 支持时间函数 clock_gettime(). #include #include #include #include #include #include #include using namespace std; unsigned my_inverse(unsigned c)/c is an odd number with 32 bit
9、s asm( /step1, cal x0 “movl %eax,%ebx;“ “incl %eax;“ “movl $4,%ecx;“ “andl %ecx,%eax;“ “shll $1,%eax;“ “addl %ebx,%eax;“ / “jmp in2_end_loop;“ /setp2 the process of Newton iteraton . 3 times “movl $3,%ecx;“ “pushq %rcx;“ “in2_start_loop:;“ “popq %rcx;“ “cmpl $0,%ecx;“ “jz in2_end_loop;“ “decl %ecx;“
10、 “pushq %rcx;“ “pushq %rax;“ “mull %ebx;“ “movl $2,%ecx;“ “subl %eax,%ecx;“ “movl %ecx,%eax;“ “popq %rcx;“ “mull %ecx;“ “jmp in2_start_loop;“ 12 “in2_end_loop:;“ “movl %eax,%edx;“ :“+d“(c) :“a“(c); return c; double my_inverse(double c) /step one:set c in 1,2 and s is the multipicity of c; int s; asm
11、( “movq %2,%rcx;“ “movq $0x7ff0000000000000,%rax;“ “andq %rcx, %rax;“ “shrq $52,%rax;“ “subq $1023,%rax;“ “movq %rax,%rdx;“ “movq $0x800fffffffffffff,%rax;“ “andq %rax,%rcx;“ “movq $0x3ff0000000000000,%rax;“ “orq %rax,%rcx;“ :“+c“(c),“=d“(s) :“m“(c); double p=-0.4705882352941176,q=1.411764705882353;
12、 double result; /step2 the process of Newton iteration; 3 times asm( “movsd %3,%xmm0;“ “mulsd %1,%xmm0;“ “addsd %2,%xmm0;“ “movsd %3,%xmm2;“/store c “movabsq $4611686018427387904,%rax;“ /the first operan is 2. thats the format of double float “movl $4,%ecx;“ “in1_start_loop:;“ “cmpl $0,%ecx;“ “jz in
13、1_end_loop;“ “decl %ecx;“ “movsd %xmm0,%xmm3;“ “mulsd %xmm2,%xmm0;“ “movq %rax,%xmm1;“ “subsd %xmm0,%xmm1;“ “mulsd %xmm3,%xmm1;“ “movsd %xmm1,%xmm0;“ “jmp in1_start_loop;“ 13 “in1_end_loop:;“ “movq %xmm0,%rdx;“ :“+d“(result) :“m“(p),“m“(q),“m“(c); /cout“r “resultendl; /step3 : add the multiplicity -
14、s; asm( “movq %1,%rdx;“ “movq $0x7ff0000000000000,%rax;“ “andq %rdx,%rax;“ “shrq $52,%rax;“ “movq $0x800fffffffffffff,%rcx;“ “andq %rcx,%rdx;“ “subq %rbx,%rax;“ “shlq $52,%rax;“ “orq %rax,%rdx;“ :“+d“(result) :“m“(result),“b“(s); /* cout“p= “setprecision(16)pendl; cout“q= “setprecision(16)qendl; cou
15、t“result= “setprecision(16)resultendl; cout“c= “setprecision(16)cendl; cout“s= “sendl; */ return result; double my_sqrt(double c) /step one:set c in 1,2 and s is the multipicity of c; coutcendl; int s; double result; double p=0.4144142135623730951,q=0.5946699141100893; asm( “movq %2,%rcx;“ “movq $0x7ff0000000000000,%rax;“ “andq %rcx, %rax;“ “shrq $52,%rax;“ “subq $1023,%rax;“ “movq %rax,%rdx;“ “movq $0x800fffffffffffff,%rax;“ “andq %rax,%rcx;“ “movq $0x3ff0000000000000,%rax;“ 14 “orq %rax,%rcx;“ :“+c“(c),“=d“(s) :“m“(c); /step 2 /the process o
