1 /*
2 Implements the sequential vectors.
3 */
4
5 #include <../src/vec/vec/impls/dvecimpl.h> /*I "petscvec.h" I*/
6 /*MC
7 VECSEQ - VECSEQ = "seq" - The basic sequential vector
8
9 Options Database Keys:
10 . -vec_type seq - sets the vector type to VECSEQ during a call to VecSetFromOptions()
11
12 Level: beginner
13
14 .seealso: `VecCreate()`, `VecSetType()`, `VecSetFromOptions()`, `VecCreateSeqWithArray()`, `VECMPI`, `VecType`, `VecCreateMPI()`, `VecCreateSeq()`
15 M*/
16
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
21
VecCreate_Seq(Vec V)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;
28
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(PetscShmgetAllocateArray(n, sizeof(PetscScalar), (void **)&array));
34 PetscCall(PetscArrayzero(array, n));
35 PetscCall(VecCreate_Seq_Private(V, array));
36
37 s = (Vec_Seq *)V->data;
38 s->array_allocated = array;
39 PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_2], 0));
40 PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_1], 0));
41 PetscCall(PetscObjectComposedDataSetReal((PetscObject)V, NormIds[NORM_INFINITY], 0));
42
43 #else
44 switch (((PetscObject)V)->precision) {
45 case PETSC_PRECISION_SINGLE: {
46 float *aarray;
47
48 PetscCall(PetscCalloc1(n, &aarray));
49 PetscCall(VecCreate_Seq_Private(V, aarray));
50
51 s = (Vec_Seq *)V->data;
52 s->array_allocated = (PetscScalar *)aarray;
53 } break;
54 case PETSC_PRECISION_DOUBLE: {
55 double *aarray;
56
57 PetscCall(PetscCalloc1(n, &aarray));
58 PetscCall(VecCreate_Seq_Private(V, aarray));
59
60 s = (Vec_Seq *)V->data;
61 s->array_allocated = (PetscScalar *)aarray;
62 } break;
63 default:
64 SETERRQ(PetscObjectComm((PetscObject)V), PETSC_ERR_SUP, "No support for mixed precision %d", (int)(((PetscObject)V)->precision));
65 }
66 #endif
67 PetscFunctionReturn(PETSC_SUCCESS);
68 }
69