C#による 非線形計画問題（スパース）

e04ug Example Program Data
 4   6              :Values of N and M
 3   4   2          :Values of NCNLN, NONLN and NJNLN
14   6  'C'  10     :Values of NNZ, IOBJ, START and NNAME
'Varble 1'  'Varble 2'  'Varble 3'  'Varble 4'  'NlnCon 1'
'NlnCon 2'  'NlnCon 3'  'LinCon 1'  'LinCon 2'  'Free Row' :End of NAMES
 1.0E+25   1   1
 1.0E+25   2   1
 1.0E+25   3   1
     1.0   5   1
    -1.0   4   1
 1.0E+25   1   2
 1.0E+25   2   2
 1.0E+25   3   2
    -1.0   5   2
     1.0   4   2
     3.0   6   3
    -1.0   1   3
    -1.0   2   4
     2.0   6   4    :End of matrix A
-0.55     -0.55      0.0       0.0     -894.8   -894.8   -1294.8  -0.55
-0.55     -1.0E+25  :End of BL
 0.55      0.55      1200.0    1200.0  -894.8   -894.8   -1294.8   1.0E+25
 1.0E+25   1.0E+25  :End of BU
 0    0    0    0   :End of ISTATE 
 0.0  0.0  0.0  0.0 :End of XS
 0.0  0.0  0.0      :End of CLAMDA 

e04ug Example Program Results
Computed values of the size of z and iz are 758, 628 respectively.


 Calls to E04UJF
 ---------------

      List
      Print level = 10

 *** E04UGF

 Parameters
 ----------

 Frequencies.
 Check frequency.........        60       Expand frequency.......     10000
 Factorization frequency.        50

 QP subproblems.
 Scale tolerance.........  9.00E-01       Minor feasibility tol..  1.05E-08
 Scale option............         1       Minor optimality tol...  1.05E-08
 Partial price...........         1       Crash tolerance........  1.00E-01
 Pivot tolerance.........  2.04E-11       Minor print level......         0
 Crash option............         0       Elastic weight.........  1.00E+02

 The SQP method.
 Minimize................
 Nonlinear objective vars         4       Major optimality tol...  1.05E-08
 Function precision......  1.72E-13       Unbounded step size....  1.00E+20
 Superbasics limit.......         4       Forward difference int.  4.15E-07
 Unbounded objective.....  1.00E+15       Central difference int.  5.56E-05
 Major step limit........  2.00E+00       Derivative linesearch..
 Derivative level........         3       Major iteration limit..      1000
 Linesearch tolerance....  9.00E-01       Verify level...........         0
 Minor iteration limit...       500       Major print level......        10
 Infinite bound size.....  1.00E+20       Iteration limit........     10000

 Hessian approximation.
 Hessian full memory.....                 Hessian updates........  99999999
 Hessian frequency.......  99999999

 Nonlinear constraints.
 Nonlinear constraints...         3       Major feasibility tol..  1.05E-08
 Nonlinear Jacobian vars.         2       Violation limit........  1.00E+01

 Miscellaneous.
 Variables...............         4       Linear constraints.....         3
 Nonlinear variables.....         4       Linear variables.......         0
 LU factor tolerance.....  5.00E+00       LU singularity tol.....  2.04E-11
 LU update tolerance.....  5.00E+00       LU density tolerance...  6.00E-01
 eps (machine precision).  1.11E-16       Monitoring file........        -1
 COLD start..............                 Infeasible exit........

 Workspace provided is                 IZ(    1000),  Z(    1000).
 To start solving the problem we need  IZ(     628),  Z(     758).

 Itn      0 -- Scale option reduced from                1 to                0.

 Itn      0 -- Feasible linear rows.

 Itn      0 -- Norm(x-x0) minimized. Sum of infeasibilities = 0.00E+00.

 confun  sets       6   out of       6   constraint gradients.
 objfun  sets       4   out of       4   objective  gradients.

 Cheap test on confun...

 The Jacobian seems to be OK.

 The largest discrepancy was    4.41E-08  in constraint     2.

 Cheap test on objfun...

 The objective gradients seem to be OK.
 Gradient projected in two directions   0.00000000000E+00   0.00000000000E+00
 Difference approximations              1.74111992322E-19   4.48742248252E-21

 Itn      0 -- All-slack basis B = I selected.

 Itn      7 -- Large multipliers.
               Elastic mode started with weight =  2.0E+02.


  Maj  Mnr    Step Merit Function Feasibl Optimal Cond Hz PD
    0   12 0.0E+00   3.199952E+05 1.7E+00 8.0E-01 2.1E+06 FF  R   i
    1    2 1.0E+00   2.463016E+05 1.2E+00 3.2E+03 4.5E+00 FF s
    2    1 1.0E+00   1.001802E+04 3.3E-02 9.2E+01 4.5E+00 FF
    3    1 1.0E+00   5.253418E+03 6.6E-04 2.5E+01 4.8E+00 FF
    4    1 1.0E+00   5.239444E+03 2.0E-06 2.8E+01 1.0E+02 FF
    5    1 1.0E+00   5.126208E+03 6.0E-04 5.9E-01 1.1E+02 FF
    6    1 1.0E+00   5.126498E+03 4.7E-07 2.9E-02 1.0E+02 FF
    7    1 1.0E+00   5.126498E+03 5.9E-10 1.5E-03 1.1E+02 TF
    8    1 1.0E+00   5.126498E+03 1.2E-12 7.6E-09 1.1E+02 TT

 Exit from NP problem after     8 major iterations,
                               21 minor iterations.

 Variable State     Value       Lower Bound  Upper Bound   Lagr Mult  Residual

 Varble 1    BS   0.118876     -0.55000      0.55000     -1.2529E-07  0.4311
 Varble 2    BS  -0.396234     -0.55000      0.55000      1.9246E-08  0.1538
 Varble 3    BS    679.945           .        1200.0      1.7001E-10   520.1
 Varble 4   SBS    1026.07           .        1200.0     -2.1918E-10   173.9

 Constrnt State     Value       Lower Bound  Upper Bound   Lagr Mult  Residual

 NlnCon 1    EQ   -894.800      -894.80      -894.80      -4.387      3.3644E-09
 NlnCon 2    EQ   -894.800      -894.80      -894.80      -4.106      6.0061E-10
 NlnCon 3    EQ   -1294.80      -1294.8      -1294.8      -5.463      3.3554E-09
 LinCon 1    BS  -0.515110     -0.55000          None          .      3.4890E-02
 LinCon 2    BS   0.515110     -0.55000          None          .       1.065
 Free Row    BS    4091.97          None         None     -1.000       4092.

 Exit E04UGF - Optimal solution found.

 Final objective value =     5126.498

 ** e04ug returned with ifail =   0

  Variable    Istate          Value                  Lagr Mult

  Varble   1    3              0.118876      -1.253e-07
  Varble   2    3             -0.396234       1.925e-08
  Varble   3    3               679.945         1.7e-10
  Varble   4    2               1026.07      -2.192e-10


  Constrnt    Istate          Value                  Lagr Mult

  NlnCon   1    0                -894.8          -4.387
  NlnCon   2    0                -894.8          -4.106
  NlnCon   3    0               -1294.8          -5.463
  LinCon   1    3              -0.51511               0
  LinCon   2    3               0.51511               0
  Free Row      3               4091.97              -1


  Final objective value =        5126.498

//      e04ug Example Program Text
//      C# version, nAG Copyright 2008
using System;
using NagLibrary;
using System.IO;
namespace NagDotNetExamples
{
  public class E04UGE
  {
    static bool defaultdata = true;
    static string datafile = "";
    // as a command line argument. It defaults to the named file specified below otherwise.
    //
    static void Main(String[] args)
    {
      if (args.Length == 1)
      {
        defaultdata = false;
        datafile = args[0];
      }
      StartExample();
    }
    public static void StartExample()
    {
      E04.E04UG_CONFUN confunE04UG = new E04.E04UG_CONFUN(confun);
      E04.E04UG_OBJFUN objfunE04UG = new E04.E04UG_OBJFUN(objfun);
      try 
      {
        DataReader sr = null;
        if(defaultdata)
        { sr = new DataReader("exampledata/e04uge.d");
      }
      else
      {
        sr = new DataReader(datafile);
      }
      double obj,   sinf; int i,  icol=0,  iobj,  j,  jcol,  m,  miniz,  minz,  n,  ncnln,  ninf,  njnln,  nname,  nnz,  nonln,  ns=0;
      string start="";
      int idummy = -11111;
      string[] names = new string[200];
      int ifail;
      Console.WriteLine("e04ug Example Program Results");
      //      Skip heading in data file.
      sr.Reset();
      sr.Reset();
      n = int.Parse(sr.Next());
      m = int.Parse(sr.Next());
      //
      //         Read ncnln, nonln and njnln from data file.
      //
      sr.Reset();
      ncnln = int.Parse(sr.Next());
      nonln = int.Parse(sr.Next());
      njnln = int.Parse(sr.Next());
      //
      //         Read nnz, iobj, start, nname and names from data file.
      //
      sr.Reset();
      nnz = int.Parse(sr.Next());
      iobj = int.Parse(sr.Next());
      start = sr.Next();
      nname = int.Parse(sr.Next());
      double[] a = new double[nnz];
      double[] bl = new double[n+m];
      double[] bu = new double[n+m];
      double[] clamda = new double[n+m];
      double[] user = new double[1];
      double[] xs = new double[n+m];
      int[] ha = new int[300];
      int[] istate = new int[n+m];
      int[] iuser = new int[1];
      int[] ka = new int[n+1];
      if (nname == n + m)
      {
        sr.Reset();
        for (i = 1 ; i <= n + m ; i++)
        {
          names[i - 1] = sr.Next();
        }
      }
      //
      //         Initialize ka.
      //
      for (i = 1 ; i <= n + 1 ; i++)
      {
        ka[i - 1] = idummy;
      }
      //
      //         Read the matrix a from data file. Set up ka.
      //
      jcol = 1;
      ka[jcol - 1] = 1;
      for (i = 1 ; i <= nnz ; i++)
      {
        //
        //            Element ( ha( i-1 ), icol-1 ) is stored in a[ i-1].
        //
        sr.Reset();
        a[i - 1] = double.Parse(sr.Next());
        ha[i - 1] = int.Parse(sr.Next());
        icol = int.Parse(sr.Next());
        //
        if (icol < jcol)
        {
          //
          //               Elements not ordered by increasing column index.
          //
          Console.WriteLine("\n  {0}{1,5}{2}{3,5}{4}{5}","Element in column",icol," found after element in column",jcol,". Problem"," abandoned.");
          goto L160;
        }
        else if (icol == jcol + 1)
        {
          //
          //               Index in a of the start of the icol-th column equals i.
          //
          ka[icol - 1] = i;
          jcol = (int)icol;
        }
        else if (icol > jcol + 1)
        {
          //
          //               Index in a of the start of the icol-th column equals i,
          //               but columns jcol+1,jcol+2,...,icol-1 are empty. Set the
          //               corresponding elements of ka to i.
          //
          for (j = jcol + 1 ; j <= icol - 1 ; j++)
          {
            ka[j - 1] = i;
          }
          ka[icol - 1] = i;
          jcol = (int)icol;
        }
      }
      //
      ka[n + 1 - 1] = nnz + 1;
      //
      if (n > icol)
      {
        //
        //            Columns n,n-1,...,icol+1 are empty. Set the corresponding
        //            elements of ka accordingly.
        //
        for (i = n ; i <= icol + 1 ; i +=  -1)
        {
          if (ka[i - 1] == idummy)
          {
            ka[i - 1] = ka[i + 1 - 1];
          }
        }
      }
      //
      //         Read bl, bu, istate, xs and clamda from data file.
      //
      sr.Reset();
      for (i = 1 ; i <= n + m ; i++)
      {
        bl[i - 1] = double.Parse(sr.Next());
      }
      sr.Reset();
      for (i = 1 ; i <= n + m ; i++)
      {
        bu[i - 1] = double.Parse(sr.Next());
      }
      if (start == "C")
      {
        sr.Reset();
        for (i = 1 ; i <= n ; i++)
        {
          istate[i - 1] = int.Parse(sr.Next());
        }
      }
      else if (start == "W")
      {
        sr.Reset();
        for (i = 1 ; i <= n + m ; i++)
        {
          istate[i - 1] = int.Parse(sr.Next());
        }
      }
      sr.Reset();
      for (i = 1 ; i <= n ; i++)
      {
        xs[i - 1] = double.Parse(sr.Next());
      }
      if (ncnln > 0)
      {
        sr.Reset();
        for (i = n + 1 ; i <= n + ncnln ; i++)
        {
          clamda[i - 1] = double.Parse(sr.Next());
        }
      }
      // Switch off all output.
      PrintManager.Warning = new PrintManager.MessageLogger(discardmessage);
      PrintManager.Message = new PrintManager.MessageLogger(discardmessage);
      //
      //         Initialise E04.e04ug 
      //
      E04.e04ugOptions options = new E04.e04ugOptions();
      //
      //            Solve the problem.
      //
      start="C";
      // First get appropriate minimum lengths of arrays, z and iz
      int leniz = 500;
      int lenz = 500;
      double[] z = new double[lenz];
      int[] iz = new int[leniz];
      E04.e04ug(confunE04UG, objfunE04UG, n, m, ncnln, nonln, njnln, iobj, nnz, a, ha, ka, bl,
      bu, start, nname, names, ref ns, xs, istate, clamda, out miniz, out minz, out ninf, out sinf,
      out obj, iz, z, options, out ifail);
      Console.WriteLine("Computed values of the size of z and iz are {0}, {1} respectively.", minz, miniz);
      //         Initialise E04.e04ug again and check for error exits
      PrintManager.Warning = new PrintManager.MessageLogger(printmessage);
      PrintManager.Message = new PrintManager.MessageLogger(printmessage);
      options = new E04.e04ugOptions();
      options.Set("List");
      options.Set("Print level = 10");
      //
      //            Solve the problem.
      //
      start="C";
      // The minimum values of lenz and leniz found above are only a preliminary estimate.
      // In practice, larger values need to be chosen frequently, we have therefore specified 
      // each of leniz and lenz as 1000.
      leniz = 1000;
      lenz = 1000;
      z = new double[lenz];
      iz = new int[leniz];
      E04.e04ug(confunE04UG, objfunE04UG, n, m, ncnln, nonln, njnln, iobj, nnz, a, ha, ka, bl,
      bu, start, nname, names, ref ns, xs, istate, clamda, out miniz, out minz, out ninf, out sinf,
      out obj, iz, z, options, out ifail);
      //
      //            Check for error exits
      //
      Console.WriteLine("");
      if (ifail >= (7))
      {
        Console.WriteLine(" ** e04ug returned with ifail = {0, 3}", ifail);
      }
      else if (ifail >= 0)
      {
        Console.WriteLine(" ** e04ug returned with ifail = {0, 3}", ifail);
        Console.WriteLine("");
        Console.WriteLine("  Variable    Istate          Value                  Lagr Mult");
        Console.WriteLine("");
        for (i = 1 ; i <= n ; i++)
        {
          Console.WriteLine("  Varble  {0,2}  {1,3}        {2,14:g6}    {3,12:g4}",i,istate[i - 1],xs[i - 1],clamda[i - 1]);
        }
        Console.WriteLine("");
        Console.WriteLine("");
        Console.WriteLine("  Constrnt    Istate          Value                  Lagr Mult");
        Console.WriteLine("");
        if (ncnln > 0)
        {
          for (i = n + 1 ; i <= n + ncnln ; i++)
          {
            j = i - n;
            Console.WriteLine("  NlnCon  {0,2}  {1,3}        {2,14:g6}    {3,12:g4}",j,istate[i - 1],xs[i - 1],clamda[i - 1]);
          }
        }
        if ((((ncnln == 0)) && ((m == 1))) && ((a[0] == 0.00e0)))
        {
          Console.WriteLine("  DummyRow    {0,3}        {1,14:g6}    {2,12:g4}",istate[n + 1 - 1],xs[n + 1 - 1],clamda[n + 1 - 1]);
        }
        else if (m > ncnln)
        {
          for (i = n + ncnln + 1 ; i <= n + m ; i++)
          {
            j = i - n - ncnln;
            if (i - n == iobj)
            {
              Console.WriteLine("  Free Row    {0,3}        {1,14:g6}    {2,12:g4}",istate[i - 1],xs[i - 1],clamda[i - 1]);
            }
            else
            {
              Console.WriteLine("  LinCon  {0,2}  {1,3}        {2,14:g6}    {3,12:g4}",j,istate[i - 1],xs[i - 1],clamda[i - 1]);
            }
          }
        }
        Console.WriteLine("");
        Console.WriteLine("");
        Console.WriteLine("  Final objective value = {0,15:g7}",obj);
      }
      else
      {
        Console.WriteLine(" ** e04ug returned with ifail = {0, 3}", ifail);
      }
      //
      L160: ;
      //
      }catch (Exception e)
      {
        Console.WriteLine(e.Message);
        Console.WriteLine( "Exception Raised");
      }
    }
    //
    public static void confun (ref int mode, int ncnln, int njnln, int nnzjac,
    double[] x, double[] f, double[] fjac, int nstate)
    {
      //      Computes the nonlinear constraint functions and their Jacobian.
      if ((mode == 0) || (mode == 2))
      {
        f[0] = 1000.00e+0 * Math.Sin( -x[0] - 0.250e+0) + 1000.00e+0 * Math.Sin( -x[1] - 0.250e+0);
        f[1] = 1000.00e+0 * Math.Sin(x[0] - 0.250e+0) + 1000.00e+0 * Math.Sin(x[0] - x[1] - 0.250e+0);
        f[2] = 1000.00e+0 * Math.Sin(x[1] - x[0] - 0.250e+0) + 1000.00e+0 * Math.Sin(x[1] - 0.250e+0);
      }
      //
      if ((mode == 1) || (mode == 2))
      {
        //
        //         Nonlinear Jacobian elements for column 1.
        //
        fjac[0] =  -1000.00e+0 * Math.Cos( -x[0] - 0.250e+0);
        fjac[1] = 1000.00e+0 * Math.Cos(x[0] - 0.250e+0) + 1000.00e+0 * Math.Cos(x[0] - x[1] - 0.250e+0);
        fjac[2] =  -1000.00e+0 * Math.Cos(x[1] - x[0] - 0.250e+0);
        //
        //         Nonlinear Jacobian elements for column 2.
        //
        fjac[3] =  -1000.00e+0 * Math.Cos( -x[1] - 0.250e+0);
        fjac[4] =  -1000.00e+0 * Math.Cos(x[0] - x[1] - 0.250e+0);
        fjac[5] = 1000.00e+0 * Math.Cos(x[1] - x[0] - 0.250e+0) + 1000.00e+0 * Math.Cos(x[1] - 0.250e+0);
      }
      //
    }
    //
    public static void objfun (ref int mode, int nonln, double[] x, ref double objf,
    double[] objgrd, int nstate)
    {
      //      Computes the nonlinear part of the objective function and its
      //      gradient
      if ((mode == 0) || (mode == 2))
      {
        objf = 1.00e-6 * Math.Pow(x[2] ,3) + 2.00e-6 * Math.Pow(x[3] ,3) / 3.00e+0;
      }
      //
      if ((mode == 1) || (mode == 2))
      {
        objgrd[0] = 0.00e+0;
        objgrd[1] = 0.00e+0;
        objgrd[2] = 3.00e-6 * (x[2])*(x[2]);
        objgrd[3] = 2.00e-6 * (x[3])*(x[3]);
      }
      //
    }
    static void discardmessage(String message)
    {
    }
    static void printmessage(String message)
    {
      Console.WriteLine(message);
    }
  }
}