螺杆转子刀具设计（2） matlab 解法

shouce

螺杆转子刀具设计 数学模型
  m: }3 E- Q* Z4 p, }) G) c
解如下方程：1679999.7*sin(x)-1679997.5*cos(y)*sin(x)-1499997.5*sin(y)*cos(x)-75000*sin(x)^3-179999.4000005*y*cos(y)*cos(x)+179999.4000005*y*sin(y)*sin(x)+75000*cos(y)*cos(x)*sin(x)=0
6 @/ p8 o" H' k求当x=-1.361时  y=?  
>> syms x y
5 x9 x' V1 O& u>> f=1679999.7*sin(x)-1679997.5*cos(y)*sin(x)-1499997.5*sin(y)*cos(x)-75000*sin(x)^3-179999.4000005*y*cos(y)*cos

(x)+179999.4000005*y*sin(y)*sin(x)+75000*cos(y)*cos(x)*sin(x)

5 l5 S$ q6 _* ]4 G' W$ j4 Of =

(7215543768789811*sin(x))/4294967296 - (2999995*cos(x)*sin(y))/2 - (3359995*cos(y)*sin(x))/2 - 75000*sin(x)^3 + 75000*cos(x)*cos(y)

- a1 X: y0 u# i1 ^# g*sin(x) - (6184732290414159*y*cos(x)*cos(y))/34359738368 + (6184732290414159*y*sin(x)*sin(y))/34359738368

/ \0 L3 h+ ~0 T>>  subs(f, x, -1.361)
+ d7 ~8 F7 s1 B. V% A$ z1 ]6 N* J% a- }
ans =
% O. B$ v( Q+ u2 Q' _
(3359995*sin(1361/1000)*cos(y))/2 - (2999995*cos(1361/1000)*sin(y))/2 - (7215543768789811*sin(1361/1000))/4294967296 + 75000*sin
+ h5 B# L  i+ _+ d
  e0 L0 V) a  j% L) r. a6 a(1361/1000)^3 - 75000*cos(1361/1000)*sin(1361/1000)*cos(y) - (6184732290414159*y*cos(1361/1000)*cos(y))/34359738368 -

: e$ F2 B5 m( B4 m  S0 ^& @/ _(6184732290414159*y*sin(1361/1000)*sin(y))/34359738368
: I4 L$ R9 R5 {  {- V

>>  a=0; b=1;
( p+ k" p& ~3 I2 xeps1=1e-8;eps2=1e-8;
N=300;
8 a" `0 V* V( ^1 M/ z; pf=@ (x) ((3359995*sin(1361/1000)*cos(x))/2 - (2999995*cos(1361/1000)*sin(x))/2 - (7215543768789811*sin(1361/1000))/4294967296 + 75000*sin(1361/1000)^3- 75000*cos(1361/1000)*sin(1361/1000)*cos(x) - (6184732290414159*x*cos(1361/1000)*cos(x))/34359738368 -(6184732290414159*x*sin(1361/1000)*sin(x))/34359738368);
Hfun=@Bisection;
[k,x,f_value]= feval(Hfun, f,a,b,eps1,eps2,N);

+ _; f8 t/ V: t# E6 e- h+ {( E1 e% `运行结果
/ |3 n7 N8 p8 D3 xk ,a ,b ,x, f
  1, 0.000000000, 1.000000000, 0.500000000,-352805.622314164,
$ Q+ f" M: }5 d9 S3 i0 K7 u3 F  2, 0.000000000, 0.500000000, 0.250000000,-92968.331400711,
  3, 0.000000000, 0.250000000, 0.125000000,-4146.780462183,
  4, 0.000000000, 0.125000000, 0.062500000,29178.936016433,
  5, 0.062500000, 0.125000000, 0.093750000,13458.364371693,
# u9 f9 Y/ K# o& K% Y  ?$ G  6, 0.093750000, 0.125000000, 0.109375000,4890.124291139,
  7, 0.109375000, 0.125000000, 0.117187500,430.092693351,
  8, 0.117187500, 0.125000000, 0.121093750,-1843.759437944,
/ |  |. o. K; N  9, 0.117187500, 0.121093750, 0.119140625,-703.184646645,
10, 0.117187500, 0.119140625, 0.118164063,-135.633470273,
11, 0.117187500, 0.118164063, 0.117675781,147.457778641,
12, 0.117675781, 0.118164063, 0.117919922,5.969190901,
13, 0.117919922, 0.118164063, 0.118041992,-64.817881139,
14, 0.117919922, 0.118041992, 0.117980957,-29.420780403,
15, 0.117919922, 0.117980957, 0.117950439,-11.724903562,
16, 0.117919922, 0.117950439, 0.117935181,-2.877633532,
2 S, c& }4 N* k1 C: [7 K1 I 17, 0.117919922, 0.117935181, 0.117927551,1.545834384,
18, 0.117927551, 0.117935181, 0.117931366,-0.665885649,
( w) _8 G) ?# q  k  j 19, 0.117927551, 0.117931366, 0.117929459,0.439977849,
2 O/ P6 @! G$ z7 `5 n& G. N0 j 20, 0.117929459, 0.117931366, 0.117930412,-0.112953030,
21, 0.117929459, 0.117930412, 0.117929935,0.163512627,
- l% n% Y7 k. n3 P+ [9 }. R) s: l 22, 0.117929935, 0.117930412, 0.117930174,0.025279853,
23, 0.117930174, 0.117930412, 0.117930293,-0.043836575,
8 Z. {# N  R3 f/ w$ J' j. I 24, 0.117930174, 0.117930293, 0.117930233,-0.009278357,
2 @( s. J7 e5 V5 R6 n 25, 0.117930174, 0.117930233, 0.117930204,0.008000749,
26, 0.117930204, 0.117930233, 0.117930219,-0.000638804,
27, 0.117930204, 0.117930219, 0.117930211,0.003680972,
>>

结论 当x=-1.361时  y=0.117930219
. O8 u+ f$ v" e' O. d8 @1 f按此方法   可依次  得出150对x与y的值

单车居士

西交大学生？？问邢老师

明月山河

你这个是代数方程？

cosxuan

算法？

shouce

cosxuan 发表于 2015-6-4 19:45
" e2 T+ v: y6 o算法？

6 b7 O: s! e1 O: K8 \3 T: q- o二分法        
% H2 v# M. u3 t还有一个程序
% Bisection.m
: L, H5 S' {0 zfunction[k,x,f_value]=Bisection(f,a,b,eps1,eps2,N)
fprintf('k ,a ,b ,x, f\n')
for k=1:N
    x=(a+b)/2;
5 p, C, }8 F+ i& ~3 Y9 D    f_value=f(x);
    fprintf( '%3d, %10.9f, %10.9f, %10.9f,%10.9f,\n'...
: E. q4 }/ |; s      ,k ,a ,b ,x, f_value)
    if abs(f_value)< eps1||0.5*(b-a)<eps2
    return
    else
        if f(x)*f(a)<0
' B5 g& B9 H% }+ f7 s          b=x;
        else
+ }' `, q0 u' }          a=x;
        end
        if k== N
7 z: L) z# |0 V; D          warning ('算法超出最大迭代数!')   
9 I) c/ a" l# t" K- q1 ~) N  `! V: D& K        end
    end

shouce

单车居士 发表于 2015-6-4 12:41
西交大学生？？问邢老师

4 n; t; E, T0 _7 ?  C3 g7 s这点问题   自己就可以解决
( R. h2 d$ m) J+ p$ }

shouce

明月山河 发表于 2015-6-4 19:12
0 }( K; B' m) P2 P( h7 i你这个是代数方程？

是的呀      不过matlab写出来的方程  括号很多  不好看   
$ V; A4 e8 ]' H0 W) F/ r# ~

shouce

明月山河 发表于 2015-6-4 19:12
8 S" f6 C3 H+ G你这个是代数方程？

; A9 y+ f/ T* R. d1 T. s
- X' Y8 N2 A" m2 L, F! z# X' A是超越方程      
# k+ ?  k+ R3 k9 d  n; v# x& [  C

[chāo yuè fāng chéng