1c4762a1bSJed Brown static char help[] = "Reads a PETSc matrix and vector from a file and reorders it.\n\
2c4762a1bSJed Brown -f0 <input_file> : first file to load (small system)\n\
3c4762a1bSJed Brown -f1 <input_file> : second file to load (larger system)\n\n";
4c4762a1bSJed Brown
5c4762a1bSJed Brown /*
6c4762a1bSJed Brown Include "petscmat.h" so that we can use matrices.
7c4762a1bSJed Brown automatically includes:
8c4762a1bSJed Brown petscsys.h - base PETSc routines petscvec.h - vectors
9c4762a1bSJed Brown petscmat.h - matrices
10c4762a1bSJed Brown petscis.h - index sets petscviewer.h - viewers
11c4762a1bSJed Brown */
12c4762a1bSJed Brown #include <petscmat.h>
13c4762a1bSJed Brown
main(int argc,char ** args)14d71ae5a4SJacob Faibussowitsch int main(int argc, char **args)
15d71ae5a4SJacob Faibussowitsch {
16c4762a1bSJed Brown Mat A; /* matrix */
17c4762a1bSJed Brown PetscViewer fd; /* viewer */
18c4762a1bSJed Brown char file[PETSC_MAX_PATH_LEN]; /* input file name */
19c4762a1bSJed Brown IS isrow, iscol; /* row and column permutations */
20c4762a1bSJed Brown MatOrderingType rtype = MATORDERINGRCM;
21c4762a1bSJed Brown PetscBool flg;
22c4762a1bSJed Brown
23327415f7SBarry Smith PetscFunctionBeginUser;
24*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help));
25c4762a1bSJed Brown /*
26c4762a1bSJed Brown Determine files from which we read the two linear systems
27c4762a1bSJed Brown (matrix and right-hand-side vector).
28c4762a1bSJed Brown */
299566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
3028b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate binary file with the -f option");
31c4762a1bSJed Brown
32c4762a1bSJed Brown /* -----------------------------------------------------------
33c4762a1bSJed Brown Beginning of loop
34c4762a1bSJed Brown ----------------------------------------------------------- */
35c4762a1bSJed Brown /*
36c4762a1bSJed Brown Loop through the reordering 2 times.
37c4762a1bSJed Brown - The intention here is to preload and solve a small system;
38c4762a1bSJed Brown then load another (larger) system and solve it as well.
39c4762a1bSJed Brown This process preloads the instructions with the smaller
40c4762a1bSJed Brown system so that more accurate performance monitoring (via
41c4762a1bSJed Brown -log_view) can be done with the larger one (that actually
42c4762a1bSJed Brown is the system of interest).
43c4762a1bSJed Brown */
44c4762a1bSJed Brown
45c4762a1bSJed Brown /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
46c4762a1bSJed Brown Load system i
47c4762a1bSJed Brown - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
48c4762a1bSJed Brown
49c4762a1bSJed Brown /*
50c4762a1bSJed Brown Open binary file. Note that we use FILE_MODE_READ to indicate
51c4762a1bSJed Brown reading from this file.
52c4762a1bSJed Brown */
539566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
54c4762a1bSJed Brown
55c4762a1bSJed Brown /*
56c4762a1bSJed Brown Load the matrix; then destroy the viewer.
57c4762a1bSJed Brown */
589566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
599566063dSJacob Faibussowitsch PetscCall(MatSetType(A, MATSEQAIJ));
609566063dSJacob Faibussowitsch PetscCall(MatLoad(A, fd));
619566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&fd));
62c4762a1bSJed Brown
639566063dSJacob Faibussowitsch PetscCall(MatGetOrdering(A, rtype, &isrow, &iscol));
649566063dSJacob Faibussowitsch PetscCall(ISView(isrow, PETSC_VIEWER_STDOUT_WORLD));
65c4762a1bSJed Brown
66c4762a1bSJed Brown /*
67c4762a1bSJed Brown Free work space. All PETSc objects should be destroyed when they
68c4762a1bSJed Brown are no longer needed.
69c4762a1bSJed Brown */
709566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A));
719566063dSJacob Faibussowitsch PetscCall(ISDestroy(&isrow));
729566063dSJacob Faibussowitsch PetscCall(ISDestroy(&iscol));
73c4762a1bSJed Brown
749566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
75b122ec5aSJacob Faibussowitsch return 0;
76c4762a1bSJed Brown }
77c4762a1bSJed Brown
78c4762a1bSJed Brown /*TEST
79c4762a1bSJed Brown
80c4762a1bSJed Brown test:
81dfd57a17SPierre Jolivet requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
82c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/medium -viewer_binary_skip_info
83c4762a1bSJed Brown
84c4762a1bSJed Brown TEST*/
85