Matlab之Gauss-Newton法

haidaowang

这一题就是已知初始点，然后求某一点到三个圆的距离的和的最小值。使用的就是大名鼎鼎的Gauss-Newton方法。这里要求雅克比行列式，由函数 jacobian可以求出。所以这就减少了许多不必要的麻烦。

代码如下：         

• 
    Q' L: Z; @; _0 }4 p' E/ H
  % using Gauss-Newton method
  
  ; o! F6 {# ~3 e, B
• 
  
  ) ], a/ h1 c+ y% initial point(0,0)
  6 v  J2 m4 P  }  |* m) n, A
  
• 
  3 _) z6 Y0 g# T1 y# t* n
  * g6 G% C3 v' D/ s% 1、with circul point (0,1),(1,1),(0,-1)
  
  2 |" C2 C3 j" G0 y$ T. q
• 
  
  ! D( e" B% O9 z5 L9 M, E% 2、with circul point (-1,0),(1,1),(1,-1)
  
  : W. d( q' [2 H6 k1 v$ D4 P, x
• 
  
  8 m3 ]6 x- D6 v5 R- \% its radius is 1
  
  
• 
  - Z& J( z/ a: P- `9 p! E
  . M6 J4 Y8 a( ?- E" i+ o% output :the sum of distance is the minist
  
  ) A  M$ T# ~# H  `2 w# U6 B1 X% v
• 
  ! }( q- M" \  b: \, d0 q6 G
  & C' G5 P' X" S8 o3 Ifunction page_237_1
  ( q6 ]/ {4 O* ~8 @. _% A4 B0 ~1 _
  7 |, [2 y5 Z3 r
• 
  
  8 d; k/ U; f5 e( `( P) |format long
  7 A7 A0 H+ ]  u6 ^& @  L- R+ j2 y7 [9 `
  3 H0 J5 }  r* T( `; u6 T+ |
• 
  1 n( D2 b7 Q) W% B* s
  ) k' M$ x! G* `syms x y
  
  
• 
  
  ) U5 x( ?1 r: ~( |8 ir1      = sqrt(x^2+(y-1)^2) - 1;
  
  
• 
  2 y) B3 j$ B8 P
  1 f  p# s7 C$ cr2      = sqrt((x-1)^2+(y-1)^2) - 1;
  
  1 k6 k. {. m1 p, {# G# ^4 |7 T
• 
  4 a: `7 v  K( M7 B8 V
  5 p$ X2 U; `2 s9 W5 Pr3      = sqrt(x^2+(y+1)^2) - 1;
  
  * D6 L7 j9 K8 Q+ _- `9 l+ I4 U1 C# O
• 
  ' @+ E4 a7 b. O4 l" N
  rxk     = [r1;r2;r3];
  - O3 \8 Y9 M+ V
  
• 
  ' Q3 |% n$ w1 q5 r/ p
  ; Z  n3 j$ T/ @6 WDrx     = jacobian([r1;r2;r3]);
  
  8 Z6 g4 O3 L5 ^, }) c1 b6 i
• 
  
  DrxT    = Drx';
  # ^. I8 m7 g" G: b/ U; V6 [7 j  H
    P+ g! z9 R; {6 W
• 
  
  DrxTdrx = DrxT*Drx;                 %Dr(xk)TDr(xk)
  
  
• 
  + T- N2 K4 F- \" }/ g) i
  0 ~$ ^- a. E2 v0 R# }DrxTrxk = -1*DrxT*rxk;              %Dr(xk)Tr(xk)
  * m  x" Z( a& n
  
• 
  
  : u4 F0 [! f( `/ c%vk     = DrxTrxk\DrxTdrx;
  
  
• 
  
  3 j: y' \4 c) J8 ODrxTdrx_= subs(DrxTdrx,{x,y},{0,0});
  ' f% \+ L& \+ {' [
  
• 
  
  DrxTrxk_= subs(DrxTrxk,{x,y},{0,0});
  7 E2 B9 d& ]. P" n3 X# V* M: x
  
• 
  
  / W8 \# l4 `9 ?1 v3 Dxy      = zeros(2,1);
  
  
• 
  7 n1 |+ C5 Q0 q% C/ s
  9 Q3 Y' S! c9 j) E! {for k = 1:15
  
  . e2 ]3 S* F0 ~+ |8 B) F& ^
• 
  
      v       = DrxTdrx_\DrxTrxk_;
  ( n. g  I* T1 K8 P1 V1 ~0 E
  
• 
  
      xy      = xy+v;      
  ! `5 N6 J4 s" D
  , \( ?# s( x7 C( Y. i: q
• 
  
  1 V9 q0 u5 O9 W# H& |    DrxTdrx_= subs(DrxTdrx,{x,y},{xy(1),xy(2)});
  
  
• 
  
      DrxTrxk_= subs(DrxTrxk,{x,y},{xy(1),xy(2)});
  
  3 c: \$ j! i) L/ |% V1 s' w
• 
  0 k* `6 S' m* k0 ?5 H$ Y+ u# ?
  end
  
  
• 
  
  ( ?2 d3 g7 j5 ^. A9 ^% Q, z9 ydisp('a with circular initial vector(0,0),Centers(0,1), (1,1), (0,?1)and all radii 1,the point is');
  " U$ ~1 u8 m& @9 P
  ( n6 L9 l  D# x) t3 \6 D3 p
• 
  * P: W& b  t9 c( ^) O
  ' k* n' o) N1 t5 w  D; ^vpa(xy,6)
  * z6 T' b  I9 e7 [4 c8 e
  
• 
  + g2 z' b) u) c2 [6 g5 b; ~
  

NingW

Matlab之Gauss-Newton法