DG fem for 1D advection equation - TypeError with dolfin-dg

dolfinx
1

code :

import dolfinx
import numpy as np
import ufl
from mpi4py import MPI
from petsc4py import PETSc

import dolfin_dg.dolfinx
import dolfin_dg.hdg_form

comm = MPI.COMM_WORLD
poly_o = 2
n_eles = [8, 16, 32]
l2errors_u = np.zeros_like(n_eles, dtype=np.double)
l2errors_p = np.zeros_like(n_eles, dtype=np.double)
hs = np.zeros_like(n_eles, dtype=np.double)


for run_no, n_ele in enumerate(n_eles):
    mesh = dolfinx.mesh.create_unit_square(
        comm, n_ele, n_ele, cell_type=dolfinx.mesh.CellType.triangle,
        ghost_mode=dolfinx.cpp.mesh.GhostMode.shared_facet)
    V = dolfinx.fem.FunctionSpace(mesh, ("DG", poly_o))
    u = dolfinx.fem.Function(V)
    v = ufl.TestFunction(V)
    x = ufl.SpatialCoordinate(mesh)
    u_soln = ufl.exp(x[0] - x[1])
    b = dolfinx.fem.Constant(mesh, np.array((1, 1), dtype=np.double))
    n = ufl.FacetNormal(mesh)
    # First order terms
    def F_c(u):
        return b * u
    flux_function = dolfin_dg.LocalLaxFriedrichs(
        flux_jacobian_eigenvalues=lambda u, n: ufl.dot(b, n))
    ho = dolfin_dg.operators.HyperbolicOperator(
        mesh, V, dolfin_dg.DGDirichletBC(ufl.ds, u_soln), F_c, flux_function)
    F = ho.generate_fem_formulation(u, v)
    # Volume source
    f = ufl.div(F_c(u_soln))
    F += - f * v * ufl.dx
    J = ufl.derivative(F, u)
    F, J = dolfinx.fem.form(F), dolfinx.fem.form(J)
    problem = dolfin_dg.dolfinx.nls.NonlinearPDE_SNESProblem(F, J, u, [])
    snes = PETSc.SNES().create(mesh.comm)
    opts = PETSc.Options()
    opts["snes_monitor"] = None
    snes.setFromOptions()
    snes.setFunction(problem.F_mono, dolfinx.fem.petsc.create_vector(F))
    snes.setJacobian(problem.J_mono, J=dolfinx.fem.petsc.create_matrix(J))
    snes.getKSP().getPC().setType("lu")
    snes.getKSP().getPC().setFactorSolverType("mumps")
    snes.solve(None, u.vector)
    l2error_u = comm.allreduce(
        dolfinx.fem.assemble.assemble_scalar(dolfinx.fem.form(
            (u - u_soln) ** 2 * ufl.dx)) ** 0.5,
        op=MPI.SUM)
    h_measure = dolfinx.cpp.mesh.h(
        mesh._cpp_object, 2, np.arange(mesh.topology.index_map(2).size_local,
                                       dtype=np.int32))
    hmin = mesh.comm.allreduce(h_measure.min(), op=MPI.MIN)
    hs[run_no] = hmin
    l2errors_u[run_no] = l2error_u

print(l2errors_u)
rates_u = np.log(l2errors_u[:-1] / l2errors_u[1:]) / np.log(hs[:-1] / hs[1:])
print("rates u: %s" % str(rates_u))

ERROR :

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[3], line 32
     28 flux_function = dolfin_dg.LocalLaxFriedrichs(
     29     flux_jacobian_eigenvalues=lambda u, n: ufl.dot(b, n))
     30 ho = dolfin_dg.operators.HyperbolicOperator(
     31     mesh, V, dolfin_dg.DGDirichletBC(ufl.ds, u_soln), F_c, flux_function)
---> 32 F = ho.generate_fem_formulation(u, v)
     34 # Volume source
     35 f = ufl.div(F_c(u_soln))

File ~/dolfin_dg/dolfin_dg/operators.py:275, in HyperbolicOperator.generate_fem_formulation(self, u, v, dx, dS)
    272 F_vec = [F_0, self.F_c]
    273 L_vec = [ufl.div]
--> 275 fos = dolfin_dg.primal.FirstOrderSystem(F_vec, L_vec, u, v)
    277 F = fos.domain(dx)
    279 n = ufl.FacetNormal(u.ufl_domain())

File ~/dolfin_dg/dolfin_dg/primal/__init__.py:82, in FirstOrderSystem.__init__(self, F_vec, L_vec, u, v, n_ibps)
     79 self.L_vec = L_vec
     80 self.F_vec = F_vec
     81 self.G_vec = [
---> 82     dolfin_dg.math.homogenize(F_vec[i], u, L_vec[i](F_vec[i + 1](u)))
     83     for i in range(len(F_vec) - 1)]
     84 self.G_vec.append(
     85     dolfin_dg.math.homogenize(F_vec[-1], u, u)
     86 )
     88 v_vec = [v]

File ~/dolfin_dg/dolfin_dg/primal/__init__.py:60, in first_order_flux.<locals>.flux_wrapper.<locals>.flux_guard(u, flux)
     58 if flux is None:
     59     flux = flux_func(u)
---> 60 return func(u, flux)

TypeError: F_c() takes 1 positional argument but 2 were given
3

Please read Read before posting: How do I get my question answered? and, in particular, format your code properly with “```”.

@nate may be able to help, but first of all you should edit your post in a way that makes it legible.

4

sorry for the inconvinence

5

No worries: just go to the first post and edit it :wink:

6

please help me to remove the error

7

I helped you formatting properly the code (but please next time do this yourself), but can’t help you any further with dolfin_dg because I’ve never used it :slight_smile:

8

I’ll take a look when I get time. I suspect there’s a lack of coverage of that class after I migrated away from that technique recently.

9

This is indeed an embarrassing error on my part. I’ll keep you updated after I’ve composed a fix. Sadly this project is not a priority for me at the moment.

10

@nate we are planning to implement dg in fenicsx if you could help us with this issue it will be very helpful for us