Design Example 1

dx1.m

global cv th rn;
cv = [ 0 0.25 -0.15];
th = [ 0  0    0];
rn = [ 1 1.5 1 ];

x = [ 0.25 -0.15]';
for (k=1:2)
   fprintf('step %d\n',k);
   [A fz] = calculate_derivatives(@sing1,x)
   fprintf('merit function: %g\n',norm(fz));
   s = -A\fz
   x = x + s
end

Output

step 1

A =

    0.5000   -0.5000
   -0.0292    0.1958


fz =

    0.1500
   -0.0183

merit function: 0.151116

s =

   -0.2425
    0.0575


x =

    0.0075
   -0.0925

step 2

A =

    0.5000   -0.5000
   -0.0302    0.0719


fz =

    0.0000
   -0.0034

merit function: 0.0034375

s =

    0.0825
    0.0825


x =

    0.0900
   -0.0100

sing1.m

function f = sing1(v);
global cv th rn;
cv(2)=v(1); % first surface curvature
cv(3)=v(2); % second surface curvature
yap = 5.0;
uco = 0.1;
scl = yap^2/2;
ya = parax([yap 0],cv,th,rn);
yc = parax([0 uco],cv,th,rn);
f(1) = -ya(4,2)/ya(1,1) - 1/20; % power error
w = ford(ya,yc,cv,th,rn,1);
f(2) = w(2)/scl;
f = f';

calculate_derivatives.m

function [A, fz] = calculate_derivatives(f,v)

fz = feval(f,v);
nv = length(v);
delta = 1e-5;
for j=1:nv
    vold = v(j);
    v(j) = vold+delta;
    fp = feval(f,v);
    v(j) = vold-delta;
    fn = feval(f,v);
    A(:,j)=(fp-fn)/(2*delta);
    v(j)=vold;
end

Maintained by John Loomis, last updated 6 Feb 2003