annotate main/optim/inst/cdiff.m @ 9930:d30cfca46e8a octave-forge

optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
author carandraug
date Fri, 30 Mar 2012 15:14:48 +0000
parents 24d6a5cdedfe
children f82660cf7dee
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
9930
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
1 ## Copyright (C) 2002 Etienne Grossmann <etienne@isr.ist.utl.pt>
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
2 ##
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
3 ## This program is free software; you can redistribute it and/or modify it under
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
4 ## the terms of the GNU General Public License as published by the Free Software
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
5 ## Foundation; either version 3 of the License, or (at your option) any later
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
6 ## version.
2370
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
7 ##
9930
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
8 ## This program is distributed in the hope that it will be useful, but WITHOUT
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
9 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
10 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
11 ## details.
2370
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
12 ##
9930
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
13 ## You should have received a copy of the GNU General Public License along with
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
14 ## this program; if not, see <http://www.gnu.org/licenses/>.
2370
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
15
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
16 ## c = cdiff (func,wrt,N,dfunc,stack,dx) - Code for num. differentiation
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
17 ## = "function df = dfunc (var1,..,dvar,..,varN) .. endfunction
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
18 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
19 ## Returns a string of octave code that defines a function 'dfunc' that
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
20 ## returns the derivative of 'func' with respect to it's 'wrt'th
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
21 ## argument.
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
22 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
23 ## The derivatives are obtained by symmetric finite difference.
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
24 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
25 ## dfunc()'s return value is in the same format as that of ndiff()
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
26 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
27 ## func : string : name of the function to differentiate
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
28 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
29 ## wrt : int : position, in argument list, of the differentiation
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
30 ## variable. Default:1
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
31 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
32 ## N : int : total number of arguments taken by 'func'.
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
33 ## If N=inf, dfunc will take variable argument list.
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
34 ## Default:wrt
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
35 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
36 ## dfunc : string : Name of the octave function that returns the
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
37 ## derivatives. Default:['d',func]
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
38 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
39 ## stack : string : Indicates whether 'func' accepts vertically
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
40 ## (stack="rstack") or horizontally (stack="cstack")
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
41 ## arguments. Any other string indicates that 'func'
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
42 ## does not allow stacking. Default:''
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
43 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
44 ## dx : real : Step used in the symmetric difference scheme.
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
45 ## Default:10*sqrt(eps)
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
46 ##
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
47 ## See also : ndiff, eval, todisk
9930
d30cfca46e8a optim: upgrade license to GPLv3+ and mention on DESCRIPTION the other package licenses
carandraug
parents: 2370
diff changeset
48
2370
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
49 function c = cdiff (func,wrt,nargs,dfunc,stack,dx)
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
50
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
51 if nargin<2,
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
52 wrt = 1 ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
53 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
54 if nargin<3,
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
55 nargs = wrt ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
56 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
57 if nargin<4 || strcmp(dfunc,""),
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
58 dfunc = ["d",func] ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
59 if exist(dfunc)>=2,
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
60 printf(["cdiff : Warning : name of derivative not specified\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
61 " and canonic name '%s' is already taken\n"],\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
62 dfunc);
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
63 ## keyboard
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
64 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
65 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
66 if nargin<5, stack = "" ; end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
67 if nargin<6, dx = 10*sqrt(eps) ; end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
68
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
69 ## verbose = 0 ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
70 ## build argstr = "var1,..,dvar,...var_nargs"
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
71 if isfinite (nargs)
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
72 argstr = sprintf("var%i,",1:nargs);
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
73 else
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
74 argstr = [sprintf("var%i,",1:wrt),"...,"];
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
75 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
76
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
77 argstr = strrep(argstr,sprintf("var%i",wrt),"dvar") ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
78 argstr = argstr(1:length(argstr)-1) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
79
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
80 if strcmp("cstack",stack) , # Horizontal stacking ################
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
81
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
82 calstr = "reshape (kron(ones(1,2*ps), dvar(:))+[-dx*eye(ps),dx*eye(ps)], sz.*[1,2*ps])";
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
83 calstr = strrep(argstr,"dvar",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
84 calstr = sprintf("%s(%s)",func,calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
85
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
86 calstr = sprintf(strcat(" res = %s;\n",
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
87 " pr = prod (size(res)) / (2*ps);\n",
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
88 " res = reshape (res,pr,2*ps);\n",
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
89 " df = (res(:,ps+1:2*ps)-res(:,1:ps)) / (2*dx);\n"),
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
90 calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
91
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
92
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
93 elseif strcmp("rstack",stack), # Vertical stacking ##################
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
94
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
95 calstr = "kron(ones(2*ps,1),dvar)+dx*[-dv;dv]" ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
96 calstr = strrep(argstr,"dvar",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
97 calstr = sprintf("%s(%s)",func,calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
98
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
99 calstr = sprintf(strcat(" dv = kron (eye(sz(2)), eye(sz(1))(:));\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
100 " res = %s;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
101 " sr = size(res)./[2*ps,1];\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
102 " pr = prod (sr);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
103 " df = (res(sr(1)*ps+1:2*sr(1)*ps,:)-res(1:sr(1)*ps,:))/(2*dx);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
104 " scramble = reshape (1:pr,sr(2),sr(1))';\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
105 " df = reshape (df',pr,ps)(scramble(:),:);\n"),\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
106 calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
107 ## sayif(verbose,"cdiff : calstr='%s'\n",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
108 else # No stacking ########################
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
109 calstr = sprintf("%s (%s)",func,argstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
110 ## "func(var1,dvar%sdv(:,%d:%d),...,varN),"
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
111 ## calstr = strrep(calstr,"dvar","dvar%sdv(:,(i-1)*sz(2)+1:i*sz(2))")(:)';
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
112
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
113 calstr = strrep(calstr,"dvar","dvar%sdv")(:)';
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
114
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
115 ## func(..,dvar+dv(:,1:sz(2)),..) - func(..)
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
116 calstr = strcat(calstr,"-",calstr) ; ## strcat(calstr,"-",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
117 calstr = sprintf(calstr,"+","-") ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
118 tmp = calstr ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
119 ## sayif(verbose,"cdiff : calstr='%s'\n",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
120 calstr = sprintf(strcat(" dv = zeros (sz); dv(1) = dx;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
121 " df0 = %s;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
122 " sr = size (df0);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
123 " df = zeros(prod (sr),ps); df(:,1) = df0(:);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
124 " for i = 2:ps,\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
125 " dv(i) = dx; dv(i-1) = 0;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
126 " df(:,i) = (%s)(:);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
127 " end;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
128 " df ./= 2*dx;\n"
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
129 ),
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
130 calstr, tmp) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
131
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
132
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
133 ## sayif(verbose,"cdiff : calstr='%s'\n",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
134
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
135 ## "func(var1,reshape(dvar(1:NV,1),SZ1,SZ2),...,varN)," ,
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
136 ## "func(var1,reshape(dvar(1:NV,2),SZ1,SZ2),...,varN)," , ...
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
137 ## "func(var1,reshape(dvar(1:NV,NP),SZ1,SZ2),...,varN)"
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
138 ## sayif(verbose,"cdiff : calstr='%s'\n",calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
139 end
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
140 argstr = strrep (argstr, "...", "varargin");
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
141 calstr = strrep (calstr, "...", "varargin{:}");
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
142
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
143 c = sprintf(strcat("function df = %s (%s)\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
144 " ## Numerical differentiation of '%s' wrt to it's %d'th argument\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
145 " ## This function has been written by 'cdiff()'\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
146 " dx = %e;\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
147 " sz = size (dvar);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
148 " ps = prod (sz);\n",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
149 "%s",\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
150 "endfunction\n"),\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
151 dfunc,argstr,\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
152 func,wrt,\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
153 dx,\
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
154 calstr) ;
24d6a5cdedfe Changed the directory structure to match the package system
hauberg
parents:
diff changeset
155