Dear all, This request stems from the "How to add a source term for PETSCFV ?" thread. Basically, I would like to output in a VTK file all the components of the gradient computed using the following piece of code (given sol is a global representation of the solution over the whole mesh described by the manager dm with a unique field fvm).
ierr = PetscFVSetComputeGradients(fvm,PETSC_TRUE);CHKERRQ(ierr); ierr = DMConvert(dm, DMPLEX, &plex);CHKERRQ(ierr); ierr = DMPlexGetDataFVM(plex, fvm, &cellGeom, &faceGeom, &gradDM);CHKERRQ(ierr); ierr = DMCreateLocalVector(plex,&locX);CHKERRQ(ierr); ierr = DMPlexInsertBoundaryValues(plex, PETSC_TRUE, locX, 0.0, faceGeom, cellGeom, NULL);CHKERRQ(ierr); ierr = DMGlobalToLocalBegin(plex, sol, INSERT_VALUES, locX);CHKERRQ(ierr); ierr = DMGlobalToLocalEnd (plex, sol, INSERT_VALUES, locX);CHKERRQ(ierr); ierr = DMCreateGlobalVector(gradDM, &grad);CHKERRQ(ierr); ierr = DMPlexReconstructGradientsFVM(plex, locX, grad);CHKERRQ(ierr); ierr = DMCreateLocalVector(gradDM, &locGrad);CHKERRQ(ierr); ierr = DMGlobalToLocalBegin(gradDM, grad, INSERT_VALUES, locGrad);CHKERRQ(ierr); ierr = DMGlobalToLocalEnd(gradDM, grad, INSERT_VALUES, locGrad);CHKERRQ(ierr); ierr = VecDestroy(&grad);CHKERRQ(ierr); ierr = PetscViewerCreate(PetscObjectComm((PetscObject)gradDM), &viewer);CHKERRQ(ierr); ierr = PetscSNPrintf(filename,sizeof filename,"%s-%03D-gradient.vtu",user->outputBasename,stepnum);CHKERRQ(ierr); ierr = PetscViewerSetType(viewer, PETSCVIEWERVTK);CHKERRQ(ierr); ierr = PetscViewerFileSetName(viewer, filename);CHKERRQ(ierr); ierr = VecView(locGrad,viewer);CHKERRQ(ierr); I am having troubles with a few things I think but I cannot figure out how to put it together in the right order ... (i) I am not sure whether VecView should be called with a local or a global vector (ii) I can see from the error log that gradDM does not have fields and thus when it goes through dmplexvtu.c it crashes and cannot write the file. I tried to add the fvm as a field, but it does not work either, it crashes on a memory corruption. Could somebody please advise ? Thank you very much in advance, Thibault Bridel-Bertomeu — Eng, MSc, PhD Research Engineer CEA/CESTA 33114 LE BARP Tel.: (+33)557046924 Mob.: (+33)611025322 Mail: [email protected] >> >>
