Actual source code: bvec3.c
petsc-3.11.4 2019-09-28
2: /*
3: Implements the sequential vectors.
4: */
6: #include <../src/vec/vec/impls/dvecimpl.h>
7: /*MC
8: VECSEQ - VECSEQ = "seq" - The basic sequential vector
10: Options Database Keys:
11: . -vec_type seq - sets the vector type to VECSEQ during a call to VecSetFromOptions()
13: Level: beginner
15: .seealso: VecCreate(), VecSetType(), VecSetFromOptions(), VecCreateSeqWithArray(), VECMPI, VecType, VecCreateMPI(), VecCreateSeq()
16: M*/
18: #if defined(PETSC_USE_MIXED_PRECISION)
19: extern PetscErrorCode VecCreate_Seq_Private(Vec,const float*);
20: extern PetscErrorCode VecCreate_Seq_Private(Vec,const double*);
21: #endif
23: PETSC_EXTERN PetscErrorCode VecCreate_Seq(Vec V)
24: {
25: Vec_Seq *s;
26: PetscScalar *array;
28: PetscInt n = PetscMax(V->map->n,V->map->N);
29: PetscMPIInt size;
32: MPI_Comm_size(PetscObjectComm((PetscObject)V),&size);
33: if (size > 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Cannot create VECSEQ on more than one process");
34: #if !defined(PETSC_USE_MIXED_PRECISION)
35: PetscMalloc1(n,&array);
36: PetscLogObjectMemory((PetscObject)V, n*sizeof(PetscScalar));
37: VecCreate_Seq_Private(V,array);
39: s = (Vec_Seq*)V->data;
40: s->array_allocated = array;
42: VecSet(V,0.0);
43: #else
44: switch (((PetscObject)V)->precision) {
45: case PETSC_PRECISION_SINGLE: {
46: float *aarray;
48: PetscMalloc1(n,&aarray);
49: PetscLogObjectMemory((PetscObject)V, n*sizeof(float));
50: PetscMemzero(aarray,n*sizeof(float));
51: VecCreate_Seq_Private(V,aarray);
53: s = (Vec_Seq*)V->data;
54: s->array_allocated = (PetscScalar*)aarray;
55: } break;
56: case PETSC_PRECISION_DOUBLE: {
57: double *aarray;
59: PetscMalloc1(n,&aarray);
60: PetscLogObjectMemory((PetscObject)V, n*sizeof(double));
61: PetscMemzero(aarray,n*sizeof(double));
62: VecCreate_Seq_Private(V,aarray);
64: s = (Vec_Seq*)V->data;
65: s->array_allocated = (PetscScalar*)aarray;
66: } break;
67: default: SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"No support for mixed precision %d",(int)(((PetscObject)V)->precision));
68: }
69: #endif
70: return(0);
71: }