@@ -61,8 +61,9 @@ std::pair<int, std::vector<int>> DiagoBlas::pdsygvx_once(const int *const desc,
6161 const char jobz = ' V' , range = ' I' , uplo = ' U' ;
6262 const int itype = 1 , il = 1 , iu = GlobalV::NBANDS, one = 1 ;
6363 int M = 0 , NZ = 0 , lwork = -1 , liwork = -1 , info = 0 ;
64+ double vl = 0 , vu = 0 ;
6465 const double abstol = 0 , orfac = -1 ;
65- std::vector<double > work (1 , 0 );
66+ std::vector<double > work (3 , 0 );
6667 std::vector<int > iwork (1 , 0 );
6768 std::vector<int > ifail (GlobalV::NLOCAL, 0 );
6869 std::vector<int > iclustr (2 * GlobalV::DSIZE);
@@ -81,8 +82,8 @@ std::pair<int, std::vector<int>> DiagoBlas::pdsygvx_once(const int *const desc,
8182 &one,
8283 &one,
8384 desc,
84- NULL ,
85- NULL ,
85+ &vl ,
86+ &vu ,
8687 &il,
8788 &iu,
8889 &abstol,
@@ -109,7 +110,7 @@ std::pair<int, std::vector<int>> DiagoBlas::pdsygvx_once(const int *const desc,
109110
110111 // GlobalV::ofs_running<<"lwork="<<work[0]<<"\t"<<"liwork="<<iwork[0]<<std::endl;
111112 lwork = work[0 ];
112- work.resize (lwork, 0 );
113+ work.resize (std::max ( lwork, 3 ) , 0 );
113114 liwork = iwork[0 ];
114115 iwork.resize (liwork, 0 );
115116
@@ -126,8 +127,8 @@ std::pair<int, std::vector<int>> DiagoBlas::pdsygvx_once(const int *const desc,
126127 &one,
127128 &one,
128129 desc,
129- NULL ,
130- NULL ,
130+ &vl ,
131+ &vu ,
131132 &il,
132133 &iu,
133134 &abstol,
@@ -184,8 +185,12 @@ std::pair<int, std::vector<int>> DiagoBlas::pzhegvx_once(const int *const desc,
184185 const int itype = 1 , il = 1 , iu = GlobalV::NBANDS, one = 1 ;
185186 int M = 0 , NZ = 0 , lwork = -1 , lrwork = -1 , liwork = -1 , info = 0 ;
186187 const double abstol = 0 , orfac = -1 ;
188+ // Note: pzhegvx_ has a bug
189+ // We must give vl,vu a value, although we do not use range 'V'
190+ // We must give rwork at least a memory of sizeof(double) * 3
191+ const double vl = 0 , vu = 0 ;
187192 std::vector<std::complex <double >> work (1 , 0 );
188- std::vector<double > rwork (1 , 0 );
193+ std::vector<double > rwork (3 , 0 );
189194 std::vector<int > iwork (1 , 0 );
190195 std::vector<int > ifail (GlobalV::NLOCAL, 0 );
191196 std::vector<int > iclustr (2 * GlobalV::DSIZE);
@@ -204,8 +209,8 @@ std::pair<int, std::vector<int>> DiagoBlas::pzhegvx_once(const int *const desc,
204209 &one,
205210 &one,
206211 desc,
207- NULL ,
208- NULL ,
212+ &vl ,
213+ &vu ,
209214 &il,
210215 &iu,
211216 &abstol,
@@ -236,7 +241,8 @@ std::pair<int, std::vector<int>> DiagoBlas::pzhegvx_once(const int *const desc,
236241 lwork = work[0 ].real ();
237242 work.resize (lwork, 0 );
238243 lrwork = rwork[0 ] + this ->degeneracy_max * GlobalV::NLOCAL;
239- rwork.resize (lrwork, 0 );
244+ int maxlrwork = std::max (lrwork,3 );
245+ rwork.resize (maxlrwork, 0 );
240246 liwork = iwork[0 ];
241247 iwork.resize (liwork, 0 );
242248
@@ -253,8 +259,8 @@ std::pair<int, std::vector<int>> DiagoBlas::pzhegvx_once(const int *const desc,
253259 &one,
254260 &one,
255261 desc,
256- NULL ,
257- NULL ,
262+ &vl ,
263+ &vu ,
258264 &il,
259265 &iu,
260266 &abstol,
0 commit comments