11690c2aeSBarry Smith static char help[] = "Tests that MatView() and MatLoad() work for MPIAIJ matrix with total nz > PETSC_INT_MAX\n\n"; 238b83642SBarry Smith 338b83642SBarry Smith #include <petscmat.h> 438b83642SBarry Smith 538b83642SBarry Smith int main(int argc, char **args) 638b83642SBarry Smith { 71690c2aeSBarry Smith PetscInt n = 1000000, nzr = (PetscInt)((double)PETSC_INT_MAX) / (3.8 * n); 838b83642SBarry Smith Mat A; 938b83642SBarry Smith PetscScalar *a; 1038b83642SBarry Smith PetscInt *ii, *jd, *jo; 1138b83642SBarry Smith PetscMPIInt rank, size; 1238b83642SBarry Smith PetscViewer viewer; 1338b83642SBarry Smith 1438b83642SBarry Smith PetscFunctionBeginUser; 15*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help)); 1638b83642SBarry Smith 1738b83642SBarry Smith PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 1838b83642SBarry Smith PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run with two MPI ranks"); 1938b83642SBarry Smith PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 2038b83642SBarry Smith PetscCall(PetscMalloc4(n + 1, &ii, n * nzr, &jd, n * nzr, &jo, n * nzr, &a)); 2138b83642SBarry Smith ii[0] = 0; 2238b83642SBarry Smith for (PetscInt i = 0; i < n; i++) { 2338b83642SBarry Smith ii[i + 1] = ii[i] + nzr; 2438b83642SBarry Smith for (PetscInt j = 0; j < nzr; j++) jd[ii[i] + j] = j; 2538b83642SBarry Smith if (rank == 0) { 2638b83642SBarry Smith for (PetscInt j = 0; j < nzr; j++) jo[ii[i] + j] = n + j - 1; 2738b83642SBarry Smith } else { 2838b83642SBarry Smith for (PetscInt j = 0; j < nzr; j++) jo[ii[i] + j] = j; 2938b83642SBarry Smith } 3038b83642SBarry Smith } 3138b83642SBarry Smith PetscCall(MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD, n, n, PETSC_DETERMINE, PETSC_DETERMINE, ii, jd, a, ii, jo, a, &A)); 3238b83642SBarry Smith PetscCall(MatView(A, PETSC_VIEWER_BINARY_WORLD)); 3338b83642SBarry Smith PetscCall(MatDestroy(&A)); 3438b83642SBarry Smith PetscCall(PetscFree4(ii, jd, jo, a)); 3538b83642SBarry Smith 3638b83642SBarry Smith PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "binaryoutput", FILE_MODE_READ, &viewer)); 3738b83642SBarry Smith PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 3838b83642SBarry Smith PetscCall(MatLoad(A, viewer)); 3938b83642SBarry Smith PetscCall(PetscViewerDestroy(&viewer)); 4038b83642SBarry Smith PetscCall(MatDestroy(&A)); 4138b83642SBarry Smith PetscCall(PetscFinalize()); 4238b83642SBarry Smith return 0; 4338b83642SBarry Smith } 4438b83642SBarry Smith 4538b83642SBarry Smith /*TEST 4638b83642SBarry Smith 4738b83642SBarry Smith test: 4838b83642SBarry Smith TODO: requires to much memory to run in the CI 4938b83642SBarry Smith nsize: 2 5038b83642SBarry Smith 5138b83642SBarry Smith TEST*/ 52