[petsc-dev] hdf5 viewer question

[Dominic Meiser]

Hi Matt,


On 08/10/2015 08:39 AM, Matthew Knepley wrote:
> On Mon, Aug 10, 2015 at 9:55 AM, Dominic Meiser <dmeiser at txcorp.com
> <mailto:dmeiser at txcorp.com>> wrote:
>
>     Hi Matt,
>
>
>     On 06/23/2015 11:24 AM, Matthew Knepley wrote:
>
>         There could be a bug with the calculation of chunksize. I will
>         run your
>         example as soon as I can.
>
>
>     Have you had a chance to run the example? Thanks.
>
>
> Unfortunately, no, and now my mail has lost testcase.c. Can you resend
> to me?

testcase.c is attached. Thanks very much for having a look at this.

Cheers,
Dominic

>
>    Thanks,
>
>       Matt
>
>
>     Dominic
>
>
>     --
>     Dominic Meiser
>     Tech-X Corporation
>     5621 Arapahoe Avenue
>     Boulder, CO 80303
>     USA
>     Telephone: 303-996-2036 <tel:303-996-2036>
>     Fax: 303-448-7756 <tel:303-448-7756>
>     www.txcorp.com <http://www.txcorp.com>
>
>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which
> their experiments lead.
> -- Norbert Wiener

-- 
Dominic Meiser
Tech-X Corporation
5621 Arapahoe Avenue
Boulder, CO 80303
USA
Telephone: 303-996-2036
Fax: 303-448-7756
www.txcorp.com
-------------- next part --------------
#include <petscvec.h>
#include <petscdmda.h>
#include <petscdm.h>
#include <petscviewer.h>
#include <petscviewerhdf5.h>

struct Field {
  PetscScalar E;
  PetscScalar H;
};
static const char *fieldNames[] = {"E", "H"};

PetscErrorCode setInitialState(Vec *field);
PetscErrorCode dump(Vec field, const char *fieldname, const char *filename);

int main(int argn, char **argv) {
  PetscErrorCode ierr;
  DM             da;
  Vec            field;

  PetscFunctionBegin;
  ierr = PetscInitialize(&argn, &argv, NULL, NULL);CHKERRQ(ierr);
  ierr = DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_PERIODIC, -1000, 2, 1, NULL, &da);CHKERRQ(ierr);
  ierr = DMDASetUniformCoordinates(da, -1.0, 1.0, 0, 0, 0, 0);CHKERRQ(ierr);
  ierr = DMDASetCoordinateName(da, 0, "x");CHKERRQ(ierr);
  ierr = DMDASetFieldNames(da, fieldNames);CHKERRQ(ierr);
  ierr = DMSetFromOptions(da);CHKERRQ(ierr);
  ierr = DMCreateLocalVector(da, &field);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject)field, "field");CHKERRQ(ierr);
  ierr = setInitialState(&field);CHKERRQ(ierr);
  ierr = dump(field, "field", "field_initial.h5");CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

PetscErrorCode setInitialState(Vec *field) {
  PetscErrorCode ierr;
  struct Field   *f;
  DM             da;
  PetscInt       i, x, m;

  PetscFunctionBegin;
  ierr = VecGetDM(*field, &da);CHKERRQ(ierr);
  ierr = DMDAVecGetArray(da, *field, &f);CHKERRQ(ierr);
  ierr = DMDAGetCorners(da, &x, NULL, NULL, &m, NULL, NULL);CHKERRQ(ierr);
  for (i = x; i < x + m; ++i) {
    f[i].E = 0.0;
    f[i].H = 0.0;
  }
  ierr = DMDAVecRestoreArray(da, *field, &f);CHKERRQ(ierr);
  DMLocalToLocalBegin(da, *field, INSERT_VALUES, *field);CHKERRQ(ierr);
  DMLocalToLocalEnd(da, *field, INSERT_VALUES, *field);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

PetscErrorCode dump(Vec field, const char *fieldname, const char *filename) {
  PetscErrorCode ierr;
  PetscViewer    v;
  DM             dm;
  Vec            tmp;
  PetscInt       numProcs;

  PetscFunctionBegin;
  ierr = VecGetDM(field, &dm);CHKERRQ(ierr);
  ierr = DMDAGetInfo(dm, 0, 0, 0, 0, &numProcs, 0, 0, 0, 0, 0, 0, 0, 0);
  if (numProcs > 1) {
    ierr = DMGetGlobalVector(dm, &tmp);CHKERRQ(ierr);
    ierr = DMLocalToGlobalBegin(dm, field, INSERT_VALUES, tmp);CHKERRQ(ierr);
    ierr = DMLocalToGlobalEnd(dm, field, INSERT_VALUES, tmp);CHKERRQ(ierr);
  } else {
    tmp = field;
  }
  ierr = PetscObjectSetName((PetscObject)tmp, fieldname);CHKERRQ(ierr);

  ierr = PetscViewerHDF5Open(MPI_COMM_WORLD, filename, FILE_MODE_WRITE, &v);CHKERRQ(ierr);
  ierr = VecView(tmp, v);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&v);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}