Sequential Quantum Gate Decomposer  v1.9.6
Powerful decomposition of general unitarias into one- and two-qubit gates gates
optimization_engines/Adam.cpp
Go to the documentation of this file.
1 /*
2 Created on Fri Jun 26 14:13:26 2020
3 Copyright 2020 Peter Rakyta, Ph.D.
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 
17 @author: Peter Rakyta, Ph.D.
18 */
24 #include "Optimization_Interface.h"
26 #include "Adam.h"
27 
28 
29 #include <fstream>
30 
31 
32 #ifdef __DFE__
33 #include "common_DFE.h"
34 #endif
35 
36 
37 
38 
45 
46 #ifdef __DFE__
47  if ( qbit_num >= 2 && get_accelerator_num() > 0 ) {
48  upload_Umtx_to_DFE();
49  }
50 #endif
51 
52 
53 
54  if (gates.size() == 0 ) {
55  return;
56  }
57 
58 
59 
60  if (solution_guess.size() == 0 ) {
61  solution_guess = Matrix_real(num_of_parameters,1);
62  }
63 
64 
65  if (optimized_parameters_mtx.size() == 0) {
66  optimized_parameters_mtx = Matrix_real(1, num_of_parameters);
67  memcpy(optimized_parameters_mtx.get_data(), solution_guess.get_data(), num_of_parameters*sizeof(double) );
68  }
69 
70  int random_shift_count = 0;
71  long long sub_iter_idx = 0;
72  double current_minimum_hold = current_minimum;
73 
74 
75  tbb::tick_count adam_start = tbb::tick_count::now();
76  CPU_time = 0.0;
77 
78  Adam optimizer;
79  optimizer.initialize_moment_and_variance( num_of_parameters );
80 
81 
82 
83  // the array storing the optimized parameters
84  Matrix_real solution_guess_tmp = Matrix_real( num_of_parameters, 1 );
85  memcpy(solution_guess_tmp.get_data(), solution_guess.get_data(), num_of_parameters*sizeof(double) );
86 
87  Matrix_real grad_mtx = Matrix_real( num_of_parameters, 1 );
88 
89 
90 
91 
92 
93 
94  int ADAM_status = 0;
95 
96 
97  long long max_inner_iterations_loc;
98  if ( config.count("max_inner_iterations_adam") > 0 ) {
99  config["max_inner_iterations_adam"].get_property( max_inner_iterations_loc );
100  }
101  else if ( config.count("max_inner_iterations") > 0 ) {
102  config["max_inner_iterations"].get_property( max_inner_iterations_loc );
103  }
104  else {
105  max_inner_iterations_loc =max_inner_iterations;
106  }
107 
108  long long iteration_threshold_of_randomization_loc;
109  if ( config.count("randomization_threshold_adam") > 0 ) {
110  config["randomization_threshold_adam"].get_property( iteration_threshold_of_randomization_loc );
111  }
112  else if ( config.count("randomization_threshold") > 0 ) {
113  config["randomization_threshold"].get_property( iteration_threshold_of_randomization_loc );
114  }
115  else {
116  iteration_threshold_of_randomization_loc = 2500000;
117  }
118 
119  long long export_circuit_2_binary_loc;
120  if ( config.count("export_circuit_2_binary_adam") > 0 ) {
121  config["export_circuit_2_binary_adam"].get_property( export_circuit_2_binary_loc );
122  }
123  else if ( config.count("export_circuit_2_binary") > 0 ) {
124  config["export_circuit_2_binary"].get_property( export_circuit_2_binary_loc );
125  }
126  else {
127  export_circuit_2_binary_loc = 0;
128  }
129 
130 
131  double optimization_tolerance_loc;
132  if ( config.count("optimization_tolerance_adam") > 0 ) {
133  config["optimization_tolerance_adam"].get_property( optimization_tolerance_loc );
134  }
135  else if ( config.count("optimization_tolerance") > 0 ) {
136  config["optimization_tolerance"].get_property( optimization_tolerance_loc );
137  }
138  else {
139  optimization_tolerance_loc = optimization_tolerance;
140  }
141 
142 
143  bool adaptive_eta_loc;
144  if ( config.count("adaptive_eta_adam") > 0 ) {
145  long long tmp;
146  config["adaptive_eta_adam"].get_property( tmp );
147  adaptive_eta_loc = (bool)tmp;
148  }
149  if ( config.count("adaptive_eta") > 0 ) {
150  long long tmp;
151  config["adaptive_eta"].get_property( tmp );
152  adaptive_eta_loc = (bool)tmp;
153  }
154  else {
155  adaptive_eta_loc = adaptive_eta;
156  }
157 
158 
159  double eta_loc;
160  if ( config.count("eta_adam") > 0 ) {
161  config["eta_adam"].get_property( eta_loc );
162  }
163  if ( config.count("eta") > 0 ) {
164  config["eta"].get_property( eta_loc );
165  }
166  else {
167  eta_loc = 1e-3;
168  }
169  optimizer.eta = eta_loc;
170 
171 
172 
173  // The number if iterations after which the current results are displed/exported
174  int output_periodicity;
175  if ( config.count("output_periodicity_cosine") > 0 ) {
176  long long value = 1;
177  config["output_periodicity_cosine"].get_property( value );
178  output_periodicity = (int) value;
179  }
180  if ( config.count("output_periodicity") > 0 ) {
181  long long value = 1;
182  config["output_periodicity"].get_property( value );
183  output_periodicity = (int) value;
184  }
185  else {
186  output_periodicity = 0;
187  }
188 
189 
190  double f0 = DBL_MAX;
191  std::stringstream sstream;
192  sstream << "max_inner_iterations: " << max_inner_iterations_loc << ", randomization threshold: " << iteration_threshold_of_randomization_loc << std::endl;
193  print(sstream, 2);
194 
195 
196  for ( long long iter_idx=0; iter_idx<max_inner_iterations_loc; iter_idx++ ) {
197 
198 
199  optimization_problem_combined( solution_guess_tmp, &f0, grad_mtx );
200 
201  prev_cost_fnv_val = f0;
202 
203  if (sub_iter_idx == 1 ) {
204  current_minimum_hold = f0;
205 
206  if ( adaptive_eta_loc ) {
207  optimizer.eta = optimizer.eta > 1e-3 ? optimizer.eta : 1e-3;
208  //std::cout << "reset learning rate to " << optimizer.eta << std::endl;
209  }
210 
211  }
212 
213 
214  if ((cost_fnc != VQE) && (current_minimum_hold*0.95 > f0 || (current_minimum_hold*0.97 > f0 && f0 < 1e-3) || (current_minimum_hold*0.99 > f0 && f0 < 1e-4) )) {
215  sub_iter_idx = 0;
216  current_minimum_hold = f0;
217  }
218 
219  if (current_minimum > f0 ) {
220  current_minimum = f0;
221  memcpy( optimized_parameters_mtx.get_data(), solution_guess_tmp.get_data(), num_of_parameters*sizeof(double) );
222  //double new_eta = 1e-3 * f0 * f0;
223 
224  if ( adaptive_eta_loc ) {
225  double new_eta = 1e-3 * f0;
226  optimizer.eta = new_eta > 1e-6 ? new_eta : 1e-6;
227  optimizer.eta = new_eta < 1e-1 ? new_eta : 1e-1;
228  }
229 
230  }
231 
232  if ( output_periodicity>0 && iter_idx % output_periodicity == 0 ) {
234  }
235 
236  if ( iter_idx % 5000 == 0 ) {
237  if (cost_fnc != VQE){
238 
239  std::stringstream sstream;
240  sstream << "ADAM: processed iterations " << (double)iter_idx/max_inner_iterations_loc*100 << "%, current minimum:" << current_minimum << ", current cost function:" << optimization_problem(solution_guess_tmp) << ", sub_iter_idx:" << sub_iter_idx <<std::endl;
241  print(sstream, 1);
242  }
243  else{
244  std::stringstream sstream;
245  sstream << "ADAM: processed iterations " << (double)iter_idx/max_inner_iterations_loc*100 << "%, current minimum:" << current_minimum <<", sub_iter_idx:" << sub_iter_idx <<std::endl;
246  print(sstream, 1);
247  }
248  if ( export_circuit_2_binary_loc > 0 ) {
249  std::string filename("initial_circuit_iteration.binary");
250  if (project_name != "") {
251  filename=project_name+ "_" +filename;
252  }
254  }
255  }
256 
257 //std::cout << grad_norm << std::endl;
258  if (f0 < optimization_tolerance_loc || random_shift_count > random_shift_count_max ) {
259  break;
260  }
261 
262 
263 
264  // calculate the gradient norm
265  double norm = 0.0;
266  for ( int grad_idx=0; grad_idx<num_of_parameters; grad_idx++ ) {
267  norm += grad_mtx[grad_idx]*grad_mtx[grad_idx];
268  }
269  norm = std::sqrt(norm);
270 
271 //grad_mtx.print_matrix();
272 /*
273  if ( ADAM_status == 0 && norm > 0.01 && optimizer.eta < 1e-4) {
274 
275  std::uniform_real_distribution<> distrib_prob(0.0, 1.0);
276  if ( distrib_prob(gen) < 0.05 ) {
277  optimizer.eta = optimizer.eta*10;
278  std::cout << "Increasing learning rate at " << f0 << " to " << optimizer.eta << std::endl;
279  }
280 
281  }
282 */
283 /*
284 
285  if ( ADAM_status == 1 && norm > 0.01 ) {
286  optimizer.eta = optimizer.eta > 1e-5 ? optimizer.eta/10 : 1e-6;
287  std::cout << "Decreasing learning rate at " << f0 << " to " << optimizer.eta << std::endl;
288  ADAM_status = 0;
289  }
290 
291  */
292 
293  if ( sub_iter_idx> iteration_threshold_of_randomization_loc || ADAM_status != 0 ) {
294 
295  //random_shift_count++;
296  sub_iter_idx = 0;
297  random_shift_count++;
298  current_minimum_hold = current_minimum;
299 
300 
301 
302  std::stringstream sstream;
303  if ( ADAM_status == 0 ) {
304  sstream << "ADAM: initiate randomization at " << f0 << ", gradient norm " << norm << std::endl;
305  }
306  else {
307  sstream << "ADAM: leaving local minimum " << f0 << ", gradient norm " << norm << " eta: " << optimizer.eta << std::endl;
308  }
309  print(sstream, 1);
310 
311  randomize_parameters(optimized_parameters_mtx, solution_guess_tmp, f0 );
312 
313  optimizer.reset();
314  optimizer.initialize_moment_and_variance( num_of_parameters );
315 
316  ADAM_status = 0;
317 
318  //optimizer.eta = 1e-3;
319 
320  }
321 
322  else {
323  ADAM_status = optimizer.update(solution_guess_tmp, grad_mtx, f0);
324  }
325 
326  sub_iter_idx++;
327 
328  }
329  sstream.str("");
330  sstream << "obtained minimum: " << current_minimum << std::endl;
331 
332 
333  tbb::tick_count adam_end = tbb::tick_count::now();
334  CPU_time = CPU_time + (adam_end-adam_start).seconds();
335  sstream << "adam time: " << CPU_time << " " << f0 << std::endl;
336 
337  print(sstream, 1);
338 
339 }
340 
341 
int update(Matrix_real &parameters, Matrix_real &grad, const double &f0)
Call to set the number of gate blocks to be optimized in one shot.
bool adaptive_eta
logical variable indicating whether adaptive learning reate is used in the ADAM algorithm ...
Header file for a class containing basic methods for the decomposition process.
void export_current_cost_fnc(double current_minimum)
Call to print out into a file the current cost function and the second Rényi entropy on the subsyste...
void print(const std::stringstream &sstream, int verbose_level=1) const
Call to print output messages in the function of the verbosity level.
Definition: logging.cpp:55
cost_function_type cost_fnc
The chosen variant of the cost function.
int get_accelerator_num()
Get the number of accelerators to be reserved on DFEs on users demand.
double optimization_problem(double *parameters)
Evaluate the optimization problem of the optimization.
double prev_cost_fnv_val
the previous value of the cost funtion to be used to evaluate bitflip errors in the cost funtion (see...
scalar * get_data() const
Call to get the pointer to the stored data.
void initialize_moment_and_variance(int parameter_num)
?????????????
std::vector< Gate * > gates
The list of stored gates.
Definition: Gates_block.h:49
std::string project_name
the name of the project
double optimization_tolerance
The maximal allowed error of the optimization problem (The error of the decomposition would scale wit...
void reset()
?????????????
Definition: common/Adam.cpp:90
double eta
Definition: Adam.h:42
void randomize_parameters(Matrix_real &input, Matrix_real &output, const double f0)
Call to randomize the parameter.
double CPU_time
time spent on optimization
int verbose
Set the verbosity level of the output messages.
Definition: logging.h:50
int size() const
Call to get the number of the allocated elements.
static void optimization_problem_combined(Matrix_real parameters, void *void_instance, double *f0, Matrix_real &grad)
Call to calculate both the cost function and the its gradient components.
std::map< std::string, Config_Element > config
config metadata utilized during the optimization
void solve_layer_optimization_problem_ADAM(int num_of_parameters, Matrix_real &solution_guess)
Call to solve layer by layer the optimization problem via ADAM algorithm.
A class for Adam optimization according to https://towardsdatascience.com/how-to-implement-an-adam-op...
Definition: Adam.h:36
Header file for the paralleized calculation of the cost function of the final optimization problem (s...
volatile double current_minimum
The current minimum of the optimization problem.
void export_gate_list_to_binary(Matrix_real &parameters, Gates_block *gates_block, const std::string &filename, int verbosity)
Use to export a quantum circuit into binary format.
int qbit_num
number of qubits spanning the matrix of the operation
Definition: Gate.h:94
Header file for DFE support in unitary simulation.
int max_inner_iterations
the maximal number of iterations for which an optimization engine tries to solve the optimization pro...
int random_shift_count_max
the maximal number of parameter randomization tries to escape a local minimum.
Matrix_real optimized_parameters_mtx
The optimized parameters for the gates.
Class to store data of complex arrays and its properties.
Definition: matrix_real.h:41