From 40b4062b835c0fc246067dd26093e77b4ac5c150 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 6 Jul 2023 14:39:36 +0200 Subject: [PATCH] reshape to expected shape --- src/grid/grid_mech_spectral_basic.f90 | 12 +++++++----- src/grid/grid_mech_spectral_polarisation.f90 | 12 +++++++----- 2 files changed, 14 insertions(+), 10 deletions(-) diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index 81d00db3f..f8b317b77 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -375,16 +375,17 @@ end subroutine grid_mechanical_spectral_basic_forward !-------------------------------------------------------------------------------------------------- -!> @brief Update coordinates +!> @brief Update coordinates. !-------------------------------------------------------------------------------------------------- -subroutine grid_mechanical_spectral_basic_updateCoords +subroutine grid_mechanical_spectral_basic_updateCoords() PetscErrorCode :: err_PETSc real(pREAL), dimension(:,:,:,:), pointer :: F + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) CHKERRQ(err_PETSc) - call utilities_updateCoords(F) + call utilities_updateCoords(reshape(F,[3,3,size(F,2),size(F,3),size(F,4)])) call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) CHKERRQ(err_PETSc) @@ -392,14 +393,15 @@ end subroutine grid_mechanical_spectral_basic_updateCoords !-------------------------------------------------------------------------------------------------- -!> @brief Write current solver and constitutive data for restart to file +!> @brief Write current solver and constitutive data for restart to file. !-------------------------------------------------------------------------------------------------- -subroutine grid_mechanical_spectral_basic_restartWrite +subroutine grid_mechanical_spectral_basic_restartWrite() PetscErrorCode :: err_PETSc integer(HID_T) :: fileHandle, groupHandle real(pREAL), dimension(:,:,:,:), pointer :: F + call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) CHKERRQ(err_PETSc) diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 73d2b89f2..3db67a7b1 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -432,16 +432,17 @@ end subroutine grid_mechanical_spectral_polarisation_forward !-------------------------------------------------------------------------------------------------- -!> @brief Update coordinates +!> @brief Update coordinates. !-------------------------------------------------------------------------------------------------- -subroutine grid_mechanical_spectral_polarisation_updateCoords +subroutine grid_mechanical_spectral_polarisation_updateCoords() PetscErrorCode :: err_PETSc real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau + call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) CHKERRQ(err_PETSc) - call utilities_updateCoords(FandF_tau(0:8,:,:,:)) + call utilities_updateCoords(reshape(FandF_tau(0:8,:,:,:),[3,3,size(FandF_tau,2),size(FandF_tau,3),size(FandF_tau,4)])) call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,err_PETSc) CHKERRQ(err_PETSc) @@ -449,14 +450,15 @@ end subroutine grid_mechanical_spectral_polarisation_updateCoords !-------------------------------------------------------------------------------------------------- -!> @brief Write current solver and constitutive data for restart to file +!> @brief Write current solver and constitutive data for restart to file. !-------------------------------------------------------------------------------------------------- -subroutine grid_mechanical_spectral_polarisation_restartWrite +subroutine grid_mechanical_spectral_polarisation_restartWrite() PetscErrorCode :: err_PETSc integer(HID_T) :: fileHandle, groupHandle real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau, F, F_tau + call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) CHKERRQ(err_PETSc) F => FandF_tau(0: 8,:,:,:)