« Previous « Start » Next »
8 NLLS Problem
The
constrained nonlinear least squares (
cls) problem
is defined as
|
|
|
|
| |
|
| s/t |
| xL |
≤ |
x |
≤ |
xU, |
| bL |
≤ |
A x |
≤ |
bU |
| cL |
≤ |
c(x) |
≤ |
cU |
|
|
(10) |
where
x,
xL,
xU 
Rn,
r(
x)
RM,
A
Rm1 × n,
bL,
bU Rm1 and
cL,
c(
x),
cU Rm2. The following file defines and
solves a problem in TOMNET.
The following code defines a problem in TOMNET.
using System;
using TOMNET;
namespace TOMNET
{
/// <summary>
/// Quick guide class for nonlinear least squares problem.
/// </summary>
static public class nllsQG
{
/// <summary>
/// Class for the residual function for nllsQG.
/// </summary>
private class nllsQG_r : IDFunction
{
/// <summary>
/// The observations.
/// </summary>
private double[] y;
/// <summary>
/// The observation times.
/// </summary>
private double[] t;
/// <summary>
/// Constant set by user.
/// </summary>
private double uP;
/// <summary>
/// Indicates if the problem has y and t set.
/// </summary>
bool hasYT = false;
/// <summary>
/// Constructor
/// </summary>
/// <param name="uP">
/// Parameter passed to function from problem.
/// </param>
public nllsQG_r(double uP)
{
this.uP = uP;
}
/// <summary>
/// Copies information from a problem to the function.
/// </summary>
/// <param name="Prob">Problem to copy information from.</param>
public void initFunc(NonlinearLeastSquaresProblem Prob)
{
this.y = Prob.y;
this.t = Prob.t;
hasYT = true;
}
/// <summary>
/// Evaluates the residual.
/// </summary>
/// <param name="r">Values of the residual.</param>
/// <param name="x">The decision variables.</param>
public void Evaluate(double[] r, double[] x)
{
if (!hasYT)
{
throw new Exception("Function not initiated yet, " +
"please run initFunc first.");
}
for (int i = 0; i < y.Length; i++)
{
r[i] = this.uP * x[0] * (Math.Exp(-x[1] * t[i]) -
Math.Exp(-x[0] * t[i])) / (x[2] * (x[0] - x[1]));
}
}
/// <summary>
/// Evaluates the gradient of the residual Jacobain.
/// </summary>
/// <param name="gx">Values of the gradient.</param>
/// <param name="x">The decision variables.</param>
public void Grad(double[] gx, double[] x)
{
if (gx.Length != y.Length * x.Length)
{
throw new TOMNETErrors.MismatchException("gx has illegal " +
"length " + gx.Length.ToString() + " (should have " +
(y.Length * x.Length).ToString() + ")");
}
double a = this.uP * x[0] / (x[2] * (x[0] - x[1]));
double b = x[0] - x[1];
double e1, e2;
for (int i = 0; i < y.Length; i++)
{
e1 = Math.Exp(-x[0] * t[i]);
e2 = Math.Exp(-x[1] * t[i]);
gx[i*3+0] = a*(t[i]*e1+(e2-e1)*(1-1/b));
gx[i*3+1] = a*(-t[i]*e2+(e2-e1)/b);
gx[i*3+2] =-a * (e2 - e1) / x[2];
}
}
}
/// <summary>
/// Class for blank Jacobian matrix for nnlsQG.
/// Needed for the class NonlinearLeastSquaresProblem.
/// </summary>
public class nllsQG_J : IDConstraints
{
/// <summary>
/// Evaluates the constraints by
/// returning an array of zeros.
/// </summary>
/// <param name="x">The decision variables.</param>
/// <returns>An array of zeros.</returns>
public double[] Evaluate(double[] x)
{
return new double[28];
}
/// <summary>
/// Evaluates the constraints.
/// </summary>
/// <param name="cx">
/// Already allocated array for the values of the constraints
/// </param>
/// <param name="x">The decision variables.</param>
public void Evaluate(double[] cx, double[] x)
{
cx.Initialize();
}
/// <summary>
/// Computes the Jacobian of the nonlinear constraints.
/// </summary>
/// <param name="Jac">Values of the Jacobian.</param>
/// <param name="x">The decision variables.</param>
public void dc(Jacobian Jac, double[] x)
{
;
}
}
/// <summary>
/// Testprogram for nllsQG.
/// </summary>
static void Main()
{
//
// The observations and the observation times.
//
double[] y = new double[28] { 30.5, 44, 43, 41.5, 38.6, 38.6,
39, 41, 37, 37, 24, 32, 29, 23, 21, 19, 17, 14, 9.5, 8.5, 7,
6, 6, 4.5, 3.6, 3, 2.2, 1.6 };
double[] t = new double[28] { 0.25, 0.5, 0.75, 1, 1.5, 2, 3, 4,
6, 8, 12, 24, 32, 48, 54, 72, 80, 96, 121, 144, 168, 192, 216,
246, 276, 324, 348, 386 };
//
// Starting point, lower and upper bounds for decision variables
//
double[] x_0 = new double[3] {6.8729, 0.0108, 0.1248};
double[] x_L = ArrayUtils.getx_L(3);
double[] x_U = ArrayUtils.getx_U(3);
//
// Nonlinear constraint bounds
//
double[] c_L = new double[0];
double[] c_U = new double[0];
//
// Linear constraints, and bounds.
//
double[] A = new double[3] { 1, 0, 0 };
double[] b_L = new double[1] { double.NegativeInfinity };
double[] b_U = new double[1] { double.PositiveInfinity };
//
// The residuals and Jacobian
//
nllsQG_r R = new nllsQG_r(5);
nllsQG_J C = new nllsQG_J();
//
// Create an instance of the problem
//
NonlinearLeastSquaresProblem Prob = new
NonlinearLeastSquaresProblem(R, C, A, b_L, b_U, x_L, x_U, x_0,
c_L, c_U, y, t, "NllsProb", double.NegativeInfinity, null,
null);
R.initFunc(Prob);
//
// Create a solver and result instance
//
SNOPT solver = new SNOPT();
Result result;
//
// Define a print file for SNOPT.
//
string printfilename = "SnoptNllsQP.txt";
solver.Options.PrintFile = printfilename;
//
// Solve the problem
//
solver.Solve(Prob, out result);
//
// Get the solution results
//
double[] solution = result.x_k;
double objective = result.f_k;
//
// Display the solution on the console.
//
Console.WriteLine(" * * * * * * * *");
Console.WriteLine(" * Solved " + Prob.Name);
Console.WriteLine(" * Printfile: " + printfilename);
Console.WriteLine(" * Objective: {0}", objective);
Console.WriteLine(" * Solution");
Console.WriteLine(" x_L[i] <= x_k[i] <= x_U[i]");
for (int i = 0; i < solution.Length; i++)
{
Console.WriteLine(" {0,10:g3} <= {1,10:g3} <= {2,10:g3} ",
Prob.x_L[i], solution[i], Prob.x_U[i]);
}
}
}
}
There are two callback function defined in the code:
r: Residual function
J: Residual gradient matrix (Jacobian)
It is also possible to define nonlinear constraints as needed.
The results should match the following command prompt outputs:
« Previous « Start » Next »
TOMNET
