20 #ifndef __mast_gcmma_optimization_interface_h__ 21 #define __mast_gcmma_optimization_interface_h__ 33 #include <libmesh/parallel.h> 37 extern void raasta_(
int *M,
int *N,
38 double *RAA0,
double *RAA,
39 double *XMIN,
double *XMAX,
40 double *DF0DX,
double *DFDX);
41 extern void asympg_(
int *ITER,
int *M,
int *N,
42 double *ALBEFA,
double *GHINIT,
43 double *GHDECR,
double *GHINCR,
44 double *XVAL,
double *XMIN,
double *XMAX,
45 double *XOLD1,
double *XOLD2,
46 double *XLOW,
double *XUPP,
47 double *ALFA,
double *BETA);
48 extern void mmasug_(
int *ITER,
int *M,
int *N,
double *GEPS,
int *IYFREE,
49 double *XVAL,
double *XMMA,
50 double *XMIN,
double *XMAX,
51 double *XLOW,
double *XUPP,
52 double *ALFA,
double *BETA,
53 double *A,
double *B,
double *C,
double *Y,
double *Z,
54 double *RAA0,
double *RAA,
double *ULAM,
55 double *F0VAL,
double *FVAL,
56 double *F0APP,
double *FAPP,
57 double *FMAX,
double *DF0DX,
double *DFDX,
58 double *P,
double *Q,
double *P0,
double *Q0,
59 double *UU,
double *GRADF,
double *DSRCH,
double *HESSF);
60 extern void conser_(
int *M,
int *ICONSE,
61 double *GEPS,
double *F0NEW,
double *F0APP,
62 double *FNEW,
double *FAPP);
63 extern void raaupd_(
int *M,
int *N,
double *GEPS,
64 double *XMMA,
double *XVAL,
65 double *XMIN,
double *XMAX,
66 double *XLOW,
double *XUPP,
67 double *F0NEW,
double *FNEW,
68 double *F0APP,
double *FAPP,
69 double *RAA0,
double *RAA);
70 extern void xupdat_(
int *N,
int *ITER,
71 double *XMMA,
double *XVAL,
72 double *XOLD1,
double *XOLD2);
73 extern void fupdat_(
int *M,
double *F0NEW,
double *FNEW,
74 double *F0VAL,
double *FVAL);
79 namespace Optimization {
84 template <
typename FunctionEvaluationType>
147 #if MAST_ENABLE_GCMMA == 1 155 Assert1(N > 0, N,
"Design variables must be greater than 0");
156 Assert2(
_XVAL.size() == N,
_XVAL.size(), N,
"Design variables must be initialized");
157 Assert2(
_XMIN.size() == N,
_XMIN.size(), N,
"Design variables must be initialized");
158 Assert2(
_XMAX.size() == N,
_XMAX.size(), N,
"Design variables must be initialized");
160 std::vector<real_t> XOLD1(N, 0.), XOLD2(N, 0.),
161 XMMA(N, 0.), XLOW(N, 0.), XUPP(N, 0.),
162 ALFA(N, 0.), BETA(N, 0.), DF0DX(N, 0.),
163 A(M, 0.), B(M, 0.), C(M, 0.), Y(M, 0.), RAA(M, 0.), ULAM(M, 0.),
164 FVAL(M, 0.), FAPP(M, 0.), FNEW(M, 0.), FMAX(M, 0.),
165 DFDX(M*N, 0.), P(M*N, 0.), Q(M*N, 0.), P0(N, 0.), Q0(N, 0.),
166 UU(M, 0.), GRADF(M, 0.), DSRCH(M, 0.), HESSF(M*(M+1)/2, 0.),
167 f0_iters(n_rel_change_iters);
169 std::vector<int> IYFREE(M, 0);
170 std::vector<bool> eval_grads(M,
false);
265 for (
uint_t i=0; i<N; i++)
266 if (max_x < fabs(
_XVAL[i]))
267 max_x = fabs(
_XVAL[i]);
273 bool terminate =
false, inner_terminate=
false;
276 std::fill(C.begin(), C.end(), std::max(1.e0*max_x, constr_penalty));
288 std::fill(eval_grads.begin(), eval_grads.end(),
true);
291 FVAL, eval_grads, DFDX);
300 asympg_(&ITER, &M, &N, &ALBEFA, &GHINIT, &GHDECR, &GHINCR,
302 &XLOW[0], &XUPP[0], &ALFA[0], &BETA[0]);
308 if (write_internal_iteration_data)
312 inner_terminate =
false;
313 while (!inner_terminate) {
318 mmasug_(&ITER, &M, &N, &GEPS, &IYFREE[0], &
_XVAL[0], &XMMA[0],
319 &
_XMIN[0], &
_XMAX[0], &XLOW[0], &XUPP[0], &ALFA[0], &BETA[0],
320 &A[0], &B[0], &C[0], &Y[0], &Z, &RAA0, &RAA[0], &ULAM[0],
321 &F0VAL, &FVAL[0], &F0APP, &FAPP[0], &FMAX[0], &DF0DX[0], &DFDX[0],
322 &P[0], &Q[0], &P0[0], &Q0[0], &UU[0], &GRADF[0], &DSRCH[0], &HESSF[0]);
327 std::fill(eval_grads.begin(), eval_grads.end(),
false);
330 FNEW, eval_grads, DFDX);
334 std::vector<real_t> XMMA_new(XMMA);
336 while (M && FNEW[0] > 1.e2) {
337 std::cout <<
"*** Backtracking: frac = " 339 <<
" constr: " << FNEW[0]
341 for (
uint_t i=0; i<XMMA.size(); i++)
342 XMMA_new[i] = XOLD1[i] + frac*(XMMA[i]-XOLD1[i]);
346 FNEW, eval_grads, DFDX);
349 for (
uint_t i=0; i<XMMA.size(); i++)
350 XMMA[i] = XMMA_new[i];
353 if (INNER >= INNMAX) {
355 <<
"** Max Inner Iter Reached: Terminating! Inner Iter = " 356 << INNER << std::endl;
357 inner_terminate =
true;
363 conser_( &M, &ICONSE, &GEPS, &F0NEW, &F0APP, &FNEW[0], &FAPP[0]);
366 <<
"** Conservative Solution: Terminating! Inner Iter = " 367 << INNER << std::endl;
368 inner_terminate =
true;
378 &F0NEW, &FNEW[0], &F0APP, &FAPP[0], &RAA0, &RAA[0]);
389 xupdat_( &N, &ITER, &XMMA[0], &
_XVAL[0], &XOLD1[0], &XOLD2[0]);
390 fupdat_( &M, &F0NEW, &FNEW[0], &F0VAL, &FVAL[0]);
395 f0_iters[(ITE-1)%n_rel_change_iters] = F0VAL;
401 if (ITE == max_iters) {
403 <<
"GCMMA: Reached maximum iterations, terminating! " 409 bool rel_change_conv =
true;
410 real_t f0_curr = f0_iters[n_rel_change_iters-1];
412 for (
uint_t i=0; i<n_rel_change_iters-1; i++) {
413 if (f0_curr > sqrt(GEPS))
414 rel_change_conv = (rel_change_conv &&
415 fabs(f0_iters[i]-f0_curr)/fabs(f0_curr) < GEPS);
417 rel_change_conv = (rel_change_conv &&
418 fabs(f0_iters[i]-f0_curr) < GEPS);
420 if (rel_change_conv) {
422 <<
"GCMMA: Converged relative change tolerance, terminating! " 427 _comm.broadcast(terminate, 0);
430 #endif //MAST_ENABLE_GCMMA == 1 436 libMesh::Parallel::Communicator &
_comm;
439 const std::vector<real_t>& XVAL,
440 const std::vector<real_t>& XMIN,
441 const std::vector<real_t>& XMAX,
442 const std::vector<real_t>& XLOW,
443 const std::vector<real_t>& XUPP,
444 const std::vector<real_t>& ALFA,
445 const std::vector<real_t>& BETA) {
449 XVAL.size(),
_feval->n_vars(),
450 "Incorrect vector size");
452 XMIN.size(),
_feval->n_vars(),
453 "Incorrect vector size");
455 XMAX.size(),
_feval->n_vars(),
456 "Incorrect vector size");
458 XLOW.size(),
_feval->n_vars(),
459 "Incorrect vector size");
461 XUPP.size(),
_feval->n_vars(),
462 "Incorrect vector size");
464 ALFA.size(),
_feval->n_vars(),
465 "Incorrect vector size");
467 BETA.size(),
_feval->n_vars(),
468 "Incorrect vector size");
471 <<
"****************************************************\n" 472 <<
" GCMMA: ASYMPTOTE DATA \n" 473 <<
"****************************************************\n" 474 << std::setw(5) <<
"Iter: " << i << std::endl
475 << std::setw(5) <<
"DV" 476 << std::setw(20) <<
"XMIN" 477 << std::setw(20) <<
"XLOW" 478 << std::setw(20) <<
"ALFA" 479 << std::setw(20) <<
"X" 480 << std::setw(20) <<
"BETA" 481 << std::setw(20) <<
"XUP" 482 << std::setw(20) <<
"XMAX" << std::endl;
487 << std::setw(20) << XMIN[j]
488 << std::setw(20) << XLOW[j]
489 << std::setw(20) << ALFA[j]
490 << std::setw(20) << XVAL[j]
491 << std::setw(20) << BETA[j]
492 << std::setw(20) << XUPP[j]
493 << std::setw(20) << XMAX[j] << std::endl;
500 std::vector<real_t> &obj_grad,
501 std::vector<real_t> &fvals,
502 std::vector<bool> &eval_grads,
503 std::vector<real_t> &grads) {
506 _comm.broadcast(x, 0);
507 _comm.broadcast(eval_obj_grad, 0);
521 "Objective function has different values on ranks");
522 Assert0(_comm.verify(obj_grad),
523 "Objective function gradient has different values on ranks");
525 "Constraint functions has different values on ranks");
527 "Constraint function gradient has different values on ranks");
538 #endif // __mast_gcmma_optimization_interface_h__
void _evaluate_wrapper(std::vector< real_t > &x, real_t &obj, bool eval_obj_grad, std::vector< real_t > &obj_grad, std::vector< real_t > &fvals, std::vector< bool > &eval_grads, std::vector< real_t > &grads)
void raasta_(int *M, int *N, double *RAA0, double *RAA, double *XMIN, double *XMAX, double *DF0DX, double *DFDX)
bool write_internal_iteration_data
void fupdat_(int *M, double *F0NEW, double *FNEW, double *F0VAL, double *FVAL)
std::vector< real_t > _XMIN
void mmasug_(int *ITER, int *M, int *N, double *GEPS, int *IYFREE, double *XVAL, double *XMMA, double *XMIN, double *XMAX, double *XLOW, double *XUPP, double *ALFA, double *BETA, double *A, double *B, double *C, double *Y, double *Z, double *RAA0, double *RAA, double *ULAM, double *F0VAL, double *FVAL, double *F0APP, double *FAPP, double *FMAX, double *DF0DX, double *DFDX, double *P, double *Q, double *P0, double *Q0, double *UU, double *GRADF, double *DSRCH, double *HESSF)
libMesh::Parallel::Communicator & _comm
virtual ~GCMMAInterface()
#define Assert1(cond, v1, msg)
FunctionEvaluationType * _feval
std::vector< real_t > _XMAX
void conser_(int *M, int *ICONSE, double *GEPS, double *F0NEW, double *F0APP, double *FNEW, double *FAPP)
void xupdat_(int *N, int *ITER, double *XMMA, double *XVAL, double *XOLD1, double *XOLD2)
void init()
initializes the design variable vector from the function object.
real_t asymptote_reduction
std::vector< real_t > _XVAL
void raaupd_(int *M, int *N, double *GEPS, double *XMMA, double *XVAL, double *XMIN, double *XMAX, double *XLOW, double *XUPP, double *F0NEW, double *FNEW, double *F0APP, double *FAPP, double *RAA0, double *RAA)
#define Assert0(cond, msg)
#define Assert2(cond, v1, v2, msg)
real_t asymptote_expansion
void asympg_(int *ITER, int *M, int *N, double *ALBEFA, double *GHINIT, double *GHDECR, double *GHINCR, double *XVAL, double *XMIN, double *XMAX, double *XOLD1, double *XOLD2, double *XLOW, double *XUPP, double *ALFA, double *BETA)
void set_function_evaluation(FunctionEvaluationType &feval)
void _output_iteration_data(uint_t i, const std::vector< real_t > &XVAL, const std::vector< real_t > &XMIN, const std::vector< real_t > &XMAX, const std::vector< real_t > &XLOW, const std::vector< real_t > &XUPP, const std::vector< real_t > &ALFA, const std::vector< real_t > &BETA)
GCMMAInterface(libMesh::Parallel::Communicator &comm)
uint_t n_rel_change_iters