PK 使人进步：highgear 向网友刘前辈提出很简单的切磋请求

只看该作者 · 2011-10-7 20:48

每个月都来那么一下？

只看该作者 · 2011-10-7 21:18

站位

只看该作者 · 2011-10-7 21:48

谢谢各位捧场。

多任务的实现说穿了很简单，就是利用定时中断打断当前运行的任务，把当前所使用的寄存器 push 到栈中，然后把切换到下一个任务的栈点，pop 出这个即将运行的任务现场，最后把运行权交给这个任务。可以看到，多任务的实现的关键点是中断，在 windows 下要实现打断运行的任务得到控制权，不是做不到，但非常复杂，所以在 PC 中应在 DOS 下实现，DOS 支持 bios 的0x10定时中断，稍加改造就可以是支持 1ms 的中断。
编译器方面， Turbo C 2.0 或是 Borland C++ 5.02 都直接支持中断响应，底层 register 的存取以及汇编的嵌入等等。

等待刘前辈回应的同时，响应 123姐姐的号召，讲解一下高精度的整数运算。整数运算在嵌入式，尤其是低端的 mcu 应用中意义重大。实际应用中使用浮点，固然可以简化程序，但代价是效率。浮点与定点运算的效率差别可达几十至几百倍，所以使用定点运算可以极大的提高运算速度。

只看该作者 · 2011-10-7 21:53

本帖最后由 highgear 于 2011-10-7 22:00 编辑

下面以sin 函数为例，浮点的sin 函数在16bit 的cpu 上（带有16x16=32bit 的乘法器）需要约 2000 周期左右。如何用定点来处理这个函数呢？关键的一点要记住就是：映射。必须把浮点数转换成相应范围内的16或32位的整数，比如 PI ---- 3.1416. 角度的范围是 [0, 2PI], 一个技巧就是映射到 [0, 0xFFFF], 对于16bit的整数来说，[0, 0xFFFF]也是周期性的，甚至不用考虑符号问题，因为[0x8000, 0xFFFF]对于signed int 来说，就是负整数。

实现sin 函数的方法很多，最常用是泰勒多项式法，这是用泰勒展开式计算 sin函数，快速的有查表法，是用空间换时间的方法。为了简单起见，这里介绍一下查表法。

既然角度的范围是[0, 0xFFFF], 我们不可能把 65536 个整数放入表中，折中的办法是用128，256个整数，用插值来提高精度。

首先定义 PI:
#define _PI 0x8000

创造一个256整数的 sine 函数表，计算公式为 32768 * sin(2* 3.1415927 * i / 256):

const short SineTable[] =
{
      0,       804,       1608,       2411  , 3212,       4011,       4808,       5602  , 6393,

创剩下的就是插值：
short  Sin(short value)
{
      unsigned short vh = ((unsigned short) value) >> 8; //得到表的索引值
      short vl = value & 0xFF; //用余数做插值
      short r = SineTable[vh];
      return r + (vl * (SineTable[(vh + 1) & 0xFF] - r) >> 8); // 插值
}
而结果为[-1, 1] 映射到[0x8000, 0x7FFF], 即实数 1 对应整数 2^15, 精度可以达到 0.0001。

以此可以方便的计算出 cos 函数：
short Cos(short value)
{
      return Sin(value + _PI / 2);
}

这个例子很简单，但还是可以看出一些整数运算的技巧，就是把实数映射到 2^N 中，然后可以利用加减乘以及移位来完成计算，尽量避免除法。

只看该作者 · 2011-10-7 22:06

24# highgear

楼主，我拜你为师吧

只看该作者 · 2011-10-7 22:33

膜拜ｉｎｇ！

1万 · 2011-10-8 04:37

如何解决解决哲学家吃饭问题?

下面连接就有，highgear老师有详细讲解，随手抄点就能拼凑了一个最基本的多任务操作系统，并在此 os 上解决哲学家就餐问题。

如何解决解决哲学家吃饭问题，见链接中的讨论：
https://bbs.21ic.com/viewthread.p ... 2%E5%AD%A6%E5%AE%B6

1万 · 2011-10-8 05:02

下面以sin 函数为例，浮点的sin 函数在16bit 的cpu 上（带有16x16=32bit 的乘法器）需要约 2000 周期左右。如何用定点来处理这个函数呢？关键的一点要记住就是：映射。必须把浮点数转换成相应范围内的16或32位的整 ...
highgear 发表于 2011-10-7 21:53

highgear老师介绍的方法很实用，实现sin 函数的方法很多，highgear老师介绍了最常用的有泰勒多项式法，但泰勒多项式法收敛太慢，快速的有查表法，但是用空间换时间的方法，比较耗内存。俺比较喜欢用切比雪夫插值法，运算简单收敛快~~~

为了简单起见，这里介绍一下用刘前辈最喜欢的C51, 编写的实现sin 函数的方法。

由于俺水平和highgear老师不是一个挡次------差很多很多，程序写的比较烂，见谅。大伙就凑合着看看吧~~~

此程序难点：
由于内部要用到16位X16位=32位乘法，乘积的小数点位置不能简单的确定和定位，需要随机调整，否则精度难于达到 0.0001。

标准C51没有16位X16位=32位乘法，为简单起见，将16位数据扩展到32位，执行长整型32位X32位=32位乘法。

1万 · 2011-10-8 05:32

本帖最后由 Cortex-M0 于 2011-10-8 09:28 编辑

为简单起见，本程序只作输入弧度 0-PI之运算，超过上述范围的需先自行转换，如范围在 PI-2PI之间，先减去一个 PI, 再将sin计算结果取负即可。

由于去掉了符号，这样的好处是变相增加了1bit位，提高了运算精度。

先定义一个输入值，范围在 0-PI弧度之间
#define INPUT 1.2345f //  输入变量

用Keil C51编译后模拟运行，显示结果如下：
FX --->  输入值，范围在 0-PI弧度之间
FY0  --->  整数求sin函数，返回值转换成浮点数供观察
FY1  --->  调用浮点数求sin函数(验证切比雪夫插值公式正确性)
FY2  --->  调用C51浮点函数库求sin(标准结论参照)

复制

#include <DEFINE51.H>

#include <MATH.H>



#define INPUT 1.2345f // 输入变量



#define K6 43120 // -1.2850635e-3 /256/256/256/2 ( >>25 ) 

#define K5 50800 // +0.012111701 /256/256/64 ( >>22 ) 

#define K4 50536 // -6.02441338e-3 /256/256/128 ( >>23 ) 

#define K3 42304 // -0.161379884 /256/256/4 ( >>18 ) 

#define K2 39551 // -2.357414e-3 /256/256/256 ( >>24 ) 

#define K1 32782 // +1.00042182 /256/128 ( >>15 ) 

#define K0 222 // -1.32048e-5 /256/256/256 ( >>24 ) 

#define KDP -8 // 小数点修正



#define FK6 -1.2850635e-3 

#define FK5 +0.012111701 

#define FK4 -6.02441338e-3 

#define FK3 -0.161379884 

#define FK2 -2.357414e-3 

#define FK1 +1.00042182 

#define FK0 -1.32048e-5 



int32 isin(uint16 Ix)

{ uint16 X1, X2, X3, X4, X5, X6;

int32 Ly;

X1 = Ix; // 2

X2 = (uint16)(((uint32)X1 * (uint32)X1) >> 16); // 4

X3 = (uint16)(((uint32)X1 * (uint32)X2) >> 15); // 5

X4 = (uint16)(((uint32)X2 * (uint32)X2) >> 15); // 7

X5 = (uint16)(((uint32)X2 * (uint32)X3) >> 16); // 9

X6 = (uint16)(((uint32)X3 * (uint32)X3) >> 16); // 10

Ly = 0; // 

Ly -= (int32)(((uint32)K6 * (uint32)X6) >>(25-10+KDP));

Ly += (int32)(((uint32)K5 * (uint32)X5) >>(22-9+KDP));

Ly -= (int32)(((uint32)K4 * (uint32)X4) >>(23-7+KDP));

Ly -= (int32)(((uint32)K3 * (uint32)X3) >>(18-5+KDP));

Ly -= (int32)(((uint32)K2 * (uint32)X2) >>(24-4+KDP));

Ly += (int32)(((uint32)K1 * (uint32)X1) >>(15-2+KDP));

Ly -= (int32)K0;

if(Ly>=0) 

{ return((uint16)(Ly>>9)); }

else 

{ return(0); }

}



float fsin(float fx)

{ float fX1, fX2, fX3, fX4, fX5, fX6;

float fY;

fX1 = fx; 

fX2 = fX1 * fX1;

fX3 = fX1 * fX2;

fX4 = fX2 * fX2;

fX5 = fX2 * fX3;

fX6 = fX3 * fX3;

fY = 0;

fY += FK6 * fX6;

fY += FK5 * fX5;

fY += FK4 * fX4;

fY += FK3 * fX3;

fY += FK2 * fX2;

fY += FK1 * fX1;

fY += FK0;

if(fY>=0) 

{ return(fY); }

else 

{ return(0); }

}



void main(void)

{ float FX, FY0, FY1, FY2;

uint16 Ia, Ib, X, Y; 

FX = INPUT;

FY1 = modf(FX, &FY0); // 分离输入变量整数部分和小数部分

Ia = (uint16)FY0; // 输入变量整数部分 --> Ia

Ib = (uint16)(FY1*256*64); // 输入变量小数部分 --> Ib

X = (Ia << 14) + Ib; // 两部分合并后 --> X(前2bit为整数部分，后14bit为小数部分)

Y = isin(X); // 调用仅使用整数变量求sin函数(返回值前1bit为整数部分，后15bit为小数部分)

FY0 = ((float)Y)/256/128; // 返回值转换成浮点数供观察

FY1 = fsin(FX); // 调用浮点数求sin函数(验证公式正确性) 

FY2 = sin(FX); // 调用C51浮点函数库求sin(标准结论参照) 

while(1);

}

1万 · 2011-10-8 05:34

1万 · 2011-10-8 05:37

俺随机测试了 0-PI之间的若干个点，由上图随机测试结果可见，整数求sin的精度达到且优于 0.0001。

1万 · 2011-10-8 05:51

上述sin(x)的运算公式：
sin(x) = k6 * x⁶ + k5 * X⁵ + k4 * x⁴ + k3 * x³ + k2 * x² + k1 * x + k0
还可以简化成：
sin(x) =（（（（（（（（（（k6 * x）+k5）* x）+k4）* x）+k3）* x）+k2）* x）+k1）+ k0
以减少乘法次数，加速运算速度。

只看该作者 · 2011-10-8 09:19

zhi chi

只看该作者 · 2011-10-8 11:20

围观……

3万 · 2011-10-8 11:38

pk。。。。又见pk。。。。

只看该作者 · 2011-10-8 11:38

支持一下

只看该作者 · 2011-10-8 14:31

此贴有技术含量。

只看该作者 · 2011-10-8 21:40

赶紧抢座围观

1万 · 2011-10-9 00:01

刷屏党到此一游:lol

蛋疼期作品：求SinX

int_sin: ; sin(x) = x - x^3/3! + x^5/5! - x^7/7! + x^9/9! - x^11/11! ..., x 单位为弧度(rad)
;为加快运算速度及减少代码量, 特经预处理使 x = (0 ~ Pi/4)
;原理是当 x = (Pi/4 ~ Pi/2) 时 sin(x) = cos(Pi/2 - x) 成立
;其他象限各区间也都有类似等式成立
;x = [r4~7] ; 205 bytes; <=3006T
;已验证, 准确, 结果有效位31位
mov a, r4
mov mult + 0, a
mov mult + 4, a
mov sinx + 0, a
mov root + 0, a
mov a, r5
mov mult + 1, a
mov mult + 5, a
mov sinx + 1, a
mov root + 1, a
mov a, r6
mov mult + 2, a
mov mult + 6, a
mov sinx + 2, a
mov root + 2, a
mov a, r7
mov mult + 3, a
mov mult + 7, a
mov sinx + 3, a
mov root + 3, a ; 种子 x --> sinx(0~3), mult(0~3), mult(4~7), root(0~3)
lcall mul_32 ; x^2
mov cosx + 0, r4
mov cosx + 1, r5
mov cosx + 2, r6
mov cosx + 3, r7 ; x^2 --> cosx(0~3)
mov mult + 10, #0

int_sin_loop:
mov mult + 0, root + 0
mov mult + 1, root + 1
mov mult + 2, root + 2
mov mult + 3, root + 3 ;x^n
mov mult + 4, cosx + 0
mov mult + 5, cosx + 1
mov mult + 6, cosx + 2
mov mult + 7, cosx + 3 ;x^2
lcall mul_32 ; x^n X x^2
mov a, r4
mov mult + 0, a
mov root + 0, a
mov a, r5
mov mult + 1, a
mov root + 1, a
mov a, r6
mov mult + 2, a
mov root + 2, a
mov a, r7
mov mult + 3, a
mov root + 3, a ; x^n --> root(0~3), mult(0~3)

mov dptr, #fact_table
mov a, mult + 10
mov r0, a
movc a, @a + dptr
mov mult + 4, a
inc r0
mov a, r0
movc a, @a + dptr
mov mult + 5, a
inc r0
mov a, r0
movc a, @a + dptr
mov mult + 6, a
inc r0
mov a, r0
movc a, @a + dptr
mov mult + 7, a ; 1/n! --> mult(4~7)
inc r0
mov a, r0
mov mult + 10, a
jnb acc.4, int_sin_xxxx
mov r0, mult + 4
mov r1, mult + 5
lcall mul_1632
mov r3, 5
mov r4, 6
mov r5, 7
clr a
mov r6, a
mov r7, a
sjmp $ + 5
int_sin_xxxx:
lcall mul_32 ; x^n/n!
mov a, r4
orl a, r5
orl a, r6
orl a, r7
jz int_sin_exit

mov a, r3
rlc a
mov a, mult + 10
jnb acc.2, add_xn
sub_xn: mov a, sinx + 0
subb a, r4
mov sinx + 0, a
mov a, sinx + 1
subb a, r5
mov sinx + 1, a
mov a, sinx + 2
subb a, r6
mov sinx + 2, a
mov a, sinx + 3
subb a, r7
mov sinx + 3, a ; x - x^n/n!
ljmp int_sin_loop
add_xn: mov a, sinx + 0
addc a, r4
mov sinx + 0, a
mov a, sinx + 1
addc a, r5
mov sinx + 1, a
mov a, sinx + 2
addc a, r6
mov sinx + 2, a
mov a, sinx + 3
addc a, r7
mov sinx + 3, a ; x + x^n/n!
ljmp int_sin_loop
int_sin_exit:
ret

sin_cos:

1万 · 2011-10-9 00:12

math_code segment code
rseg math_code
fact_table: db 0xab, 0xaa, 0xaa, 0x2a ;1/3!  = 1/6
db 0x22, 0x22, 0x22, 0x02 ;1/5!  = 1/120
db 0xd0, 0x00, 0x0d, 0x00 ;1/7!  = 1/5040
db 0x3c, 0x2e, 0x00, 0x00 ;1/9!  = 1/362880
db 0x6c, 0x00, 0x00, 0x00 ;1/11! = 1/39916800 ;(pi/4)^11/11! = 0x00000008
db 0x00, 0x00, 0x00, 0x00
pi4: db 0x69, 0x21, 0xa2, 0xda, 0x0f, 0xc9 ;pi/4  = 3.141592653589793.../4
pi2: db 0xd2, 0x42, 0x44, 0xb5, 0x1f, 0x92 ;pi/2  = 3.141592653589793.../2, 最高字节为0x01, 不参加运算

.
.
.
.

mul_1632: ;r0r1 X mult[0~3] --> r2~r7; ;101 bytes; <= 108t
;未验证, 冒然使用后果自负!
;破坏 Acc, B, F0, PSW
mov a, mult ;a0 X b0
mov b, r0
mul ab
mov r2, a
mov r3, b
mov a, mult + 2 ;a2 X b0
mov b, r0
mul ab
mov r4, a
mov r5, b
mov a, mult + 3 ;a3 X b1 ;a0b0, a2b0, a3b1 --> r2~r7
mov b, r1
mul ab
mov r6, a
mov r7, b

mov a, mult ;a0 X b1
mov b, r1
mul ab
add a, r3
mov r3, a
mov a, b
addc a, r4
mov r4, a
mov f0, c
mov a, mult + 2 ;a2 X b1
mov b, r1
mul ab
mov c, f0
addc a, r5
mov r5, a
mov a, b
addc a, r6
mov r6, a
jnc $ + 3
inc r7 ;a0b0, a0b1, a2b0, a2b1, a3b1

mov a, mult + 1 ;a1 X b0
mov b, r0
mul ab
add a, r3
mov r3, a
mov a, b
addc a, r4
mov r4, a
mov f0, c
mov a, mult + 3 ;a3 X b0
mov b, r0
mul ab
mov c, f0
addc a, r5
mov r5, a
mov a, b
addc a, r6
mov r6, a
jnc $ + 3
inc r7 ;a0b0, a0b1, a1b0, a2b0, a2b1, a3b0, a3b1

mov a, mult + 1 ;a1 X b1 ;a0b0, a0b1, a1b0, a1b1, a2b0, a2b1, a3b0, a3b1
mov b, r1
mul ab
add a, r4
mov r4, a
mov a, b
mul_160:addc a, r5
mov r5, a
jnc $ + 7
inc r6
mov a, r6
jnz $ + 3
inc r7
ret

mul_32: ;mul[t0~3] X mult[4~7] --> r0~7;  ;46 bytes, <= 267T
;不需要验证了吧 :-P
;破坏 [mult(4~9)], Acc, B, F0, PSW
mov r0, mult + 4
mov r1, mult + 5
lcall mul_1632
mov r0, mult + 6
mov r1, mult + 7
mov mult + 4, r2
mov mult + 5, r3
mov mult + 6, r4
mov mult + 7, r5
mov mult + 8, r6
mov mult + 9, r7
lcall mul_1632
mov r0, mult + 4
mov r1, mult + 5
mov a, mult + 6
add a, r2
mov r2, a
mov a, mult + 7
addc a, r3
mov r3, a
mov a, mult + 8
addc a, r4
mov r4, a
mov a, mult + 9
sjmp mul_160

PK 使人进步：highgear 向网友刘前辈提出很简单的切磋请求

社区建设奖章

欢快之小溪

伴坛终老

技术奇才奖章

技术新星奖章

经常做客

核心会员奖章

坚毅之洋流

七世轮回

技术高手奖章

无冕之王奖章

沉静之湖泊

缘定三生

技术领袖奖章

精英会员奖章

奔腾之江水

十世金身

技术导师奖章

甘甜之泉水

希望之星奖章