ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
unitary.h
1 /*
2  * This source code is part of
3  *
4  * E R K A L E
5  * -
6  * HF/DFT from Hel
7  *
8  * Written by Susi Lehtola, 2010-2013
9  * Copyright (c) 2010-2013, Susi Lehtola
10  *
11  * This program is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU General Public License
13  * as published by the Free Software Foundation; either version 2
14  * of the License, or (at your option) any later version.
15  */
16 
17 #ifndef ERKALE_UNITARY
18 #define ERKALE_UNITARY
19 
20 #include "global.h"
21 class Timer;
22 #include <armadillo>
23 
24 
38 enum unitmethod {
40  POLY_F,
42  POLY_DF,
44  FOURIER_DF,
46  ARMIJO
47 };
48 
49 enum unitacc {
51  SDSA,
53  CGPR,
55  CGFR,
57  CGHS
58 };
59 
60 
63  protected:
65  arma::cx_mat W;
67  double f;
69  int q;
71  int sign;
72 
73  public:
75  UnitaryFunction(int q, bool max);
77  virtual ~UnitaryFunction();
78 
80  virtual void setW(const arma::cx_mat & W);
82  arma::cx_mat getW() const;
83 
85  int getq() const;
87  double getf() const;
89  int getsign() const;
90 
92  virtual UnitaryFunction *copy() const=0;
94  virtual double cost_func(const arma::cx_mat & W)=0;
96  virtual arma::cx_mat cost_der(const arma::cx_mat & W)=0;
98  virtual void cost_func_der(const arma::cx_mat & W, double & f, arma::cx_mat & der)=0;
99 
101  virtual std::string legend() const;
103  virtual std::string status(bool lfmt=false);
105  virtual bool converged();
106 };
107 
110  private:
112  arma::cx_mat G;
114  arma::cx_mat H;
116  arma::cx_mat Hvec;
118  arma::vec Hval;
120  double Tmu;
121 
122  protected:
124  bool verbose;
126  bool real;
127 
129  double Gthr;
131  double Fthr;
132 
135 
140 
142  bool debug;
143 
145  FILE *log;
146 
148  virtual void print_legend(const UnitaryFunction *f) const;
150  virtual void print_progress(size_t k, UnitaryFunction *f, const UnitaryFunction *fold) const;
152  virtual void print_time(const Timer & t) const;
154  virtual void print_step(enum unitmethod & met, double step) const;
155 
157  void check_unitary(const arma::cx_mat & W) const;
159  void check_derivative(const UnitaryFunction *f);
161  void classify(const arma::cx_mat & W) const;
162 
164  void update_gradient(const arma::cx_mat & W, UnitaryFunction *f);
166  void update_search_direction(int q);
167 
169  arma::cx_mat get_rotation(double step) const;
171  double step_der(const arma::cx_mat & W, const arma::cx_mat & der) const;
172 
174  void set_q(int q);
175 
177  void armijo_step(UnitaryFunction* & f);
183  void fourier_step_df(UnitaryFunction* & f);
184 
185  public:
187  UnitaryOptimizer(double Gthr, double Fthr, bool verbose=true, bool real=false);
189  virtual ~UnitaryOptimizer();
190 
192  void open_log(const std::string & fname);
194  void set_debug(bool dbg);
195 
197  void set_poly(int deg);
199  void set_fourier(int Nsamples, int Nperiods);
201  void set_thr(double Gtol, double Ftol);
202 
204  double optimize(UnitaryFunction* & f, enum unitmethod met, enum unitacc acc, size_t maxiter);
205 };
206 
207 
209 double bracket(const arma::cx_mat & X, const arma::cx_mat & Y);
210 
212 arma::cx_vec fourier_shift(const arma::cx_vec & v);
213 
215 arma::vec fit_polynomial(const arma::vec & x, const arma::vec & y, int deg=-1);
217 arma::vec fit_polynomial_fdf(const arma::vec & x, const arma::vec & y, const arma::vec & dy, int deg=-1);
218 
220 arma::vec derivative_coefficients(const arma::vec & c);
221 
223 arma::cx_vec solve_roots_cplx(const arma::cx_vec & v);
225 arma::cx_vec solve_roots_cplx(const arma::vec & v);
227 arma::vec solve_roots(const arma::vec & v);
229 double smallest_positive(const arma::vec & v);
230 
231 
233 class Brockett : public UnitaryFunction {
235  arma::cx_mat sigma;
237  arma::mat Nmat;
238 
240  std::string legend() const;
242  std::string status(bool lfmt=false);
243 
245  double diagonality() const;
247  double unitarity() const;
248 
249  public:
250  Brockett(size_t N, unsigned long int seed=0);
251  ~Brockett();
252 
254  Brockett *copy() const;
256  double cost_func(const arma::cx_mat & W);
258  arma::cx_mat cost_der(const arma::cx_mat & W);
260  void cost_func_der(const arma::cx_mat & W, double & f, arma::cx_mat & der);
261 };
262 
263 
264 #endif
int fourier_periods
Amount of quasi-periods for Fourier method (N_T = 1, 2, ...)
Definition: unitary.h:137
UnitaryFunction(int q, bool max)
Constructor.
Definition: unitary.cpp:26
arma::cx_mat Hvec
Eigenvectors of search direction.
Definition: unitary.h:116
void open_log(const std::string &fname)
Open log file.
Definition: unitary.cpp:99
void classify(const arma::cx_mat &W) const
Classify matrix.
Definition: unitary.cpp:387
int fourier_samples
Amount of samples per one period (K = 3, 4, or 5)
Definition: unitary.h:139
void update_gradient(const arma::cx_mat &W, UnitaryFunction *f)
Get new gradient direction.
Definition: unitary.cpp:149
A timer routine.
Definition: timer.h:43
double Gthr
Convergence threshold wrt norm of Riemannian derivative.
Definition: unitary.h:129
int sign
Maximization or minimization?
Definition: unitary.h:71
void check_derivative(const UnitaryFunction *f)
Check that the programmed cost function and its derivative are OK.
Definition: unitary.cpp:406
double Fthr
Convergence threshold wrt relative change in function.
Definition: unitary.h:131
double cost_func(const arma::cx_mat &W)
Evaluate cost function.
Definition: unitary.cpp:1124
arma::cx_mat get_rotation(double step) const
Get rotation matrix with wanted step size.
Definition: unitary.cpp:130
virtual std::string status(bool lfmt=false)
Print status information, possibly in a longer format.
Definition: unitary.cpp:64
void polynomial_step_f(UnitaryFunction *&f)
Polynomial step (fit function)
Definition: unitary.cpp:448
virtual arma::cx_mat cost_der(const arma::cx_mat &W)=0
Evaluate derivative of cost function.
double step_der(const arma::cx_mat &W, const arma::cx_mat &der) const
Get derivative wrt step length.
Definition: unitary.cpp:644
double unitarity() const
Compute unitarity criterion.
Definition: unitary.cpp:1174
std::string legend() const
Print legend.
Definition: unitary.cpp:1140
virtual ~UnitaryFunction()
Destructor.
Definition: unitary.cpp:31
FILE * log
Log file.
Definition: unitary.h:145
double Tmu
Maximum step size.
Definition: unitary.h:120
virtual void print_progress(size_t k, UnitaryFunction *f, const UnitaryFunction *fold) const
Print progress.
Definition: unitary.cpp:338
arma::cx_mat H
Search direction.
Definition: unitary.h:114
arma::cx_mat G
Gradient.
Definition: unitary.h:112
virtual UnitaryFunction * copy() const =0
Copy constructor.
UnitaryOptimizer(double Gthr, double Fthr, bool verbose=true, bool real=false)
Constructor.
Definition: unitary.cpp:70
void set_q(int q)
Set degree.
void set_fourier(int Nsamples, int Nperiods)
Set Fourier search options.
Definition: unitary.cpp:144
Unitary function optimizer, used to hold values during the optimization.
Definition: unitary.h:62
void set_thr(double Gtol, double Ftol)
Set convergence threshold.
Definition: unitary.cpp:94
void fourier_step_df(UnitaryFunction *&f)
Fourier step.
Definition: unitary.cpp:747
void check_unitary(const arma::cx_mat &W) const
Check that the matrix is unitary.
Definition: unitary.cpp:119
arma::cx_mat sigma
Sigma matrix.
Definition: unitary.h:235
void set_debug(bool dbg)
Set debug mode.
Definition: unitary.cpp:90
virtual void cost_func_der(const arma::cx_mat &W, double &f, arma::cx_mat &der)=0
Evaluate cost function and its derivative.
int getsign() const
Get sign.
Definition: unitary.cpp:50
arma::cx_mat getW() const
Get matrix.
Definition: unitary.cpp:38
std::string status(bool lfmt=false)
Print progress.
Definition: unitary.cpp:1146
void armijo_step(UnitaryFunction *&f)
Armijo step.
Definition: unitary.cpp:648
virtual void setW(const arma::cx_mat &W)
Set matrix.
Definition: unitary.cpp:34
void polynomial_step_df(UnitaryFunction *&f)
Polynomial step (fit only derivative)
Definition: unitary.cpp:533
virtual void print_step(enum unitmethod &met, double step) const
Print chosen step length.
Definition: unitary.cpp:365
int polynomial_degree
Degree of polynomial used for fit: a_0 + a_1*mu + ... + a_(d-1)*mu^(d-1)
Definition: unitary.h:134
bool real
Operate with real or complex matrices?
Definition: unitary.h:126
virtual std::string legend() const
Get status legend.
Definition: unitary.cpp:59
virtual bool converged()
Check convergence.
Definition: unitary.cpp:54
bool verbose
Verbose operation?
Definition: unitary.h:124
double getf() const
Get function value.
Definition: unitary.cpp:46
virtual double cost_func(const arma::cx_mat &W)=0
Evaluate cost function.
double diagonality() const
Compute diagonality criterion.
Definition: unitary.cpp:1155
virtual void print_legend(const UnitaryFunction *f) const
Print legend.
Definition: unitary.cpp:334
virtual void print_time(const Timer &t) const
Print time.
Definition: unitary.cpp:355
Brockett * copy() const
Copy constructor.
Definition: unitary.cpp:1120
Unitary optimization worker.
Definition: unitary.h:109
Brockett.
Definition: unitary.h:233
int q
Order in W.
Definition: unitary.h:69
double optimize(UnitaryFunction *&f, enum unitmethod met, enum unitacc acc, size_t maxiter)
Unitary optimization.
Definition: unitary.cpp:175
int getq() const
Get q.
Definition: unitary.cpp:42
bool debug
Debugging mode - print out line search every iteration.
Definition: unitary.h:142
arma::cx_mat cost_der(const arma::cx_mat &W)
Evaluate derivative of cost function.
Definition: unitary.cpp:1130
void set_poly(int deg)
Set polynomial search options.
Definition: unitary.cpp:140
void update_search_direction(int q)
Compute new search direction (diagonalize H) and max step length.
Definition: unitary.cpp:161
arma::mat Nmat
N matrix.
Definition: unitary.h:237
virtual ~UnitaryOptimizer()
Destructor.
Definition: unitary.cpp:85
arma::cx_mat W
Present matrix.
Definition: unitary.h:65
arma::vec Hval
Eigenvalues of search direction.
Definition: unitary.h:118
double f
Present value.
Definition: unitary.h:67
void cost_func_der(const arma::cx_mat &W, double &f, arma::cx_mat &der)
Evaluate cost function and its derivative.
Definition: unitary.cpp:1135