# Building mass and stiffness matrix for time dependent PDE with Gridap

**URL:** <https://discourse.julialang.org/t/building-mass-and-stiffness-matrix-for-time-dependent-pde-with-gridap/126164>\
**Category:** Modelling & Simulations\
**Tags:** question, gridap\
**Created:** [February 21, 2025, 5:29pm UTC](https://discourse.julialang.org/t/building-mass-and-stiffness-matrix-for-time-dependent-pde-with-gridap/126164 "2025-02-21T17:29:58Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![Cirdans-Home](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cirdans-home/32/215388_2.png) [@Cirdans-Home](https://discourse.julialang.org/u/Cirdans-Home)\
**Post date:** [February 21, 2025, 5:29pm UTC](https://discourse.julialang.org/t/building-mass-and-stiffness-matrix-for-time-dependent-pde-with-gridap/126164/1 "2025-02-21T17:29:58Z")

</div>

Hi,

I have implemented an algorithm for the solution of system of ODEs that I want to test against a method of line discretization of a PDE.

My routine needs to take as input the mass matrix (assumed to be time-independent) and the right-hand side. Something that should be applied to:

M y’ = K y + f(t)  
y[0] = y0

Looking at the Poisson tutorial I have tried doing something of the form:

```julia
u0f(x) = 2.0
ftrue(x,t) = 1.0 + 0.5*sin(2*pi*t)
f(x) = ftrue(x,0.0)
h(x) = 3.0
m(u,v) = ∫( v*u )*dΩ
a(u,v) = ∫( ∇(v)⋅∇(u) )*dΩ
b(v) = ∫( v*f )*dΩ + ∫( v*h )*dΓ
u0(v) = ∫( v*u0f )*dΩ
massop = AffineFEOperator(m,u0,Ug,V0)
op = AffineFEOperator(a,b,Ug,V0)

Mass = get_matrix(massop)
y0 = get_vector(massop)
L = -get_matrix(op)

function fx(t,ftrue,a,Ug,V0,h,dΩ,dΓ)
    f(x) = ftrue(x,t)
    b(v) = ∫( v*f )*dΩ + ∫( v*h )*dΓ
    op = AffineFEOperator(a,b,Ug,V0)
    return get_vector(op)
end

ff(t) = fx(t,ftrue,a,Ug,V0,h,dΩ,dΓ)

```

so that I can pass to the ODE solver the different pieces:

```julia
tspan = [0,1]
n = 100
t,y = odesolve(Mass,L,ff,tspan,n,y0)

```

The solution seems to numerically diverge. I have tried the same code against a finite difference discretization of an analogous problem, and things seem to work fine. I suppose I’m doing something wrong with the mass matrix/boundary assembly.

I have also tried looking at the

```julia
TransientLinearFEOperator( )

```

but those object do not seem to have matrix assembly routines.

Any pointers would be very helpful.

Thank you.
