Actual source code: vecs.c


  2: #include <petscvec.h>

  4: PetscErrorCode VecsDestroy(Vecs x)
  5: {
  6:   PetscFunctionBegin;
  7:   PetscCall(VecDestroy(&(x)->v));
  8:   PetscCall(PetscFree(x));
  9:   PetscFunctionReturn(PETSC_SUCCESS);
 10: }

 12: PetscErrorCode VecsCreateSeq(MPI_Comm comm, PetscInt p, PetscInt m, Vecs *x)
 13: {
 14:   PetscFunctionBegin;
 15:   PetscCall(PetscNew(x));
 16:   PetscCall(VecCreateSeq(comm, p * m, &(*x)->v));
 17:   (*x)->n = m;
 18:   PetscFunctionReturn(PETSC_SUCCESS);
 19: }

 21: PetscErrorCode VecsCreateSeqWithArray(MPI_Comm comm, PetscInt p, PetscInt m, PetscScalar *a, Vecs *x)
 22: {
 23:   PetscFunctionBegin;
 24:   PetscCall(PetscNew(x));
 25:   PetscCall(VecCreateSeqWithArray(comm, 1, p * m, a, &(*x)->v));
 26:   (*x)->n = m;
 27:   PetscFunctionReturn(PETSC_SUCCESS);
 28: }

 30: PetscErrorCode VecsDuplicate(Vecs x, Vecs *y)
 31: {
 32:   PetscFunctionBegin;
 33:   PetscCall(PetscNew(y));
 34:   PetscCall(VecDuplicate(x->v, &(*y)->v));
 35:   (*y)->n = x->n;
 36:   PetscFunctionReturn(PETSC_SUCCESS);
 37: }