dolfin team mailing list archive

[Kent-Andre Mardal <question143195@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Open => Answered

Kent-Andre Mardal proposed the following answer:
I did something like this a while ago and attach the code.

Here, I use the Trilinos backend to handle the LA on a global level.
I guess you may perform the static condensation on the element level
and this may be more efficient.

Anyway, the code is attach.

Kent

On 28 January 2011 16:55, David Ham
<question143195@xxxxxxxxxxxxxxxxxxxxx>wrote:

> New question #143195 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/143195
>
>
> Suppose I'm solving the linear shallow water equations with discontinuous
> velocity and continuous height:
>
> V = VectorFunctionSpace(mesh, 'DG', 1, dim=2)# Velocity space
>
> H = FunctionSpace(mesh, 'CG', 2)             # Height space
>
> W=V*H                                        # Mixed space of both.
>
> (v, q)=TestFunctions(W)
> (u, h)=TrialFunctions(W)
>
> M=inner(v,q)*dx+inner(u,h)*dx
>
> # Divergence term.
> Ct=-inner(u,grad(q))*dx+inner(avg(u),jump(q,n))*dS
> # Pressure gradient operator
> C=-g*Depth*adjoint(Ct)
>
> # Coriolis term
> F=params["f"]*inner(v,as_vector([-u[1],u[0]]))*dx
>
> If I then apply backward Euler timesteping, I get a system whose block
> structure is:
>
> /   M_u-F  C    \    /u_new \         / rhs_u\
> |                       |  |             |    =   |          |
> \   Ct         M_h/  \ h_new /         \ rhs_h/
>
> This system is (a) big  and (b) non-symmetric. However, the standard trick
> at this stage is to use static condensation to turn this into a discrete
> wave equation which is much smaller and is SPD:
>
> (-Ct(M_u-F)^(-1)C + M_h) h_new = Ct(M_u-F)^(-1)rhs_u + rhs_h
>
> The reason this is affordable is that M_u-F contains no derivatives and
> therefore no facet integrals. This means that the inversion commutes with
> assembly so one can directly assemble the inverse matrix by first inverting
> the local tensor  and then conducting global assembly. (M_u-F)^(-1) has the
> same sparsity as M_u-F.
>
> So, the question is, how do I get Dolfin to assemble (M_u-F)^(-1) rather
> than M_u-F?
>
> (The next question will be, can I multiply sparse matrices as above in
> Dolfin, but one thing at a time!)
>
> Many thanks,
>
> David
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin<https://launchpad.net/%7Edolfin>
> Post to     : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin<https://launchpad.net/%7Edolfin>
> More help   : https://help.launchpad.net/ListHelp
>

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