1 /* _______________________________________________________________________
2
3 DAKOTA: Design Analysis Kit for Optimization and Terascale Applications
4 Copyright 2014-2020 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
5 This software is distributed under the GNU Lesser General Public License.
6 For more information, see the README file in the top Dakota directory.
7 _______________________________________________________________________ */
8
9 #include <cstdlib>
10 #include <iostream>
11 #include <fstream>
12 #include <vector>
13 #include <string>
14 #include <cmath>
15 //#include <thread> // for sleep_for
16
17
18 //KRD modified this starting from text_book.cpp
main(int argc,char ** argv)19 int main(int argc, char** argv)
20 {
21
22 std::ifstream fin(argv[1]);
23 if (!fin) {
24 std::cerr << "\nError: failure opening " << argv[1] << std::endl;
25 exit(-1);
26 }
27 size_t i, j, k, num_vars, num_fns, num_deriv_vars;
28 std::string vars_text, fns_text, dvv_text;
29
30 // Get the parameter std::vector and ignore the labels
31 fin >> num_vars >> vars_text;
32 std::vector<double> x(num_vars);
33 for (i=0; i<num_vars; i++) {
34 fin >> x[i];
35 fin.ignore(256, '\n');
36 }
37
38 // Get the ASV std::vector and ignore the labels
39 fin >> num_fns >> fns_text;
40 if(num_fns!=1) {
41 std::cerr << "\nError: Shubert doesn't support constraints but you said #constraints=" << num_fns-1 << std::endl;
42 exit(-1);
43 }
44 std::vector<int> ASV(num_fns);
45 fin >> ASV[0];
46 fin.ignore(256, '\n');
47
48 // Get the DVV std::vector and ignore the labels
49 fin >> num_deriv_vars >> dvv_text;
50 std::vector<int> DVV(num_deriv_vars);
51 for (i=0; i<num_deriv_vars; i++) {
52 fin >> DVV[i];
53 fin.ignore(256, '\n');
54 }
55
56 //srand ( (unsigned int) (time(NULL)/x[0]) );
57 //std::this_thread::sleep_for
58 // (std::chrono::seconds((int)(3.0*((double)rand()/RAND_MAX))));
59
60 //std::this_thread::sleep_for(std::chrono::seconds(5));
61
62 // Compute the results and output them directly to argv[2] (the NO_FILTER
63 // option is used). Response tags are now optional; output them for ease
64 // of results readability.
65 std::ofstream fout(argv[2]);
66 if (!fout) {
67 std::cerr << "\nError: failure creating " << argv[2] << std::endl;
68 exit(-1);
69 }
70 fout.precision(15); // 16 total digits
71 fout.setf(std::ios::scientific);
72 fout.setf(std::ios::right);
73
74 // **** f:
75 // f(\underline{x})=\prod_{i=1}^num_vars w(x_i) //latex notation
76 // w(x_i) = \sum_{j=1}^5 j*cos(x_i*(j+1)+j) //w for shubert
77 std::vector<double> w(num_vars);
78 double jdbl;
79 if (ASV[0] >=1)
80 for (i=0; i<num_vars; i++) {
81 w[i]=0.0;
82 for (j=1; j<=5; ++j) {
83 jdbl=static_cast<double>(j);
84 w[i]+=jdbl*std::cos(x[i]*(jdbl+1.0)+jdbl);
85 }
86 }
87
88 if (ASV[0] & 1) {
89 double val = 1.0;
90 for (i=0; i<num_vars; i++)
91 val*=w[i];
92 fout << " " << val << " f\n";
93 }
94
95 // **** df/dx:
96 // df/dx_i=(\prod_{j=1}^{i-1} w(x_j)) * dw(x_i)/dx_i * (\prod_{j=i+1}^{num_vars} w(x_j))
97 // dw/dx_i= \sum_{j=1}^5 j*(j+1)*-sin(x_i*(j+1)+j) //for shubert
98 std::vector<double> d1w(num_vars);
99 if (ASV[0] >=2)
100 for (i=0; i<num_deriv_vars; i++) {
101 int ii=DVV[i]-1;
102 d1w[ii]=0.0;
103 for (j=1; j<=5; ++j) {
104 jdbl=static_cast<double>(j);
105 d1w[ii]+=jdbl*(jdbl+1.0)*-std::sin(x[ii]*(jdbl+1.0)+jdbl);
106 }
107 }
108 double d1;
109 if (ASV[0] & 2) {
110 fout << "[ ";
111 for (i=0; i<num_deriv_vars; i++) {
112 int ii=DVV[i]-1;
113 d1=d1w[ii];
114 for (j=0; j<ii; ++j)
115 d1*=w[j];
116 for (j=ii+1; j<num_vars; ++j)
117 d1*=w[j];
118 fout << d1 << ' ';
119 }
120 fout << "]\n";
121 }
122
123 // **** d^2f/dx^2:
124 // if i<k d^2f/(dx_i*dx_k)=(\prod_{j=1)^{i-1} w(x_j)) * dw(x_i)/dx_i ...
125 // *(\prod_{j=i+1)^{k-1} w(x_j)) * dw(x_k)/dx_k ...
126 // *(\prod_{j=k+1)^{num_vars} w(x_j))
127 // if i==k d^2f/dx_i^2 =(\prod_{j=1)^{i-1} w(x_j)) * d^2w(x_i)/dx_i^2 * (\prod_{j=i+1)^{num_vars} w(x_j))
128 // d^2w(x_i)/dx_i^2 = \sum{j=1}^5 j*(j+1)*(j+1)*-cos(x_i*(j+1)+j) //for shubert
129 if (ASV[0] & 4) {
130 fout << "[[ ";
131 double d2;
132 for (i=0; i<num_deriv_vars; ++i) {
133 for (k=0; k<num_deriv_vars; ++k) {
134 int ii=DVV[i]-1;
135 int kk=DVV[k]-1;
136 if (ii==kk) {
137 d2=0.0;
138 for (j=1; j<=5; ++j) {
139 jdbl=static_cast<double>(j);
140 d2+=jdbl*(jdbl+1.0)*(jdbl+1.0)*-std::cos(x[ii]*(jdbl+1.0)+jdbl);
141 }
142 for (j=0; j<ii; ++j)
143 d2*=w[j];
144 for (j=ii+1; j<num_vars; ++j)
145 d2*=w[j];
146 }
147 else {
148 d2=d1w[ii]*d1w[kk];
149 if(kk<ii) {
150 j=ii;
151 ii=kk;
152 kk=j;
153 }
154 for (j=0; j<ii; ++j)
155 d2*=w[j];
156 for (j=ii+1; j<kk; ++j)
157 d2*=w[j];
158 for (j=kk+1; j<num_vars; ++j)
159 d2*=w[j];
160 }
161 // for (j=0; j<num_vars; ++j)
162 // if((j!=ii)&&(j!=kk))
163 // d2*=w[j];
164 fout << d2 << ' ';
165 }
166 }
167 fout << "]]\n";
168 }
169
170 fout.flush();
171 fout.close();
172 return 0;
173 }
174