31 SymMatrix(
const char* fname):
LinOp(0,0,SYMMETRIC,2),value() { this->load(fname); }
39 size_t size()
const {
return nlin()*(nlin()+1)/2; };
50 double*
data()
const {
return value.get(); }
55 return data()[(i<=j) ? i+j*(j+1)/2 : j+i*(i+1)/2];
61 return data()[(i<=j) ? i+j*(j+1)/2 : j+i*(i+1)/2];
70 void solveLin(
Vector* B,
const int nbvect);
85 void operator *=(
const double x);
86 void operator /=(
const double x) { (*this)*=(1/x); }
94 void save(
const char* filename)
const;
95 void load(
const char* filename);
97 void save(
const std::string& s)
const {
save(s.c_str()); }
98 void load(
const std::string& s) {
load(s.c_str()); }
111 BLAS_INT* pivots=
new BLAS_INT[
nlin()];
120 std::cout <<
"solveLin not defined" << std::endl;
132 BLAS_INT *pivots=
new BLAS_INT[
nlin()];
137 for(
int i = 0; i < nbvect; i++)
143 std::cout <<
"solveLin not defined" << std::endl;
152 const size_t sz =
size();
153 for (
size_t i=0; i<sz; ++i)
163 const size_t sz =
size();
164 for (
size_t i=0; i<sz; ++i)
179 std::cerr <<
"Positive definite inverse not defined" << std::endl;
189 BLAS_INT *pivots=
new BLAS_INT[
nlin()];
194 std::cout <<
"Big problem in det (DSPTRF)" << std::endl;
195 for (
size_t i = 0; i<
nlin(); i++){
196 if (pivots[i] >= 0) {
199 if (i <
nlin()-1 && pivots[i] == pivots[i+1]) {
200 d *= (invA(i,i)*invA(i+1,i+1)-invA(i,i+1)*invA(i+1,i));
203 std::cout <<
"Big problem in det" << std::endl;
209 throw OpenMEEG::maths::LinearAlgebraError(
"Determinant not defined without LAPACK");
258 BLAS_INT* pivots =
new BLAS_INT[
nlin()];
261 DSPTRF(
'U',M,invA.
data(),pivots,Info);
263 double* work =
new double[
nlin()*64];
264 DSPTRI(
'U',M,invA.
data(),pivots,work,Info);
271 throw OpenMEEG::maths::LinearAlgebraError(
"Inverse not implemented, requires LAPACK");
278 BLAS_INT* pivots =
new BLAS_INT[
nlin()];
281 DSPTRF(
'U',M,
data(),pivots,Info);
283 double* work =
new double[
nlin()*64];
284 DSPTRI(
'U',M,
data(),pivots,work,Info);
291 throw OpenMEEG::maths::LinearAlgebraError(
"Inverse not implemented, requires LAPACK");
300 DSPMV(CblasUpper,M,1.0,
data(),v.
data(),1,0.0,y.
data(),1);
305 y(i)+=(*this)(i,j)*v(j);
Matrix class Matrix class.
void save(const std::string &s) const
Matrix submat(const Index istart, const Index isize, const Index jstart, const Index jsize) const
SymMatrix posdefinverse() const
Matrix operator*(const Matrix &B) const
SymMatrix operator/(const double x) const
SymMatrix(const char *fname)
Vector solveLin(const Vector &B) const
Matrix solveLin(Matrix &B) const
const SymMatrix & operator=(const double d)
Matrix operator*(const SymMatrix &B) const
SymMatrix operator*(const double x) const
void operator+=(const SymMatrix &B)
SymMatrix submat(const Index istart, const Index iend) const
SymMatrix operator+(const SymMatrix &B) const
Matrix operator()(const Index i_start, const Index i_end, const Index j_start, const Index j_end) const
void setlin(const Index i, const Vector &v)
double & operator()(const Index i, const Index j)
Vector getlin(const Index i) const
void save(const char *filename) const
SymMatrix inverse() const
void reference_data(const double *array)
SymMatrix(const SymMatrix &S, const DeepCopy)
double operator()(const Index i, const Index j) const
void load(const char *filename)
SymMatrix(const Matrix &A)
SymMatrix(const Vector &v)
void load(const std::string &s)
void operator-=(const SymMatrix &B)
SymMatrix(Dimension M, Dimension N)
SymMatrix operator-(const SymMatrix &B) const
Vect3 operator*(const double d, const Vect3 &V)
double det(const Vect3 &V1, const Vect3 &V2, const Vect3 &V3)
BLAS_INT sizet_to_int(const unsigned &num)