FindMinimum
- To: mathgroup at smc.vnet.net
- Subject: [mg65686] FindMinimum
- From: "Chris Chiasson" <chris.chiasson at gmail.com>
- Date: Sun, 16 Apr 2006 01:44:43 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
'lo Gurus,
If I wanted to make a routine for finding the minimum of a
multidimensional function, should I use FindMinimum for a 1 dimensional
minimization after picking a search direction? I am doing so now and
have run into problems where the precision of the position of the last
step is "polluting" the present step's 1 dimensional minimization call.
BTW, I know Mathematica already has functionality for minimizing
multidimensional functions.
Of course, I could set the precision of the last position to be the
same as the working precision (or higher), but that makes some of the
reported digits in the answer incorrect.
What should I be going if I want to code a certain type of minimization
algorithm into Mathematica? Is there an API for adding new methods a la
differential equations?
Here is some example code -- sorry for the density.
Thank you,
XVectorQ[xvec_] := MatchQ[xvec, {{_Symbol}..}]
XValueVectorQ[xvec_] := MatchQ[xvec, {{_?NumericQ}..}]
XRepValueVectorQ[xvec_] := MatchQ[xvec, {{_}..}]
vectorrules[(vec__)?XVectorQ, (valvec_)?XRepValueVectorQ] :=
MapThread[#1[[1]] -> #2[[1]] & , {vec, valvec}];
grad[f_, (xvec_)?XVectorQ] := (D[f, {#1, 1}] & ) /@ xvec;
hess[f_, (xvec_)?XVectorQ] := First[Outer[D, {f}, Flatten[xvec],
Reverse[Flatten[xvec]]]];
Options[FindMinimumSteepestDescent] = {MaxIterations -> 100,
AccuracyGoal -> Automatic, PrecisionGoal -> Automatic,
WorkingPrecision -> MachinePrecision, StepMonitor -> None};
FindMinimumSteepestDescent[f_, {(xvec_)?XVectorQ,
(xvalvecinput_)?XValueVectorQ}, (opts___)?OptionQ] :=
Module[{accgoal, sdcount = 0, fstar = f /. vectorrules[xvec,
xvalvecinput],
fval, grd = grad[f, xvec], maxiter = MaxIterations /. {opts} /.
Options[FindMinimumSteepestDescent], precgoal, S, stepminsol,
wcp = WorkingPrecision /. {opts} /.
Options[FindMinimumSteepestDescent] /. MachinePrecision ->
$MachinePrecision, xvalvec, xvalvecstar = xvalvecinput,
\[Alpha], \[Alpha]star},
{accgoal, precgoal} = {AccuracyGoal, PrecisionGoal} /. {opts} /.
Options[FindMinimumSteepestDescent] /. Automatic -> wcp/2;
While[(sdcount < maxiter && !Abs[fstar - fval] < 10^(-accgoal) +
Abs[fval]/10^precgoal) || sdcount == 0,
{fval, xvalvec} = {fstar, xvalvecstar}; ++sdcount;
S = -grd /. vectorrules[xvec, xvalvec]; stepminsol =
FindMinimum[f /. vectorrules[xvec, xvalvec + \[Alpha]*S],
{\[Alpha], 0},
WorkingPrecision -> wcp, AccuracyGoal -> accgoal,
PrecisionGoal -> precgoal]; {fstar, \[Alpha]star} =
{stepminsol[[1]], \[Alpha] /. stepminsol[[2]]};
xvalvecstar = xvalvec + \[Alpha]star*S; ReleaseHold[
ReplacePart[Hold[Block[Null, StepMonitor /. {opts} /.
Options[FindMinimumSteepestDescent]]],
vectorrules[xvec, xvalvecstar] /. Rule -> ReplaceMe, {1, 1}]
/.
ReplaceMe -> Set]]; If[sdcount >= maxiter,
Print["maximum number of iterations reached"]];
{fstar, vectorrules[xvec, xvalvecstar]}]
K[i_]=500+200*(5/3-i)^2;
W[j_]=50*j;
L[0]=10;
\[CapitalDelta]L[i_]=Sqrt[(X[i+1]-X[i])^2+(Y[i+1]-Y[i])^2]-L[0];
X[1]=0;
X[7]=60;
Y[1|7]=0;
varz=Most[Most[Rest[Rest[Flatten[({X[#1],Y[#1]}&)/@Range[7]]]]]]
PE=Sum[(1/2)*K[i]*\[CapitalDelta]L[i]^2,{i,6}]+Sum[W[i]*Y[i],{i,5}]
initrepz={X[i_]\[RuleDelayed]10*(i-1),Y[i_]\[Rule]0}
symbrepz={X[blah_]\[RuleDelayed]
ToExpression[StringJoin["X\[UnderBracket]",ToString[blah]]],
Y[blah_]\[RuleDelayed]ToExpression[
StringJoin["Y\[UnderBracket]",ToString[blah]]]}
FindMinimumSteepestDescent[
PE/.symbrepz,{List/@varz/.symbrepz,List/@varz/.initrepz},
WorkingPrecision\[Rule]64]