Matlab之Gauss-Newton法

haidaowang

这一题就是已知初始点，然后求某一点到三个圆的距离的和的最小值。使用的就是大名鼎鼎的Gauss-Newton方法。这里要求雅克比行列式，由函数 jacobian可以求出。所以这就减少了许多不必要的麻烦。
6 S, Z4 s7 k# E3 v. o4 |
5 W, k' a: g( [4 r' U+ u

代码如下：         
! ^1 `# U, z/ z- l: _' ^

( U3 V* a5 h/ }; k

• 
  # x9 v8 J6 ]. C( B' B* b9 w
  , R0 F# M5 G$ P4 G/ k% using Gauss-Newton method
  
  
• 
  % O. G. q# A: C7 r5 ?
  % initial point(0,0)
  
  ' s( a2 B! L- g7 ^2 ~
• 
  3 m/ m8 G! e/ g+ |+ R( W  u
  % 1、with circul point (0,1),(1,1),(0,-1)
  * l  u4 p- f$ u' |& t' g# p, L
  / z% U& ]1 f& {9 u
• 
  8 i  @* R4 D& ]6 D, `
  % 2、with circul point (-1,0),(1,1),(1,-1)
  
  1 P, R9 f) n7 t. g; f
• 
  & n% i8 L7 Z4 F8 W( b" @' [
  % its radius is 1
  3 P, ~. Q& G. B* L$ B1 W
  2 ^  Q* x' Z& Y8 k& \
• 
  
  % output :the sum of distance is the minist
  
  
• 
  - T  I( @3 `+ r' n( B& n( w
  function page_237_1
  
  
• 
  4 t1 P- p* b# v$ k
  format long
  
  ! `; O( |+ |! k6 a- G9 v, j
• 
    q9 c2 J0 O% d( o- e
  syms x y
  , b; N( x: X( M4 e+ {
  
• 
  ! m/ s, x% n" y7 b9 P( [+ n
  * M( V# b$ Q: D9 W7 i3 p" Ir1      = sqrt(x^2+(y-1)^2) - 1;
  
  + x0 o$ V+ l3 ~' d- b0 e
• 
  $ ^# s1 Q( n* z2 D( B
  9 d- n# u+ Y& u) t% Jr2      = sqrt((x-1)^2+(y-1)^2) - 1;
  
  3 \$ g7 `  P- I9 w0 i
• 
  ' ~* p" G- ]/ p. b
  r3      = sqrt(x^2+(y+1)^2) - 1;
  
  6 P+ I1 Y$ N0 u  i4 V2 Q, T
• 
  * \  S" Y' V7 L, M3 \
  ( `+ U, Y8 f: L' xrxk     = [r1;r2;r3];
  7 f% |; p4 p1 A& k/ Z: i
  
• 
  / m7 T- }$ O, A- Q$ ]6 S, W5 E
  # R( n0 \8 z0 x5 \* wDrx     = jacobian([r1;r2;r3]);
  5 ~& E1 y0 p" T% ^; I
  
• 
  
  DrxT    = Drx';
  # C: r/ X$ j' x! L. k8 Q) N$ u
  5 _; ~& Z+ o% R1 z! d
• 
  
  DrxTdrx = DrxT*Drx;                 %Dr(xk)TDr(xk)
  . c, v* ~6 `9 m! W1 B6 w
  6 A: Q$ |. t0 p% n: x, Q6 }
• 
  
  5 a4 X4 `. c& T* a9 {  @0 Q6 kDrxTrxk = -1*DrxT*rxk;              %Dr(xk)Tr(xk)
  , t- O3 m! D# O# v) R0 T
  
• 
  7 B+ r" a0 v2 H; ~5 t$ Z8 j
  4 {$ f0 Q$ z( [- m! {2 n) u8 F% S2 Y%vk     = DrxTrxk\DrxTdrx;
  * \) q! X% X! N1 o: i3 z& e
  
• 
  0 ?7 `* c$ A+ D% _' B( |
  0 ]8 x+ a. V% F1 Y, Q: j/ {- EDrxTdrx_= subs(DrxTdrx,{x,y},{0,0});
  # ]. P3 l: W6 M6 k# r3 |
  
• 
  
  DrxTrxk_= subs(DrxTrxk,{x,y},{0,0});
  
  / S7 H6 V$ T9 y8 K# k
• 
  , H8 A2 K( L) a; h
  xy      = zeros(2,1);
  
  ) a$ I7 A: R5 Z' Y0 W
• 
  3 d3 g8 J. V: [) N: O: W* R: j; K# ]
  5 g$ V" ?; s9 E% p  t% }9 m! Ofor k = 1:15
  
  8 Y3 l1 V' y" Y1 i; Z2 G
• 
  
      v       = DrxTdrx_\DrxTrxk_;
  
  , }. W: p8 H  x$ E: C& q0 o
• 
  
  ) v8 B9 H  _. |9 M  s    xy      = xy+v;      
  
  / \3 y7 P3 C4 p6 A; |# q- O
• 
  . P4 u( @: ?4 E  {2 q! R) \9 A+ Z
      DrxTdrx_= subs(DrxTdrx,{x,y},{xy(1),xy(2)});
  " N1 D6 E+ i- |( N( h) [5 ~* h1 H
  + I& `; Z2 d6 [- e1 y/ R
• 
  
      DrxTrxk_= subs(DrxTrxk,{x,y},{xy(1),xy(2)});
  ; e0 |. A/ M8 e  Y
  & T( q6 E) `9 c0 G1 K
• 
  
  . B! V$ n* s0 [/ p" T; nend
  7 z2 Z0 K; j2 ?5 n" S, s
  
• 
  ! }& \5 a2 |6 i* i
  disp('a with circular initial vector(0,0),Centers(0,1), (1,1), (0,?1)and all radii 1,the point is');
  
  ! x* P) y1 q( l7 o' m% S7 [6 w
• 
  
  2 D% }, O5 R6 [* Zvpa(xy,6)
  
  
• 
  + ?2 g5 h3 j/ |( J$ \7 {" f3 K& h
  

1 Z& n/ W# |. ]

NingW

Matlab之Gauss-Newton法