Nrthugu

I am new to fenics and I really need help to debug my code. Basically I am trying to solve a time dependent partial deferential equation and I would like to compute the convergence rate of the error norms that I have. however I keep get this error which I don't understand what does it mean or even how to get over it!

rates[degree][error_type] = 



KeyError: 1

and this is my full code

T = 2.0       # final time

beta =  1.2   # parameter beta

tol = 3E-16

bound_a=0.0

bound_b=1.0

r=1.0

m=0.1





num_levels=3

max_degree=3

h = {}  # discretization parameter: h[degree][level]

E = {}  # error measure(s): E[degree][level][error_type]



# Iterate over degrees and mesh refinement levels

degrees = range(1, max_degree + 1)

for degree in degrees:

    h[degree] = 

    E[degree] = 

    for i in range(num_levels):

        num_steps = pow(4,i) # number of time steps

        dt = T / num_steps   # time step size

        # Create mesh and define function space

        mesh_division= int(sqrt(num_steps))

        mesh = IntervalMesh(mesh_division,bound_a,bound_b)

        h[degree].append(1.0 / mesh_division)

        V = FunctionSpace(mesh, 'P', 1)



        # Define boundary condition

        u_D = Expression('(b-a)/2 + beta*t +(x[0]-a)*(x[0]-b)', degree=1, beta=beta,a=bound_a,b=bound_b, t=0.0)



        def boundary(x, on_boundary):

            return on_boundary and (abs((x[0]-bound_a)*(x[0]-bound_b))<tol)





        bc = DirichletBC(V, u_D, boundary)



        # Define initial value

        u_n = interpolate(u_D, V)



        # Define variational problem

        u = TrialFunction(V)

        v = TestFunction(V)

        f = Constant(0.0)



        F = u*v*dx + dt*inner(grad(u), grad(v))*dx + dt*inner(grad(u_n), grad(u))*dx -dt*r*(m-u_n)*v*dx- (u_n + dt*f)*v*dx

        a, L = lhs(F), rhs(F)



        # Time-stepping

        u = Function(V)

        time = 0.0



        for n in range(num_steps):



            # Update current time

            time += dt

            u_D.t = time



            # Compute solution

            solve(a == L, u, bc)

            u_e= interpolate(u_D, V)



            def compute_errors(u_e, u):

                # Infinity norm based on nodal values

                E1 = abs(u_e.vector().get_local() - u.vector().get_local()).max()



                # L2 norm

                E2 = errornorm(u_e, u, norm_type='L2')



                # H1 seminorm

                E3 = errornorm(u_e, u, norm_type='H10')



                # Collect error measures in a dictionary with self-explanatory keys

                errors = {'infinity norm': E1,

                          'L2 norm': E2,

                          'H10 seminorm': E3}



                return errors



            errors = compute_errors(u_e, u)

            E[degree].append(errors)

            #Compute convergence rates

            etypes = list(E[1][0].keys())

            rates = {}

            for error_type in sorted(etypes):

                 rates[degree][error_type] = 

                 for i in range(1, num_levels):

                     Ei = E[degree][i][error_type]

                     r = ln(Ei / Eim1) / ln(h[degree][i] / h[degree][i - 1])

                     rates[degree][error_type].append(round(r, 2))



            # Update previous solution

            u_n.assign(u)

any help will be appreciated. Many thanks