dolfin team mailing list archive

[Evan Lezar <question125019@xxxxxxxxxxxxxxxxxxxxx>]

Question #125019 on DOLFIN changed:
https://answers.launchpad.net/dolfin/+question/125019

Evan Lezar posted a new comment:
At the moment I am using a bit of a hack - converting to numpy arrays.
Which suits me fine for testing purposes.

The code to do this is given below. There S_petc and T_petc are two PETc
matrices, and boundary_indicators is an integer array of 0's and 1's
from DirichletBC.get_bc().

I have not had a chance to look into how this would be done in PETC, as
I have some deadlines that are quite pressing at the moment.


def convert_dolfin_matrices_to_numpy_and_remove_dirichlet_dofs ( S_petc, T_petc, boundary_indicators ):
    """
    Convert the dolfin mass and stiffness matrices to numpy arrays and remove the rows and columns
    associated with the dirichlet boundary conditions
    """
    M = S_petc.size(0)
    dof_map = np.where(boundary_indicators == 0)[0]
    
    free_dofs = len(dof_map)
    S = np.zeros((free_dofs, free_dofs), np.float32, order='F')
    T = np.zeros((free_dofs, free_dofs), np.float32, order='F')
    
    one_row = np.zeros(M, np.float32)
    for i in range(free_dofs):
        petc_data = S_petc.getrow(dof_map[i])
        one_row[petc_data[0]] = petc_data[1].astype(np.float32)
        S[i,:] = one_row[dof_map]
        one_row[:] = 0;
        petc_data = T_petc.getrow(dof_map[i])
        one_row[petc_data[0]] = petc_data[1].astype(np.float32)
        T[i,:] = one_row[dof_map] 
        one_row[:] = 0;
        
    return S, T, dof_map, free_dofs, M

-- 
You received this question notification because you are a member of
DOLFIN Team, which is an answer contact for DOLFIN.