/****************************************************
*  Examples #1 to #3 of nonlinear systems for NNES  *
*                                                   *
* Release1.1: Added example of size 4 (03/30/07).   *
****************************************************/
#include <basis.h>

int nfetot;

REAL Sqr(REAL);

void fcn1(bool *overfl, int n, REAL *fvec, REAL *xc) {

  *overfl = FALSE;
  nfetot++;

/* Nonlinear systems to solve 
   Example #1:

   X1**2 + X1 + X2**2 - 2           = 0
   X1**2 + X2 - X2**2 - 1 + log(X1) = 0

   Example #2  (for comparison with Brown's method):
   10.0*x + x2 + x3  + x4 - 20.0 + Sqr(sin(x1)) + Sqr(cos(x2)) = 0
   x1 + 20.0*x2 + x3 + x4 - 48.0 + one/pow^6                   = 0
   Sqr(x1 + x2) + 30.0*x3 + x4 - 97.0 + log(x1) + log(x2+x3)   = 0
   x1     + x2  + x3 + 40.0*x4 - 166.0 + Sqr(x1)               = 0

   Example #3:  
     Hiebert's 2nd Chemical Engineering Problem

     source: Hiebert; Sandia Technical Report #SAND80-0181
             Sandia National Laboratories, Albuquerque, NM (1980)

   X1 + X2 + X4 - .001 = 0
   X5 + X6 -55         = 0
   X1 + X2 + X3 + 2X5 + X6 - 110.001 = 0
   X1 - 0.1X2        = 0
   X1 - 10000 X3 X4  = 0
   X5 - 5.5e15 X3 X6 = 0

   solution: (8.264e-5, 8.264e-4, 9.091e-5, 9.091e-5, 55, 1.1e-10) */

if (n==2) {
  fvec[1] = xc[1]*xc[1] + xc[1] + xc[2]*xc[2] -2.0;
  fvec[2] = xc[1]*xc[1] + xc[2] - xc[2]*xc[2] - 1.0 + log(xc[1]);
}  
else if (n==4) {
  fvec[1] = 10.0*xc[1] + xc[2] + xc[3] + xc[4] - 20.0 + Sqr(sin(xc[1])) + Sqr(cos(xc[2]));
  fvec[2] = xc[1] + 20.0*xc[2] + xc[3] + xc[4] - 48.0 + 1.0/pow(xc[1],6);
  fvec[3] = Sqr(xc[1] + xc[2]) + 30.0*xc[3] + xc[4] - 97.0 + log(xc[1]) + log(xc[2]+xc[3]);
  fvec[4] = xc[1] + xc[2] + xc[3] + 40.0*xc[4] - 166.0 + Sqr(xc[1]);
}
else {
  fvec[1] = xc[1] + xc[2] + xc[4] - 0.001;
  fvec[2] = xc[5] + xc[6] - 55.0;
  fvec[3] = xc[1] + xc[2] + xc[3] + TWO * xc[5] + xc[6] - 110.001;
  fvec[4] = xc[1] - 0.1 * xc[2];
  fvec[5] = xc[1] - 1.0e04 * xc[3] * xc[4];
  fvec[6] = xc[5] - 5.5e15 * xc[3] * xc[6];
}

} // fcn1()


//Dummt function
void jacob(bool *overfl, int n, REAL **jac, REAL *xc) {
  int i, j;
  *overfl = FALSE;
  for (i=0; i<=n; i++)
    for (j=0; j<=n; j++)
	  jac[i][j] = ZERO;
}

//end of file fcn.cpp