Problem for the Real space

[zhyxue0]

I tried this example The (R) space by using different parameters of the solver（gmres and direct solver), but all failed.

Here is the code:

from firedrake import *

m = UnitSquareMesh(25, 25)
V = FunctionSpace(m, 'CG', 1)
R = FunctionSpace(m, 'R', 0)
W = V * R
u, r = TrialFunctions(W)
v, s = TestFunctions(W)

a = inner(grad(u), grad(v))*dx + u*s*dx + v*r*dx
L = -v*ds(3) + v*ds(4)

solver_parameters={
                   "ksp_type": "preonly",
                   "pc_type": "lu",
                #    "pc_factor_mat_solver_type": "mumps",
                   "pc_factor_mat_solver_type": "umfpack",
                #    "pc_factor_makt_solver_type": "superlu_dist",
                   'ksp_converged_reason': None,
                   'ksp_monitor_true_residual': None,
                #    'ksp_view': None
                  }

w = Function(W)

solve(a == L, w, solver_parameters = solver_parameters)

u, s = split(w)
exact = Function(V)
x, y = SpatialCoordinate(m)
exact.interpolate(y - 0.5)
print(sqrt(assemble((u - exact)*(u - exact)*dx)))

Here is the feedback：

Traceback (most recent call last):
  File "/home/firedrake/src/PyOP2/pyop2/caching.py", line 163, in __new__
    return cache[key]
KeyError: ((MixedDataSet((DataSet(Set((676, 676, 676), 'set_#x7fec80ff2c10'), (1,), 'None_nodes_dset'), GlobalDataSet(Global((1,), array([0.]), dtype('float64'), 'global_#x7feb83dd64c0')))), MixedDataSet((DataSet(Set((676, 676, 676), 'set_#x7fec80ff2c10'), (1,), 'None_nodes_dset'), GlobalDataSet(Global((1,), array([0.]), dtype('float64'), 'global_#x7feb83dd64c0'))))), frozenset({((MixedMap((Map(Set((1250, 1250, 1250), 'Cells'), Set((676, 676, 676), 'set_#x7fec80ff2c10'), 3, None, 'None_cell_node'), None)), MixedMap((Map(Set((1250, 1250, 1250), 'Cells'), Set((676, 676, 676), 'set_#x7fec80ff2c10'), 3, None, 'None_cell_node'), None))), (<IterationRegion.ALL: 4>,))}), False, False)

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/firedrake/src/firedrake/firedrake/assemble.py", line 174, in allocate_matrix
    sparsity = op2.Sparsity((test.function_space().dof_dset,
  File "/home/firedrake/src/PyOP2/pyop2/caching.py", line 165, in __new__
    obj = make_obj()
  File "/home/firedrake/src/PyOP2/pyop2/caching.py", line 146, in make_obj
    obj.__init__(*args, **kwargs)
  File "/home/firedrake/src/PyOP2/pyop2/types/mat.py", line 125, in __init__
    raise ex.SparsityFormatError("Mixed monolithic matrices with Global rows or columns are not supported.")
pyop2.exceptions.SparsityFormatError: Mixed monolithic matrices with Global rows or columns are not supported.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "Test_R.py", line 28, in <module>
    solve(a == L, w, solver_parameters = solver_parameters)
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/firedrake/src/firedrake/firedrake/adjoint/solving.py", line 53, in wrapper
    output = solve(*args, **kwargs)
  File "/home/firedrake/src/firedrake/firedrake/solving.py", line 130, in solve
    _solve_varproblem(*args, **kwargs)
  File "/home/firedrake/src/firedrake/firedrake/solving.py", line 156, in _solve_varproblem
    solver = vs.LinearVariationalSolver(problem, solver_parameters=solver_parameters,
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/firedrake/src/pyadjoint/pyadjoint/tape.py", line 46, in wrapper
    return function(*args, **kwargs)
  File "/home/firedrake/src/firedrake/firedrake/adjoint/variational_solver.py", line 45, in wrapper
    init(self, problem, *args, **kwargs)
  File "/home/firedrake/src/firedrake/firedrake/variational_solver.py", line 209, in __init__
    ctx.set_jacobian(self.snes)
  File "/home/firedrake/src/firedrake/firedrake/solving_utils.py", line 302, in set_jacobian
    snes.setJacobian(self.form_jacobian, J=self._jac.petscmat,
  File "/home/firedrake/src/PyOP2/pyop2/utils.py", line 62, in __get__
    obj.__dict__[self.__name__] = result = self.fget(obj)
  File "/home/firedrake/src/firedrake/firedrake/solving_utils.py", line 504, in _jac
    return allocate_matrix(self.J,
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/firedrake/src/firedrake/firedrake/assemble.py", line 181, in allocate_matrix
    raise ValueError("Monolithic matrix assembly not supported for systems "
ValueError: Monolithic matrix assembly not supported for systems with R-space blocks

I noticed that the feedback mentioned that Mixed Monolithic Matrices with Global Rows or Columns Are Not Supported and Monolithic matrix assembly not supported for systems with R-space blocks.

Can I only use the default solver Schur complement elimination PC？ Or there are other choices？

[dham]

The constraint is that only assembling a nested matrix is supported. (Equivalent to passing mat_type="nest" to solve). This means that direct solvers won't work ( though you can use them inside a Schur compliment) however iterative solvers such as GMRes should work.

[zhyxue0]

Thanks, I change the solver_parameters, but it still does not work

solver_parameters={
                   "ksp_type": "gmres",
                   'ksp_converged_reason': None,
                   'ksp_monitor_true_residual': None,
                   'ksp_view': None
                  }

Do I need to change other parameters?

The constraint is that only assembling a nested matrix is supported. (Equivalent to passing mat_type="nest" to solve). This means that direct solvers won't work ( though you can use them inside a Schur compliment) however iterative solvers such as GMRes should work.

[dham]

try passing mat_type="nest" to solve.

[zhyxue0]

try passing mat_type="nest" to solve.

Thanks, it works.

solver_parameters={
                   "ksp_type": "gmres",
                   'mat_type': 'nest',          
                   'ksp_rtol': 1e-15,
                   'ksp_atol': 1e-50,
                   'ksp_divtol': 1e4,
                   'ksp_max_it': 10000,
                   'ksp_converged_reason': None,
                   'ksp_monitor_true_residual': None,
                   'ksp_view': None
                  }