[petsc-users] Problems with PetscObjectViewFromOptions in Fortran
Barry Smith
bsmith at mcs.anl.gov
Sat May 20 15:09:55 CDT 2017
> On May 18, 2017, at 5:11 AM, Fabian.Jakub <Fabian.Jakub at physik.uni-muenchen.de> wrote:
>
> Dear Petsc Team,
>
> I have a problem with object viewing through PetscObjectViewFromOptions
>
> The C Version works fine, e.g.
>
> static char help[] = "Testing multiple PetscObjectViewFromOptions";
> #include <petsc.h>
>
> int main(int argc,char **argv) {
> DM dmA, dmB;
> PetscInitialize(&argc,&argv,(char*)0,help);
>
> PetscErrorCode ierr;
>
> ierr = DMPlexCreate(PETSC_COMM_WORLD, &dmA); CHKERRQ(ierr);
> ierr = PetscObjectSetName((PetscObject) dmA, "DMPlex_A"); CHKERRQ(ierr);
>
> ierr = DMPlexCreate(PETSC_COMM_WORLD, &dmB); CHKERRQ(ierr);
> ierr = PetscObjectSetName((PetscObject) dmB, "DMPlex_B"); CHKERRQ(ierr);
>
> PetscObjectViewFromOptions((PetscObject) dmA, NULL, "-dmA");
> PetscObjectViewFromOptions((PetscObject) dmB, NULL, "-dmB");
>
> ierr = DMDestroy(&dmA); CHKERRQ(ierr);
> ierr = DMDestroy(&dmB); CHKERRQ(ierr);
>
> PetscFinalize();
> }
>
> and running it with -help, correctly produces the options and views as:
>
> -dmA
>
> -dmB
>
>
> but the equivalent in Fortran, e.g.:
>
> program main
> #include "petsc/finclude/petsc.h"
> use petsc
> implicit none
>
> PetscErrorCode :: ierr
>
> DM :: dmA, dmB
>
> call PetscInitialize(PETSC_NULL_CHARACTER, ierr); CHKERRQ(ierr)
>
> call DMPlexCreate(PETSC_COMM_WORLD, dmA, ierr);CHKERRQ(ierr)
> call PetscObjectSetName(dmA, 'DMPlex_A', ierr);CHKERRQ(ierr)
>
> call DMPlexCreate(PETSC_COMM_WORLD, dmB, ierr);CHKERRQ(ierr)
> call PetscObjectSetName(dmB, 'DMPlex_B', ierr);CHKERRQ(ierr)
>
> call PetscObjectViewFromOptions(dmA, PETSC_NULL_CHARACTER, "-dmA",
> ierr); CHKERRQ(ierr)
> call PetscObjectViewFromOptions(dmB, PETSC_NULL_CHARACTER, "-dmB",
The second argument is a PETScObject, not a character string. This is what is causing the error. You should replace the PETSC_NULL_CHARACTER with PETSC_NULL_OBJECT in PETSc 3.7.x or earlier or with PETSC_NULL_VEC with the master branch development version of PETSc.
I have added more error checking to the branch barry/errorcheck-fortran-petscobjectviewfromoptions that will detect this error in the future.
Thanks for reporting the problem,
Barry
> ierr); CHKERRQ(ierr)
>
> call DMDestroy(dmA, ierr);CHKERRQ(ierr)
> call DMDestroy(dmB, ierr);CHKERRQ(ierr)
>
> call PetscFinalize(ierr)
>
> end program
>
> produces the options to be:
>
> -dmA-dmB
>
> -dmB
>
>
>
> While this works as expected when running with:
> ./example -dmA-dmB -dmB
>
> This is not intuitive.
>
> Is the hickup on my side or is it somewhere in the Fortran stubs?
>
> Please, let me know if you need more details on the build or if you
> cannot reproduce this.
>
> Many thanks,
>
> Fabian
>
>
>
> Petsc Development GIT revision: v3.7.6-3910-gd04c6f6 GIT Date:
> 2017-05-15 17:09:20 -0500
> ./configure \
> --with-cc=$(which mpicc) \
> --with-fc=$(which mpif90) \
> --with-cxx=$(which mpicxx) \
> --with-fortran \
> --with-fortran-interfaces \
> --with-shared-libraries=1 \
> --download-hdf5 \
> --download-triangle \
> --download-ctetgen \
> --with-cmake=$(which cmake) \
> --with-debugging=1 \
> COPTFLAGS='-O2 ' \
> FOPTFLAGS='-O2 ' \
> \
> && make all test
>
> GNU Fortran (Ubuntu 5.4.0-6ubuntu1~16.04.4) 5.4.0 20160609
> (Open MPI) 1.10.2
>
>
> Complete output of -help -info (Fortran Version):
> [0] petscinitialize_internal(): (Fortran):PETSc successfully started:
> procs 1
> [0] PetscGetHostName(): Rejecting domainname, likely is NIS
> met-ws-740m19.(none)
> [0] petscinitialize_internal(): Running on machine: met-ws-740m19
> ------Additional PETSc component options--------
> -log_exclude: <vec,mat,pc.ksp,snes>
> -info_exclude: <null,vec,mat,pc,ksp,snes,ts>
> -----------------------------------------------
> [0] PetscCommDuplicate(): Duplicating a communicator 47693199447680
> 11260976 max tags = 2147483647
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199447680 11260976
>
> -dmA-dmB ascii[:[filename][:[format][:append]]]: Prints object to
> stdout or ASCII file (PetscOptionsGetViewer)
> -dmA-dmB binary[:[filename][:[format][:append]]]: Saves object to a
> binary file (PetscOptionsGetViewer)
> -dmA-dmB draw[:drawtype[:filename]] Draws object (PetscOptionsGetViewer)
> -dmA-dmB socket[:port]: Pushes object to a Unix socket
> (PetscOptionsGetViewer)
> -dmA-dmB saws[:communicatorname]: Publishes object to SAWs
> (PetscOptionsGetViewer)
>
> DM Object: DMPlex_A 1 MPI processes
> type: plex
> [0] PetscCommDuplicate(): Duplicating a communicator 47693199449728
> 13799472 max tags = 2147483647
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199449728 13799472
> DMPlex_A in 0 dimensions:
> 0-cells: 0
>
> -dmB ascii[:[filename][:[format][:append]]]: Prints object to stdout
> or ASCII file (PetscOptionsGetViewer)
> -dmB binary[:[filename][:[format][:append]]]: Saves object to a
> binary file (PetscOptionsGetViewer)
> -dmB draw[:drawtype[:filename]] Draws object (PetscOptionsGetViewer)
> -dmB socket[:port]: Pushes object to a Unix socket
> (PetscOptionsGetViewer)
> -dmB saws[:communicatorname]: Publishes object to SAWs
> (PetscOptionsGetViewer)
>
> DM Object: DMPlex_B 1 MPI processes
> type: plex
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199449728 13799472
> [0] PetscCommDuplicate(): Using internal PETSc communicator
> 47693199449728 13799472
> DMPlex_B in 0 dimensions:
> 0-cells: 0
> [0] Petsc_DelComm_Inner(): Removing reference to PETSc communicator
> embedded in a user MPI_Comm 13799472
> [0] Petsc_DelComm_Outer(): User MPI_Comm 47693199449728 is being freed
> after removing reference from inner PETSc comm to this outer comm
> [0] PetscCommDestroy(): Deleting PETSc MPI_Comm 13799472
> [0] Petsc_DelCounter(): Deleting counter data in an MPI_Comm 13799472
> [0] PetscFinalize(): PetscFinalize() called
> [0] PetscGetHostName(): Rejecting domainname, likely is NIS
> met-ws-740m19.(none)
> [0] PetscFOpen(): Opening file Log.0
> [0] Petsc_DelViewer(): Removing viewer data attribute in an MPI_Comm
> 11260976
> [0] Petsc_DelComm_Inner(): Removing reference to PETSc communicator
> embedded in a user MPI_Comm 11260976
> [0] Petsc_DelComm_Outer(): User MPI_Comm 47693199447680 is being freed
> after removing reference from inner PETSc comm to this outer comm
> [0] PetscCommDestroy(): Deleting PETSc MPI_Comm 11260976
> [0] Petsc_DelCounter(): Deleting counter data in an MPI_Comm 11260976
>
More information about the petsc-users
mailing list