LINEBURG


<< . .

 14
( 16)



. . >>

0.2400
0.4800
0.4946
0.5092
0.5092
0.5092
0.7546
1.0000
1.0000
1.0000
1.0000
1.0000

project3_1(9,1)
f=
1.5792
g=
-1.1642e-10
xm =
3.0000 5.0000
3.0007 -3.8736
0.2249 -2.4065
0.2249 -2.4065
-3.4249 -1.0539
-3.4249 -1.0539
-3.4249 -1.0539
-3.4249 -1.0539
-3.4249 -1.0539
-3.4249 -1.0539
s=
0
0.3164
0.4882
0.4882
178 OPTIMAL CONTROL MODELS IN FINANCE

1.0000
1.0000
1.0000
1.0000
1.0000
1.0000

project3_1(9,2)
f=
0.9830
g=
2.2104e-10
xm =
3.0000 5.0000
3.7478 -2.6935
-0.1429 -3.2994
-0.4100 -3.0087
-0.6527 -2.7209
-0.6527 -2.7209
-0.6527 -2.7209
-2.8224 -0.6915
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
-2.8407 0.4891
s=
0
0.2322
0.4644
0.4814
0.4983
0.4983
0.4983
0.7492
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
1.0000
APPENDIX B: Some Computation Results 179

1.0000
1.0000
This page intentionally left blank
Appendix C
Differential Equation Solver from the SCOM
Package




function xk=nqq(pd,nx,nu,nn,fil,ma,t,it,hs,um,xm,lm,ps)
% Calling the Runge-Kutta integration

yin=fil(1,:);xk(1,:)=yin;
while it < nn+1
[y2,it2,t2]=rqq(pd,ma,t,it,hs,yin,nx,nu,nn,um,xm,lm,ps);
xk(it+1,:)=y2;
it=it2;t=t2;
yin=y2;
end



function [y1,it,t]=rqq(pd,ma,t,it,hs,yin,nx,nu,nn,um,xm,lm,ps)
% Runge-Kutta stages

fp=zeros(ma,1)™; tt=zeros(ma,1)™;
P=0; q=1;
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
p=0.5;q=2;t=t+0.5*hs;
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
t=t+0.5*hs;p=1;q=1;
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
it=it+sign(hs);
it;
y1=yin+tt/6;
182 OPTIMAL CONTROL MODELS IN FINANCE


function tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps)
% Calling the function

z=yin+p*hs*fp;
ff=feval(pd,t,it,z,yin,hs,um,xm,lm,ps);
fp=ff;
tz=[tt+q*hs*fp;fp] ;



function ix=jqq(vv,hh,xi)
% Linear interpolation

global um xm jm
nn=length(vv™)-1;
nz=-1/hh;
if xi < 0, xi = 0; end
fr=-nz*mod(xi,hh);
bb=ceil(xi*nz)+1;

if bb>1, vw=vv(bb-1,:);
else vw=vv(1);
end

ix=(1-fr)*vv(bb,:) + fr*vw;
Appendix D
SCOM Package




subs=cell(1,9);
subs={™dnx™ , ™dnj™ , ™dnf™ , ™dnc™};

par=[2, 2, 20, 0, 0, 1, 6, 0.2, 0.1, 0.1, 0.75];
% nx, nu, nn, npa, grad, c, T, r, p, d , b=k/r\
% nx = number of states, nu = number of controls,
nn = number of subdivisions
% npa = 0 (reserved), grad = 1 if gradients calculated,
otherwise 0
% The remaining parameters, specific to the problem,
are passed to the subroutines

xinit=[2 1]; nn=par(3); % Initial values for the state(s)
u0=zeros(nn,2); % Starting values for computing the control
ul=zeros(nn,2); % Lower bound (vector or matrix) for the control
uu=ones(nn,2); % Upper bound (vector or matrix) for the control

figure
Control=constr(™fqq™ ,u0, [] ,ul,uu, [] ,par,subs,xinit)
% Calls constr package

[Objective,Constraint,State,Integral]= cqq(Control,par,subs,xinit)
% Computes optimal state, etc., from the optimal control
got from constr

% If gradients are calculated, then used instead
% [Objective, Constraint, State, Integral, Costate,Gradient]
=cqq(Control,par,subs,xinit)
% The following lines plot one state component against another,
and plot two control components, ehich are step-functions,
against time
184 OPTIMAL CONTROL MODELS IN FINANCE

plot(State(:,1),State(:,2),™x-™)
xlabel(™State 1™)
ylabel(™State 2™)

figure
t2=[0:0.001:0.999];
plot(t2,Control(floor(nn*t2)+1,1),™r™)

hold on
plot(t2,Control(floor(nn*t2)+1,2),™g™)
%The following codes comprise the SCOM package;
the user does not alter them

% Get function values as required by constr



function [f,g]=fqq(uz,ps,q,xinit) % Input data to the subroutine

if ps(5)==1
[a1,a2,a3,a4,a5,a6]=cqq(uz,ps,q,xinit);
% Output data
elseif ps(5)==0
[a1,a2,a3,a4]=cqq(uz,ps,q,xinit);
end

f=a1;
g=a2;



function [df,dg]=gqq(uu,ps,q,xinit)
% Get gradient values as required by constr

[a1,a2,a3,a4,a5,a6]=cqq(uu,ps,q,xinit);
df=a6™;
nn=ps(3);
pk=char(q(5));
dg=feval(pk,uu,nn);



function [f,g,xm,f0,lm,gr]=cqq(um,ps,q,xinit)
% Solve differential equations

nx=ps(1);nu=ps(2);nn=ps(3);npa=ps(4); rec=0;npar=1;
xm=zeros(nn+1,nx); xm(1,:)=xinit; lm=zeros(nn+1,nx);
ma=nx;t=0;it=1;hs=1/nn;
px=char(q(1));
xm=nqq(px,nx,nu,nn,xm,ma,t,it,hs,um,xm,lm,ps); %compute state
ma=1;t=0;it=1;hs=1/nn;
zz=zeros(nn+1,1); zz(1,:)=0;
APPENDIX D: SCOM Package 185

pj=char(q(2));

jm=nqq(pj,nx,nu,nn,zz,ma,t,it,hs,um,xm,lm,ps);%compute integral
xf=xm(nn+1,:) ;
pf=char(q(3));
f=jm(nn+1) + feval(pf,xf,um,xm,ps); %objective
pc=char(q(4));

for ii=1:nn
g(ii)=feval(pc,ii,hs,um,xm,lm,ps) ; %control constraint
end

f0=jm(nn+1);
%final state

if ps(5)==1
%if gradients are supplied

ma=nx;t=1;it=nn;hs=-1/nn;
lm=zeros(nn+1,nx);
pa=char(q(8));
lm(nn+1,:)=feval(pa,nn,xf,um,xm,ps);
pq=char(q(6));
lm=lqq(pq,nx,nu,nn,lm,ma,t,it,hs,um,xm,lm,ps);
%costate



function xk=nqq(pd,nx,nu,nn,fil,ma,t,it,hs,um,xm,lm,ps)
%Runge-Kutta (organize steps for fourth order RK
integration of DEs)

yin=fil(1,:);xk(1,:)=yin;

while it < nn+1
[y2,it2,t2]=rqq(pd,ma,t,it,hs,yin,nx,nu,nn,um,xm,lm,ps);
xk(it+1,:)=y2;
it=it2;t=t2;
yin=y2;
end

function [y1,it,t]=rqq(pd,ma,t,it,hs,yin,nx,nu,nn,um,xm,lm,ps)
% Runge-Kutta : increments

fp=zeros(ma,1)™; tt=zeros(ma,1)™;
P=0; q=1;
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
p=0.5;q=2;t=t+0.5*hs;

tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
186 OPTIMAL CONTROL MODELS IN FINANCE

tt=tz(1,:); fp=tz(2,:);
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);
tt=tz(1,:); fp=tz(2,:);
t=t+0.5*hs;p=1;q=1;
tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps);

tt=tz(1,:); fp=tz(2,:);
it=it+sign(hs);
y1=yin+tt/6;

function tz=kqq(pd,ma,t,it,hs,p,q,fp,yin,tt,um,xm,lm,ps)
% Runge-Kutta : get function values

z=yin+p*hs*fp;
ff=feval(pd,t,it,z,yin,hs,um,xm,lm,ps);
fp=ff;
tz=[tt+q*hs*fp;fp] ;




function xk=lqq(pd,nx,nu,nn,fil,ma,t,it,hs,um,xm,lm,ps)
%Runge-Kutta for adjoint differential equation
(time reversed)

yin=fil(nn+1,:);xk(nn+1,:)=yin;

while it > 0 %< nn+1
[y2,it2,t2]=rqq(pd,ma,t,it,hs,yin,nx,nu,nn,um,xm,lm,ps);
xk(it,:)=y2;
it=it2; t=t2;
yin=y2;
end



function ix=iqq(vv,hh,xi)
% Linear interpolation (forwards)
nn=length(vv™)-1;
nz=1/hh;
fr=nz*mod(xi,hh);
bb=floor(xi*nz)+1;

if bb <= nn, vu=vv(bb+1,:);
else vu=vv(nn);
end

ix=(1-fr)*vv(bb,:) + fr*vu;



function ix=jqq(vv,hh,xi)
% Linear interpolation (backwards)
APPENDIX D: SCOM Package 187

global um xm jm
nn=length(vv™)-1;
nz=-1/hh;

if xi < 0, xi = 0;
end
fr=-nz*mod(xi,hh);
bb=ceil(xi*nz)+1;

if bb>1, vw=vv(bb-1,:);
else vw=vv(1);
end

ix=(1-fr)*vv(bb,:) + fr*vw;

Note that linear interpolation of the state is required in integrating the objective function, and
in solving the adjoint differential equation. The latter is solved backwards in time t, so needs
backwards interpolation.
This page intentionally left blank
Appendix E
Format of Problem Optimization




The control function approximates the control function by a step-function, dividing
[0,1] into nn equal subintervals. (Because the dynamic equation has a smoothing effect, set
function controls are usually a sufficient approximation.) Function values (and often also gra-
dients with respect to control) for the objective function are obtained by solving differential
equations. They are then supplied to the optimization program constr in MATLAB™s Optimiza-
tion Toolbox. The computation is considerably faster if gradients are supplied, but this is not
suitable for some problems, especially if the functions are not well behaved outside the feasible
region of the problem. If gradients are used, then the adjoint differential equation is required.
If a constraint is required, then it is added to s a positive penalty parameter. The user
must supply a calling program (defining all parameters), and user subroutines for the functions
of the given control problem. This use of MATLAB minimizes the amount of programming
required for a given control problem, since MATLAB handles matrix calculations and pass-
ing of parameters very effectively. (But another programming language may be preferred for a
large control problem, if faster computation is needed.) The acronym SCOM stands for step-
function control optimization on Macintosh, since the intention was to compute optimal control
on a desktop computer, rather than on a mainframe or workstation. Note that MATLAB on a
Windows computer does the same job.

<< . .

 14
( 16)



. . >>

Copyright Design by: Sunlight webdesign