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 <[email protected]> To: Roland Richter <[email protected]> Cc: PETSc <[email protected]> Subject: Re: [petsc-users] Efficient FFTShift-implementation for vectors/matrices Message-ID: <[email protected]> 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-co... < https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-co...
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 <[email protected]> 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