PEXSI
 All Classes Namespaces Files Functions Variables Typedefs Pages
ppexsi.hpp
Go to the documentation of this file.
1 /*
2  Copyright (c) 2012 The Regents of the University of California,
3  through Lawrence Berkeley National Laboratory.
4 
5  Author: Lin Lin
6 
7  This file is part of PEXSI. All rights reserved.
8 
9  Redistribution and use in source and binary forms, with or without
10  modification, are permitted provided that the following conditions are met:
11 
12  (1) Redistributions of source code must retain the above copyright notice, this
13  list of conditions and the following disclaimer.
14  (2) Redistributions in binary form must reproduce the above copyright notice,
15  this list of conditions and the following disclaimer in the documentation
16  and/or other materials provided with the distribution.
17  (3) Neither the name of the University of California, Lawrence Berkeley
18  National Laboratory, U.S. Dept. of Energy nor the names of its contributors may
19  be used to endorse or promote products derived from this software without
20  specific prior written permission.
21 
22  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
23  ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
24  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
25  DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
26  ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
27  (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
28  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
29  ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
30  (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
31  SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
32 
33  You are under no obligation whatsoever to provide any bug fixes, patches, or
34  upgrades to the features, functionality or performance of the source code
35  ("Enhancements") to anyone; however, if you choose to make your Enhancements
36  available either publicly, or directly to Lawrence Berkeley National
37  Laboratory, without imposing a separate written license agreement for such
38  Enhancements, then you hereby grant the following license: a non-exclusive,
39  royalty-free perpetual license to install, use, modify, prepare derivative
40  works, incorporate into other computer software, distribute, and sublicense
41  such enhancements or derivative works thereof, in binary and source code form.
42 */
47 #ifndef _PPEXSI_HPP_
48 #define _PPEXSI_HPP_
49 #include "pexsi/environment.hpp"
50 #include "pexsi/sparse_matrix.hpp"
51 #include "pexsi/NumVec.hpp"
52 #include "pexsi/utility.hpp"
53 #include "pexsi/pole.hpp"
54 #include "pexsi/mpi_interf.hpp"
55 #include "pexsi/SuperLUGrid.hpp"
57 #include "pexsi/pselinv.hpp"
58 #include "pexsi/pselinv_unsym.hpp"
59 //#include "pexsi/ngchol_interf.hpp"
60 //#include "pexsi/c_pexsi_interface.h"
61 
62 namespace PEXSI{
63 
68  class PPEXSIData{
69  private:
70  // *********************************************************************
71  // Computational variables
72  // *********************************************************************
73 
74  std::vector<Complex> zshift_; // Complex shift for the pole expansion
75  std::vector<Complex> zweightRho_; // Complex weight for the pole expansion for density
76  std::vector<Complex> zweightRhoDrvMu_; // Complex weight for the pole expansion for derivative of the Fermi-Dirac with respect to the chemical potential
77  std::vector<Complex> zweightRhoDrvT_; // Complex weight for the pole expansion for derivative of the Fermi-Dirac with respect to the temperature T (1/beta, in au)
78  std::vector<Complex> zweightHelmholtz_; // Complex shift for the pole expansion for Helmholtz free energy
79  std::vector<Complex> zweightForce_; // Complex weight for the pole expansion for force
80 
81  // Outer layer communicator. Also used for distributing the
82  // DistSparseMatrix. Each DistSparseMatrix is replicated in the row
83  // (numPoleGroup) direction of gridPole.
84  const GridType* gridPole_;
85  const GridType* gridSelInv_; // Inner layer communicator for SelInv
86 
87  // Inner layer communicator for SuperLU factorization
88  const SuperLUGrid<Real>* gridSuperLUReal_;
89  const SuperLUGrid<Complex>* gridSuperLUComplex_;
90 
91 
92  DistSparseMatrix<Real> HRealMat_;
93  DistSparseMatrix<Real> SRealMat_;
94 
95  DistSparseMatrix<Real> shiftRealMat_;
96  DistSparseMatrix<Complex> shiftComplexMat_;
97 
98  DistSparseMatrix<Real> shiftInvRealMat_;
99  DistSparseMatrix<Complex> shiftInvComplexMat_;
100 
101  DistSparseMatrix<Real> rhoRealMat_; // Density matrix
102  DistSparseMatrix<Real> rhoDrvMuRealMat_; // Derivative of the Fermi-Dirac with respect to mu
103  DistSparseMatrix<Real> rhoDrvTRealMat_; // Derivative of the Fermi-Dirac with respect to T
104  DistSparseMatrix<Real> freeEnergyDensityRealMat_; // Helmholtz free energy density matrix
105  DistSparseMatrix<Real> energyDensityRealMat_; // Energy density matrix for computing the Pulay force
106 
107  // SuperLUMatrix and PMatrix structures These structures are saved
108  // to avoid repetitive symbolic factorization process, and saved in
109  // pointer form because of the constructors.
110  SuperLUMatrix<Real>* luRealMat_;
111  SuperLUMatrix<Complex>* luComplexMat_;
112 
113  SuperLUOptions luOpt_;
114  PSelInvOptions selinvOpt_;
115 
116  PMatrix<Real>* PMRealMat_;
117  PMatrix<Complex>* PMComplexMat_;
118  PMatrixUnsym<Real>* PMRealUnsymMat_;
119  PMatrixUnsym<Complex>* PMComplexUnsymMat_;
120 
121  // Whether the matrices have been loaded into HRealMat_ and
122  // SRealMat_
123  bool isMatrixLoaded_;
124  // Whether the matrices (luMat and PMat) have obtained symbolic
125  // information
126  bool isRealSymmetricSymbolicFactorized_;
127  bool isComplexSymmetricSymbolicFactorized_;
128  bool isRealUnsymmetricSymbolicFactorized_;
129  bool isComplexUnsymmetricSymbolicFactorized_;
130  // Supernode partition for the real matrix
131  SuperNodeType superReal_;
132  // Supernode partition for the complex matrix
133  SuperNodeType superComplex_;
134 
135  // Saves all the indices of diagonal elements in H, so that
136  // H.nzvalLocal(diagIdxLocal_[j]) are diagonal elements for all j.
137  // This is manly used when S is implicitly given as an identity matrix.
138  std::vector<Int> diagIdxLocal_;
139 
140  // Energy computed from Tr[H*DM]
141  Real totalEnergyH_;
142  // Energy computed from Tr[S*EDM]
143  Real totalEnergyS_;
144  // Free energy
145  Real totalFreeEnergy_;
146 
147 
148  // *********************************************************************
149  // Saved variables for nonlinear iterations
150  // *********************************************************************
151 
152  public:
153  PPEXSIData(
154  MPI_Comm comm,
155  Int numProcRow,
156  Int numProcCol,
157  Int outputFileIndex );
158 
159  ~PPEXSIData();
160 
161  void LoadRealSymmetricMatrix(
162  Int nrows,
163  Int nnz,
164  Int nnzLocal,
165  Int numColLocal,
166  Int* colptrLocal,
167  Int* rowindLocal,
168  Real* HnzvalLocal,
169  Int isSIdentity,
170  Real* SnzvalLocal,
171  Int verbosity );
172 
173  void LoadRealUnsymmetricMatrix(
174  Int nrows,
175  Int nnz,
176  Int nnzLocal,
177  Int numColLocal,
178  Int* colptrLocal,
179  Int* rowindLocal,
180  Real* HnzvalLocal,
181  Int isSIdentity,
182  Real* SnzvalLocal,
183  Int verbosity );
184 
185 
186 
202  std::string ColPerm,
203  Int numProcSymbFact,
204  Int verbosity );
205 
224  std::string ColPerm,
225  std::string RowPerm,
226  Int numProcSymbFact,
227  Int Transpose,
228  double* AnzvalLocal,
229  Int verbosity );
230 
231 
247  std::string ColPerm,
248  Int numProcSymbFact,
249  Int verbosity );
250 
269  std::string ColPerm,
270  std::string RowPerm,
271  Int numProcSymbFact,
272  Int Transpose,
273  double* AnzvalLocal,
274  Int verbosity );
275 
276 
277 
278  void SelInvRealSymmetricMatrix(
279  double* AnzvalLocal,
280  Int verbosity,
281  double* AinvnzvalLocal );
282 
283  void SelInvRealUnsymmetricMatrix(
284  double* AnzvalLocal,
285  Int verbosity,
286  double* AinvnzvalLocal );
287 
288 
289  void SelInvComplexSymmetricMatrix(
290  double* AnzvalLocal,
291  Int verbosity,
292  double* AinvnzvalLocal );
293 
294  void SelInvComplexUnsymmetricMatrix(
295  double* AnzvalLocal,
296  Int verbosity,
297  double* AinvnzvalLocal );
298 
299 
300 
330  const std::vector<Real>& shiftVec,
331  std::vector<Real>& inertiaVec,
332  Int verbosity );
333 
334 
361  Int numPole,
362  Real temperature,
363  Real gap,
364  Real deltaE,
365  Real mu,
366  Real numElectronExact,
367  Real numElectronTolerance,
368  Int verbosity,
369  Real& numElectron,
370  Real& numElectronDrvMu );
371 
372 
374  void DFTDriver(
375  Real numElectronExact,
376  Real temperature,
377  Real gap,
378  Real deltaE,
379  Int numPole,
380  Int isInertiaCount,
381  Int maxPEXSIIter,
382  Real muMin0,
383  Real muMax0,
384  Real mu0,
385  Real muInertiaTolerance,
386  Real muInertiaExpansion,
387  Real muPEXSISafeGuard,
388  Real numElectronPEXSITolerance,
389  Int matrixType,
390  Int isSymbolicFactorize,
391  Int ordering,
392  Int numProcSymbFact,
393  Int verbosity,
394  Real& muPEXSI,
395  Real& numElectronPEXSI,
396  Real& muMinInertia,
397  Real& muMaxInertia,
398  Int& numTotalInertiaIter,
399  Int& numTotalPEXSIIter );
400 
401 
402  // *********************************************************************
403  // Access data
404  // *********************************************************************
405 
406  const GridType* GridPole() const {return gridPole_;}
407 
412  const DistSparseMatrix<Real>& RhoRealMat() const {return rhoRealMat_;}
413 
419  const DistSparseMatrix<Real>& EnergyDensityRealMat() const {return energyDensityRealMat_;}
420 
429  const DistSparseMatrix<Real>& FreeEnergyDensityRealMat() const {return freeEnergyDensityRealMat_;}
430 
431  Real TotalEnergyH() const {return totalEnergyH_;}
432 
433  Real TotalEnergyS() const {return totalEnergyS_;}
434 
435  Real TotalFreeEnergy() const {return totalFreeEnergy_;}
436 
437 
438  }; // PPEXSIData
439 
440 
441 } // namespace PEXSI
442 #endif // _PPEXSI_HPP_
Environmental variables.
A thin interface for passing parameters to set the SuperLU options.
Definition: superlu_dist_internal.hpp:62
const DistSparseMatrix< Real > & FreeEnergyDensityRealMat() const
Total Helmholtz free energy matrix (band energy part only).
Definition: ppexsi.hpp:429
SuperNodeType describes mapping between supernode and column, the permutation information, and potentially the elimination tree (not implemented here).
Definition: pselinv.hpp:165
void SymbolicFactorizeComplexSymmetricMatrix(std::string ColPerm, Int numProcSymbFact, Int verbosity)
Symbolically factorize the loaded matrices for complex arithmetic factorization and selected inversio...
Definition: ppexsi.cpp:725
Interface with SuperLU_Dist (version 3.0 and later)
void DFTDriver(Real numElectronExact, Real temperature, Real gap, Real deltaE, Int numPole, Int isInertiaCount, Int maxPEXSIIter, Real muMin0, Real muMax0, Real mu0, Real muInertiaTolerance, Real muInertiaExpansion, Real muPEXSISafeGuard, Real numElectronPEXSITolerance, Int matrixType, Int isSymbolicFactorize, Int ordering, Int numProcSymbFact, Int verbosity, Real &muPEXSI, Real &numElectronPEXSI, Real &muMinInertia, Real &muMaxInertia, Int &numTotalInertiaIter, Int &numTotalPEXSIIter)
Main driver for solving KSDFT.
Definition: ppexsi.cpp:2197
const DistSparseMatrix< Real > & EnergyDensityRealMat() const
Energy density matrix.
Definition: ppexsi.hpp:419
Main class for parallel PEXSI.
Definition: ppexsi.hpp:68
Main file for parallel selected inversion on unsymmetric matrices.
Main file for parallel selected inversion.
SuperLU processor grid.
Definition: SuperLUMatrix.hpp:206
void CalculateFermiOperatorReal(Int numPole, Real temperature, Real gap, Real deltaE, Real mu, Real numElectronExact, Real numElectronTolerance, Int verbosity, Real &numElectron, Real &numElectronDrvMu)
Compute the Fermi operator for a given chemical potential for real symmetric matrices.
Definition: ppexsi.cpp:1595
Definition: SuperLUMatrix.hpp:337
Sparse matrix and Distributed sparse matrix in compressed column format.
GridType is the PSelInv way of defining the grid.
Definition: pselinv.hpp:128
Interface with MPI to facilitate communication.
void CalculateNegativeInertiaReal(const std::vector< Real > &shiftVec, std::vector< Real > &inertiaVec, Int verbosity)
Compute the negative inertia (the number of eigenvalues below a shift) for real symmetric matrices...
Definition: ppexsi.cpp:1435
const DistSparseMatrix< Real > & RhoRealMat() const
Density matrix.
Definition: ppexsi.hpp:412
Various utility subroutines.
void SymbolicFactorizeComplexUnsymmetricMatrix(std::string ColPerm, std::string RowPerm, Int numProcSymbFact, Int Transpose, double *AnzvalLocal, Int verbosity)
Symbolically factorize the loaded matrices for complex arithmetic factorization and selected inversio...
Definition: ppexsi.cpp:834
void SymbolicFactorizeRealSymmetricMatrix(std::string ColPerm, Int numProcSymbFact, Int verbosity)
Symbolically factorize the loaded matrices for real arithmetic factorization and selected inversion...
Definition: ppexsi.cpp:473
Definition: SuperLUGrid.hpp:124
Definition: SuperLUGrid.hpp:111
Numerical vector.
void SymbolicFactorizeRealUnsymmetricMatrix(std::string ColPerm, std::string RowPerm, Int numProcSymbFact, Int Transpose, double *AnzvalLocal, Int verbosity)
Symbolically factorize the loaded matrices for real arithmetic factorization and selected inversion...
Definition: ppexsi.cpp:586
A thin interface for passing parameters to set the PSelInv options.
Definition: pselinv.hpp:101
Pole expansion subroutines.