Commit 0597e537 authored by Davis King's avatar Davis King

polished LAPACK bindings

--HG--
extra : convert_revision : svn%3Afdd8eb12-d10e-0410-9acb-85c331704f74/trunk%403820
parent 93dce7ca
......@@ -205,13 +205,14 @@ namespace dlib
typename T,
long NR1, long NR2, long NR3, long NR4, long NR5,
long NC1, long NC2, long NC3, long NC4, long NC5,
typename MM
typename MM,
typename layout
>
int gees (
const char jobz,
matrix<T,NR1,NC1,MM,column_major_layout>& a,
matrix<T,NR2,NC2,MM,column_major_layout>& wr,
matrix<T,NR3,NC3,MM,column_major_layout>& wi,
matrix<T,NR2,NC2,MM,layout>& wr,
matrix<T,NR3,NC3,MM,layout>& wi,
matrix<T,NR4,NC4,MM,column_major_layout>& vs,
matrix<T,NR5,NC5,MM,column_major_layout>& work
)
......
......@@ -166,14 +166,15 @@ namespace dlib
typename T,
long NR1, long NR2, long NR3, long NR4, long NR5, long NR6,
long NC1, long NC2, long NC3, long NC4, long NC5, long NC6,
typename MM
typename MM,
typename layout
>
int geev (
const char jobvl,
const char jobvr,
matrix<T,NR1,NC1,MM,column_major_layout>& a,
matrix<T,NR2,NC2,MM,column_major_layout>& wr,
matrix<T,NR3,NC3,MM,column_major_layout>& wi,
matrix<T,NR2,NC2,MM,layout>& wr,
matrix<T,NR3,NC3,MM,layout>& wi,
matrix<T,NR4,NC4,MM,column_major_layout>& vl,
matrix<T,NR5,NC5,MM,column_major_layout>& vr,
matrix<T,NR6,NC6,MM,column_major_layout>& work
......
......@@ -259,13 +259,16 @@ namespace dlib
const char jobz,
matrix<T,NR1,NC1,MM,row_major_layout>& a,
matrix<T,NR2,NC2,MM,row_major_layout>& s,
matrix<T,NR3,NC3,MM,row_major_layout>& vt,
matrix<T,NR4,NC4,MM,row_major_layout>& u,
matrix<T,NR3,NC3,MM,row_major_layout>& u_,
matrix<T,NR4,NC4,MM,row_major_layout>& vt_,
matrix<T,NR5,NC5,MM,row_major_layout>& work,
matrix<integer,NR6,NC6,MM,row_major_layout>& iwork
)
{
// Row major order matrices are transposed from LAPACK's point of view.
matrix<T,NR3,NC3,MM,row_major_layout>& u = vt_;
matrix<T,NR4,NC4,MM,row_major_layout>& vt = u_;
const long m = a.nc();
const long n = a.nr();
......
......@@ -246,12 +246,14 @@ namespace dlib
const char jobvt,
matrix<T,NR1,NC1,MM,row_major_layout>& a,
matrix<T,NR2,NC2,MM,row_major_layout>& s,
matrix<T,NR3,NC3,MM,row_major_layout>& vt,
matrix<T,NR4,NC4,MM,row_major_layout>& u,
matrix<T,NR3,NC3,MM,row_major_layout>& u_,
matrix<T,NR4,NC4,MM,row_major_layout>& vt_,
matrix<T,NR5,NC5,MM,row_major_layout>& work
)
{
// Row major order matrices are transposed from LAPACK's point of view.
matrix<T,NR3,NC3,MM,row_major_layout>& u = vt_;
matrix<T,NR4,NC4,MM,row_major_layout>& vt = u_;
const long m = a.nc();
const long n = a.nr();
......
......@@ -152,6 +152,51 @@ namespace dlib
return info;
}
// ------------------------------------------------------------------------------------
template <
typename T,
long NR1, long NR2, long NR3,
long NC1, long NC2, long NC3,
typename MM
>
int syev (
char jobz,
char uplo,
matrix<T,NR1,NC1,MM,row_major_layout>& a,
matrix<T,NR2,NC2,MM,row_major_layout>& w,
matrix<T,NR3,NC3,MM,row_major_layout>& work
)
{
if (uplo == 'L')
uplo = 'U';
else
uplo = 'L';
const long n = a.nr();
w.set_size(n,1);
// figure out how big the workspace needs to be.
T work_size = 1;
int info = binding::syev(jobz, uplo, n, &a(0,0),
a.nc(), &w(0,0), &work_size, -1);
if (info != 0)
return info;
if (work.size() < work_size)
work.set_size(work_size, 1);
// compute the actual decomposition
info = binding::syev(jobz, uplo, n, &a(0,0),
a.nc(), &w(0,0), &work(0,0), work.size());
a = trans(a);
}
// ------------------------------------------------------------------------------------
}
......
......@@ -350,6 +350,80 @@ namespace dlib
&w(0,0), &z(0,0), z.nr(), &isuppz(0,0), &work(0,0), work.size(),
&iwork(0,0), iwork.size());
return info;
}
// ------------------------------------------------------------------------------------
template <
typename T,
long NR1, long NR2, long NR3, long NR4, long NR5, long NR6,
long NC1, long NC2, long NC3, long NC4, long NC5, long NC6,
typename MM
>
int syevr (
const char jobz,
const char range,
char uplo,
matrix<T,NR1,NC1,MM,row_major_layout>& a,
const double vl,
const double vu,
const integer il,
const integer iu,
const double abstol,
integer& num_eigenvalues_found,
matrix<T,NR2,NC2,MM,row_major_layout>& w,
matrix<T,NR3,NC3,MM,row_major_layout>& z,
matrix<integer,NR4,NC4,MM,row_major_layout>& isuppz,
matrix<T,NR5,NC5,MM,row_major_layout>& work,
matrix<integer,NR6,NC6,MM,row_major_layout>& iwork
)
{
if (uplo == 'L')
uplo = 'U';
else
uplo = 'L';
const long n = a.nr();
w.set_size(n,1);
isuppz.set_size(2*n, 1);
if (jobz == 'V')
{
z.set_size(n,n);
}
else
{
z.set_size(1,1);
}
// figure out how big the workspace needs to be.
T work_size = 1;
integer iwork_size = 1;
int info = binding::syevr(jobz, range, uplo, n, &a(0,0),
a.nc(), vl, vu, il, iu, abstol, &num_eigenvalues_found,
&w(0,0), &z(0,0), z.nc(), &isuppz(0,0), &work_size, -1,
&iwork_size, -1);
if (info != 0)
return info;
if (work.size() < work_size)
work.set_size(work_size, 1);
if (iwork.size() < iwork_size)
iwork.set_size(iwork_size, 1);
// compute the actual decomposition
info = binding::syevr(jobz, range, uplo, n, &a(0,0),
a.nc(), vl, vu, il, iu, abstol, &num_eigenvalues_found,
&w(0,0), &z(0,0), z.nc(), &isuppz(0,0), &work(0,0), work.size(),
&iwork(0,0), iwork.size());
z = trans(z);
return info;
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment