[petsc-users] petsc-users Digest, Vol 149, Issue 47

[Sajid Ali]

 You could use VecPermulte and MatPermute from the PETSc API to permute the
vectors/matrices. However, MatPermute creates a new matrix and even though
VecPermute permutes the vector (locally) in-place, it allocates a temporary
array and frees the original array.

Since you are working on dense matrices and vectors and desire to avoid
temporary allocation, you could use BLAS level 1 swap function (used by
VecSwap to swap two vectors) which will probably be the most optimized
version for the hardware you're using (since it's implemented by platform
specific intrinsics and assembly).



----------------------------------------------------------------------
>
> Message: 1
> Date: Tue, 18 May 2021 11:44:36 -0500
> From: Barry Smith <bsmith at petsc.dev>
> To: Roland Richter <roland.richter at ntnu.no>
> Cc: PETSc <petsc-users at mcs.anl.gov>
> Subject: Re: [petsc-users] Efficient FFTShift-implementation for
>         vectors/matrices
> Message-ID: <83012843-A364-4CE1-A92B-1768EB11591E at petsc.dev>
> Content-Type: text/plain; charset="us-ascii"
>
>
>   I found a variety of things on the web, below. I don't understand this
> but for the even case it seems one simply modifies the input matrix before
> the FFT http://www.fftw.org/faq/section3.html#centerorigin
>
>
> https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-code
> <
> https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-code
> >
>
>   https://www.dsprelated.com/showthread/comp.dsp/20790-1.php <
> https://www.dsprelated.com/showthread/comp.dsp/20790-1.php>
>
>
>
> > On May 18, 2021, at 9:48 AM, Roland Richter <roland.richter at ntnu.no>
> wrote:
> >
> > Dear all,
> >
> > I tried to implement the function fftshift from numpy (i.e. swap the
> > half-spaces of all axis) for row vectors in a matrix by using the
> > following code
> >
> > void fft_shift(Mat &fft_matrix) {
> >         PetscScalar *mat_ptr;
> >         MatDenseGetArray (fft_matrix, &mat_ptr);
> >         PetscInt r_0, r_1;
> >         MatGetOwnershipRange(fft_matrix, &r_0, &r_1);
> >         PetscInt local_row_num = r_1 - r_0;
> >         arma::cx_mat temp_mat(local_row_num, Ntime, arma::fill::zeros);
> >         for(int i = 0; i < Ntime; ++i) {
> >             const PetscInt row_shift = i * local_row_num;
> >             for(int j = 0; j < local_row_num; ++j) {
> >                 const PetscInt cur_pos = j + row_shift;
> >                 if(i < (int)(Ntime / 2))
> >                     temp_mat(j, i + int(Ntime / 2)) = *(mat_ptr +
> cur_pos);
> >                 else
> >                     temp_mat(j, i - int(Ntime / 2)) = *(mat_ptr +
> cur_pos);
> >             }
> >         }
> >         for(int i = 0; i < Ntime; ++i) {
> >             const PetscInt row_shift = i * local_row_num;
> >             for(int j = 0; j < local_row_num; ++j) {
> >                 const PetscInt cur_pos = j + row_shift;
> >                 *(mat_ptr + cur_pos) = temp_mat(j, i);
> >             }
> >         }
> >         MatDenseRestoreArray (fft_matrix, &mat_ptr);
> > }
> >
> > but I do not like the approach of having a second matrix as temporary
> > storage space. Are there more efficient approaches possible using
> > PETSc-functions?
> >
> > Thanks!
> >
> > Regards,
> >
> > Roland Richter
> >
>
> -------------- next part --------------
> An HTML attachment was scrubbed...
> URL: <
> http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210518/b8710455/attachment-0001.html
> >
>
>

-- 
Sajid Ali (he/him) | PhD Candidate
Applied Physics
Northwestern University
s-sajid-ali.github.io
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210518/37903fdd/attachment.html>