1 #ifndef __TAO_MFQNLS_H 2 #define __TAO_MFQNLS_H 3 #include <petsc-private/taoimpl.h> 4 #include <petsc.h> 5 #include <petscblaslapack.h> 6 7 typedef struct { 8 PetscInt npmax; /* Max number of interpolation points (>n+1) (def: 2n+1) */ 9 PetscInt nmax; /* Max(n*(n+1)/2, 5*npmax) */ 10 PetscInt m,n; 11 Vec *Xhist; 12 Vec *Fhist; 13 PetscReal *Fres; /* (nfmax) */ 14 PetscReal *RES; /* npxm */ 15 PetscReal *work; /* (n) */ 16 PetscReal *work2; /* (n) */ 17 PetscReal *work3; /* (n) */ 18 PetscReal *xmin; /* (n) */ 19 PetscReal *mwork; /* (m) */ 20 PetscReal *Disp; /* nxn */ 21 PetscReal *Fdiff;/* nxm */ 22 PetscReal *H; /* model hessians (mxnxn) */ 23 PetscReal *Hres; /* nxn */ 24 PetscReal *Gres; /* n */ 25 PetscReal *Gdel; /* mxn */ 26 PetscReal *Hdel; /* mxnxn */ 27 PetscReal *Gpoints; /* nxn */ 28 PetscReal *C; /* m */ 29 PetscReal *Xsubproblem; /* n */ 30 PetscInt *indices; /* 1,2,3...m */ 31 PetscInt minindex; 32 PetscInt nmodelpoints; 33 PetscInt *model_indices; /* n */ 34 PetscInt *interp_indices; /* n */ 35 PetscBLASInt *iwork; /* n */ 36 PetscInt nHist; 37 VecScatter scatterf,scatterx; 38 Vec localf, localx, localfmin, localxmin; 39 Vec workxvec; 40 PetscMPIInt size; 41 42 43 PetscReal delta; /* Trust region radius (>0) */ 44 PetscBool usegqt; 45 Mat Hs; 46 Vec b; 47 48 PetscReal deltamax; 49 PetscReal deltamin; 50 PetscReal c1; /* Factor for checking validity */ 51 PetscReal c2; /* Factor for linear poisedness */ 52 PetscReal theta1; /* Pivot threshold for validity */ 53 PetscReal theta2; /* Pivot threshold for additional points */ 54 PetscReal gamma0; /* parameter for shrinking trust region (<1) */ 55 PetscReal gamma1; /* parameter for enlarging trust region (>2) */ 56 PetscReal eta0; /* parameter 1 for accepting point (0 <= eta0 < eta1)*/ 57 PetscReal eta1; /* parameter 2 for accepting point (eta0 < eta1 < 1)*/ 58 PetscReal gqt_rtol; /* parameter used by gqt */ 59 PetscInt gqt_maxits; /* parameter used by gqt */ 60 /* QR factorization data */ 61 PetscInt q_is_I; 62 PetscReal *Q; /* npmax x npmax */ 63 PetscReal *Q_tmp; /* npmax x npmax */ 64 PetscReal *tau; /* scalar factors of H(i) */ 65 PetscReal *tau_tmp; /* scalar factors of H(i) */ 66 PetscReal *npmaxwork; /* work vector of length npmax */ 67 PetscBLASInt *npmaxiwork; /* integer work vector of length npmax */ 68 /* morepoints and getquadnlsmfq */ 69 PetscReal *L; /* n*(n+1)/2 x npmax */ 70 PetscReal *L_tmp; /* n*(n+1)/2 x npmax */ 71 PetscReal *L_save; /* n*(n+1)/2 x npmax */ 72 PetscReal *Z; /* npmax x npmax-(n+1) */ 73 PetscReal *M; /* npmax x n+1 */ 74 PetscReal *N; /* npmax x n*(n+1)/2 */ 75 PetscReal *alpha; /* n+1 */ 76 PetscReal *beta; /* r(n+1)/2 */ 77 PetscReal *omega; /* npmax - np - 1 */ 78 79 Tao subtao; 80 Vec subxl,subxu,subx,subpdel,subndel,subb; 81 Mat subH; 82 83 } TAO_POUNDERS; 84 85 86 PetscErrorCode gqt(PetscInt n, PetscReal *a, PetscInt lda, PetscReal *b, PetscReal delta, PetscReal rtol, PetscReal atol, PetscInt itmax, PetscReal *par, PetscReal *f, PetscReal *x, PetscInt *info, PetscInt *its, PetscReal *z, PetscReal *wa1, PetscReal *wa2); 87 88 PetscErrorCode gqtwrap(Tao tao,PetscReal *gnorm, PetscReal *qmin); 89 PetscErrorCode phi2eval(PetscReal *x, PetscInt n, PetscReal *phi); 90 PetscErrorCode getquadpounders(TAO_POUNDERS *mfqP); 91 PetscErrorCode morepoints(TAO_POUNDERS *mfqP); 92 PetscErrorCode addpoint(Tao tao, TAO_POUNDERS *mfqP, PetscInt index); 93 PetscErrorCode modelimprove(Tao tao, TAO_POUNDERS *mfqP, PetscInt addallpoints); 94 PetscErrorCode affpoints(TAO_POUNDERS *mfqP, PetscReal *xmin, PetscReal c); 95 96 EXTERN_C_BEGIN 97 void dgqt_(PetscInt *n, PetscReal *a, PetscInt *lda, PetscReal *b, PetscReal *delta, PetscReal *rtol, PetscReal *atol, PetscInt *itmax, PetscReal *par, PetscReal *f, PetscReal *x, PetscInt *info, int *its, PetscReal *z, PetscReal *wa1, PetscReal *wa2); 98 EXTERN_C_END 99 #endif /* ifndef __TAO_MFQNLS */ 100