Download Assignment 2 Solutions | Optimization | MATH 5750 and more Assignments Optimization Techniques in Engineering in PDF only on Docsity!
Optimization
Math 5750/
Assignment 2
1. Derivative approximation using finite difference approximation as well as central difference approxi-
mation.
(a) For the derivative approximation code see listing (1). See figure (1) for the error plots.
10 −20^10 −15^10 −10^10 −5^100
10 −
10 −
10 −
10 −
10 −
10 −
100
10 2
Forward Difference Central Difference
Figure 1: Forward Finite Difference log log plot
Listing 1: Derivative Approximation
r a n g e = 1 0. ˆ [ 0 : − 1 : − 2 0 ] ; f o r w a r d D i f f = ( exp ( repmat ( 1 , 1 , 2 1 ) + r a n g e ) −... repmat ( exp ( 1 ) , 1 , 2 1 ) ). / r a n g e ;
f E r r o r = abs ( f o r w a r d D i f f − exp ( 1 ) ) ;
c e n t r a l D i f f = ( exp ( repmat ( 1 , 1 , 2 1 ) + r a n g e ) −... exp ( repmat ( 1 , 1 , 2 1 ) − r a n g e ) ). /... ( 2. ∗ r a n g e ) ;
c E r r o r = abs ( c e n t r a l D i f f − exp ( 1 ) ) ;
l o g l o g ( range , f E r r o r , ’−o r ’ ,... range , c E r r o r , ’−sb ’ , ’ LineWidth ’ , 2 ) ;
h = legend ( ’ Forward D i f f e r e n c e ’ , ’ C e n t r a l D i f f e r e n c e ’ , 2 ) ; s e t ( h , ’ I n t e r p r e t e r ’ , ’ none ’ , ’ L o c a t i o n ’ , ’ SouthWest ’ ) ;
(b) It appears that both the forward difference approximation as well as the central difference approx-
imation adhere to the floting point error ≤ = O(10−^16 ) with the forward difference h∗^ = O(≤^1 /^2 )
and the central difference h∗^ = O(≤^1 /^3 ). As expected, the error for the forward difference again
increases at 10−^8 where the central difference error increases at ≈ 10 −^5.
2. Curve Fitting, Non-Linear Least Squares. Given data (t, y), find the non-linear least square solution
for:
y ≈ β 1 e
−λ 1 t
+ β 2 e
−λ 2 t
The new objective function φ would be:
min
β 1 ,β 2 ,λ 1 ,λ 2
∑^ m
j=
‖ β 1 e
−λ 1 tj
+ β 2 e
−λ 2 tj
− yj ‖
or in vector notation:
min
β 1 ,β 2 ,λ 1 ,λ 2
‖ β 1 e
−λ 1 T
+ β 2 e
−λ 2 T
− Y ‖
The gradient φ would be:
[
∂xi
]
β 1 , β 2 λ 1 , λ 2
)T
e
−λ 1 T
(β 1 e
−λ 1 T
+ β 2 e
−λ 2 T
− Y )
e−λ^2 T^ (β 1 e−λ^1 T^ + β 2 e−λ^2 T^ − Y )
−T e−λ^1 T^ (β 1 e−λ^1 T^ + β 2 e−λ^2 T^ − Y )
−T e−λ^2 T^ (β 1 e−λ^1 T^ + β 2 e−λ^2 T^ − Y )
When using BFGS (as required in the problem) it is only necessary to approximate the Hessian for
the first iteration. The subsequent iterations approximate the Hessian. As an apporixmation to the
Hessian we will use the identity matrix. Listing (2) shows the actual computation of the function as
well as its gradient.
Listing 2: Data Fitting Example
% f u n c t i o n [ f , g ,H] = r o s e n b r o c k ( x , p a r s ) % % e x p f i t f u n c t i o n f ( x ) = % \ d f r a c { 1 }{ 2 } \ p a r a l l e l \ b e t a 1 eˆ{−\lambda 1 T} + \ b e t a 2 eˆ{−\lambda 2 % T} − Y \ p a r a l l e l ˆ % % I n p u t s % x v e c t o r o f v a r i a b l e s ( \ b e t a 1 , \ b e t a 2 , \ lambda 1 , \ lambda 2 )
% p a r s s t r u c t u r e c o n t a i n i n g parameters % T, t e s t p o i n t s % Y, r e s u l t i n g p o i n t s % Outputs % f f ( x ) % % g r a d i e n t %eˆ{−\ lambda 1 T}(\ b e t a 1 eˆ{−\ lambda 1 T} + \ b e t a 2 eˆ{−\lambda 2 T}−Y) %eˆ{−\ lambda 2 T}(\ b e t a 1 eˆ{−\ lambda 1 T} + \ b e t a 2 eˆ{−\lambda 2 T}−Y) %−T eˆ{−\lambda 1 T}(\ b e t a 1 eˆ{−\lambda 1 T} + \ b e t a 2 eˆ{−\ lambda 2 T}−Y) %−T eˆ{−\lambda 2 T}(\ b e t a 1 eˆ{−\lambda 1 T} + \ b e t a 2 eˆ{−\ lambda 2 T}−Y) function [ f , g ,H] = e x p f i t ( x , p a r s ) i f ˜ i s f i e l d ( p a r s , ’T ’ ) error ( ’ M i s s i n g T v a l u e s ’ ) ; end ; i f ˜ i s f i e l d ( p a r s , ’Y ’ ) error ( ’ M i s s i n g Y v a l u e s ’ ) ; end ;
% pre−computation t o r e d u c e ops exp1 = exp(−x ( 3 ). ∗ p a r s. T ) ; exp2 = exp(−x ( 4 ). ∗ p a r s. T ) ; f u n c = ( x ( 1 ). ∗ exp1 ) + ( x ( 2 ). ∗ exp2 ) − p a r s .Y;
f =. 5 ∗ norm( f u n c ) ;
% g r a d i e n t computation i f ( nargout >= 2 ) g = [ exp1 ’ ∗ f u n c ;... exp2 ’ ∗ f u n c ;... exp1 ’ ∗ ( f u n c. ∗ −p a r s. T ) ;... exp2 ’ ∗ ( f u n c. ∗ −p a r s. T ) ] ;
i f ( nargout >= 3 ) % Hessian a p p r o x i m a t i o n H = eye ( s i z e ( x , 1 ) ) ; end ; end ;
Additionally, listing (3) displays the approximation code.
% S e t d e f a u l t kappa. i f ( nargin ==1), kappa = 1 ; end
% I n i t i a l i z a t i o n. h = 1/ n ; t = h / 2 : h : 1 ; c = h / ( 2 ∗ kappa ∗ sqrt ( pi ) ) ; d = 1 / ( 4 ∗ kappa ˆ 2 ) ;
% Compute t h e m a t r i x A. k = c ∗ t. ˆ ( − 1. 5 ). ∗ exp(−d. / t ) ; r = zeros ( 1 , length ( t ) ) ; r ( 1 ) = k ( 1 ) ; A = t o e p l i t z ( k , r ) ;
% Compute t h e v e c t o r s x and b. i f ( nargout >1) x = zeros ( n , 1 ) ; f o r i =1:n/ t i = i ∗20/ n ; i f ( t i < 2 ) x ( i ) = 0. 7 5 ∗ t i ˆ 2 / 4 ; e l s e i f ( t i < 3 ) x ( i ) = 0. 7 5 + ( t i −2)∗(3− t i ) ; e l s e x ( i ) = 0. 7 5 ∗ exp( −( t i − 3 ) ∗ 2 ) ; end end x ( n /2+1: n ) = zeros ( 1 , n / 2 ) ; b = A∗x ; end
Listing 5: Heat Demo
% Demo f o r h e a t .m from % R e g u l a r i z a t i o n Tools , Per C h r i s t i a n Hansen % h t t p : / /www2 .imm. d t u. dk /˜ pch / R e g u t o o l s / % t h i s f i l e models t h e n o t o r i o u s l y i l l −posed i n v e r s e h e a t e q u a t i o n % as a l i n e a r l e a s t s q u a r e s problem.
% use t h e s e parameters f o r t h e homework ( or h i g h e r ones ) n = 1 0 0 ; m = 5 0 ; kappa =1; [ A, b , x t r u e ] = h e a t ( n , kappa ) ; A = f u l l (A ) ; % c u t A and x i n o r d e r t o g e t an n x m r e c t a n g u l a r system A = A ( : , 1 :m) ; x t r u e=x t r u e ( 1 :m) ;
% p u t some n o i s e i n t h e d a t a b n o i s y = b + 1 e − 4 ∗randn ( n , 1 ) ;
f i g u r e ( 1 ) ; c l f ; plot ( x t r u e ) ; t i t l e ( ’ t r u e s o l u t i o n ’ ) ;
f i g u r e ( 2 ) ; c l f ; plot ( b ) ; t i t l e ( ’ n o i s e l e s s da t a ’ ) ;
f i g u r e ( 3 ) ; c l f ; plot ( b n o i s y ) t i t l e ( ’ n o i s y d a t a ’ ) ;
f i g u r e ( 4 ) ; c l f ; % i n matlab x=A\ b <=> x = p i n v (A)∗ b % s o l v e t h e problem w i t h t h e n o i s e l e s s data ‘ x1=A\b ; % s o l v e t h e problem w i t h t h e n o i s y d a t a x2=A\ b n o i s y ; plot ( 1 :m, x1 , 1 : m, x2 ) ; t i t l e ( ’ c o m p a r i s o n o f s o l u t i o n s ’ ) ;
% compute t h e SVD o f A f i g u r e ( 5 ) ; c l f ; [ U, S ,V] = svd (A ) ; plot ( diag ( S ) ) ; t i t l e ( ’ svd (A) ’ ) ;
