// e04he Example Program Text // C# version, NAG Copyright 2008 using System; using NagLibrary; using System.Globalization; using System.IO; namespace NagDotNetExamples { public class E04HEE { static int m = 15; static int nt = 3; static double[,] t = new double[m, nt]; static double[] y = new double[m]; static string datafile = "ExampleData/e04hee.d"; // as a command line argument. It defaults to the named file specified below otherwise. // static void Main(String[] args) { if (args.Length == 1) { datafile = args[0]; } StartExample(); } public static void StartExample() { E04.E04YA_LSQFUN lsqfunE04YA = new E04.E04YA_LSQFUN(lsqfun); E04.E04YB_LSQFUN lsqfunE04YB = new E04.E04YB_LSQFUN(lsqfun); E04.E04YB_LSQHES lsqhesE04YB = new E04.E04YB_LSQHES(lsqhes); E04.E04HE_LSQFUN lsqfunE04HE = new E04.E04HE_LSQFUN(lsqfun); E04.E04HE_LSQHES lsqhesE04HE = new E04.E04HE_LSQHES(lsqhes); E04.E04HE_LSQMON lsqmonE04HE = new E04.E04HE_LSQMON(lsqmon); try { DataReader sr = new DataReader(datafile); double eta, fsumsq, stepmx, xtol; int i, ifail, iprint, j, maxcal, nf, niter; int n = 3; int ldfjac = m; int ldv = n; int lb = n * (n + 1) / 2; double[] b = new double[lb]; double[,] fjac = new double[ldfjac, n]; double[] fvec = new double[m]; double[] g = new double[n]; double[] s = new double[n]; double[,] v = new double[ldv, n]; double[] x = new double[n]; Console.WriteLine("e04he Example Program Results"); // Skip heading in data file sr.Reset(); // Observations of tj (j = 1, 2, 3) are held in t[i-1, j-1] // (i = 1, 2, . . . , 15) for (i = 1; i <= m; i++) { sr.Reset(); y[i - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture); for (j = 1; j <= nt; j++) { t[i - 1, j - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture); } } // Set up an arbitrary point at which to check the derivatives x[0] = 0.190e0; x[1] = -1.340e0; x[2] = 0.880e0; // Check the 1st derivatives // E04.e04ya(m, n, lsqfunE04YA, x, fvec, fjac, out ifail); // if (ifail == 0) { // // Check the evaluation of B E04.e04yb(m, n, lsqfunE04YB, lsqhesE04YB, x, fvec, fjac, b, out ifail); // } if (ifail != 0) { Console.WriteLine(""); Console.WriteLine(" ** e04ya returned with ifail = {0, 3}", ifail); goto L40; } // // Continue setting parameters for e04he // Set iprint to 1 to obtain output from lsqmonE04HE at each iteration * iprint = -1; maxcal = 50 * n; eta = 0.90e0; xtol = 10.00e0 * Math.Sqrt(X02.x02aj()); // We estimate that the minimum will be within 10 units of the // starting point stepmx = 10.00e0; // Set up the starting point x[0] = 0.50e0; x[1] = 1.00e0; x[2] = 1.50e0; // E04.e04he(m, n, lsqfunE04HE, lsqhesE04HE, lsqmonE04HE, iprint, maxcal, eta, xtol, stepmx, x, out fsumsq, fvec, fjac, s, v, out niter, out nf, out ifail); // if (ifail >= 0) { if (ifail != 0) { Console.WriteLine(""); Console.WriteLine(" {0}{1,3}{2}", "Error exit type", ifail, " - see method document"); } if (ifail != 1) { Console.WriteLine(""); Console.WriteLine(" {0}{1,12:f4}", "On exit, the sum of squares is", fsumsq); Console.Write(" " + " {0}", "at the point"); for (j = 1; j <= n; j++) { Console.Write(" " + " {0, 12:f4}", x[j - 1]); } Console.WriteLine(" "); lsqgrd(m, n, fvec, fjac, g); Console.Write(" " + " {0}", "The corresponding gradient is"); for (j = 1; j <= n; j++) { Console.Write(" " + " {0, 12:e3}", g[j - 1]); } Console.WriteLine(" "); Console.WriteLine(" {0}", " (machine dependent)"); Console.WriteLine(" {0}\n", "and the residuals are"); for (i = 1; i <= m; i++) { Console.WriteLine(" " + " {0, 9:e1}", fvec[i - 1]); } Console.WriteLine(" "); } } else { Console.WriteLine(""); Console.WriteLine(" ** e04he returned with ifail = {0, 3}", ifail); } L40: ; // } catch (Exception e) { Console.WriteLine(e.Message); Console.WriteLine( "Exception Raised"); } } // public static void lsqfun(ref int iflag, int m, int n, double[] xc, double[] fvec, double[,] fjac) { // Routine to evaluate the residuals and their 1st derivatives double denom, dummy; int i = 0; for (i = 1; i <= m; i++) { denom = xc[1] * t[i - 1, 1] + xc[2] * t[i - 1, 2]; fvec[i - 1] = xc[0] + t[i - 1, 0] / denom - y[i - 1]; fjac[i - 1, 0] = 1.00e0; dummy = -1.00e0 / (denom * denom); fjac[i - 1, 1] = t[i - 1, 0] * t[i - 1, 1] * dummy; fjac[i - 1, 2] = t[i - 1, 0] * t[i - 1, 2] * dummy; } } // public static void lsqhes(ref int iflag, int m, int n, double[] fvec, double[] xc, double[] b) { // Routine to compute the lower triangle of the matrix B // (stored by rows in the array B) double dummy, sum22, sum32, sum33; int i = 0; b[0] = 0.00e0; b[1] = 0.00e0; sum22 = 0.00e0; sum32 = 0.00e0; sum33 = 0.00e0; for (i = 1; i <= m; i++) { dummy = 2.00e0 * t[i - 1, 0] / Math.Pow((xc[1] * t[i - 1, 1] + xc[2] * t[i - 1, 2]), 3); sum22 = sum22 + fvec[i - 1] * dummy * (t[i - 1, 1]) * (t[i - 1, 1]); sum32 = sum32 + fvec[i - 1] * dummy * t[i - 1, 1] * t[i - 1, 2]; sum33 = sum33 + fvec[i - 1] * dummy * (t[i - 1, 2]) * (t[i - 1, 2]); } b[2] = sum22; b[3] = 0.00e0; b[4] = sum32; b[5] = sum33; } // public static void lsqmon(int m, int n, double[] xc, double[] fvec, double[,] fjac, double[] s, int igrade, int niter, int nf) { // Monitoring method double fsumsq, gtg; int j = 0; int ifail; double[] g = new double[3]; fsumsq = F06.f06ea(m, fvec, 1, fvec, 1, out ifail); lsqgrd(m, n, fvec, fjac, g); gtg = F06.f06ea(n, g, 1, g, 1, out ifail); Console.WriteLine(""); Console.WriteLine(" {0}", " Itns F evals SUMSQ GTG grade"); Console.WriteLine(" {0,4} {1,5} {2,13:e5} {3,9:e1} {4,3}", niter, nf, fsumsq, gtg, igrade); Console.WriteLine(""); Console.WriteLine(" {0}", " x g Singular values"); for (j = 1; j <= n; j++) { Console.WriteLine(" {0,13:e5} {1,9:e1} {2,9:e1}", xc[j - 1], g[j - 1], s[j - 1]); } // } // public static void lsqgrd(int m, int n, double[] fvec, double[,] fjac, double[] g) { // Routine to evaluate gradient of the sum of squares double sum = 0.0; int i, j; for (j = 1; j <= n; j++) { sum = 0.00e0; for (i = 1; i <= m; i++) { sum = sum + fjac[i - 1, j - 1] * fvec[i - 1]; } g[j - 1] = sum + sum; } } } }