petsc4py.PETSc.Vec#
- class petsc4py.PETSc.Vec#
Bases:
Object
A vector object.
See also
Enumerations
Vector assembly option.
The vector type.
Methods Summary
abs
()Replace each entry (xₙ) in the vector by abs|xₙ|.
appendOptionsPrefix
(prefix)Append to the prefix used for searching for options in the database.
assemble
()Assemble the vector.
Begin an assembling stage of the vector.
Finish the assembling stage initiated with
assemblyBegin
.attachDLPackInfo
([vec, dltensor])Attach tensor information from another vector or DLPack tensor.
axpby
(alpha, beta, x)Compute and store y = ɑ·x + β·y.
axpy
(alpha, x)Compute and store y = ɑ·x + y.
aypx
(alpha, x)Compute and store y = x + ɑ·y.
bindToCPU
(flg)Bind vector operations execution on the CPU.
Return whether the vector has been bound to the CPU.
chop
(tol)Set all vector entries less than some absolute tolerance to zero.
Clear tensor information.
concatenate
(vecs)Concatenate vectors into a single vector.
Conjugate the vector.
copy
([result])Return a copy of the vector.
create
([comm])Create a vector object.
createCUDAWithArrays
([cpuarray, cudahandle, ...])Create a
Type.CUDA
vector with optional arrays.createGhost
(ghosts, size[, bsize, comm])Create a parallel vector with ghost padding on each processor.
createGhostWithArray
(ghosts, array[, size, ...])Create a parallel vector with ghost padding and provided arrays.
createHIPWithArrays
([cpuarray, hiphandle, ...])Create a
Type.HIP
vector with optional arrays.Create a local vector.
createMPI
(size[, bsize, comm])Create a parallel
Type.MPI
vector.createNest
(vecs[, isets, comm])Create a
Type.NEST
vector containing multiple nested subvectors.createSeq
(size[, bsize, comm])Create a sequential
Type.SEQ
vector.createShared
(size[, bsize, comm])Create a
Type.SHARED
vector that uses shared memory.createViennaCLWithArrays
([cpuarray, ...])Create a
Type.VIENNACL
vector with optional arrays.createWithArray
(array[, size, bsize, comm])Create a vector using a provided array.
createWithDLPack
(dltensor[, size, bsize, comm])Create a vector wrapping a DLPack object, sharing the same memory.
destroy
()Destroy the vector.
dot
(vec)Return the dot product with
vec
.dotBegin
(vec)Begin computing the dot product.
dotEnd
(vec)Finish computing the dot product initiated with
dotBegin
.dotNorm2
(vec)Return the dot product with
vec
and its squared norm.duplicate
([array])Create a new vector with the same type, optionally with data.
equal
(vec)Return whether the vector is equal to another.
exp
()Replace each entry (xₙ) in the vector by exp(xₙ).
getArray
([readonly])Return local portion of the vector as an
ndarray
.Return the block size of the vector.
getBuffer
([readonly])Return a buffered view of the local portion of the vector.
Return the OpenCL context associated with the vector.
getCLMemHandle
([mode])Return the OpenCL buffer associated with the vector.
Return the OpenCL command queue associated with the vector.
getCUDAHandle
([mode])Return a pointer to the device buffer.
getDM
()Return the
DM
associated to the vector.Return ghosting indices of a ghost vector.
getHIPHandle
([mode])Return a pointer to the device buffer.
getLGMap
()Return the local-to-global mapping.
Return the local size of the vector.
getLocalVector
(lvec[, readonly])Maps the local portion of the vector into a local vector.
Return all the vectors contained in the nested vector.
Return the offloading status of the vector.
Return the prefix used for searching for options in the database.
Return the locally owned range of indices
(start, end)
.Return the range of indices owned by each process.
getSize
()Return the global size of the vector.
getSizes
()Return the vector sizes.
getSubVector
(iset[, subvec])Return a subvector from given indices.
getType
()Return the type of the vector.
getValue
(index)Return a single value from the vector.
getValues
(indices[, values])Return values from certain locations in the vector.
getValuesStagStencil
(indices[, values])Not implemented.
ghostUpdate
([addv, mode])Update ghosted vector entries.
ghostUpdateBegin
([addv, mode])Begin updating ghosted vector entries.
ghostUpdateEnd
([addv, mode])Finish updating ghosted vector entries initiated with
ghostUpdateBegin
.isaxpy
(idx, alpha, x)Add a scaled reduced-space vector to a subset of the vector.
isset
(idx, alpha)Set specific elements of the vector to the same value.
load
(viewer)Load a vector.
Return a context manager for viewing ghost vectors in local form.
log
()Replace each entry in the vector by its natural logarithm.
mDot
(vecs[, out])Compute Xᴴ·y with X an array of vectors.
mDotBegin
(vecs, out)Starts a split phase multiple dot product computation.
mDotEnd
(vecs, out)Ends a split phase multiple dot product computation.
max
()Return the vector entry with maximum real part and its location.
maxPointwiseDivide
(vec)Return the maximum of the component-wise absolute value division.
maxpy
(alphas, vecs)Compute and store y = Σₙ(ɑₙ·Xₙ) + y with X an array of vectors.
min
()Return the vector entry with minimum real part and its location.
mtDot
(vecs[, out])Compute Xᵀ·y with X an array of vectors.
mtDotBegin
(vecs, out)Starts a split phase transpose multiple dot product computation.
mtDotEnd
(vecs, out)Ends a split phase transpose multiple dot product computation.
norm
([norm_type])Compute the vector norm.
normBegin
([norm_type])Begin computing the vector norm.
normEnd
([norm_type])Finish computations initiated with
normBegin
.Normalize the vector by its 2-norm.
permute
(order[, invert])Permute the vector in-place with a provided ordering.
placeArray
(array)Set the local portion of the vector to a provided array.
pointwiseDivide
(x, y)Compute and store the component-wise division of two vectors.
pointwiseMax
(x, y)Compute and store the component-wise maximum of two vectors.
pointwiseMaxAbs
(x, y)Compute and store the component-wise maximum absolute values.
pointwiseMin
(x, y)Compute and store the component-wise minimum of two vectors.
pointwiseMult
(x, y)Compute and store the component-wise multiplication of two vectors.
Replace each entry in the vector by its reciprocal.
resetArray
([force])Reset the vector to use its default array.
Restore a pointer to the OpenCL buffer obtained with
getCLMemHandle
.restoreCUDAHandle
(handle[, mode])Restore a pointer to the device buffer obtained with
getCUDAHandle
.restoreHIPHandle
(handle[, mode])Restore a pointer to the device buffer obtained with
getHIPHandle
.restoreLocalVector
(lvec[, readonly])Unmap a local access obtained with
getLocalVector
.restoreSubVector
(iset, subvec)Restore a subvector extracted using
getSubVector
.scale
(alpha)Scale all entries of the vector.
set
(alpha)Set all components of the vector to the same value.
setArray
(array)Set values for the local portion of the vector.
setBlockSize
(bsize)Set the block size of the vector.
setDM
(dm)Associate a
DM
to the vector.Configure the vector from the options database.
setLGMap
(lgmap)Set the local-to-global mapping.
setMPIGhost
(ghosts)Set the ghost points for a ghosted vector.
setNestSubVecs
(sx[, idxm])Set the component vectors at specified indices in the nested vector.
setOption
(option, flag)Set option.
setOptionsPrefix
(prefix)Set the prefix used for searching for options in the database.
setRandom
([random])Set all components of the vector to random numbers.
setSizes
(size[, bsize])Set the local and global sizes of the vector.
setType
(vec_type)Set the vector type.
setUp
()Set up the internal data structures for using the vector.
setValue
(index, value[, addv])Insert or add a single value in the vector.
setValueLocal
(index, value[, addv])Insert or add a single value in the vector using a local numbering.
setValues
(indices, values[, addv])Insert or add multiple values in the vector.
setValuesBlocked
(indices, values[, addv])Insert or add blocks of values in the vector.
setValuesBlockedLocal
(indices, values[, addv])Insert or add blocks of values in the vector with a local numbering.
setValuesLocal
(indices, values[, addv])Insert or add multiple values in the vector with a local numbering.
setValuesStagStencil
(indices, values[, addv])Not implemented.
shift
(alpha)Shift all entries in the vector.
sqrtabs
()Replace each entry (xₙ) in the vector by √|xₙ|.
strideGather
(field, vec[, addv])Insert component values into a single-component vector.
strideMax
(field)Return the maximum of entries in a subvector.
strideMin
(field)Return the minimum of entries in a subvector.
strideNorm
(field[, norm_type])Return the norm of entries in a subvector.
strideScale
(field, alpha)Scale a component of the vector.
strideScatter
(field, vec[, addv])Scatter entries into a component of another vector.
strideSum
(field)Sum subvector entries.
sum
()Return the sum of all the entries of the vector.
swap
(vec)Swap the content of two vectors.
tDot
(vec)Return the indefinite dot product with
vec
.tDotBegin
(vec)Begin computing the indefinite dot product.
tDotEnd
(vec)Finish computing the indefinite dot product initiated with
tDotBegin
.toDLPack
([mode])Return a DLPack
PyCapsule
wrapping the vector data.view
([viewer])Display the vector.
waxpy
(alpha, x, y)Compute and store w = ɑ·x + y.
Set all entries in the vector to zero.
Attributes Summary
Alias for
array_w
.Read-only
ndarray
containing the local portion of the vector.Writeable
ndarray
containing the local portion of the vector.The block size.
Alias for
buffer_w
.Read-only buffered view of the local portion of the vector.
Writeable buffered view of the local portion of the vector.
The local vector size.
The locally owned range of indices in the form
[low, high)
.The range of indices owned by each process.
The global vector size.
The local and global vector sizes.
Methods Documentation
- abs()#
Replace each entry (xₙ) in the vector by abs|xₙ|.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2293
- Return type:
- appendOptionsPrefix(prefix)#
Append to the prefix used for searching for options in the database.
Logically collective.
- assemble()#
Assemble the vector.
Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:3047
- Return type:
- assemblyBegin()#
Begin an assembling stage of the vector.
Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:3023
- Return type:
- assemblyEnd()#
Finish the assembling stage initiated with
assemblyBegin
.Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:3035
- Return type:
- attachDLPackInfo(vec=None, dltensor=None)#
Attach tensor information from another vector or DLPack tensor.
Logically collective.
This tensor information is required when converting a
Vec
to a DLPack object.- Parameters:
vec (Vec | None) – Vector with attached tensor information. This is typically created by calling
createWithDLPack
.dltensor – DLPack tensor. This will only be used if
vec
isNone
.
- Return type:
Notes
This operation does not copy any data from
vec
ordltensor
.See also
- axpby(alpha, beta, x)#
Compute and store y = ɑ·x + β·y.
Logically collective.
- Parameters:
- Return type:
- axpy(alpha, x)#
Compute and store y = ɑ·x + y.
Logically collective.
- aypx(alpha, x)#
Compute and store y = x + ɑ·y.
Logically collective.
- Parameters:
- Return type:
- bindToCPU(flg)#
Bind vector operations execution on the CPU.
Logically collective.
See also
- boundToCPU()#
Return whether the vector has been bound to the CPU.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1409
- Return type:
- chop(tol)#
Set all vector entries less than some absolute tolerance to zero.
Collective.
- Parameters:
tol (float) – The absolute tolerance below which entries are set to zero.
- Return type:
See also
- clearDLPackInfo()#
Clear tensor information.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:705
- Return type:
- classmethod concatenate(vecs)#
Concatenate vectors into a single vector.
Collective.
- Parameters:
- Returns:
- Return type:
See also
- conjugate()#
Conjugate the vector.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2305
- Return type:
- copy(result=None)#
Return a copy of the vector.
Logically collective.
This operation copies vector entries to the new vector.
- create(comm=None)#
Create a vector object.
Collective.
After creation the vector type can then be set with
setType
.- Parameters:
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.- Return type:
- createCUDAWithArrays(cpuarray=None, cudahandle=None, size=None, bsize=None, comm=None)#
Create a
Type.CUDA
vector with optional arrays.Collective.
- Parameters:
cpuarray (Sequence[Scalar] | None) – Host array. Will be lazily allocated if not provided.
cudahandle (Any | None) – Address of the array on the GPU. Will be lazily allocated if not provided.
size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
- createGhost(ghosts, size, bsize=None, comm=None)#
Create a parallel vector with ghost padding on each processor.
Collective.
- Parameters:
size (LayoutSizeSpec) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
- createGhostWithArray(ghosts, array, size=None, bsize=None, comm=None)#
Create a parallel vector with ghost padding and provided arrays.
Collective.
- Parameters:
ghosts (Sequence[int]) – Global indices of ghost points.
array (Sequence[Scalar]) – Array to store the vector values. Must be at least as large as the local size of the vector (including ghost points).
size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
- createHIPWithArrays(cpuarray=None, hiphandle=None, size=None, bsize=None, comm=None)#
Create a
Type.HIP
vector with optional arrays.Collective.
- Parameters:
cpuarray (Sequence[Scalar] | None) – Host array. Will be lazily allocated if not provided.
hiphandle (Any | None) – Address of the array on the GPU. Will be lazily allocated if not provided.
size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
- createLocalVector()#
Create a local vector.
Not collective.
- Returns:
The local vector.
- Return type:
See also
- createMPI(size, bsize=None, comm=None)#
Create a parallel
Type.MPI
vector.Collective.
- Parameters:
size (LayoutSizeSpec) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
See also
- createNest(vecs, isets=None, comm=None)#
Create a
Type.NEST
vector containing multiple nested subvectors.Collective.
- Parameters:
- Return type:
See also
- createSeq(size, bsize=None, comm=None)#
Create a sequential
Type.SEQ
vector.Collective.
- Parameters:
- Return type:
See also
Create a
Type.SHARED
vector that uses shared memory.Collective.
- Parameters:
size (LayoutSizeSpec) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
See also
- createViennaCLWithArrays(cpuarray=None, viennaclvechandle=None, size=None, bsize=None, comm=None)#
Create a
Type.VIENNACL
vector with optional arrays.Collective.
- Parameters:
cpuarray (Sequence[Scalar] | None) – Host array. Will be lazily allocated if not provided.
viennaclvechandle (Any | None) – Address of the array on the GPU. Will be lazily allocated if not provided.
size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
- createWithArray(array, size=None, bsize=None, comm=None)#
Create a vector using a provided array.
Collective.
This method will create either a
Type.SEQ
orType.MPI
depending on the size of the communicator.- Parameters:
array (Sequence[Scalar]) – Array to store the vector values. Must be at least as large as the local size of the vector.
size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
See also
- createWithDLPack(dltensor, size=None, bsize=None, comm=None)#
Create a vector wrapping a DLPack object, sharing the same memory.
Collective.
This operation does not modify the storage of the original tensor and should be used with contiguous tensors only. If the tensor is stored in row-major order (e.g. PyTorch tensors), the resulting vector will look like an unrolled tensor using row-major order.
The resulting vector type will be one of
Type.SEQ
,Type.MPI
,Type.SEQCUDA
,Type.MPICUDA
,Type.SEQHIP
orType.MPIHIP
depending on the type ofdltensor
and the number of processes in the communicator.- Parameters:
dltensor – Either an object with a
__dlpack__
method or a DLPack tensor object.size (LayoutSizeSpec | None) – Vector size.
comm (Comm | None) – MPI communicator, defaults to
Sys.getDefaultComm
.
- Return type:
Self
- destroy()#
Destroy the vector.
Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:163
- Return type:
- dot(vec)#
Return the dot product with
vec
.Collective.
For complex numbers this computes yᴴ·x with
self
as x,vec
as y and where yᴴ denotes the conjugate transpose of y.Use
tDot
for the indefinite form yᵀ·x where yᵀ denotes the transpose of y.
- dotBegin(vec)#
Begin computing the dot product.
Collective.
This should be paired with a call to
dotEnd
.See also
- dotNorm2(vec)#
Return the dot product with
vec
and its squared norm.Collective.
See also
- duplicate(array=None)#
Create a new vector with the same type, optionally with data.
Collective.
- Parameters:
array (Sequence[Scalar] | None) – Optional values to store in the new vector.
- Return type:
See also
- equal(vec)#
Return whether the vector is equal to another.
Collective.
See also
- exp()#
Replace each entry (xₙ) in the vector by exp(xₙ).
Logically collective.
Source code at petsc4py/PETSc/Vec.pyx:2257
- Return type:
- getArray(readonly=False)#
Return local portion of the vector as an
ndarray
.Logically collective.
- Parameters:
readonly (bool) – Request read-only access.
- Return type:
- getBlockSize()#
Return the block size of the vector.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1149
- Return type:
- getBuffer(readonly=False)#
Return a buffered view of the local portion of the vector.
Logically collective.
- Parameters:
readonly (bool) – Request read-only access.
- Returns:
Buffer object
wrapping the local portion of the vector data. This can be used either as a context manager providing access as a numpy array or can be passed to array constructors accepting buffered objects such asnumpy.asarray
.- Return type:
Examples
Accessing the data with a context manager:
>>> vec = PETSc.Vec().createWithArray([1, 2, 3]) >>> with vec.getBuffer() as arr: ... arr array([1., 2., 3.])
Converting the buffer to an
ndarray
:>>> buf = PETSc.Vec().createWithArray([1, 2, 3]).getBuffer() >>> np.asarray(buf) array([1., 2., 3.])
See also
- getCLContextHandle()#
Return the OpenCL context associated with the vector.
Not collective.
- Returns:
Pointer to underlying CL context. This can be used with
pyopencl
throughpyopencl.Context.from_int_ptr
.- Return type:
See also
- getCLMemHandle(mode='rw')#
Return the OpenCL buffer associated with the vector.
Not collective.
- Returns:
Pointer to the device buffer. This can be used with
pyopencl
throughpyopencl.Context.from_int_ptr
.- Return type:
- Parameters:
mode (AccessModeSpec)
Notes
This method may incur a host-to-device copy if the device data is out of date and
mode
is"r"
or"rw"
.
- getCLQueueHandle()#
Return the OpenCL command queue associated with the vector.
Not collective.
- Returns:
Pointer to underlying CL command queue. This can be used with
pyopencl
throughpyopencl.Context.from_int_ptr
.- Return type:
See also
- getCUDAHandle(mode='rw')#
Return a pointer to the device buffer.
Not collective.
The returned pointer should be released using
restoreCUDAHandle
with the same access mode.- Returns:
CUDA device pointer.
- Return type:
- Parameters:
mode (AccessModeSpec)
Notes
This method may incur a host-to-device copy if the device data is out of date and
mode
is"r"
or"rw"
.
- getDM()#
Return the
DM
associated to the vector.Not collective.
Source code at petsc4py/PETSc/Vec.pyx:3519
- Return type:
- getGhostIS()#
Return ghosting indices of a ghost vector.
Collective.
- Returns:
Indices of ghosts.
- Return type:
See also
- getHIPHandle(mode='rw')#
Return a pointer to the device buffer.
Not collective.
The returned pointer should be released using
restoreHIPHandle
with the same access mode.- Returns:
HIP device pointer.
- Return type:
- Parameters:
mode (AccessModeSpec)
Notes
This method may incur a host-to-device copy if the device data is out of date and
mode
is"r"
or"rw"
.
- getLGMap()#
Return the local-to-global mapping.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2899
- Return type:
- getLocalSize()#
Return the local size of the vector.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1107
- Return type:
- getLocalVector(lvec, readonly=False)#
Maps the local portion of the vector into a local vector.
Logically collective.
- Parameters:
lvec (Vec) – The local vector obtained from
createLocalVector
.readonly (bool) – Request read-only access.
- Return type:
- getNestSubVecs()#
Return all the vectors contained in the nested vector.
Not collective.
See also
- getOffloadMask()#
Return the offloading status of the vector.
Not collective.
Common return values include:
1:
PETSC_OFFLOAD_CPU
- CPU has valid entries2:
PETSC_OFFLOAD_GPU
- GPU has valid entries3:
PETSC_OFFLOAD_BOTH
- CPU and GPU are in sync
- Returns:
Enum value from
PetscOffloadMask
describing the offloading status.- Return type:
See also
- getOptionsPrefix()#
Return the prefix used for searching for options in the database.
Not collective.
Source code at petsc4py/PETSc/Vec.pyx:1014
- Return type:
- getOwnershipRange()#
Return the locally owned range of indices
(start, end)
.Not collective.
- Returns:
- Return type:
See also
- getOwnershipRanges()#
Return the range of indices owned by each process.
Not collective.
The returned array is the result of exclusive scan of the local sizes.
See also
Source code at petsc4py/PETSc/Vec.pyx:1184
- Return type:
- getSize()#
Return the global size of the vector.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1093
- Return type:
- getSizes()#
Return the vector sizes.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1121
- Return type:
- getSubVector(iset, subvec=None)#
Return a subvector from given indices.
Collective.
Once finished with the subvector it should be returned with
restoreSubVector
.- Parameters:
- Return type:
See also
- getType()#
Return the type of the vector.
Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1079
- Return type:
- getValue(index)#
Return a single value from the vector.
Not collective.
Only values locally stored may be accessed.
See also
- getValues(indices, values=None)#
Return values from certain locations in the vector.
Not collective.
Only values locally stored may be accessed.
- Parameters:
- Return type:
See also
- getValuesStagStencil(indices, values=None)#
Not implemented.
Source code at petsc4py/PETSc/Vec.pyx:2767
- Return type:
- ghostUpdate(addv=None, mode=None)#
Update ghosted vector entries.
Neighborwise collective.
- Parameters:
addv (InsertModeSpec) – Insertion mode.
mode (ScatterModeSpec) – Scatter mode.
- Return type:
Examples
To accumulate ghost region values onto owning processes:
>>> vec.ghostUpdate(InsertMode.ADD_VALUES, ScatterMode.REVERSE)
Update ghost regions:
>>> vec.ghostUpdate(InsertMode.INSERT_VALUES, ScatterMode.FORWARD)
See also
- ghostUpdateBegin(addv=None, mode=None)#
Begin updating ghosted vector entries.
Neighborwise collective.
See also
ghostUpdateEnd
,ghostUpdate
,createGhost
,VecGhostUpdateBegin
Source code at petsc4py/PETSc/Vec.pyx:3297
- Parameters:
addv (InsertModeSpec)
mode (ScatterModeSpec)
- Return type:
- ghostUpdateEnd(addv=None, mode=None)#
Finish updating ghosted vector entries initiated with
ghostUpdateBegin
.Neighborwise collective.
See also
ghostUpdateBegin
,ghostUpdate
,createGhost
,VecGhostUpdateEnd
Source code at petsc4py/PETSc/Vec.pyx:3314
- Parameters:
addv (InsertModeSpec)
mode (ScatterModeSpec)
- Return type:
- isaxpy(idx, alpha, x)#
Add a scaled reduced-space vector to a subset of the vector.
Logically collective.
This is equivalent to
y[idx[i]] += alpha*x[i]
.- Parameters:
- Return type:
- isset(idx, alpha)#
Set specific elements of the vector to the same value.
Not collective.
- Parameters:
- Return type:
See also
- load(viewer)#
Load a vector.
Collective.
- localForm()#
Return a context manager for viewing ghost vectors in local form.
Logically collective.
- Returns:
Context manager yielding the vector in local (ghosted) form.
- Return type:
Notes
This operation does not perform a copy. To obtain up-to-date ghost values
ghostUpdateBegin
andghostUpdateEnd
must be called first.Non-ghost values can be found at
values[0:nlocal]
and ghost values atvalues[nlocal:nlocal+nghost]
.Examples
>>> with vec.localForm() as lf: ... # compute with lf
- log()#
Replace each entry in the vector by its natural logarithm.
Logically collective.
Source code at petsc4py/PETSc/Vec.pyx:2269
- Return type:
- mDot(vecs, out=None)#
Compute Xᴴ·y with X an array of vectors.
Collective.
- Parameters:
vecs (Sequence[Vec]) – Array of vectors.
out (ArrayScalar | None) – Optional placeholder for the result.
- Return type:
- mDotBegin(vecs, out)#
Starts a split phase multiple dot product computation.
Collective.
- Parameters:
out (ArrayScalar) – Placeholder for the result.
- Return type:
See also
- mDotEnd(vecs, out)#
Ends a split phase multiple dot product computation.
Collective.
- Parameters:
out (ArrayScalar) – Placeholder for the result.
- Return type:
See also
- max()#
Return the vector entry with maximum real part and its location.
Collective.
- Returns:
- Return type:
- maxPointwiseDivide(vec)#
Return the maximum of the component-wise absolute value division.
Logically collective.
Equivalent to
result = max_i abs(x[i] / y[i])
.See also
pointwiseMin
,pointwiseMax
,pointwiseMaxAbs
,VecMaxPointwiseDivide
- maxpy(alphas, vecs)#
Compute and store y = Σₙ(ɑₙ·Xₙ) + y with X an array of vectors.
Logically collective.
Equivalent to
y[:] = alphas[i]*vecs[i, :] + y[:]
.- Parameters:
- Return type:
- min()#
Return the vector entry with minimum real part and its location.
Collective.
- Returns:
- Return type:
- mtDot(vecs, out=None)#
Compute Xᵀ·y with X an array of vectors.
Collective.
- Parameters:
vecs (Sequence[Vec]) – Array of vectors.
out (ArrayScalar | None) – Optional placeholder for the result.
- Return type:
See also
- mtDotBegin(vecs, out)#
Starts a split phase transpose multiple dot product computation.
Collective.
- Parameters:
out (ArrayScalar) – Placeholder for the result.
- Return type:
See also
- mtDotEnd(vecs, out)#
Ends a split phase transpose multiple dot product computation.
Collective.
- Parameters:
out (ArrayScalar) – Placeholder for the result.
- Return type:
See also
- norm(norm_type=None)#
Compute the vector norm.
Collective.
A 2-tuple is returned if
NormType.NORM_1_AND_2
is specified.Source code at petsc4py/PETSc/Vec.pyx:2090
- Parameters:
norm_type (NormTypeSpec)
- Return type:
- normBegin(norm_type=None)#
Begin computing the vector norm.
Collective.
This should be paired with a call to
normEnd
.See also
Source code at petsc4py/PETSc/Vec.pyx:2112
- Parameters:
norm_type (NormTypeSpec)
- Return type:
- normEnd(norm_type=None)#
Finish computations initiated with
normBegin
.Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2131
- Parameters:
norm_type (NormTypeSpec)
- Return type:
- normalize()#
Normalize the vector by its 2-norm.
Collective.
- Returns:
The vector norm before normalization.
- Return type:
See also
- permute(order, invert=False)#
Permute the vector in-place with a provided ordering.
Collective.
- Parameters:
- Return type:
See also
- placeArray(array)#
Set the local portion of the vector to a provided array.
Not collective.
See also
- pointwiseDivide(x, y)#
Compute and store the component-wise division of two vectors.
Logically collective.
Equivalent to
w[i] = x[i] / y[i]
.See also
- pointwiseMax(x, y)#
Compute and store the component-wise maximum of two vectors.
Logically collective.
Equivalent to
w[i] = max(x[i], y[i])
.- Parameters:
- Return type:
See also
- pointwiseMaxAbs(x, y)#
Compute and store the component-wise maximum absolute values.
Logically collective.
Equivalent to
w[i] = max(abs(x[i]), abs(y[i]))
.- Parameters:
- Return type:
See also
- pointwiseMin(x, y)#
Compute and store the component-wise minimum of two vectors.
Logically collective.
Equivalent to
w[i] = min(x[i], y[i])
.- Parameters:
- Return type:
See also
- pointwiseMult(x, y)#
Compute and store the component-wise multiplication of two vectors.
Logically collective.
Equivalent to
w[i] = x[i] * y[i]
.- Parameters:
- Return type:
See also
- reciprocal()#
Replace each entry in the vector by its reciprocal.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2245
- Return type:
- resetArray(force=False)#
Reset the vector to use its default array.
Not collective.
- Parameters:
force (bool) – Force the calling of
VecResetArray
even if no user array has been placed withplaceArray
.- Returns:
The array previously provided by the user with
placeArray
. Can beNone
ifforce
isTrue
and no array was placed before.- Return type:
See also
- restoreCLMemHandle()#
Restore a pointer to the OpenCL buffer obtained with
getCLMemHandle
.Not collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1670
- Return type:
- restoreCUDAHandle(handle, mode='rw')#
Restore a pointer to the device buffer obtained with
getCUDAHandle
.Not collective.
- Parameters:
handle (Any) – CUDA device pointer.
mode (AccessModeSpec) – Access mode.
- Return type:
- restoreHIPHandle(handle, mode='rw')#
Restore a pointer to the device buffer obtained with
getHIPHandle
.Not collective.
- Parameters:
handle (Any) – HIP device pointer.
mode (AccessModeSpec) – Access mode.
- Return type:
- restoreLocalVector(lvec, readonly=False)#
Unmap a local access obtained with
getLocalVector
.Logically collective.
- Parameters:
- Return type:
- restoreSubVector(iset, subvec)#
Restore a subvector extracted using
getSubVector
.Collective.
- Parameters:
- Return type:
See also
- scale(alpha)#
Scale all entries of the vector.
Collective.
This method sets each entry (xₙ) in the vector to ɑ·xₙ.
- set(alpha)#
Set all components of the vector to the same value.
Collective.
See also
- setArray(array)#
Set values for the local portion of the vector.
Logically collective.
See also
- setBlockSize(bsize)#
Set the block size of the vector.
Logically collective.
See also
- setFromOptions()#
Configure the vector from the options database.
Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:1042
- Return type:
- setLGMap(lgmap)#
Set the local-to-global mapping.
Logically collective.
This allows users to insert vector entries using a local numbering with
setValuesLocal
.See also
setValues
,setValuesLocal
,getLGMap
,VecSetLocalToGlobalMapping
- setMPIGhost(ghosts)#
Set the ghost points for a ghosted vector.
Collective.
See also
- setNestSubVecs(sx, idxm=None)#
Set the component vectors at specified indices in the nested vector.
Not collective.
- Parameters:
- Return type:
See also
- setOption(option, flag)#
Set option.
Collective.
See also
- setOptionsPrefix(prefix)#
Set the prefix used for searching for options in the database.
Logically collective.
- setRandom(random=None)#
Set all components of the vector to random numbers.
Collective.
- Parameters:
random (Random | None) – Random number generator. If
None
then one will be created internally.- Return type:
See also
- setSizes(size, bsize=None)#
Set the local and global sizes of the vector.
Collective.
- Parameters:
size (LayoutSizeSpec) – Vector size.
- Return type:
See also
- setType(vec_type)#
Set the vector type.
Collective.
See also
- setUp()#
Set up the internal data structures for using the vector.
Collective.
Source code at petsc4py/PETSc/Vec.pyx:1054
- Return type:
- setValue(index, value, addv=None)#
Insert or add a single value in the vector.
Not collective.
- Parameters:
index (int) – Location to write to. Negative indices are ignored.
value (Scalar) – Value to insert at
index
.addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValue
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.See also
- setValueLocal(index, value, addv=None)#
Insert or add a single value in the vector using a local numbering.
Not collective.
- Parameters:
index (int) – Location to write to.
value (Scalar) – Value to insert at
index
.addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValueLocal
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.See also
- setValues(indices, values, addv=None)#
Insert or add multiple values in the vector.
Not collective.
- Parameters:
indices (Sequence[int]) – Locations to write to. Negative indices are ignored.
addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValues
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.See also
- setValuesBlocked(indices, values, addv=None)#
Insert or add blocks of values in the vector.
Not collective.
Equivalent to
x[bs*indices[i]+j] = y[bs*i+j]
for0 <= i < len(indices)
,0 <= j < bs
andbs
block_size
.- Parameters:
indices (Sequence[int]) – Block indices to write to. Negative indices are ignored.
values (Sequence[Scalar]) – Values to insert at
indices
. Should have lengthlen(indices) * vec.block_size
.addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValuesBlocked
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.See also
- setValuesBlockedLocal(indices, values, addv=None)#
Insert or add blocks of values in the vector with a local numbering.
Not collective.
Equivalent to
x[bs*indices[i]+j] = y[bs*i+j]
for0 <= i < len(indices)
,0 <= j < bs
andbs
block_size
.- Parameters:
values (Sequence[Scalar]) – Values to insert at
indices
. Should have lengthlen(indices) * vec.block_size
.addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValuesBlockedLocal
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.
- setValuesLocal(indices, values, addv=None)#
Insert or add multiple values in the vector with a local numbering.
Not collective.
- Parameters:
addv (InsertModeSpec) – Insertion mode.
- Return type:
Notes
The values may be cached so
assemblyBegin
andassemblyEnd
must be called after all calls of this method are completed.Multiple calls to
setValuesLocal
cannot be made with different values foraddv
without intermediate calls toassemblyBegin
andassemblyEnd
.See also
- setValuesStagStencil(indices, values, addv=None)#
Not implemented.
Source code at petsc4py/PETSc/Vec.pyx:2880
- Return type:
- shift(alpha)#
Shift all entries in the vector.
Collective.
This method sets each entry (xₙ) in the vector to xₙ + ɑ.
- sqrtabs()#
Replace each entry (xₙ) in the vector by √|xₙ|.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2281
- Return type:
- strideGather(field, vec, addv=None)#
Insert component values into a single-component vector.
Collective.
The current vector is expected to be multi-component (
block_size
greater than1
) and the target vector is expected to be single-component.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.vec (Vec) – Single-component vector to be inserted into.
addv (InsertModeSpec) – Insertion mode.
- Return type:
See also
- strideMax(field)#
Return the maximum of entries in a subvector.
Collective.
Equivalent to
max(x[field], x[field+bs], x[field+2*bs], ...)
wherebs
isblock_size
.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.- Returns:
- Return type:
See also
- strideMin(field)#
Return the minimum of entries in a subvector.
Collective.
Equivalent to
min(x[field], x[field+bs], x[field+2*bs], ...)
wherebs
isblock_size
.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.- Returns:
- Return type:
See also
- strideNorm(field, norm_type=None)#
Return the norm of entries in a subvector.
Collective.
Equivalent to
norm(x[field], x[field+bs], x[field+2*bs], ...)
wherebs
isblock_size
.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.norm_type (NormTypeSpec) – The norm type.
- Return type:
See also
- strideScale(field, alpha)#
Scale a component of the vector.
Logically collective.
- Parameters:
- Return type:
See also
- strideScatter(field, vec, addv=None)#
Scatter entries into a component of another vector.
Collective.
The current vector is expected to be single-component (
block_size
of1
) and the target vector is expected to be multi-component.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.vec (Vec) – Multi-component vector to be scattered into.
addv (InsertModeSpec) – Insertion mode.
- Return type:
See also
- strideSum(field)#
Sum subvector entries.
Collective.
Equivalent to
sum(x[field], x[field+bs], x[field+2*bs], ...)
wherebs
isblock_size
.- Parameters:
field (int) – Component index. Must be between
0
andvec.block_size
.- Return type:
See also
- sum()#
Return the sum of all the entries of the vector.
Collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2166
- Return type:
- swap(vec)#
Swap the content of two vectors.
Logically collective.
See also
- tDot(vec)#
Return the indefinite dot product with
vec
.Collective.
This computes yᵀ·x with
self
as x,vec
as y and where yᵀ denotes the transpose of y.
- tDotBegin(vec)#
Begin computing the indefinite dot product.
Collective.
This should be paired with a call to
tDotEnd
.See also
- tDotEnd(vec)#
Finish computing the indefinite dot product initiated with
tDotBegin
.Collective.
See also
- toDLPack(mode='rw')#
Return a DLPack
PyCapsule
wrapping the vector data.Collective.
- Parameters:
mode (AccessModeSpec) – Access mode for the vector.
- Returns:
Capsule of a DLPack tensor wrapping a
Vec
.- Return type:
Notes
It is important that the access mode is respected by the consumer as this is not enforced internally.
See also
- view(viewer=None)#
Display the vector.
Collective.
- waxpy(alpha, x, y)#
Compute and store w = ɑ·x + y.
Logically collective.
- Parameters:
- Return type:
- zeroEntries()#
Set all entries in the vector to zero.
Logically collective.
See also
Source code at petsc4py/PETSc/Vec.pyx:2358
- Return type:
Attributes Documentation
- block_size#
The block size.
- buffer_r#
Read-only buffered view of the local portion of the vector.
- buffer_w#
Writeable buffered view of the local portion of the vector.
- local_size#
The local vector size.
- owner_range#
The locally owned range of indices in the form
[low, high)
.
- owner_ranges#
The range of indices owned by each process.
- size#
The global vector size.
- sizes#
The local and global vector sizes.