From 720d76a4baaec4872eb02ae2fa8ffc8f1fc34204 Mon Sep 17 00:00:00 2001
From: Jukka Aho <ahojukka5@gmail.com>
Date: Sat, 19 Aug 2017 14:15:11 +0300
Subject: [PATCH] Eliminate Dirichlet boundary conditions differently

Slicing is very expensive operation on sparse matrices. Avoid slicing
by calculating matrix product `K_red = C1*K*C1`, where C1 is diagonal
matrix. This should be much cheaper.
---
 src/solvers_modal.jl | 13 +++++++------
 1 file changed, 7 insertions(+), 6 deletions(-)

diff --git a/src/solvers_modal.jl b/src/solvers_modal.jl
index 6ab4c48..df57048 100644
--- a/src/solvers_modal.jl
+++ b/src/solvers_modal.jl
@@ -44,12 +44,13 @@ function eliminate_boundary_conditions!(K_red::SparseMatrixCSC,
     @assert nnz(g) == 0
     @assert C1 == C2
     @assert isdiag(C1)
-    nz = get_nonzero_rows(C1)
-    info("bc $(problem.name): $(length(nz)) nonzeros, $nz")
-    K_red[nz,:] = 0.0
-    K_red[:,nz] = 0.0
-    M_red[nz,:] = 0.0
-    M_red[:,nz] = 0.0
+    fixed_dofs = get_nonzero_rows(C1)
+    P = ones(ndim)
+    P[fixed_dofs] = 0.0
+    Q = spdiagm(P)
+    K_red[:,:] = Q' * K_red * Q
+    M_red[:,:] = Q' * M_red * Q
+    return
 end
 
 """ Given data vector, return slave displacements. """