From: janssen Date: Mon, 22 Mar 2010 14:41:39 +0000 (+0000) Subject: modified to be able to treat non quadratic matrices with distribute_local_to_global X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=92555e70e84cd08161d6997ac1f31c66b385f563;p=dealii-svn.git modified to be able to treat non quadratic matrices with distribute_local_to_global git-svn-id: https://svn.dealii.org/trunk@20875 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/constraint_matrix.h b/deal.II/lac/include/lac/constraint_matrix.h index 79ffcfe2b3..d66b846ce5 100644 --- a/deal.II/lac/include/lac/constraint_matrix.h +++ b/deal.II/lac/include/lac/constraint_matrix.h @@ -1305,6 +1305,18 @@ class ConstraintMatrix : public Subscriptor void distribute_local_to_global (const FullMatrix &local_matrix, const std::vector &local_dof_indices, + MatrixType &global_matrix) const; + + /** + * Does the same as the function above + * but can treat + * non quadratic matrices. + */ + template + void + distribute_local_to_global (const FullMatrix &local_matrix, + const std::vector &row_indices, + const std::vector &col_indices, MatrixType &global_matrix) const; /** diff --git a/deal.II/lac/include/lac/constraint_matrix.templates.h b/deal.II/lac/include/lac/constraint_matrix.templates.h index 62c5937314..ff84b0660d 100644 --- a/deal.II/lac/include/lac/constraint_matrix.templates.h +++ b/deal.II/lac/include/lac/constraint_matrix.templates.h @@ -1136,8 +1136,11 @@ namespace internals { GlobalRowsFromLocal (const unsigned int n_local_rows) : - total_row_indices (n_local_rows), n_active_rows (n_local_rows), - n_inhomogeneous_rows (0){}; + total_row_indices (n_local_rows), + n_active_rows (n_local_rows), + n_inhomogeneous_rows (0) + {} + // implemented below void insert_index (const unsigned int global_row, @@ -1322,6 +1325,7 @@ namespace internals // constraints are really empty. const unsigned int length = size(); #ifdef DEBUG + AssertIndexRange (length, total_row_indices.size()+1); for (unsigned int i=0; i &local_matrix) { - const unsigned int loc_col = global_rows.local_row(j); + const unsigned int loc_col = global_cols.local_row(j); double col_val; // case 1: row has direct contribution in @@ -1469,9 +1474,9 @@ namespace internals // account for indirect contributions by // constraints in column - for (unsigned int p=0; p= column_end, - ExcIndexRange (column_end, 0, global_rows.size())); + Assert (global_cols.size() >= column_end, + ExcIndexRange (column_end, 0, global_cols.size())); const unsigned int loc_row = global_rows.local_row(i); // fast function if there are no indirect @@ -1527,17 +1533,18 @@ namespace internals // element is zero. if (global_rows.have_indirect_rows() == false) { - Assert(loc_row < local_matrix.m(), ExcInternalError()); + AssertIndexRange(loc_row, local_matrix.m()); const double * matrix_ptr = &local_matrix(loc_row, 0); for (unsigned int j=column_start; j (col_val); - *col_ptr++ = global_rows.global_row(j); + *col_ptr++ = global_cols.global_row(j); } } } @@ -1549,7 +1556,7 @@ namespace internals { for (unsigned int j=column_start; j (col_val); - *col_ptr++ = global_rows.global_row(j); + *col_ptr++ = global_cols.global_row(j); } } } @@ -1666,7 +1673,7 @@ namespace internals { for (unsigned int j=column_start; j &local_dof_indices, internals::GlobalRowsFromLocal &global_rows) const { const unsigned int n_local_dofs = local_dof_indices.size(); + AssertDimension (n_local_dofs, global_rows.size()); // when distributing the local data to // the global matrix, we can quite // cheaply sort the indices (obviously, @@ -2293,7 +2301,7 @@ ConstraintMatrix::distribute_local_to_global ( { unsigned int * col_ptr = &cols[0]; number * val_ptr = &vals[0]; - resolve_matrix_row (global_rows, i, 0, n_actual_dofs, + resolve_matrix_row (global_rows, global_rows, i, 0, n_actual_dofs, local_matrix, col_ptr, val_ptr); const unsigned int n_values = col_ptr - &cols[0]; Assert (n_values == (unsigned int)(val_ptr - &vals[0]), @@ -2328,6 +2336,57 @@ ConstraintMatrix::distribute_local_to_global ( +template +void +ConstraintMatrix::distribute_local_to_global ( + const FullMatrix &local_matrix, + const std::vector &row_indices, + const std::vector &col_indices, + MatrixType &global_matrix) const +{ + typedef double number; + + AssertDimension (local_matrix.m(), row_indices.size()); + AssertDimension (local_matrix.n(), col_indices.size()); + //Assert (sorted == true, ExcMatrixNotClosed()); + + const unsigned int n_local_row_dofs = row_indices.size(); + const unsigned int n_local_col_dofs = col_indices.size(); + internals::GlobalRowsFromLocal global_rows (n_local_row_dofs); + internals::GlobalRowsFromLocal global_cols (n_local_col_dofs); + make_sorted_row_list (row_indices, global_rows); + make_sorted_row_list (col_indices, global_cols); + + const unsigned int n_actual_row_dofs = global_rows.size(); + const unsigned int n_actual_col_dofs = global_cols.size(); + + // create arrays for the column data + // (indices and values) that will then be + // written into the matrix. Shortcut for + // deal.II sparse matrix + std::vector cols (n_actual_col_dofs); + std::vector vals (n_actual_col_dofs); + + // now do the actual job. + for (unsigned int i=0; i 0) + global_matrix.add(row, n_values, &cols[0], &vals[0], false, true); + } +} + + // similar function as above, but now // specialized for block matrices. See // the other function for additional @@ -2412,7 +2471,7 @@ distribute_local_to_global (const FullMatrix &local_matrix, { unsigned int * col_ptr = &cols[0]; number * val_ptr = &vals[0]; - resolve_matrix_row (global_rows, i, start_block, + resolve_matrix_row (global_rows, global_rows, i, start_block, end_block, local_matrix, col_ptr, val_ptr); const unsigned int n_values = col_ptr - &cols[0]; Assert (n_values == (unsigned int)(val_ptr - &vals[0]),