State space matrices including dirichlet BC

Question

Hello
I want to compute linearized state space matrices of a non-linear system, like incompressible Navier-Stokes discretized with say Taylor-Hood elements. I know how to generate the matrices to obtain a linearized system like

Mdu/dt = Au + B*p
0 = B^T *u

but I want to completely eliminate the dirichlet degrees of freedom for the velocity. If I write velocity as

u = uf + ud

where

uf = free dofs
ud = dirichlet dofs

then I want a system like this

Mfd(uf)/dt = Afuf + Adud + Bp

where I simply neglect the term d(ud)/dt. So here Mf, Af are smaller matrices from which dirichlet dofs have been completely removed.

Is there a simple way to construct these matrices ?

Thanks
praveen

Answer 1

I am not sure about a simple way. I have done it like this:

def condense_sysmatsbybcs(Ma,Aa,BTa,Ba,fv,fp,velbcs):
        """resolve the Dirichlet BCs and condense the sysmats

        to the inner nodes"""

        #auxiliary vector containing the dirichlet data
        auxu = np.zeros((len(fv),1))
        bcinds = []
        for bc in velbcs:
                bcdict = bc.get_boundary_values()
                auxu[bcdict.keys(),0] = bcdict.values()
                bcinds.extend(bcdict.keys())

        # putting the bcs into the right hand sides
        fvbc = - Aa*auxu # '*' is np.dot for csr matrices
        fpbc = - Ba*auxu

        # indices of the innernodes
        innerinds = np.setdiff1d(range(len(fv)),bcinds).astype(np.int32)

        # extract the inner nodes equation coefficients
        Mc = Ma[innerinds,:][:,innerinds]
        Ac = Aa[innerinds,:][:,innerinds]
        fvbc= fvbc[innerinds,:]
        Bc = Ba[:,innerinds]
        BTc = BTa[innerinds,:]

        bcvals = auxu[bcinds]

        return Mc, Ac, BTc, Bc, fvbc, fpbc, bcinds, bcvals, innerinds

where Ma,Aa,BTa,Ba are the unreduced scipy matrices representing mass, stiffness, gradient, divergence with respect to mixed finite elements, where fv,fp are the rhs of the momentum and the continuity equation, and where velbcs is a list of the dolfin's DirichletBCs. See this python module and the other files in this github repository for related functions.

Comments

Thanks. This is very useful. I have to first convert all matrices to scipy matrices, is it so ?If I understand correctly the assemble function creates a PetscMatrix.

---

In my approach, yes! Cf. this function and note the command

parameters.linear_algebra_backend = "uBLAS"

which sets the LA backend. I am not so deep in FEniCS, to tell you what is the internal representation.

---

Thanks. Thats is very clear. I am using mac binary so I dont have scipy. I will have to compile fenics myself to do this.

---

Ok, you can try this. Although, there might be other options with different packages. You may want to remove the 'answered' tag from my answer to attract a different solution.