HΦ  3.2.0
CalcSpectrumByLanczos.c File Reference

Get the spectrum function by continued fraction expansions.
Ref. E.Dagotto, Rev. Mod. Phys. 66 (1994), 763. More...

#include "CalcSpectrumByLanczos.h"
#include "Lanczos_EigenValue.h"
#include "FileIO.h"
#include "wrapperMPI.h"
#include "common/setmemory.h"
#include "CalcTime.h"
+ Include dependency graph for CalcSpectrumByLanczos.c:

Go to the source code of this file.

Functions

int CalcSpectrumByLanczos (struct EDMainCalStruct *X, double complex *tmp_v1, double dnorm, int Nomega, double complex *dcSpectrum, double complex *dcomega)
 A main function to calculate spectrum by continued fraction expansions. More...
 
int GetSpectrumByTridiagonalMatrixComponents (double *tmp_alpha, double *tmp_beta, double dnorm, double complex dcomega, double complex *dcSpectrum, unsigned long int ilLanczosStp)
 Calculate the spectrum by using tridiagonal matrix components obtained by the Lanczos_GetTridiagonalMatrixComponents function. More...
 

Detailed Description

Get the spectrum function by continued fraction expansions.
Ref. E.Dagotto, Rev. Mod. Phys. 66 (1994), 763.

Version
1.1
Author
Kazuyoshi Yoshimi (The University of Tokyo)

Definition in file CalcSpectrumByLanczos.c.

Function Documentation

◆ CalcSpectrumByLanczos()

int CalcSpectrumByLanczos ( struct EDMainCalStruct X,
double complex *  tmp_v1,
double  dnorm,
int  Nomega,
double complex *  dcSpectrum,
double complex *  dcomega 
)

A main function to calculate spectrum by continued fraction expansions.

Parameters
X[in,out] Struct for getting and giving calculation information
tmp_v1[in] Normalized excited state.
dnorm[in] Norm of the excited state before normalization.
Nomega[in] Total number of frequencies.
dcSpectrum[out] Calculated spectrum.
dcomega[in] Target frequencies.
Return values
0normally finished
-1unnormally finished
Version
1.1
Author
Kazuyoshi Yoshimi (The University of Tokyo)

Definition at line 43 of file CalcSpectrumByLanczos.c.

Referenced by CalcSpectrum().

51 {
52  unsigned long int i;
53  int iret;
54  unsigned long int liLanczosStp = X->Bind.Def.Lanczos_max;
55  unsigned long int liLanczosStp_vec=0;
56 
57  if(X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents_VEC ||
58  X->Bind.Def.iFlgCalcSpec == RECALC_INOUT_TMComponents_VEC) {
59 
60  StartTimer(6201);
61  if(ReadInitialVector( &(X->Bind), v0, v1, &liLanczosStp_vec)!=0){
62  StopTimer(6201);
63  exitMPI(-1);
64  }
65  StopTimer(6201);
66  }
67 
68  //Read diagonal components
69  if(X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents ||
70  X->Bind.Def.iFlgCalcSpec ==RECALC_FROM_TMComponents_VEC||
71  X->Bind.Def.iFlgCalcSpec == RECALC_INOUT_TMComponents_VEC)
72  {
73  StartTimer(6202);
74  int iFlgTMComp=0;
75  if(X->Bind.Def.iFlgCalcSpec == RECALC_INOUT_TMComponents_VEC ||
76  X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents_VEC)
77  {
78  iFlgTMComp=0;
79  }
80  else{
81  iFlgTMComp=1;
82  }
83  iret=ReadTMComponents(&(X->Bind), &dnorm, &liLanczosStp, iFlgTMComp);
84  if(iret !=TRUE){
85  fprintf(stdoutMPI, " Error: Fail to read TMcomponents\n");
86  return FALSE;
87  }
88 
89  if(X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents){
91  }
92  else if(X->Bind.Def.iFlgCalcSpec == RECALC_INOUT_TMComponents_VEC||
93  X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents_VEC){
94  if(liLanczosStp_vec !=liLanczosStp){
95  fprintf(stdoutMPI, " Error: Input files for vector and TMcomponents are incoorect.\n");
96  fprintf(stdoutMPI, " Error: Input vector %ld th stps, TMcomponents %ld th stps.\n", liLanczosStp_vec, liLanczosStp);
97  return FALSE;
98  }
99  X->Bind.Def.Lanczos_restart=liLanczosStp;
100  liLanczosStp = liLanczosStp+X->Bind.Def.Lanczos_max;
101  }
102  StopTimer(6202);
103  }
104 
105  // calculate ai, bi
106  if (X->Bind.Def.iFlgCalcSpec == RECALC_NOT ||
107  X->Bind.Def.iFlgCalcSpec == RECALC_OUTPUT_TMComponents_VEC ||
108  X->Bind.Def.iFlgCalcSpec == RECALC_FROM_TMComponents_VEC ||
109  X->Bind.Def.iFlgCalcSpec == RECALC_INOUT_TMComponents_VEC
110  )
111  {
112  fprintf(stdoutMPI, " Start: Calculate tridiagonal matrix components.\n");
114  // Functions in Lanczos_EigenValue
115  StartTimer(6203);
116  iret = Lanczos_GetTridiagonalMatrixComponents(&(X->Bind), alpha, beta, tmp_v1, &(liLanczosStp));
117  StopTimer(6203);
118  if (iret != TRUE) {
119  //Error Message will be added.
120  return FALSE;
121  }
122  fprintf(stdoutMPI, " End: Calculate tridiagonal matrix components.\n\n");
124  StartTimer(6204);
125  OutputTMComponents(&(X->Bind), alpha,beta, dnorm, liLanczosStp);
126  StopTimer(6204);
127  }//X->Bind.Def.iFlgCalcSpec == RECALC_NOT || RECALC_FROM_TMComponents_VEC
128 
129  fprintf(stdoutMPI, " Start: Caclulate spectrum from tridiagonal matrix components.\n");
131  StartTimer(6205);
132  for( i = 0 ; i < Nomega; i++) {
133  iret = GetSpectrumByTridiagonalMatrixComponents(alpha, beta, dnorm, dcomega[i], &dcSpectrum[i], liLanczosStp);
134  if (iret != TRUE) {
135  //ToDo: Error Message will be added.
136  //ReAlloc alpha, beta and Set alpha_start and beta_start in Lanczos_EigenValue
137  return FALSE;
138  }
139  dcSpectrum[i] = dnorm * dcSpectrum[i];
140  }
141  StopTimer(6205);
142  fprintf(stdoutMPI, " End: Caclulate spectrum from tridiagonal matrix components.\n\n");
144 
145  //output vectors for recalculation
146  if(X->Bind.Def.iFlgCalcSpec==RECALC_OUTPUT_TMComponents_VEC ||
147  X->Bind.Def.iFlgCalcSpec==RECALC_INOUT_TMComponents_VEC){
148  StartTimer(6206);
149  if(OutputLanczosVector(&(X->Bind), v0, v1, liLanczosStp)!=0){
150  StopTimer(6206);
151  exitMPI(-1);
152  }
153  StopTimer(6206);
154  }
155 
156  return TRUE;
157 }
void exitMPI(int errorcode)
MPI Abortation wrapper.
Definition: wrapperMPI.c:86
struct DefineList Def
Definision of system (Hamiltonian) etc.
Definition: struct.h:410
void StartTimer(int n)
function for initializing elapse time [start]
Definition: time.c:71
void StopTimer(int n)
function for calculating elapse time [elapse time=StartTimer-StopTimer]
Definition: time.c:83
double complex * v1
Definition: global.h:35
int Lanczos_GetTridiagonalMatrixComponents(struct BindStruct *X, double *_alpha, double *_beta, double complex *tmp_v1, unsigned long int *liLanczos_step)
Calculate tridiagonal matrix components by Lanczos method.
#define TRUE
Definition: global.h:26
int OutputLanczosVector(struct BindStruct *X, double complex *tmp_v0, double complex *tmp_v1, unsigned long int liLanczosStp_vec)
Output initial vectors for the restart calculation.
int Lanczos_restart
Number of iterations performed in the restart computation.
Definition: struct.h:75
const char * c_GetTridiagonalStart
Definition: LogMessage.c:66
const char * c_CalcSpectrumFromTridiagonalEnd
Definition: LogMessage.c:69
int ReadInitialVector(struct BindStruct *X, double complex *_v0, double complex *_v1, unsigned long int *liLanczosStp_vec)
Read initial vectors for the restart calculation.
double * alpha
Definition: global.h:44
#define FALSE
Definition: global.h:25
const char * cFileNameTimeKeep
Definition: global.c:23
double * beta
Definition: global.h:44
int OutputTMComponents(struct BindStruct *X, double *_alpha, double *_beta, double _dnorm, int liLanczosStp)
Output tridiagonal matrix components obtained by the Lanczos method.
int GetSpectrumByTridiagonalMatrixComponents(double *tmp_alpha, double *tmp_beta, double dnorm, double complex dcomega, double complex *dcSpectrum, unsigned long int ilLanczosStp)
Calculate the spectrum by using tridiagonal matrix components obtained by the Lanczos_GetTridiagonalM...
int iFlgCalcSpec
Input parameter CalcSpec in teh CalcMod file.
Definition: struct.h:216
const char * c_CalcSpectrumFromTridiagonalStart
Definition: LogMessage.c:68
const char * c_GetTridiagonalEnd
Definition: LogMessage.c:67
int ReadTMComponents(struct BindStruct *X, double *_dnorm, unsigned long int *_i_max, int iFlg)
Read tridiagonal matrix components obtained by the Lanczos method. .
double complex * v0
Definition: global.h:34
struct BindStruct Bind
Binded struct.
Definition: struct.h:420
int TimeKeeper(struct BindStruct *X, const char *cFileName, const char *cTimeKeeper_Message, const char *cWriteType)
Functions for writing a time log.
Definition: log.c:42
unsigned int Lanczos_max
Maximum number of iterations.
Definition: struct.h:74
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165
+ Here is the caller graph for this function:

◆ GetSpectrumByTridiagonalMatrixComponents()

int GetSpectrumByTridiagonalMatrixComponents ( double *  tmp_alpha,
double *  tmp_beta,
double  dnorm,
double complex  dcomega,
double complex *  dcSpectrum,
unsigned long int  ilLanczosStp 
)

Calculate the spectrum by using tridiagonal matrix components obtained by the Lanczos_GetTridiagonalMatrixComponents function.

Parameters
tmp_alpha[in] Tridiagonal matrix components.
tmp_beta[in] Tridiagonal matrix components.
dnorm[in] Norm for the excited state.
dcomega[in] Target frequency.
dcSpectrum[out] Spectrum at dcomega.
ilLanczosStp[in] Lanczos step required to get tridiagonal matrix components.
Return values
FALSEFail to get the spectrum
TRUESuccess to get the spectrum

Definition at line 169 of file CalcSpectrumByLanczos.c.

References eps, eps_Energy, FALSE, stdoutMPI, and TRUE.

177 {
178  unsigned long int istp=2;
179  double complex dcDn;
180  double complex dcb0;
181  double complex dcbn, dcan;
182  double complex dcDeltahn;
183  double complex dch;
184 
185  if(ilLanczosStp < 1){
186  fprintf(stdoutMPI, "Error: LanczosStep must be greater than 1.\n");
187  return FALSE;
188  }
189 
190  dcb0 = dcomega-tmp_alpha[1];
191  if(ilLanczosStp ==1) {
192  if(cabs(dcb0)<eps_Energy){
193  dcb0=eps_Energy;
194  }
195  *dcSpectrum = dnorm / (dcb0);
196  return TRUE;
197  }
198 
199  dcbn = dcomega-tmp_alpha[2];
200  dcan = -pow(tmp_beta[1],2);
201  dcDn=1.0/dcbn;
202  dcDeltahn = dcan*dcDn;
203  dch=dcb0+dcDeltahn;
204 
205  for(istp=2; istp<=ilLanczosStp; istp++){
206  dcbn = dcomega-tmp_alpha[istp+1];
207  dcan =-pow(tmp_beta[istp],2);
208  dcDn = (dcbn+dcan*dcDn);
209  if(cabs(dcDn)<eps_Energy){
210  dcDn=eps_Energy;
211  }
212  dcDn =1.0/dcDn;
213  dcDeltahn = (dcbn*dcDn-1.0)*dcDeltahn;
214  dch += dcDeltahn;
215  if(cabs(dcDeltahn/dch)<cabs(dcb0)*eps) break;
216  }
217  *dcSpectrum = dnorm/(dch);
218  return TRUE;
219 }
#define TRUE
Definition: global.h:26
double eps_Energy
Definition: global.h:156
double eps
Definition: global.h:153
#define FALSE
Definition: global.h:25
FILE * stdoutMPI
File pointer to the standard output defined in InitializeMPI()
Definition: global.h:165