Actual source code: bvec3.c

  1: /*
  2:    Implements the sequential vectors.
  3: */

  5: #include <../src/vec/vec/impls/dvecimpl.h>
  6: /*MC
  7:    VECSEQ - VECSEQ = "seq" - The basic sequential vector

  9:    Options Database Keys:
 10: . -vec_type seq - sets the vector type to VECSEQ during a call to VecSetFromOptions()

 12:   Level: beginner

 14: .seealso: `VecCreate()`, `VecSetType()`, `VecSetFromOptions()`, `VecCreateSeqWithArray()`, `VECMPI`, `VecType`, `VecCreateMPI()`, `VecCreateSeq()`
 15: M*/

 17: #if defined(PETSC_USE_MIXED_PRECISION)
 18: extern PetscErrorCode VecCreate_Seq_Private(Vec, const float *);
 19: extern PetscErrorCode VecCreate_Seq_Private(Vec, const double *);
 20: #endif

 22: PetscErrorCode VecCreate_Seq(Vec V)
 23: {
 24:   Vec_Seq     *s;
 25:   PetscScalar *array;
 26:   PetscInt     n = PetscMax(V->map->n, V->map->N);
 27:   PetscMPIInt  size;

 29:   PetscFunctionBegin;
 30:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)V), &size));
 31:   PetscCheck(size <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cannot create VECSEQ on more than one process");
 32: #if !defined(PETSC_USE_MIXED_PRECISION)
 33:   PetscCall(PetscCalloc1(n, &array));
 34:   PetscCall(VecCreate_Seq_Private(V, array));

 36:   s                  = (Vec_Seq *)V->data;
 37:   s->array_allocated = array;
 38:   PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_2], 0));
 39:   PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_1], 0));
 40:   PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_INFINITY], 0));

 42: #else
 43:   switch (((PetscObject)V)->precision) {
 44:   case PETSC_PRECISION_SINGLE: {
 45:     float *aarray;

 47:     PetscCall(PetscCalloc1(n, &aarray));
 48:     PetscCall(VecCreate_Seq_Private(V, aarray));

 50:     s                  = (Vec_Seq *)V->data;
 51:     s->array_allocated = (PetscScalar *)aarray;
 52:   } break;
 53:   case PETSC_PRECISION_DOUBLE: {
 54:     double *aarray;

 56:     PetscCall(PetscCalloc1(n, &aarray));
 57:     PetscCall(VecCreate_Seq_Private(V, aarray));

 59:     s                  = (Vec_Seq *)V->data;
 60:     s->array_allocated = (PetscScalar *)aarray;
 61:   } break;
 62:   default:
 63:     SETERRQ(PetscObjectComm((PetscObject)V), PETSC_ERR_SUP, "No support for mixed precision %d", (int)(((PetscObject)V)->precision));
 64:   }
 65: #endif
 66:   PetscFunctionReturn(PETSC_SUCCESS);
 67: }