fewpjfiewjipfwe
parent
e89224a20e
commit
efdae14fae
@ -1,82 +0,0 @@
|
||||
#include <petscdm.h>
|
||||
#include <petscdmda.h>
|
||||
#include <petscsys.h>
|
||||
#include <petscviewerhdf5.h>
|
||||
|
||||
static char help[] = "Test to write HDF5 file from PETSc DMDA Vec.\n\n";
|
||||
|
||||
int main(int argc, char **argv) {
|
||||
DM da2D;
|
||||
PetscInt i, j, ixs, ixm, iys, iym;
|
||||
PetscViewer H5viewer;
|
||||
PetscScalar xm = -1.0, xp = 1.0;
|
||||
PetscScalar ym = -1.0, yp = 1.0;
|
||||
PetscScalar value = 1.0, dx, dy;
|
||||
PetscInt Nx = 40, Ny = 40;
|
||||
Vec gauss, input;
|
||||
PetscScalar **gauss_ptr;
|
||||
PetscReal norm;
|
||||
const char *vecname;
|
||||
|
||||
dx = (xp - xm) / (Nx - 1);
|
||||
dy = (yp - ym) / (Ny - 1);
|
||||
|
||||
/* Initialize the Petsc context */
|
||||
PetscFunctionBeginUser;
|
||||
PetscCall(PetscInitialize(&argc, &argv, NULL, help));
|
||||
PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, Nx, Ny,
|
||||
PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL, &da2D));
|
||||
PetscCall(DMSetFromOptions(da2D));
|
||||
PetscCall(DMSetUp(da2D));
|
||||
|
||||
/* Set the coordinates */
|
||||
PetscCall(DMDASetUniformCoordinates(da2D, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0));
|
||||
|
||||
/* Declare gauss as a DMDA component */
|
||||
PetscCall(DMCreateGlobalVector(da2D, &gauss));
|
||||
PetscCall(PetscObjectSetName((PetscObject)gauss, "pressure"));
|
||||
|
||||
/* Initialize vector gauss with a constant value (=1) */
|
||||
PetscCall(VecSet(gauss, value));
|
||||
|
||||
/* Get the coordinates of the corners for each process */
|
||||
PetscCall(DMDAGetCorners(da2D, &ixs, &iys, 0, &ixm, &iym, 0));
|
||||
|
||||
/* Build the gaussian profile (exp(-x^2-y^2)) */
|
||||
PetscCall(DMDAVecGetArray(da2D, gauss, &gauss_ptr));
|
||||
for (j = iys; j < iys + iym; j++) {
|
||||
for (i = ixs; i < ixs + ixm; i++)
|
||||
gauss_ptr[j][i] = PetscExpScalar(-(xm + i * dx) * (xm + i * dx) - (ym + j * dy) * (ym + j * dy));
|
||||
}
|
||||
PetscCall(DMDAVecRestoreArray(da2D, gauss, &gauss_ptr));
|
||||
|
||||
/* Create the HDF5 viewer */
|
||||
PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "gauss.h5", FILE_MODE_WRITE, &H5viewer));
|
||||
PetscCall(PetscViewerSetFromOptions(H5viewer));
|
||||
|
||||
/* Write the H5 file */
|
||||
PetscCall(VecView(gauss, H5viewer));
|
||||
|
||||
/* Close the viewer */
|
||||
PetscCall(PetscViewerDestroy(&H5viewer));
|
||||
|
||||
PetscCall(VecDuplicate(gauss, &input));
|
||||
PetscCall(PetscObjectGetName((PetscObject)gauss, &vecname));
|
||||
PetscCall(PetscObjectSetName((PetscObject)input, vecname));
|
||||
|
||||
/* Create the HDF5 viewer for reading */
|
||||
PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "gauss.h5", FILE_MODE_READ, &H5viewer));
|
||||
PetscCall(PetscViewerSetFromOptions(H5viewer));
|
||||
PetscCall(VecLoad(input, H5viewer));
|
||||
PetscCall(PetscViewerDestroy(&H5viewer));
|
||||
|
||||
PetscCall(VecAXPY(input, -1.0, gauss));
|
||||
PetscCall(VecNorm(input, NORM_2, &norm));
|
||||
PetscCheck(norm <= 1.e-6, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Vec read in does not match vector written out");
|
||||
|
||||
PetscCall(VecDestroy(&input));
|
||||
PetscCall(VecDestroy(&gauss));
|
||||
PetscCall(DMDestroy(&da2D));
|
||||
PetscCall(PetscFinalize());
|
||||
return 0;
|
||||
}
|
||||
@ -0,0 +1,12 @@
|
||||
#!/bin/bash
|
||||
#SBATCH --job-name=experiment-1
|
||||
#SBATCH --nodes=20
|
||||
#SBATCH --output=%x_%j.log
|
||||
|
||||
for i in {1..30}
|
||||
do
|
||||
echo "Iteration $i"
|
||||
mpirun -np $i ./build/arnoldi -l 100
|
||||
sleep 1
|
||||
done
|
||||
|
||||
@ -1,276 +0,0 @@
|
||||
static char help[] = "Basic vector routines.\n\n";
|
||||
|
||||
/*
|
||||
Include "petscvec.h" so that we can use vectors. Note that this file
|
||||
automatically includes:
|
||||
petscsys.h - base PETSc routines petscis.h - index sets
|
||||
petscviewer.h - viewers
|
||||
*/
|
||||
|
||||
#include <lapacke.h>
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
|
||||
#include <petscvec.h>
|
||||
|
||||
int example_lapack() {
|
||||
int n = 2; // Dimension of the matrix A
|
||||
int nrhs = 1; // Number of right-hand sides
|
||||
int lda = 2; // Leading dimension of A
|
||||
int ldb = 2; // Leading dimension of B
|
||||
int info; // Status output
|
||||
int ipiv[2]; // Pivot indices for the LU decomposition
|
||||
|
||||
// Define matrix A (row-major order)
|
||||
double A[4] = {3, 1, 1, 2};
|
||||
|
||||
// Define right-hand side vector b
|
||||
double b[2] = {9, 8};
|
||||
|
||||
// Call LAPACK dgesv to solve Ax = b
|
||||
info = LAPACKE_dgesv(LAPACK_ROW_MAJOR, n, nrhs, A, lda, ipiv, b, ldb);
|
||||
|
||||
// Check for success
|
||||
if (info == 0) {
|
||||
printf("Solution:\n");
|
||||
for (int i = 0; i < n; i++) {
|
||||
printf("x[%d] = %f\n", i, b[i]);
|
||||
}
|
||||
} else if (info < 0) {
|
||||
printf("Argument %d had an illegal value.\n", -info);
|
||||
} else {
|
||||
printf("Matrix is singular. Solution could not be computed.\n");
|
||||
}
|
||||
|
||||
return info;
|
||||
}
|
||||
|
||||
int main(int argc, char **argv) {
|
||||
printf("LAPACKE Example:\n");
|
||||
example_lapack();
|
||||
|
||||
printf("PETSC Example:\n");
|
||||
Vec x, y, w; /* vectors */
|
||||
Vec *z; /* array of vectors */
|
||||
PetscReal norm, v, v1, v2, maxval;
|
||||
PetscInt n = 20, maxind;
|
||||
PetscScalar one = 1.0, two = 2.0, three = 3.0, dots[3], dot;
|
||||
|
||||
PetscFunctionBeginUser;
|
||||
PetscCall(PetscInitialize(&argc, &argv, NULL, help));
|
||||
PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
|
||||
|
||||
/*
|
||||
Create a vector, specifying only its global dimension.
|
||||
When using VecCreate(), VecSetSizes() and VecSetFromOptions(), the vector
|
||||
format (currently parallel, shared, or sequential) is determined at
|
||||
runtime. Also, the parallel partitioning of the vector is determined by
|
||||
PETSc at runtime.
|
||||
*/
|
||||
PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
|
||||
PetscCall(VecSetSizes(x, PETSC_DECIDE, n));
|
||||
PetscCall(VecSetFromOptions(x));
|
||||
|
||||
/*
|
||||
Duplicate some work vectors (of the same format and
|
||||
partitioning as the initial vector).
|
||||
*/
|
||||
PetscCall(VecDuplicate(x, &y));
|
||||
PetscCall(VecDuplicate(x, &w));
|
||||
|
||||
/*
|
||||
Duplicate more work vectors (of the same format and
|
||||
partitioning as the initial vector). Here we duplicate
|
||||
an array of vectors, which is often more convenient than
|
||||
duplicating individual ones.
|
||||
*/
|
||||
PetscCall(VecDuplicateVecs(x, 3, &z));
|
||||
/*
|
||||
Set the vectors to entries to a constant value.
|
||||
*/
|
||||
PetscCall(VecSet(x, one));
|
||||
PetscCall(VecSet(y, two));
|
||||
PetscCall(VecSet(z[0], one));
|
||||
PetscCall(VecSet(z[1], two));
|
||||
PetscCall(VecSet(z[2], three));
|
||||
/*
|
||||
Demonstrate various basic vector routines.
|
||||
*/
|
||||
PetscCall(VecDot(x, y, &dot));
|
||||
PetscCall(VecMDot(x, 3, z, dots));
|
||||
|
||||
/*
|
||||
Note: If using a complex numbers version of PETSc, then
|
||||
PETSC_USE_COMPLEX is defined in the makefiles; otherwise,
|
||||
(when using real numbers) it is undefined.
|
||||
*/
|
||||
|
||||
PetscCall(
|
||||
PetscPrintf(PETSC_COMM_WORLD, "Vector length %" PetscInt_FMT "\n", n));
|
||||
|
||||
PetscCall(VecMax(x, &maxind, &maxval));
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD,
|
||||
"VecMax %g, VecInd %" PetscInt_FMT "\n",
|
||||
(double)maxval, maxind));
|
||||
|
||||
PetscCall(VecMin(x, &maxind, &maxval));
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD,
|
||||
"VecMin %g, VecInd %" PetscInt_FMT "\n",
|
||||
(double)maxval, maxind));
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD,
|
||||
"All other values should be near zero\n"));
|
||||
|
||||
PetscCall(VecScale(x, two));
|
||||
PetscCall(VecNorm(x, NORM_2, &norm));
|
||||
v = norm - 2.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecScale %g\n", (double)v));
|
||||
|
||||
PetscCall(VecCopy(x, w));
|
||||
PetscCall(VecNorm(w, NORM_2, &norm));
|
||||
v = norm - 2.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecCopy %g\n", (double)v));
|
||||
|
||||
PetscCall(VecAXPY(y, three, x));
|
||||
PetscCall(VecNorm(y, NORM_2, &norm));
|
||||
v = norm - 8.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecAXPY %g\n", (double)v));
|
||||
|
||||
PetscCall(VecAYPX(y, two, x));
|
||||
PetscCall(VecNorm(y, NORM_2, &norm));
|
||||
v = norm - 18.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecAYPX %g\n", (double)v));
|
||||
|
||||
PetscCall(VecSwap(x, y));
|
||||
PetscCall(VecNorm(y, NORM_2, &norm));
|
||||
v = norm - 2.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecSwap %g\n", (double)v));
|
||||
PetscCall(VecNorm(x, NORM_2, &norm));
|
||||
v = norm - 18.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecSwap %g\n", (double)v));
|
||||
|
||||
PetscCall(VecWAXPY(w, two, x, y));
|
||||
PetscCall(VecNorm(w, NORM_2, &norm));
|
||||
v = norm - 38.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecWAXPY %g\n", (double)v));
|
||||
|
||||
PetscCall(VecPointwiseMult(w, y, x));
|
||||
PetscCall(VecNorm(w, NORM_2, &norm));
|
||||
v = norm - 36.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(
|
||||
PetscPrintf(PETSC_COMM_WORLD, "VecPointwiseMult %g\n", (double)v));
|
||||
|
||||
PetscCall(VecPointwiseDivide(w, x, y));
|
||||
PetscCall(VecNorm(w, NORM_2, &norm));
|
||||
v = norm - 9.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(
|
||||
PetscPrintf(PETSC_COMM_WORLD, "VecPointwiseDivide %g\n", (double)v));
|
||||
|
||||
PetscCall(VecSetValue(y, 0, 0.0, INSERT_VALUES));
|
||||
PetscCall(VecAssemblyBegin(y));
|
||||
PetscCall(VecAssemblyEnd(y));
|
||||
PetscCall(VecPointwiseDivide(w, x, y));
|
||||
PetscCall(VecNorm(w, NORM_2, &norm));
|
||||
v = norm - 9.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(
|
||||
PetscPrintf(PETSC_COMM_WORLD, "VecPointwiseDivide %g\n", (double)v));
|
||||
|
||||
dots[0] = one;
|
||||
dots[1] = three;
|
||||
dots[2] = two;
|
||||
|
||||
PetscCall(VecSet(x, one));
|
||||
PetscCall(VecMAXPY(x, 3, dots, z));
|
||||
PetscCall(VecNorm(z[0], NORM_2, &norm));
|
||||
v = norm - PetscSqrtReal((PetscReal)n);
|
||||
if (v > -PETSC_SMALL && v < PETSC_SMALL)
|
||||
v = 0.0;
|
||||
PetscCall(VecNorm(z[1], NORM_2, &norm));
|
||||
v1 = norm - 2.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v1 > -PETSC_SMALL && v1 < PETSC_SMALL)
|
||||
v1 = 0.0;
|
||||
PetscCall(VecNorm(z[2], NORM_2, &norm));
|
||||
v2 = norm - 3.0 * PetscSqrtReal((PetscReal)n);
|
||||
if (v2 > -PETSC_SMALL && v2 < PETSC_SMALL)
|
||||
v2 = 0.0;
|
||||
PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecMAXPY %g %g %g \n", (double)v,
|
||||
(double)v1, (double)v2));
|
||||
|
||||
/*
|
||||
Free work space. All PETSc objects should be destroyed when they
|
||||
are no longer needed.
|
||||
*/
|
||||
PetscCall(VecDestroy(&x));
|
||||
PetscCall(VecDestroy(&y));
|
||||
PetscCall(VecDestroy(&w));
|
||||
PetscCall(VecDestroyVecs(3, &z));
|
||||
PetscCall(PetscFinalize());
|
||||
return 0;
|
||||
}
|
||||
|
||||
/*TEST
|
||||
|
||||
testset:
|
||||
output_file: output/ex1_1.out
|
||||
# This is a test where the exact numbers are critical
|
||||
diff_args: -j
|
||||
|
||||
test:
|
||||
|
||||
test:
|
||||
suffix: cuda
|
||||
args: -vec_type cuda
|
||||
requires: cuda
|
||||
|
||||
test:
|
||||
suffix: kokkos
|
||||
args: -vec_type kokkos
|
||||
requires: kokkos_kernels
|
||||
|
||||
test:
|
||||
suffix: hip
|
||||
args: -vec_type hip
|
||||
requires: hip
|
||||
|
||||
test:
|
||||
suffix: 2
|
||||
nsize: 2
|
||||
|
||||
test:
|
||||
suffix: 2_cuda
|
||||
nsize: 2
|
||||
args: -vec_type cuda
|
||||
requires: cuda
|
||||
|
||||
test:
|
||||
suffix: 2_kokkos
|
||||
nsize: 2
|
||||
args: -vec_type kokkos
|
||||
requires: kokkos_kernels
|
||||
|
||||
test:
|
||||
suffix: 2_hip
|
||||
nsize: 2
|
||||
args: -vec_type hip
|
||||
requires: hip
|
||||
|
||||
TEST*/
|
||||
Loading…
Reference in New Issue