1*f918ec44SMatthew G. Knepley #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2*f918ec44SMatthew G. Knepley 3*f918ec44SMatthew G. Knepley #ifdef PETSC_HAVE_LIBCEED 4*f918ec44SMatthew G. Knepley #include <petscdmplexceed.h> 5*f918ec44SMatthew G. Knepley 6*f918ec44SMatthew G. Knepley /* Define the map from the local vector (Lvector) to the cells (Evector) */ 7*f918ec44SMatthew G. Knepley PetscErrorCode DMPlexGetCeedRestriction(DM dm, CeedElemRestriction *ERestrict) 8*f918ec44SMatthew G. Knepley { 9*f918ec44SMatthew G. Knepley PetscErrorCode ierr; 10*f918ec44SMatthew G. Knepley 11*f918ec44SMatthew G. Knepley PetscFunctionBeginUser; 12*f918ec44SMatthew G. Knepley PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 13*f918ec44SMatthew G. Knepley PetscValidPointer(ERestrict, 2); 14*f918ec44SMatthew G. Knepley if (!dm->ceedERestrict) { 15*f918ec44SMatthew G. Knepley PetscDS ds; 16*f918ec44SMatthew G. Knepley PetscFE fe; 17*f918ec44SMatthew G. Knepley PetscSpace sp; 18*f918ec44SMatthew G. Knepley PetscSection s; 19*f918ec44SMatthew G. Knepley PetscInt Nf, *Nc, c, P, cStart, cEnd, Ncell, cell, lsize, *erestrict, eoffset; 20*f918ec44SMatthew G. Knepley PetscBool simplex; 21*f918ec44SMatthew G. Knepley Ceed ceed; 22*f918ec44SMatthew G. Knepley 23*f918ec44SMatthew G. Knepley ierr = DMPlexIsSimplex(dm, &simplex);CHKERRQ(ierr); 24*f918ec44SMatthew G. Knepley if (simplex) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "LibCEED does not work with simplices"); 25*f918ec44SMatthew G. Knepley ierr = DMGetDS(dm, &ds);CHKERRQ(ierr); 26*f918ec44SMatthew G. Knepley ierr = PetscDSGetNumFields(ds, &Nf); 27*f918ec44SMatthew G. Knepley if (Nf != 1) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_SUP, "LibCEED only works with one field right now"); 28*f918ec44SMatthew G. Knepley ierr = PetscDSGetComponents(ds, &Nc);CHKERRQ(ierr); 29*f918ec44SMatthew G. Knepley ierr = PetscDSGetDiscretization(ds, 0, (PetscObject *) &fe);CHKERRQ(ierr); 30*f918ec44SMatthew G. Knepley ierr = PetscFEGetBasisSpace(fe, &sp);CHKERRQ(ierr); 31*f918ec44SMatthew G. Knepley ierr = PetscSpaceGetDegree(sp, &P, NULL);CHKERRQ(ierr); 32*f918ec44SMatthew G. Knepley ++P; 33*f918ec44SMatthew G. Knepley ierr = DMGetLocalSection(dm, &s);CHKERRQ(ierr); 34*f918ec44SMatthew G. Knepley ierr = DMPlexGetHeightStratum(dm, 0, &cStart,& cEnd);CHKERRQ(ierr); 35*f918ec44SMatthew G. Knepley Ncell = cEnd - cStart; 36*f918ec44SMatthew G. Knepley 37*f918ec44SMatthew G. Knepley ierr = PetscMalloc1(Ncell*P*P, &erestrict);CHKERRQ(ierr); 38*f918ec44SMatthew G. Knepley for (cell = cStart, eoffset = 0; cell < cEnd; ++cell) { 39*f918ec44SMatthew G. Knepley PetscInt Ni, *ind, i; 40*f918ec44SMatthew G. Knepley 41*f918ec44SMatthew G. Knepley ierr = DMPlexGetClosureIndices(dm, s, s, cell, PETSC_TRUE, &Ni, &ind, NULL, NULL);CHKERRQ(ierr); 42*f918ec44SMatthew G. Knepley for (i = 0; i < Ni; i += Nc[0]) { 43*f918ec44SMatthew G. Knepley for (c = 0; c < Nc[0]; ++c) { 44*f918ec44SMatthew G. Knepley if (ind[i+c] != ind[i] + c) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Cell %D closure indices not interlaced", cell); 45*f918ec44SMatthew G. Knepley } 46*f918ec44SMatthew G. Knepley erestrict[eoffset++] = ind[i]; 47*f918ec44SMatthew G. Knepley } 48*f918ec44SMatthew G. Knepley ierr = DMPlexRestoreClosureIndices(dm, s, s, cell, PETSC_TRUE, &Ni, &ind, NULL, NULL);CHKERRQ(ierr); 49*f918ec44SMatthew G. Knepley } 50*f918ec44SMatthew G. Knepley if (eoffset != Ncell*P*P) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Size of array %D != %D restricted dofs", Ncell*P*P, eoffset); 51*f918ec44SMatthew G. Knepley 52*f918ec44SMatthew G. Knepley ierr = DMGetCeed(dm, &ceed);CHKERRQ(ierr); 53*f918ec44SMatthew G. Knepley ierr = PetscSectionGetStorageSize(s, &lsize);CHKERRQ(ierr); 54*f918ec44SMatthew G. Knepley ierr = CeedElemRestrictionCreate(ceed, Ncell, P*P, Nc[0], 1, lsize, CEED_MEM_HOST, CEED_COPY_VALUES, erestrict, &dm->ceedERestrict);CHKERRQ(ierr); 55*f918ec44SMatthew G. Knepley ierr = PetscFree(erestrict);CHKERRQ(ierr); 56*f918ec44SMatthew G. Knepley } 57*f918ec44SMatthew G. Knepley *ERestrict = dm->ceedERestrict; 58*f918ec44SMatthew G. Knepley PetscFunctionReturn(0); 59*f918ec44SMatthew G. Knepley } 60*f918ec44SMatthew G. Knepley 61*f918ec44SMatthew G. Knepley #endif 62