1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 #include <petscdmplextransform.h> 3 4 PetscErrorCode DMPlexExtrude(DM dm, PetscInt layers, PetscReal thickness, PetscBool tensor, PetscBool symmetric, const PetscReal normal[], const PetscReal thicknesses[], DM *edm) 5 { 6 DMPlexTransform tr; 7 DM cdm, ecdm; 8 const char *prefix; 9 PetscOptions options; 10 PetscErrorCode ierr; 11 12 PetscFunctionBegin; 13 ierr = DMPlexTransformCreate(PetscObjectComm((PetscObject) dm), &tr);CHKERRQ(ierr); 14 ierr = DMPlexTransformSetDM(tr, dm);CHKERRQ(ierr); 15 ierr = DMPlexTransformSetType(tr, DMPLEXEXTRUDE);CHKERRQ(ierr); 16 ierr = PetscObjectGetOptionsPrefix((PetscObject) dm, &prefix);CHKERRQ(ierr); 17 ierr = PetscObjectSetOptionsPrefix((PetscObject) tr, prefix);CHKERRQ(ierr); 18 ierr = PetscObjectGetOptions((PetscObject) dm, &options);CHKERRQ(ierr); 19 ierr = PetscObjectSetOptions((PetscObject) tr, options);CHKERRQ(ierr); 20 ierr = DMPlexTransformExtrudeSetLayers(tr, layers);CHKERRQ(ierr); 21 if (thickness > 0.) {ierr = DMPlexTransformExtrudeSetThickness(tr, thickness);CHKERRQ(ierr);} 22 ierr = DMPlexTransformExtrudeSetTensor(tr, tensor);CHKERRQ(ierr); 23 ierr = DMPlexTransformExtrudeSetSymmetric(tr, symmetric);CHKERRQ(ierr); 24 if (normal) {ierr = DMPlexTransformExtrudeSetNormal(tr, normal);CHKERRQ(ierr);} 25 if (thicknesses) {ierr = DMPlexTransformExtrudeSetThicknesses(tr, layers, thicknesses);CHKERRQ(ierr);} 26 ierr = DMPlexTransformSetFromOptions(tr);CHKERRQ(ierr); 27 ierr = PetscObjectSetOptions((PetscObject) tr, NULL);CHKERRQ(ierr); 28 ierr = DMPlexTransformSetUp(tr);CHKERRQ(ierr); 29 ierr = PetscObjectViewFromOptions((PetscObject) tr, NULL, "-dm_plex_transform_view");CHKERRQ(ierr); 30 ierr = DMPlexTransformApply(tr, dm, edm);CHKERRQ(ierr); 31 ierr = DMCopyDisc(dm, *edm);CHKERRQ(ierr); 32 ierr = DMGetCoordinateDM(dm, &cdm);CHKERRQ(ierr); 33 ierr = DMGetCoordinateDM(*edm, &ecdm);CHKERRQ(ierr); 34 ierr = DMCopyDisc(cdm, ecdm);CHKERRQ(ierr); 35 ierr = DMPlexTransformCreateDiscLabels(tr, *edm);CHKERRQ(ierr); 36 ierr = DMPlexTransformDestroy(&tr);CHKERRQ(ierr); 37 if (*edm) { 38 ((DM_Plex *) (*edm)->data)->printFEM = ((DM_Plex *) dm->data)->printFEM; 39 ((DM_Plex *) (*edm)->data)->printL2 = ((DM_Plex *) dm->data)->printL2; 40 } 41 PetscFunctionReturn(0); 42 } 43 44 PetscErrorCode DMExtrude_Plex(DM dm, PetscInt layers, DM *edm) 45 { 46 PetscErrorCode ierr; 47 48 PetscFunctionBegin; 49 ierr = DMPlexExtrude(dm, layers, PETSC_DETERMINE, PETSC_TRUE, PETSC_FALSE, NULL, NULL, edm);CHKERRQ(ierr); 50 ierr = DMViewFromOptions(*edm, NULL, "-check_extrude");CHKERRQ(ierr); 51 PetscFunctionReturn(0); 52 } 53