Add npt

c/bobyqa_c.f90

@@ -13,7 +13,7 @@ module bobyqa_c_mod
 contains
 
 
-subroutine bobyqa_c(cobj_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, iprint, info) bind(C)
+subroutine bobyqa_c(cobj_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, info) bind(C)
 use, intrinsic :: iso_c_binding, only : C_DOUBLE, C_INT, C_FUNPTR
 use, non_intrinsic :: cintrf_mod, only : COBJ
 use, non_intrinsic :: consts_mod, only : RP, IK
@@ -36,13 +36,15 @@ subroutine bobyqa_c(cobj_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxf
 real(C_DOUBLE), intent(in), value :: rhoend
 real(C_DOUBLE), intent(in), value :: ftarget
 integer(C_INT), intent(in), value :: maxfun
+integer(C_INT), intent(in), value :: npt
 integer(C_INT), intent(in), value :: iprint
 integer(C_INT), intent(out) :: info
 
 ! Local variables
 integer(IK) :: info_loc
 integer(IK) :: iprint_loc
 integer(IK) :: maxfun_loc
+integer(IK) :: npt_loc
 integer(IK) :: nf_loc
 real(RP) :: f_loc
 real(RP) :: rhobeg_loc
@@ -60,11 +62,12 @@ subroutine bobyqa_c(cobj_ptr, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxf
 rhoend_loc = real(rhoend, kind(rhoend))
 ftarget_loc = real(ftarget, kind(ftarget))
 maxfun_loc = int(maxfun, kind(maxfun_loc))
+npt_loc = int(npt, kind(npt_loc))
 iprint_loc = int(iprint, kind(iprint_loc))
 
 ! Call the Fortran code
 call bobyqa(calfun, x_loc, f_loc, xl=xl_loc, xu=xu_loc, nf=nf_loc, rhobeg=rhobeg_loc, rhoend=rhoend_loc, &
-    & ftarget=ftarget_loc, maxfun=maxfun_loc, iprint=iprint_loc, info=info_loc)
+    & ftarget=ftarget_loc, maxfun=maxfun_loc, npt=npt_loc, iprint=iprint_loc, info=info_loc)
 
 ! Write the outputs
 x = real(x_loc, kind(x))

c/examples/bobyqa/bobyqa_example.c

@@ -23,8 +23,9 @@ int main(int argc, char * argv[])
   const double ftarget = -INFINITY;
   const int iprint = PRIMA_MSG_EXIT;
   const int maxfun = 1000;
+  const int npt = 2*n+1;
   int nf = 0;
-  const int rc = prima_bobyqa(&fun, n, x, &f, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+  const int rc = prima_bobyqa(&fun, n, x, &f, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint);
   const char *msg = prima_get_rc_string(rc);
   printf("x*={%g, %g} rc=%d msg='%s' evals=%d\n", x[0], x[1], rc, msg, nf);
   return (fabs(x[0]-3)>2e-2 || fabs(x[1]-2)>2e-2);

c/examples/lincoa/lincoa_example.c

@@ -35,8 +35,9 @@ int main(int argc, char * argv[])
   const double ftarget = -INFINITY;
   const int iprint = PRIMA_MSG_EXIT;
   const int maxfun = 100000;
+  const int npt = 2*n+1;
   int nf = 0;
-  const int rc = prima_lincoa(&fun, n, x, &f, &cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+  const int rc = prima_lincoa(&fun, n, x, &f, &cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint);
   const char *msg = prima_get_rc_string(rc);
   printf("x*={%g, %g} f*=%g cstrv=%g rc=%d msg='%s' evals=%d\n", x[0], x[1], f, cstrv, rc, msg, maxfun);
   return (fabs(x[0]-3)>2e-2 || fabs(x[1]-2)>2e-2);

c/examples/newuoa/newuoa_example.c

@@ -21,8 +21,9 @@ int main(int argc, char * argv[])
   const double ftarget = -INFINITY;
   const int iprint = PRIMA_MSG_EXIT;
   const int maxfun = 1000;
+  const int npt = 2*n+1;
   int nf = 0;
-  const int rc = prima_newuoa(&fun, n, x, &f, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+  const int rc = prima_newuoa(&fun, n, x, &f, &nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint);
   const char *msg = prima_get_rc_string(rc);
   printf("x*={%g, %g} rc=%d msg='%s' evals=%d\n", x[0], x[1], rc, msg, nf);
   return (fabs(x[0]-3)>2e-2 || fabs(x[1]-2)>2e-2);

c/include/prima/prima.h

@@ -84,6 +84,7 @@ typedef void (*prima_objcon)(const double x[], double *f, double constr[]);
  * ftarget   : target function value; optimization stops when the value <= ftarget
  *             can be set to -INFINITY to disable
  * maxfun    : maximum number of function evaluations
+ * npt       : number of interpolation conditions, recommended: 2*n+1
  * iprint    : verbosity level, see the prima_message enum
  * m_nlcon   : number of non-linear constraints (>=0)
  * calcfc    : function to minimize and constraints (see prima_objcon)
@@ -103,11 +104,11 @@ typedef void (*prima_objcon)(const double x[], double *f, double constr[]);
 PRIMAC_API
 int prima_bobyqa(const prima_obj calfun, const int n, double x[n], double *f,
                  const double xl[n], const double xu[n],
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint);
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint);
 
 PRIMAC_API
 int prima_newuoa(const prima_obj calfun, const int n, double x[n], double *f,
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint);
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint);
 
 PRIMAC_API
 int prima_uobyqa(const prima_obj calfun, const int n, double x[n], double *f,
@@ -127,7 +128,7 @@ int prima_lincoa(const prima_obj calfun, const int n, double x[n], double *f,
                  const int m_ineq, const double Aineq[m_ineq*n], const double bineq[m_ineq],
                  const int m_eq, const double Aeq[m_eq*n], const double beq[m_eq],
                  const double xl[n], const double xu[n],
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint);
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint);
 
 #ifdef __cplusplus
 }

c/lincoa_c.f90

@@ -14,7 +14,7 @@ module lincoa_c_mod
 
 
 subroutine lincoa_c(cobj_ptr, n, x, f, cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, &
-    & nf, rhobeg, rhoend, ftarget, maxfun, iprint, info) bind(C)
+    & nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, info) bind(C)
 use, intrinsic :: iso_c_binding, only : C_DOUBLE, C_INT, C_FUNPTR
 use, non_intrinsic :: cintrf_mod, only : COBJ
 use, non_intrinsic :: consts_mod, only : RP, IK
@@ -44,13 +44,15 @@ subroutine lincoa_c(cobj_ptr, n, x, f, cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, b
 real(C_DOUBLE), intent(in), value :: rhoend
 real(C_DOUBLE), intent(in), value :: ftarget
 integer(C_INT), intent(in), value :: maxfun
+integer(C_INT), intent(in), value :: npt
 integer(C_INT), intent(in), value :: iprint
 integer(C_INT), intent(out) :: info
 
 ! Local variables
 integer(IK) :: info_loc
 integer(IK) :: iprint_loc
 integer(IK) :: maxfun_loc
+integer(IK) :: npt_loc
 integer(IK) :: nf_loc
 real(RP) :: Aineq_loc(m_ineq, n)
 real(RP) :: bineq_loc(m_ineq)
@@ -79,14 +81,15 @@ subroutine lincoa_c(cobj_ptr, n, x, f, cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, b
 rhoend_loc = real(rhoend, kind(rhoend))
 ftarget_loc = real(ftarget, kind(ftarget))
 maxfun_loc = int(maxfun, kind(maxfun_loc))
+npt_loc = int(npt, kind(npt_loc))
 iprint_loc = int(iprint, kind(iprint_loc))
 
 ! Call the Fortran code
 call lincoa(calfun, x_loc, f_loc, cstrv=cstrv_loc, &
     & Aineq=Aineq_loc, bineq=bineq_loc, Aeq=Aeq_loc, beq=beq_loc, &
     & xl=xl_loc, xu=xu_loc, nf=nf_loc, &
     & rhobeg=rhobeg_loc, rhoend=rhoend_loc, &
-    & ftarget=ftarget_loc, maxfun=maxfun_loc, iprint=iprint_loc, info=info_loc)
+    & ftarget=ftarget_loc, maxfun=maxfun_loc, npt=npt_loc, iprint=iprint_loc, info=info_loc)
 
 ! Write the outputs
 x = real(x_loc, kind(x))

c/newuoa_c.f90

@@ -12,7 +12,7 @@ module newuoa_c_mod
 contains
 
 
-subroutine newuoa_c(cobj_ptr, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, iprint, info) bind(C)
+subroutine newuoa_c(cobj_ptr, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, info) bind(C)
 use, intrinsic :: iso_c_binding, only : C_DOUBLE, C_INT, C_FUNPTR
 use, non_intrinsic :: cintrf_mod, only : COBJ
 use, non_intrinsic :: consts_mod, only : RP, IK
@@ -33,13 +33,15 @@ subroutine newuoa_c(cobj_ptr, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, ipri
 real(C_DOUBLE), intent(in), value :: rhoend
 real(C_DOUBLE), intent(in), value :: ftarget
 integer(C_INT), intent(in), value :: maxfun
+integer(C_INT), intent(in), value :: npt
 integer(C_INT), intent(in), value :: iprint
 integer(C_INT), intent(out) :: info
 
 ! Local variables
 integer(IK) :: info_loc
 integer(IK) :: iprint_loc
 integer(IK) :: maxfun_loc
+integer(IK) :: npt_loc
 integer(IK) :: nf_loc
 real(RP) :: f_loc
 real(RP) :: rhobeg_loc
@@ -53,11 +55,12 @@ subroutine newuoa_c(cobj_ptr, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, ipri
 rhoend_loc = real(rhoend, kind(rhoend))
 ftarget_loc = real(ftarget, kind(ftarget))
 maxfun_loc = int(maxfun, kind(maxfun_loc))
+npt_loc = int(npt, kind(npt_loc))
 iprint_loc = int(iprint, kind(iprint_loc))
 
 ! Call the Fortran code
 call newuoa(calfun, x_loc, f_loc, nf=nf_loc, rhobeg=rhobeg_loc, rhoend=rhoend_loc, ftarget=ftarget_loc, &
-    & maxfun=maxfun_loc, iprint=iprint_loc, info=info_loc)
+    & maxfun=maxfun_loc, npt=npt_loc, iprint=iprint_loc, info=info_loc)
 
 ! Write the outputs
 x = real(x_loc, kind(x))

c/prima.c

@@ -9,17 +9,17 @@ int cobyla_c(const int m_nlcon, const prima_objcon calcfc, const int n, double x
              const double xl[n], const double xu[n],
              int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint, int *info);
 int bobyqa_c(prima_obj calfun, const int n, double x[n], double *f, const double xl[n], const double xu[n],
-             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint, int *info);
+             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint, int *info);
 int newuoa_c(prima_obj calfun, const int n, double x[n], double *f,
-             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint, int *info);
+             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint, int *info);
 int uobyqa_c(prima_obj calfun, const int n, double x[n], double *f,
              int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint, int *info);
 int lincoa_c(prima_obj calfun, const int n, double x[n], double *f,
              double *cstrv,
              const int m_ineq, const double Aineq[m_ineq*n], const double bineq[m_ineq],
              const int m_eq, const double Aeq[m_eq*n], const double beq[m_eq],
              const double xl[n], const double xu[n],
-             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint, int *info);
+             int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint, int *info);
 
 /* these functions just call the fortran compatibility layer and return the status code */
 int prima_cobyla(const int m_nlcon, const prima_objcon calcfc, const int n, double x[n], double *f, double *cstrv, double nlconstr[m_nlcon],
@@ -34,18 +34,18 @@ int prima_cobyla(const int m_nlcon, const prima_objcon calcfc, const int n, doub
 }
 
 int prima_bobyqa(const prima_obj calfun, const int n, double x[n], double *f, const double xl[n], const double xu[n],
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint)
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint)
 {
   int info = 0;
-  bobyqa_c(calfun, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, iprint, &info);
+  bobyqa_c(calfun, n, x, f, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, &info);
   return info;
 }
 
 int prima_newuoa(const prima_obj calfun, const int n, double x[n], double *f,
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint)
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint)
 {
   int info = 0;
-  newuoa_c(calfun, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, iprint, &info);
+  newuoa_c(calfun, n, x, f, nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, &info);
   return info;
 }
 
@@ -61,10 +61,10 @@ int prima_lincoa(const prima_obj calfun, const int n, double x[n], double *f, do
                  const int m_ineq, const double Aineq[m_ineq*n], const double bineq[m_ineq],
                  const int m_eq, const double Aeq[m_eq*n], const double beq[m_eq],
                  const double xl[n], const double xu[n],
-                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, const int iprint)
+                 int *nf, const double rhobeg, const double rhoend, const double ftarget, const int maxfun, int npt, const int iprint)
 {
   int info = 0;
-  lincoa_c(calfun, n, x, f, cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, iprint, &info);
+  lincoa_c(calfun, n, x, f, cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, nf, rhobeg, rhoend, ftarget, maxfun, npt, iprint, &info);
   return info;
 }
 

c/tests/t_large.c

@@ -106,7 +106,7 @@ int main(int argc, char * argv[])
   if(strcmp(algo, "bobyqa") == 0)
   {
     n = 1600;
-    rc = prima_bobyqa(&fun, n, x, &f, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+    rc = prima_bobyqa(&fun, n, x, &f, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, 2*n+1, iprint);
   }
   else if(strcmp(algo, "cobyla") == 0)
   {
@@ -118,12 +118,12 @@ int main(int argc, char * argv[])
   {
     n = 1000;
     m_ineq = 1000;
-    rc = prima_lincoa(&fun, n, x, &f, &cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+    rc = prima_lincoa(&fun, n, x, &f, &cstrv, m_ineq, Aineq, bineq, m_eq, Aeq, beq, xl, xu, &nf, rhobeg, rhoend, ftarget, maxfun, 2*n+1, iprint);
   }
   else if(strcmp(algo, "newuoa") == 0)
   {
     n = 1600;
-    rc = prima_newuoa(&fun, n, x, &f, &nf, rhobeg, rhoend, ftarget, maxfun, iprint);
+    rc = prima_newuoa(&fun, n, x, &f, &nf, rhobeg, rhoend, ftarget, maxfun, 2*n+1, iprint);
   }
   else if(strcmp(algo, "uobyqa") == 0)
   {