Switch to local numbering in dof maps

Issue #119 resolved
Prof Garth Wells created an issue

Use local (to a process) numbering in dof maps. See mailing list discussion at http://fenicsproject.org/pipermail/fenics/2013-September/000512.html.

Comments (24)

  1. Prof Garth Wells reporter

    Dofs still have a unique global index, and the dofmap will need to provide a local-to-global map for all dofs (both owned and not owned). This map is passed to the linear algebra backend at construction of a vector or matrix.

  2. Prof Garth Wells reporter

    Similar to what we have now. The N local dofs have indices in [0, N), and non-local dofs have indices in dofs in [N, ...).

  3. Prof Garth Wells reporter

    Lot of questions! In principle it should all work just as it does now. The difference is that the linear algebra backends will take indices [0, N), and the dofmap needs to support reflect this. If one was to apply a map would could move back-and-forth between the proposed and current dof maps. The motivation behind the change is that local-to-global is cheap (an array look-up), but global-to-local is expensive (requiring a map). For field-split we need cheap local-to-global.

  4. Johan Hake

    Ok, so:

    sub_dofmap.cell_dofs
    

    will return a vector with numbers from 0-N, where N is the local range from the parent dofmap?

    I am working on the sub function assignment and I want to understand how the DofMap works now and tomorrow :)

  5. Prof Garth Wells reporter

    @blechta It depends somewhat on the backend, but ideally it will be possible to index via local or global indices.

  6. Jan Blechta

    @garth-wells I don't see how you can access off-diagonal entries (i.e. in column corresponding to off-process DOF) with local indices.

  7. Jan Blechta

    @garth-wells You can get anything from local row, can't you? Additionally you can set anything.

    I just say that switching to indexing by local indices introduces more complexity in GenericTensor as same entry will be indexed differently form different ranks.

  8. Jan Blechta

    @garth-wells In what sense don't we have a function to get matrix entries? Can you tell me whether following (indeed non-sensical) code would run with your changes? Basically I'd like to write a wrapper for grouping Matrices/Vectors into another (for my special needs, probably not welcomed into master) and I need to know how indexing of Tensors will work. I would be, for example, pleased to hear, everything will work the same way or it will not definitely work in a current way, from you.

    from dolfin import *
    import numpy as np
    
    comm = mpi_comm_world()
    rank = MPI.rank(comm)
    size = MPI.size(comm)
    
    mesh = UnitSquareMesh(4, 5)
    V = FunctionSpace(mesh, 'CG', 1)
    
    u, v = TrialFunction(V), TestFunction(V)
    a = inner(grad(u), grad(v))*dx
    A = assemble(a)
    print rank, A.size(0), A.size(1), A.local_range(0)
    
    # We can get any local row
    try:
        print rank, A.getrow(0)
    except RuntimeError:
        pass
    
    # We can set any nonlocal row
    j = np.array([0, 1, 2], 'uintp')
    v = np.array([42., 43., 44.])
    if rank > 0:
        A.setrow(0, j, v)
    A.apply('insert')
    
    # We can get any block on local rows, we use global indices
    i = np.array([14, 15, 21, 23], 'intc')
    j = np.array([14, 15, 21, 23], 'intc')
    v = np.ndarray((i.size, j.size))
    try:
        A.get(v, i, j)
    except RuntimeError:
        print rank, 'Cannot get nonlocal row!'
    else:
        print rank, v
    
  9. Prof Garth Wells reporter

    @blechta For matrices and vectors that can presently have entries set using global indices, this will still be possible. The parallel backends allow (in most cases) entries to be set using local and global indices. An exception is PETSc nested matrices, which require local indices and which is one of the the motivation for a change.

  10. Log in to comment