本文整理汇总了C++中SiconosMatrix类的典型用法代码示例。如果您正苦于以下问题:C++ SiconosMatrix类的具体用法?C++ SiconosMatrix怎么用?C++ SiconosMatrix使用的例子?那么, 这里精选的类代码示例或许可以为您提供帮助。
在下文中一共展示了SiconosMatrix类的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。
示例1: SimpleMatrix
// construct by copy of SiconosMatrix
OSNSMatrix::OSNSMatrix(const SiconosMatrix& MSource):
_dimRow(MSource.size(0)), _dimColumn(MSource.size(1)), _storageType(NM_DENSE)
{
_numericsMat.reset(new NumericsMatrix);
NM_null(_numericsMat.get());
_M1.reset(new SimpleMatrix(MSource));
}
示例2: setJacobianRhsx
void DynamicalSystem::setJacobianRhsx(const SiconosMatrix& newValue)
{
// check dimensions ...
if (newValue.size(0) != _n || newValue.size(1) != _n)
RuntimeException::selfThrow("DynamicalSystem::setJacobianRhsx - inconsistent sizes between jacobianRhsx input and n - Maybe you forget to set n?");
if (_jacxRhs)
*_jacxRhs = newValue;
else
_jacxRhs.reset(new SimpleMatrix(newValue));
}
示例3: SolveByLeastSquares
void SimpleMatrix::SolveByLeastSquares(SiconosMatrix &B)
{
if (B.isBlock())
SiconosMatrixException::selfThrow("SimpleMatrix::SolveByLeastSquares(Siconos Matrix &B) failed. Not yet implemented for M being a BlockMatrix.");
int info = 0;
#ifdef USE_OPTIMAL_WORKSPACE
info += lapack::gels(*mat.Dense, *(B.dense()), lapack::optimal_workspace());
#endif
#ifdef USE_MINIMAL_WORKSPACE
info += lapack::gels(*mat.Dense, *(B.dense()), lapack::minimal_workspace());
#endif
if (info != 0)
SiconosMatrixException::selfThrow("SimpleMatrix::SolveByLeastSquares failed.");
}
示例4: private_addprod
void private_addprod(const SiconosMatrix& A, unsigned int startRow, unsigned int startCol, const BlockVector& x, SiconosVector& y)
{
assert(!(A.isPLUFactorized()) && "A is PLUFactorized in prod !!");
assert(!A.isBlock() && "private_addprod(A,start,x,y) error: not yet implemented for block matrix.");
VectorOfVectors::const_iterator it;
unsigned int startColBis = startCol;
for (it = x.begin(); it != x.end(); ++it)
{
private_addprod(A, startRow, startColBis, **it, y);
startColBis += (*it)->size();
}
}
示例5: setW
void EulerMoreauOSI::setW(const SiconosMatrix& newValue, SP::DynamicalSystem ds)
{
// Check if ds is in the OSI
if (!OSIDynamicalSystems->isIn(ds))
RuntimeException::selfThrow("EulerMoreauOSI::setW(newVal,ds) - ds does not belong to this Integrator ...");
// Check dimensions consistency
unsigned int line = newValue.size(0);
unsigned int col = newValue.size(1);
if (line != col) // Check that newValue is square
RuntimeException::selfThrow("EulerMoreauOSI::setW(newVal,ds) - newVal is not square! ");
if (!ds)
RuntimeException::selfThrow("EulerMoreauOSI::setW(newVal,ds) - ds == NULL.");
unsigned int sizeW = ds->getDim(); // n for first order systems, ndof for lagrangian.
unsigned int dsN = ds->number();
if (line != sizeW) // check consistency between newValue and dynamical system size
RuntimeException::selfThrow("EulerMoreauOSI::setW(newVal,ds) - unconsistent dimension between newVal and dynamical system to be integrated ");
// Memory allocation for W, if required
if (!WMap[dsN]) // allocate a new W if required
{
WMap[dsN].reset(new SimpleMatrix(newValue));
}
else // or fill-in an existing one if dimensions are consistent.
{
if (line == WMap[dsN]->size(0) && col == WMap[dsN]->size(1))
*(WMap[dsN]) = newValue;
else
RuntimeException::selfThrow("EulerMoreauOSI - setW: inconsistent dimensions with problem size for given input matrix W");
}
}
示例6: if
SimpleMatrix operator * (double a, const SiconosMatrix & A)
{
// To compute B = a * A
unsigned int numA = A.getNum();
if (numA == 6) // if A = 0
{
//DenseMat p(zero_matrix(A.size(0),A.size(1)));
//return p;
return A;
}
else if (numA == 7)
{
return (DenseMat)(a**A.identity());
}
else if (numA == 0) // A block
{
SimpleMatrix tmp(A); // ... copy ...
tmp *= a;
return tmp;
}
else if (numA == 1) // dense)
return (DenseMat)(a** A.dense());
else if (numA == 2)
return (TriangMat)(a ** A.triang());
else if (numA == 3)
return (SymMat)(a ** A.sym());
else if (numA == 4)
return (SparseMat)(a ** A.sparse());
else //if(numA==5)
return (BandedMat)(a ** A.banded());
}
示例7: trans
void SimpleMatrix::trans(const SiconosMatrix &m)
{
if (m.isBlock())
SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, not yet implemented for m being a BlockMatrix.");
if (&m == this)
trans();//SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, m = this, use this->trans().");
else
{
unsigned int numM = m.num();
switch (numM)
{
case 1:
if (_num != 1)
SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, try to transpose a dense matrix into another type.");
noalias(*mat.Dense) = ublas::trans(*m.dense());
break;
case 2:
if (_num != 1)
SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, try to transpose a triangular matrix into a non-dense one.");
noalias(*mat.Dense) = ublas::trans(*m.triang());
break;
case 3:
*this = m;
break;
case 4:
if (_num == 1)
noalias(*mat.Dense) = ublas::trans(*m.sparse());
else if (_num == 4)
noalias(*mat.Sparse) = ublas::trans(*m.sparse());
else
SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, try to transpose a sparse matrix into a forbidden type (not dense nor sparse).");
break;
case 5:
if (_num == 1)
noalias(*mat.Dense) = ublas::trans(*m.banded());
else if (_num == 5)
noalias(*mat.Banded) = ublas::trans(*m.banded());
else
SiconosMatrixException::selfThrow("SimpleMatrix::trans(m) failed, try to transpose a banded matrix into a forbidden type (not dense nor banded).");
break;
case 6:
*this = m;
break;
case 7:
*this = m;
}
// unsigned int tmp = _dimRow;
// _dimRow = _dimCol;
// _dimCol = tmp;
resetLU();
}
}
示例8: private_prod
// x block, y siconos
void private_prod(const SiconosMatrix& A, unsigned int startRow, const SiconosVector& x, SiconosVector& y, bool init)
{
assert(!(A.isPLUFactorized()) && "A is PLUFactorized in prod !!");
// Computes y = subA *x (or += if init = false), subA being a sub-matrix of A, between el. of index (row) startRow and startRow + sizeY
if (init) // y = subA * x , else y += subA * x
y.zero();
private_addprod(A, startRow, 0, x, y);
}
示例9: subprod
void subprod(const SiconosMatrix& A, const BlockVector& x, SiconosVector& y, const Index& coord, bool init)
{
assert(!(A.isPLUFactorized()) && "A is PLUFactorized in prod !!");
// Number of the subvector of x that handles element at position coord[4]
std::size_t firstBlockNum = x.getNumVectorAtPos(coord[4]);
// Number of the subvector of x that handles element at position coord[5]
unsigned int lastBlockNum = x.getNumVectorAtPos(coord[5]);
Index subCoord = coord;
SPC::SiconosVector tmp = x[firstBlockNum];
std::size_t subSize = tmp->size(); // Size of the sub-vector
const SP::Index xTab = x.tabIndex();
if (firstBlockNum != 0)
{
subCoord[4] -= (*xTab)[firstBlockNum - 1];
subCoord[5] = std::min(coord[5] - (*xTab)[firstBlockNum - 1], subSize);
}
else
subCoord[5] = std::min(coord[5], subSize);
if (firstBlockNum == lastBlockNum)
{
subprod(A, *tmp, y, subCoord, init);
}
else
{
unsigned int xPos = 0 ; // Position in x of the current sub-vector of x
bool firstLoop = true;
subCoord[3] = coord[2] + subCoord[5] - subCoord[4];
for (VectorOfVectors::const_iterator it = x.begin(); it != x.end(); ++it)
{
if ((*it)->getNum() == 0)
SiconosMatrixException::selfThrow("subprod(A,x,y) error: not yet implemented for x block of blocks ...");
if (xPos >= firstBlockNum && xPos <= lastBlockNum)
{
tmp = x[xPos];
if (firstLoop)
{
subprod(A, *tmp, y, subCoord, init);
firstLoop = false;
}
else
{
subCoord[2] += subCoord[5] - subCoord[4]; // !! old values for 4 and 5
subSize = tmp->size();
subCoord[4] = 0;
subCoord[5] = std::min(coord[5] - (*xTab)[xPos - 1], subSize);
subCoord[3] = subCoord[2] + subCoord[5] - subCoord[4];
subprod(A, *tmp, y, subCoord, false);
}
}
xPos++;
}
}
}
示例10: isComparableTo
//=====================
// matrices comparison
//=====================
bool isComparableTo(const SiconosMatrix& m1, const SiconosMatrix& m2)
{
// return:
// - true if one of the matrices is a Simple and if they have the same dimensions.
// - true if both are block but with blocks which are facing each other of the same size.
// - false in other cases
if ((!m1.isBlock() || !m2.isBlock()) && (m1.size(0) == m2.size(0)) && (m1.size(1) == m2.size(1)))
return true;
const SP::Index I1R = m1.tabRow();
const SP::Index I2R = m2.tabRow();
const SP::Index I1C = m1.tabCol();
const SP::Index I2C = m2.tabCol();
return ((*I1R == *I2R) && (*I1C == *I2C));
}
示例11: expm
void expm(SiconosMatrix& A, SiconosMatrix& Exp, bool computeAndAdd)
{
// Implemented only for dense matrices.
// Note FP : Maybe it works for others but it has not been
// tested here --> to be done
// Do not work with sparse.
A.resetLU();
Exp.resetLU();
assert(Exp.num() == 1 || A.num() == 1);
if(computeAndAdd)
*Exp.dense() += expm_pad(*A.dense());
else
*Exp.dense() = expm_pad(*A.dense());
}
示例12: prod
void prod(const SiconosVector& x, const SiconosMatrix& A, BlockVector& y, bool init)
{
assert(!(A.isPLUFactorized()) && "A is PLUFactorized in prod !!");
unsigned int startRow = 0;
VectorOfVectors::const_iterator it;
// For Each subvector of y, y[i], private_prod computes y[i] = subA x, subA being a submatrix of A corresponding to y[i] position.
// private_prod takes into account the fact that x and y[i] may be block vectors.
for (it = y.begin(); it != y.end(); ++it)
{
private_prod(createSPtrConstSiconosVector(x), createSPtrConstSiconosMatrix(A), startRow, *it, init);
startRow += (*it)->size();
}
}
示例13: addBlock
void SimpleMatrix::addBlock(unsigned int row_min, unsigned int col_min, const SiconosMatrix& m)
{
// add m to current matrix elements, starting from row row_min and column col_min, to the values of the matrix m.
// m may be a BlockMatrix.
if (_num == 6 || _num == 7)
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(pos,..., m) forbidden for zero or identity matrix.");
if (&m == this)
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(pos,..., m): m = this.");
if (row_min >= size(0))
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(row,col): row is out of range");
if (col_min >= size(1))
SiconosMatrixException::selfThrow("SimpleMatrix::addBloc(row,col)k: col is out of range");
unsigned int row_max, col_max;
row_max = m.size(0) + row_min;
col_max = m.size(1) + col_min;
if (row_max > size(0))
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(row,col,m): m.row + row is out of range.");
if (col_max > size(1))
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(row,col,m): m.col + col is out of range.");
unsigned int numM = m.num();
if (numM == 0) // if m is a block matrix ...
{
const BlockMatrix& mB = static_cast<const BlockMatrix&>(m);
BlocksMat::const_iterator1 it;
BlocksMat::const_iterator2 it2;
unsigned int posRow = row_min;
unsigned int posCol = col_min;
for (it = mB._mat->begin1(); it != mB._mat->end1(); ++it)
{
for (it2 = it.begin(); it2 != it.end(); ++it2)
{
addBlock(posRow, posCol, **it2);
posCol += (*it2)->size(1);
}
posRow += (*it)->size(0);
posCol = 0;
}
}
else if (numM == 6) // if m = 0
{
// nothing to do !
}
else // if m is a SimpleMatrix
{
if (_num == 1)
{
switch (numM)
{
case 1:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.dense());
break;
case 2:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.triang());
break;
case 3:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.sym());
break;
case 4:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.sparse());
break;
case 5:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.banded());
break;
case 7:
noalias(ublas::subrange(*mat.Dense, row_min, row_max, col_min, col_max)) += *(m.identity());
break;
default:
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(...,m): wrong matrix type for m.");
break;
}
}
else
SiconosMatrixException::selfThrow("SimpleMatrix::addBlock(...): implemeted only for dense matrices.");
resetLU();
}
}
示例14: write
bool write(const std::string& fileName, const std::string& mode, const SiconosMatrix& m, const std::string& outputType)
{
// Open file and various checks
std::ofstream outfile;
if (mode == "ascii")
outfile.open(fileName.c_str(), std::ofstream::out);
else if (mode == "binary")
outfile.open(fileName.c_str(), std::ofstream::binary);
else
SiconosMatrixException::selfThrow("ioMatrix::write Incorrect mode for writing");
if (!outfile.good())
SiconosMatrixException::selfThrow("ioMatrix:: write error : Fail to open \"" + fileName + "\"");
if (m.isBlock())
SiconosMatrixException::selfThrow("ioMatrix:: write error : not yet implemented for BlockMatrix");
outfile.precision(15);
outfile.setf(std::ios::scientific);
// Writing
if (outputType != "noDim")
outfile << m.size(0) << " " << m.size(1) << std::endl;
if (m.getNum() == 1)
{
// DenseMat * p = m.dense();
DenseMat::iterator1 row;
DenseMat::iterator2 col;
double tmp;
for (unsigned int i = 0; i < m.size(0); i++)
{
for (unsigned int j = 0; j < m.size(1); j++)
{
tmp = m(i, j);
if (fabs(tmp) < std::numeric_limits<double>::min()) tmp = 0.0;
outfile << tmp << " " ;
assert(outfile.good());
}
outfile << std::endl;
}
}
else if (m.getNum() == 2)
{
TriangMat * p = m.triang();
TriangMat::iterator1 row;
for (row = p->begin1(); row != p->end1() ; ++row)
{
std::copy(row.begin(), row.end(), std::ostream_iterator<double>(outfile, " "));
outfile << std::endl;
}
}
else if (m.getNum() == 3)
{
SymMat * p = m.sym();
SymMat::iterator1 row;
for (row = p->begin1(); row != p->end1() ; ++row)
{
std::copy(row.begin(), row.end(), std::ostream_iterator<double>(outfile, " "));
outfile << std::endl;
}
}
else if (m.getNum() == 4)
{
SparseMat * p = m.sparse();
SparseMat::iterator1 row;
for (row = p->begin1(); row != p->end1() ; ++row)
{
std::copy(row.begin(), row.end(), std::ostream_iterator<double>(outfile, " "));
outfile << std::endl;
}
}
else
{
BandedMat * p = m.banded();
BandedMat::iterator1 row;
for (row = p->begin1(); row != p->end1() ; ++row)
{
std::copy(row.begin(), row.end(), std::ostream_iterator<double>(outfile, " "));
outfile << std::endl;
}
}
outfile.close();
return true;
}
示例15: SimpleMatrix
const SimpleMatrix operator + (const SiconosMatrix& A, const SiconosMatrix& B)
{
// To compute C = A + B
if ((A.size(0) != B.size(0)) || (A.size(1) != B.size(1)))
SiconosMatrixException::selfThrow("Matrix operator + : inconsistent sizes");
unsigned int numA = A.getNum();
unsigned int numB = B.getNum();
// == A or B equal to null ==
if (numA == 6) // A = 0
{
if (numB == 6) // B = 0
return SimpleMatrix(A.size(0), A.size(1));
else
return SimpleMatrix(B);
}
if (numB == 6)
return SimpleMatrix(A);
// == A and B different from 0 ==
if (numA == numB && numA != 0) // all matrices are of the same type and NOT block
{
if (numA == 1)
return (DenseMat)(*A.dense() + *B.dense());
else if (numA == 2)
return (TriangMat)(*A.triang() + *B.triang());
else if (numA == 3)
return (SymMat)(*A.sym() + *B.sym());
else if (numA == 4)
{
SparseMat tmp(*A.sparse());
tmp += *B.sparse();
return tmp;
// return (SparseMat)(*A.sparse() + *B.sparse());
}
else //if(numA==5)
{
BandedMat tmp(*A.banded());
tmp += *B.banded();
return tmp;
}
}
else if (numA != 0 && numB != 0 && numA != numB) // A and B of different types and none is block
{
if (numA == 1)
{
if (numB == 2)
return (DenseMat)(*A.dense() + *B.triang());
else if (numB == 3)
return (DenseMat)(*A.dense() + *B.sym());
else if (numB == 4)
return (DenseMat)(*A.dense() + *B.sparse());
else if (numB == 5)
return (DenseMat)(*A.dense() + *B.banded());
else // if(numB ==7)
return (DenseMat)(*A.dense() + *B.identity());
}
else if (numA == 2)
{
if (numB == 1)
return (DenseMat)(*A.triang() + *B.dense());
else if (numB == 3)
return (DenseMat)(*A.triang() + *B.sym());
else if (numB == 4)
return (DenseMat)(*A.triang() + *B.sparse());
else if (numB == 5)
return (DenseMat)(*A.triang() + *B.banded());
else // if(numB ==7:
return (DenseMat)(*A.triang() + *B.identity());
}
else if (numA == 3)
{
if (numB == 1)
return (DenseMat)(*A.sym() + *B.dense());
else if (numB == 2)
return (DenseMat)(*A.sym() + *B.triang());
else if (numB == 4)
return (DenseMat)(*A.sym() + *B.sparse());
else if (numB == 5)
return (DenseMat)(*A.sym() + *B.banded());
else // if(numB ==7)
return (DenseMat)(*A.sym() + *B.identity());
}
else if (numA == 4)
{
if (numB == 1)
return (DenseMat)(*A.sparse() + *B.dense());
else if (numB == 2)
return (DenseMat)(*A.sparse() + *B.triang());
else if (numB == 3)
return (DenseMat)(*A.sparse() + *B.sym());
else if (numB == 5)
return (DenseMat)(*A.sparse() + *B.banded());
else // if(numB ==7)
return (DenseMat)(*A.sparse() + *B.identity());
}
//.........这里部分代码省略.........