Hello,
I want to export and convert a matrix and load the matrix in an ublas::compressed_matrix<ScalarType> ublas_matrix;
I have a function to convert the original matrix from LDU to CSR format but I have no idea about ho to load the ublas_matrix?! I am an absolute beginner.
Function to convert the original matrix:
void ldu2csr
(
const Foam::lduMatrix & matrix,
scalar* vals,
uint* c_idx,
uint* r_idx
)
{
int n_rows = matrix.diag().size();
//
// Calculate each row size. Sizes are shifted, because in the next part
// array r_idx is modified twice.
//
for (int i=0; i < matrix.upper().size(); i++)
{
int ri1 = matrix.lduAddr().lowerAddr()[i];
int ri2 = matrix.lduAddr().upperAddr()[i];
r_idx[ri1+2] ++;
r_idx[ri2+2] ++;
}
//
// Compute row offsets. Offsets are shifted by one position,
// because they are used as START positions while filling values.
//
r_idx[0] = 0;
r_idx[1] = 0;
for (int i=1; i<n_rows; i++)
{
r_idx[i+1] += r_idx[i];
}
//
// Fill in CSR matrix.
// Order below is important to keep column indices sorted.
//
// lower triangle
for (int i=0; i < matrix.lower().size() ; i++)
{
int row = matrix.lduAddr().upperAddr()[i] +1;
int column = matrix.lduAddr().lowerAddr()[i];
int idx = r_idx[row];
vals[idx] = matrix.lower()[i];
c_idx[idx] = column;
r_idx[row]++;
}
// diagonal
for (int i=0; i<matrix.diag().size(); i++)
{
int idx = r_idx[i+1];
vals[idx] = matrix.diag()[i];
c_idx[idx] = i; // i is row and column index
r_idx[i+1]++;
}
// upper triangle
for (int i=0; i < matrix.upper().size() ; i++)
{
int row = matrix.lduAddr().lowerAddr()[i] +1;
int column = matrix.lduAddr().upperAddr()[i];
int idx = r_idx[row];
vals[idx] = matrix.upper()[i];
c_idx[idx] = column;
r_idx[row]++;
}
}
My idea of a function to load the ublas_matrix:
void load_umatrix
(
const Foam::lduMatrix & matrix, boost::numeric::ublas::compressed_matrix<scalar>& ublas_matrix
)
{
uint N_DIAG = matrix.diag().size();
uint lduMatrixSIZE = matrix.lower().size() + matrix.upper().size() + matrix.diag().size(); // nnz
// allocat mem for CSR sparse matrix
scalar * v = (scalar *)calloc(lduMatrixSIZE, sizeof(scalar));
uint * c = (uint *)calloc(lduMatrixSIZE, sizeof(uint));
uint * r = (uint *)calloc(N_DIAG + 2, sizeof(uint));
ldu2csr(matrix,v,c,r);
ublas_matrix.set(r,c,v,N_DIAG,N_DIAG,lduMatrixSIZE); //???
// free and release the matrix mem
free(v); free(r); free(c); //colloc()
}
I keep getting error messages related to "incomplete types" and that the matrix cannot be defined.
Any hint to fix the problem is highly appreciated.
Klaus