* Bug #9208 fixed - Optimization: added three optional output arguments to optim() 59/12459/6
Paul BIGNIER [Fri, 6 Sep 2013 14:31:06 +0000 (16:31 +0200)]
'iters' and 'evals', to retrieve #iterations and #evaluations,
and 'err' for a termination indicator.

Change-Id: Ief8cdceff4ca88c5cd2e2133326a5674fc1565ca

SEP/INDEX
SEP/SEP_105_optim_output_arguments.odt [new file with mode: 0644]
scilab/CHANGES_5.5.X
scilab/modules/optimization/help/en_US/optim.xml
scilab/modules/optimization/sci_gateway/fortran/sci_f_optim.f
scilab/modules/optimization/tests/nonreg_tests/bug_9208.dia.ref [new file with mode: 0644]
scilab/modules/optimization/tests/nonreg_tests/bug_9208.tst [new file with mode: 0644]

index 337538a..4aa290c 100644 (file)
--- a/SEP/INDEX
+++ b/SEP/INDEX
@@ -98,3 +98,6 @@ SEP #099: Matrices of strings as input argument for clear
 SEP #100: Extends history size
 SEP #101: File browser
 SEP #102: members function.
+SEP #103: Localization : allow multiple domains and add addlocalizationdomain function.
+SEP #104: /* RESERVED */
+SEP #105: Three new optional output arguments to optim().
diff --git a/SEP/SEP_105_optim_output_arguments.odt b/SEP/SEP_105_optim_output_arguments.odt
new file mode 100644 (file)
index 0000000..1f8006b
Binary files /dev/null and b/SEP/SEP_105_optim_output_arguments.odt differ
index bd856ef..9c880a3 100644 (file)
@@ -383,6 +383,9 @@ Bug Fixes
 
 * Bug #9109 fixed - nfreq tagged as obsolete.
 
+* Bug #9208 fixed - Added three optional output arguments to optim(),
+                    to retrieve #iterations, #evaluations and a termination indicator.
+
 * Bug #9385 fixed - The type checking in trigonometric functions has been added.
 
 * Bug #9394 fixed - is_params recognized "plist" as an existing field.
index cd4cd4e..3810fd1 100644 (file)
@@ -25,6 +25,9 @@
             [fopt, xopt] = optim(...)
             [fopt, xopt, gopt] = optim(...)
             [fopt, xopt, gopt, work] = optim(...)
+            [fopt, xopt, gopt, work, iters] = optim(...)
+            [fopt, xopt, gopt, work, iters, evals] = optim(...)
+            [fopt, xopt, gopt, work, iters, evals, err] = optim(...)
         </synopsis>
     </refsynopsisdiv>
     <refsection>
                             "ar",nap,iter
                             "ar",nap,iter,epsg
                             "ar",nap,iter,epsg,epsf
-                            "ar",nap,iter,epsg,epsf,epsx            
+                            "ar",nap,iter,epsg,epsf,epsx
                         </screen>
                     </para>
                     <para>where:</para>
                                     </para>
                                 </listitem>
                             </itemizedlist>
-                            
                         </listitem>
                         <listitem>
                             <para>
                     </para>
                 </listitem>
             </varlistentry>
+            <varlistentry>
+                <term>iters</term>
+                <listitem>
+                    <para>
+                        scalar, the number of iterations that is displayed when <literal>imp=2</literal>.
+                    </para>
+                </listitem>
+            </varlistentry>
+            <varlistentry>
+                <term>evals</term>
+                <listitem>
+                    <para>
+                        scalar, the number of <literal>cost</literal> function evaluations
+                        that is displayed when <literal>imp=2</literal>.
+                    </para>
+                </listitem>
+            </varlistentry>
+            <varlistentry>
+                <term>err</term>
+                <listitem>
+                    <para>
+                        scalar, a termination indicator.
+                        The success flag is <literal>9</literal>.
+                        <literal>err=1</literal>: Norm of projected gradient lower than...
+                        <literal>err=2</literal>: At last iteration f decreases by less than...
+                        <literal>err=3</literal>: Optimization stops because of too small variations for x.
+                        <literal>err=4</literal>: Optim stops: maximum number of calls to f is reached.
+                        <literal>err=5</literal>: Optim stops: maximum number of iterations is reached.
+                        <literal>err=6</literal>: Optim stops: too small variations in gradient direction.
+                        <literal>err=7</literal>: Stop during calculation of descent direction.
+                        <literal>err=8</literal>: Stop during calculation of estimated hessian.
+                        <literal>err=9</literal>: End of optimization, successful completion.
+                        <literal>err=10</literal>: End of optimization (linear search fails).
+                    </para>
+                </listitem>
+            </varlistentry>
         </variablelist>
     </refsection>
     <refsection>
             constrained nonlinear optimization problems:
         </para>
         <screen>min f(x)
-            binf &lt;= x &lt;= bsup      
+            binf &lt;= x &lt;= bsup
         </screen>
         <para>
             where <literal>binf</literal> is the lower bound and
                                 <para>with the following declarations:</para>
                                 <screen>integer ind,n ti(*)
                                     double precision x(n),f,g(n),td(*)
-                                    real tr(*)      
+                                    real tr(*)
                                 </screen>
                                 <para>
                                     The argument <literal>ind</literal> is described
             [fopt, xopt, gopt] = optim(cost, "b", [-1; 0; 2], [0.5; 1; 4], x0, "gc", "ar", 3, 100, 1e-6, 1e-6, [1e-3; 1e-3; 1e-3])
             
             // Additionnal messages are printed in the console.
-            [fopt, xopt] = optim(cost, x0, imp = 3)    
+            [fopt, xopt] = optim(cost, x0, imp = 3)
         </programlisting>
     </refsection>
     <refsection>
             endfunction
             
             x0 = [1; -1; 1];
-            [fopt, xopt] = optim(cost, x0, imp = -1)   
+            [fopt, xopt] = optim(cost, x0, imp = -1)
         </programlisting>
         <para>The previous script produces the following output.</para>
         <screen>--&gt;[fopt, xopt] = optim(cost, x0, imp = -1)
             2.
             3.
             fopt  =
-            0.    
+            0.
         </screen>
         <para>In the following example, we solve the Rosenbrock test case. For
             each iteration of the algorithm, we plot the current value of x into a 2D
             plot(xopt(1) , xopt(2) , "r*")
             
             // 5. Plot the optimization process, during optimization
-            [fopt, xopt] = optim ( rosenbrockPlot , x0 , imp = -1)    
+            [fopt, xopt] = optim ( rosenbrockPlot , x0 , imp = -1)
         </programlisting>
         <scilab:image>
             function [f, g, ind]=rosenbrock(x, ind)
             
             x0 = [-1.2 1.0];
             
-            [ fopt , xopt ] = optim ( rosenbrockCost , x0 )    
+            [ fopt , xopt ] = optim ( rosenbrockCost , x0 )
         </programlisting>
         <para>In the following example, we use the derivative function to solve
             Rosenbrock's problem. Given that the step computation strategy is not the
             endfunction
             
             x0 = [-1.2 1.0];
-            [fopt , xopt] = optim ( rosenbrockCost2 , x0 )    
+            [fopt , xopt] = optim ( rosenbrockCost2 , x0 )
         </programlisting>
     </refsection>
     <refsection>
             c = 3.0;
             d = 4.0;
             costf = list ( costfunction , a , b , c, d );
-            [fopt , xopt] = optim ( costf , x0 , imp = 2)    
+            [fopt , xopt] = optim ( costf , x0 , imp = 2)
         </programlisting>
         <para>In complex cases, the cost function may have so many parameters,
             that having a function which takes all arguments as inputs is not
             parameters.c = c;
             parameters.d = d;
             costf = list ( costfunction , parameters );
-            [fopt , xopt] = optim ( costf , x0 , imp = 2)    
+            [fopt , xopt] = optim ( costf , x0 , imp = 2)
         </programlisting>
         <para>In the following example, the parameters are defined before the
             optimizer is called. They are directly used in the cost function.
             b = 2.0;
             c = 3.0;
             d = 4.0;
-            [ fopt , xopt ] = optim ( costfunction , x0 , imp = 2 )   
+            [ fopt , xopt ] = optim ( costfunction , x0 , imp = 2 )
         </programlisting>
         <para>While the previous example perfectly works, there is a risk that the
             same variables are used by some internal function used by
             b = 2.0;
             c = 3.0;
             d = 4.0;
-            [ fopt , xopt ] = optim ( costfunction , x0 , imp = 2 )    
+            [ fopt , xopt ] = optim ( costfunction , x0 , imp = 2 )
         </programlisting>
         <para>While the previous example perfectly works, there is a risk that the
             same variables are used by some internal function used by
             [ f , g , index ] = rosenbrock ( x0 , index );
             gnd = derivative ( simplified , x0.' );
             mprintf("Exact derivative:[%s]\n" , strcat ( string(g) , " " ));
-            mprintf("Numerical derivative:[%s]\n" , strcat ( string(gnd) , " " ));    
+            mprintf("Numerical derivative:[%s]\n" , strcat ( string(gnd) , " " ));
         </programlisting>
         <para>The previous script produces the following output. Obviously, the
             difference between the two gradient is enormous, which shows that the
             wrong formula has been used in the gradient.
         </para>
         <programlisting role="example">      Exact derivative:[-638 -88]
-            Numerical derivative:[-215.6 -88]   
+            Numerical derivative:[-215.6 -88]
         </programlisting>
     </refsection>
     <refsection>
             //solve the problem
             x0 = [40; 10; 50];
             p = 100;
-            [f, xo, go] = optim('rosenf', x0, 'td', p)    
+            [f, xo, go] = optim('rosenf', x0, 'td', p)
         </programlisting>
     </refsection>
     <refsection>
             //
             // Solve the problem
             //
-            [f, x, g] = optim('rosenf', x0, 'in');    
+            [f, x, g] = optim('rosenf', x0, 'in');
         </programlisting>
     </refsection>
     <refsection>
             'c     if n&lt;=2 returns ind=0'
             'c     f.bonnans, oct 86'
             '      implicit double precision (a-h,o-z)'
-            [etc...]    
+            [etc...]
         </programlisting>
     </refsection>
     <refsection role="see also">
index 404c5be..e0e7540 100644 (file)
@@ -885,7 +885,7 @@ c     sauvegarde de g
       if(lhs1.eq.3) goto 320
 c     
 c     sauvegarde de tv (tableau interne a l' optimiseur - pour hot start
-      if (lhs1.eq.4) then
+      if (lhs1.ge.4) then
          istv=0
          if(ialg.eq.1) istv=1
          if (istv.eq.0) then
@@ -934,6 +934,45 @@ c     if (l+ntv+nitv+1.ge.sadr(lizs)) then
  317     continue
       endif
       lstk(top+1)=l + ntv + nitv
+      if(lhs1.eq.4) goto 320
+c     
+c     sauvegarde de niter
+      top2=top2 + 1
+      lstk(top2)=lstk(top+1)
+      il=iadr(lstk(top2))
+      istk(il)=1
+      istk(il+1)=1
+      istk(il+2)=1
+      istk(il+3)=0
+      l=sadr(il+4)
+      stk(l)=itmax
+      lstk(top+1)=l+1
+      if(lhs1.eq.5) go to 320
+c     
+c     sauvegarde de napm
+      top2=top2 + 1
+      lstk(top2)=lstk(top+1)
+      il=iadr(lstk(top2))
+      istk(il)=1
+      istk(il+1)=1
+      istk(il+2)=1
+      istk(il+3)=0
+      l=sadr(il+4)
+      stk(l)=napm
+      lstk(top+1)=l+1
+      if(lhs1.eq.6) go to 320
+c     
+c     sauvegarde de indopt
+      top2=top2 + 1
+      lstk(top2)=lstk(top+1)
+      il=iadr(lstk(top2))
+      istk(il)=1
+      istk(il+1)=1
+      istk(il+2)=1
+      istk(il+3)=0
+      l=sadr(il+4)
+      stk(l)=indopt
+      lstk(top+1)=l+1
 c     
 c     sauvegarde de izs et dzs
  320  if (lhs.eq.lhs1) go to 350
@@ -1012,4 +1051,3 @@ c     commentaires finaux
       return
       end
 c     --------------------------
-
diff --git a/scilab/modules/optimization/tests/nonreg_tests/bug_9208.dia.ref b/scilab/modules/optimization/tests/nonreg_tests/bug_9208.dia.ref
new file mode 100644 (file)
index 0000000..870c244
--- /dev/null
@@ -0,0 +1,51 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2013 - Scilab Enterprises - Paul Bignier
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+//
+// <-- CLI SHELL MODE -->
+//
+// <-- Non-regression test for bug 9208 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=9208
+//
+// <-- Short Description -->
+// Added three optional output arguments to optim(), to retrieve #iterations,
+// #evaluations and a termination flag.
+//
+// Example 1
+a  = 1.0;
+b  = 2.0;
+c  = 3.0;
+d  = 4.0;
+x0 = [1 1];
+function [f, g, ind] = costfunction(x, ind, a, b, c, d)
+    f = a * ( x(1) - c ) ^2 + b * ( x(2) - d )^2;
+    g(1) = 2 * a * ( x(1) - c );
+    g(2) = 2 * b * ( x(2) - d );
+endfunction
+costf = list ( costfunction , a , b , c, d );
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 );
+assert_checkequal([iters evals err], [10 11 1]);
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 , "ar",nap=5 );
+assert_checkequal([iters evals err], [4 5 4]);
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 , "ar",nap=100,iter=5 );
+assert_checkequal([iters evals err], [6 7 5]);
+// Example 2
+x0 = [-1.2 1.0];
+function f = rosenbrock(x)
+    f = 100.0 *(x(2)-x(1)^2)^2 + (1-x(1))^2;
+endfunction
+function [f, g, ind] = rosenbrockCost2(x, ind)
+    f = rosenbrock ( x );
+    g = derivative ( rosenbrock , x.' , order = 4 );
+endfunction
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 );
+assert_checkequal([iters evals err], [37 50 9]);
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 , "ar",nap=10 );
+assert_checkequal([iters evals err], [8 10 4]);
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 , "ar",nap=100,iter=10 );
+assert_checkequal([iters evals err], [11 14 5]);
diff --git a/scilab/modules/optimization/tests/nonreg_tests/bug_9208.tst b/scilab/modules/optimization/tests/nonreg_tests/bug_9208.tst
new file mode 100644 (file)
index 0000000..4802b25
--- /dev/null
@@ -0,0 +1,60 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2013 - Scilab Enterprises - Paul Bignier
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+//
+// <-- CLI SHELL MODE -->
+//
+// <-- Non-regression test for bug 9208 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=9208
+//
+// <-- Short Description -->
+// Added three optional output arguments to optim(), to retrieve #iterations,
+// #evaluations and a termination flag.
+//
+
+// Example 1
+a  = 1.0;
+b  = 2.0;
+c  = 3.0;
+d  = 4.0;
+x0 = [1 1];
+function [f, g, ind] = costfunction(x, ind, a, b, c, d)
+    f = a * ( x(1) - c ) ^2 + b * ( x(2) - d )^2;
+    g(1) = 2 * a * ( x(1) - c );
+    g(2) = 2 * b * ( x(2) - d );
+endfunction
+costf = list ( costfunction , a , b , c, d );
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 );
+assert_checkequal([iters evals err], [10 11 1]);
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 , "ar",nap=5 );
+assert_checkequal([iters evals err], [4 5 4]);
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( costf , x0 , "ar",nap=100,iter=5 );
+assert_checkequal([iters evals err], [6 7 5]);
+
+
+// Example 2
+x0 = [-1.2 1.0];
+function f = rosenbrock(x)
+    f = 100.0 *(x(2)-x(1)^2)^2 + (1-x(1))^2;
+endfunction
+function [f, g, ind] = rosenbrockCost2(x, ind)
+    f = rosenbrock ( x );
+    g = derivative ( rosenbrock , x.' , order = 4 );
+endfunction
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 );
+assert_checkequal([iters evals err], [37 50 9]);
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 , "ar",nap=10 );
+assert_checkequal([iters evals err], [8 10 4]);
+
+[fopt, xopt, w, g, iters, evals, err] = optim ( rosenbrockCost2 , x0 , "ar",nap=100,iter=10 );
+assert_checkequal([iters evals err], [11 14 5]);