unified naming scheme with damage/thermal

This commit is contained in:
Martin Diehl 2023-07-17 02:52:26 +02:00
parent 8682df0e86
commit ddeb218728
3 changed files with 88 additions and 88 deletions

View File

@ -52,9 +52,9 @@ module grid_mechanical_FEM
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! PETSc data ! PETSc data
DM :: mechanical_grid DM :: DM_mech
SNES :: SNES_mechanical SNES :: SNES_mech
Vec :: solution_current, solution_lastInc, solution_rate Vec :: u_PETSc, u_lastInc_PETSc, uDot_PETSc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! common pointwise data ! common pointwise data
@ -163,9 +163,9 @@ subroutine grid_mechanical_FEM_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! initialize solver specific parts of PETSc ! initialize solver specific parts of PETSc
call SNESCreate(PETSC_COMM_WORLD,SNES_mechanical,err_PETSc) call SNESCreate(PETSC_COMM_WORLD,SNES_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetOptionsPrefix(SNES_mechanical,'mechanical_',err_PETSc) call SNESSetOptionsPrefix(SNES_mech,'mechanical_',err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call MPI_Allgather(int(cells3,MPI_INTEGER_KIND),1_MPI_INTEGER_KIND,MPI_INTEGER,& call MPI_Allgather(int(cells3,MPI_INTEGER_KIND),1_MPI_INTEGER_KIND,MPI_INTEGER,&
cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI) cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI)
@ -177,44 +177,44 @@ subroutine grid_mechanical_FEM_init
1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), & 1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), &
3_pPETSCINT, 1_pPETSCINT, & ! #dof (u, vector), ghost boundary width (domain overlap) 3_pPETSCINT, 1_pPETSCINT, & ! #dof (u, vector), ghost boundary width (domain overlap)
[int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells [int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells
mechanical_grid,err_PETSc) DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetFromOptions(mechanical_grid,err_PETSc) call DMsetFromOptions(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetUp(mechanical_grid,err_PETSc) call DMsetUp(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDASetUniformCoordinates(mechanical_grid,0.0_pREAL,geomSize(1),0.0_pREAL,geomSize(2),0.0_pREAL,geomSize(3),err_PETSc) call DMDASetUniformCoordinates(DM_mech,0.0_pREAL,geomSize(1),0.0_pREAL,geomSize(2),0.0_pREAL,geomSize(3),err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMCreateGlobalVector(mechanical_grid,solution_current,err_PETSc) call DMCreateGlobalVector(DM_mech,u_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMCreateGlobalVector(mechanical_grid,solution_lastInc,err_PETSc) call DMCreateGlobalVector(DM_mech,u_lastInc_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMCreateGlobalVector(mechanical_grid,solution_rate ,err_PETSc) call DMCreateGlobalVector(DM_mech,uDot_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMSNESSetFunctionLocal(mechanical_grid,formResidual,PETSC_NULL_SNES,err_PETSc) call DMSNESSetFunctionLocal(DM_mech,formResidual,PETSC_NULL_SNES,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMSNESSetJacobianLocal(mechanical_grid,formJacobian,PETSC_NULL_SNES,err_PETSc) call DMSNESSetJacobianLocal(DM_mech,formJacobian,PETSC_NULL_SNES,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetConvergenceTest(SNES_mechanical,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "_converged" call SNESSetConvergenceTest(SNES_mech,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "_converged"
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetMaxLinearSolveFailures(SNES_mechanical, huge(1_pPETSCINT), err_PETSc) ! ignore linear solve failures call SNESSetMaxLinearSolveFailures(SNES_mech, huge(1_pPETSCINT), err_PETSc) ! ignore linear solve failures
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetDM(SNES_mechanical,mechanical_grid,err_PETSc) call SNESSetDM(SNES_mech,DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetFromOptions(SNES_mechanical,err_PETSc) ! pull it all together with additional cli arguments call SNESSetFromOptions(SNES_mech,err_PETSc) ! pull it all together with additional cli arguments
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init fields ! init fields
call VecSet(solution_current,0.0_pREAL,err_PETSc) call VecSet(u_PETSc,0.0_pREAL,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call VecSet(solution_lastInc,0.0_pREAL,err_PETSc) call VecSet(u_lastInc_PETSc,0.0_pREAL,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call VecSet(solution_rate ,0.0_pREAL,err_PETSc) call VecSet(uDot_PETSc ,0.0_pREAL,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDAVecGetArrayF90(mechanical_grid,solution_current,u,err_PETSc) call DMDAVecGetArrayF90(DM_mech,u_PETSc,u,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDAVecGetArrayF90(mechanical_grid,solution_lastInc,u_lastInc,err_PETSc) call DMDAVecGetArrayF90(DM_mech,u_lastInc_PETSc,u_lastInc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
delta = geomSize/real(cells,pREAL) ! grid spacing delta = geomSize/real(cells,pREAL) ! grid spacing
@ -271,9 +271,9 @@ subroutine grid_mechanical_FEM_init
call utilities_constitutiveResponse(P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 call utilities_constitutiveResponse(P_current,P_av,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2
F, & ! target F F, & ! target F
0.0_pREAL) ! time increment 0.0_pREAL) ! time increment
call DMDAVecRestoreArrayF90(mechanical_grid,solution_current,u,err_PETSc) call DMDAVecRestoreArrayF90(DM_mech,u_PETSc,u,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDAVecRestoreArrayF90(mechanical_grid,solution_lastInc,u_lastInc,err_PETSc) call DMDAVecRestoreArrayF90(DM_mech,u_lastInc_PETSc,u_lastInc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
restartRead2: if (CLI_restartInc > 0) then restartRead2: if (CLI_restartInc > 0) then
@ -315,9 +315,9 @@ function grid_mechanical_FEM_solution(incInfoIn) result(solution)
! update stiffness (and gamma operator) ! update stiffness (and gamma operator)
S = utilities_maskedCompliance(params%rotation_BC,params%stress_mask,C_volAvg) S = utilities_maskedCompliance(params%rotation_BC,params%stress_mask,C_volAvg)
call SNESsolve(SNES_mechanical,PETSC_NULL_VEC,solution_current,err_PETSc) call SNESsolve(SNES_mech,PETSC_NULL_VEC,u_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetConvergedReason(SNES_mechanical,reason,err_PETSc) call SNESGetConvergedReason(SNES_mech,reason,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
solution%converged = reason > 0 solution%converged = reason > 0
@ -374,15 +374,15 @@ subroutine grid_mechanical_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remai
end if end if
if (guess) then if (guess) then
call VecWAXPY(solution_rate,-1.0_pREAL,solution_lastInc,solution_current,err_PETSc) call VecWAXPY(uDot_PETSc,-1.0_pREAL,u_lastInc_PETSc,u_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call VecScale(solution_rate,1.0_pREAL/Delta_t_old,err_PETSc) call VecScale(uDot_PETSc,1.0_pREAL/Delta_t_old,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
else else
call VecSet(solution_rate,0.0_pREAL,err_PETSc) call VecSet(uDot_PETSc,0.0_pREAL,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end if end if
call VecCopy(solution_current,solution_lastInc,err_PETSc) call VecCopy(u_PETSc,u_lastInc_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
F_lastInc = F F_lastInc = F
@ -398,7 +398,7 @@ subroutine grid_mechanical_FEM_forward(cutBack,guess,Delta_t,Delta_t_old,t_remai
if (stress_BC%myType=='dot_P') P_aim = P_aim & if (stress_BC%myType=='dot_P') P_aim = P_aim &
+ merge(.0_pREAL,stress_BC%values,stress_BC%mask)*Delta_t + merge(.0_pREAL,stress_BC%values,stress_BC%mask)*Delta_t
call VecAXPY(solution_current,Delta_t,solution_rate,err_PETSc) call VecAXPY(u_PETSc,Delta_t,uDot_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -430,9 +430,9 @@ subroutine grid_mechanical_FEM_restartWrite()
PetscScalar, dimension(:,:,:,:), pointer :: u,u_lastInc PetscScalar, dimension(:,:,:,:), pointer :: u,u_lastInc
call DMDAVecGetArrayReadF90(mechanical_grid,solution_current,u,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,u_PETSc,u,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDAVecGetArrayReadF90(mechanical_grid,solution_lastInc,u_lastInc,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,u_lastInc_PETSc,u_lastInc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
print'(1x,a)', 'saving solver data required for restart'; flush(IO_STDOUT) print'(1x,a)', 'saving solver data required for restart'; flush(IO_STDOUT)
@ -459,9 +459,9 @@ subroutine grid_mechanical_FEM_restartWrite()
call HDF5_closeFile(fileHandle) call HDF5_closeFile(fileHandle)
end if end if
call DMDAVecRestoreArrayReadF90(mechanical_grid,solution_current,u,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,u_PETSc,u,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDAVecRestoreArrayReadF90(mechanical_grid,solution_lastInc,u_lastInc,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,u_lastInc_PETSc,u_lastInc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end subroutine grid_mechanical_FEM_restartWrite end subroutine grid_mechanical_FEM_restartWrite
@ -531,9 +531,9 @@ subroutine formResidual(da_local,x_local, &
integer(MPI_INTEGER_KIND) :: err_MPI integer(MPI_INTEGER_KIND) :: err_MPI
real(pREAL), dimension(3,3,3,3) :: devNull real(pREAL), dimension(3,3,3,3) :: devNull
call SNESGetNumberFunctionEvals(SNES_mechanical,nfuncs,err_PETSc) call SNESGetNumberFunctionEvals(SNES_mech,nfuncs,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetIterationNumber(SNES_mechanical,PETScIter,err_PETSc) call SNESGetIterationNumber(SNES_mech,PETScIter,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)

View File

@ -51,9 +51,9 @@ module grid_mechanical_spectral_basic
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! PETSc data ! PETSc data
DM :: da DM :: DM_mech
SNES :: SNES_mechanical SNES :: SNES_mech
Vec :: solution_vec Vec :: F_PETSc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! common pointwise data ! common pointwise data
@ -162,9 +162,9 @@ subroutine grid_mechanical_spectral_basic_init()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! initialize solver specific parts of PETSc ! initialize solver specific parts of PETSc
call SNESCreate(PETSC_COMM_WORLD,SNES_mechanical,err_PETSc) call SNESCreate(PETSC_COMM_WORLD,SNES_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetOptionsPrefix(SNES_mechanical,'mechanical_',err_PETSc) call SNESSetOptionsPrefix(SNES_mech,'mechanical_',err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call MPI_Allgather(int(cells3,MPI_INTEGER_KIND),1_MPI_INTEGER_KIND,MPI_INTEGER,& call MPI_Allgather(int(cells3,MPI_INTEGER_KIND),1_MPI_INTEGER_KIND,MPI_INTEGER,&
cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI) cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI)
@ -176,26 +176,26 @@ subroutine grid_mechanical_spectral_basic_init()
1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), & 1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), &
9_pPETSCINT, 0_pPETSCINT, & ! #dof (F, tensor), ghost boundary width (domain overlap) 9_pPETSCINT, 0_pPETSCINT, & ! #dof (F, tensor), ghost boundary width (domain overlap)
[int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells [int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells
da,err_PETSc) ! handle, error DM_mech,err_PETSc) ! handle, error
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetFromOptions(da,err_PETSc) call DMsetFromOptions(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetUp(da,err_PETSc) call DMsetUp(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMcreateGlobalVector(da,solution_vec,err_PETSc) ! global solution vector (cells x 9, i.e. every def grad tensor) call DMcreateGlobalVector(DM_mech,F_PETSc,err_PETSc) ! global solution vector (cells x 9, i.e. every def grad tensor)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector call DMDASNESsetFunctionLocal(DM_mech,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESsetConvergenceTest(SNES_mechanical,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged" call SNESsetConvergenceTest(SNES_mech,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged"
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetDM(SNES_mechanical,da,err_PETSc) call SNESSetDM(SNES_mech,DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESsetFromOptions(SNES_mechanical,err_PETSc) ! pull it all together with additional CLI arguments call SNESsetFromOptions(SNES_mech,err_PETSc) ! pull it all together with additional CLI arguments
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init fields ! init fields
call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) ! places pointer on PETSc data call DMDAVecGetArrayF90(DM_mech,F_PETSc,F,err_PETSc) ! places pointer on PETSc data
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
restartRead: if (CLI_restartInc > 0) then restartRead: if (CLI_restartInc > 0) then
@ -231,7 +231,7 @@ subroutine grid_mechanical_spectral_basic_init()
call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2
reshape(F,shape(F_lastInc)), & ! target F reshape(F,shape(F_lastInc)), & ! target F
0.0_pREAL) ! time increment 0.0_pREAL) ! time increment
call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) ! deassociate pointer call DMDAVecRestoreArrayF90(DM_mech,F_PETSc,F,err_PETSc) ! deassociate pointer
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
restartRead2: if (CLI_restartInc > 0) then restartRead2: if (CLI_restartInc > 0) then
@ -280,9 +280,9 @@ function grid_mechanical_spectral_basic_solution(incInfoIn) result(solution)
S = utilities_maskedCompliance(params%rotation_BC,params%stress_mask,C_volAvg) S = utilities_maskedCompliance(params%rotation_BC,params%stress_mask,C_volAvg)
if (num%update_gamma) call utilities_updateGamma(C_minMaxAvg) if (num%update_gamma) call utilities_updateGamma(C_minMaxAvg)
call SNESsolve(SNES_mechanical,PETSC_NULL_VEC,solution_vec,err_PETSc) call SNESsolve(SNES_mech,PETSC_NULL_VEC,F_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetConvergedReason(SNES_mechanical,reason,err_PETSc) call SNESGetConvergedReason(SNES_mech,reason,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
solution%converged = reason > 0 solution%converged = reason > 0
@ -317,7 +317,7 @@ subroutine grid_mechanical_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_
real(pREAL), pointer, dimension(:,:,:,:) :: F real(pREAL), pointer, dimension(:,:,:,:) :: F
call DMDAVecGetArrayF90(da,solution_vec,F,err_PETSc) call DMDAVecGetArrayF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
if (cutBack) then if (cutBack) then
@ -361,7 +361,7 @@ subroutine grid_mechanical_spectral_basic_forward(cutBack,guess,Delta_t,Delta_t_
F = reshape(utilities_forwardField(Delta_t,F_lastInc,Fdot, & ! estimate of F at end of time+Delta_t that matches rotated F_aim on average F = reshape(utilities_forwardField(Delta_t,F_lastInc,Fdot, & ! estimate of F at end of time+Delta_t that matches rotated F_aim on average
rotation_BC%rotate(F_aim,active=.true.)),[9,cells(1),cells(2),cells3]) rotation_BC%rotate(F_aim,active=.true.)),[9,cells(1),cells(2),cells3])
call DMDAVecRestoreArrayF90(da,solution_vec,F,err_PETSc) call DMDAVecRestoreArrayF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -381,10 +381,10 @@ subroutine grid_mechanical_spectral_basic_updateCoords()
PetscErrorCode :: err_PETSc PetscErrorCode :: err_PETSc
real(pREAL), dimension(:,:,:,:), pointer :: F real(pREAL), dimension(:,:,:,:), pointer :: F
call DMDAVecGetArrayReadF90(da,solution_vec,F,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call utilities_updateCoords(reshape(F,[3,3,size(F,2),size(F,3),size(F,4)])) call utilities_updateCoords(reshape(F,[3,3,size(F,2),size(F,3),size(F,4)]))
call DMDAVecRestoreArrayReadF90(da,solution_vec,F,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end subroutine grid_mechanical_spectral_basic_updateCoords end subroutine grid_mechanical_spectral_basic_updateCoords
@ -399,7 +399,7 @@ subroutine grid_mechanical_spectral_basic_restartWrite()
integer(HID_T) :: fileHandle, groupHandle integer(HID_T) :: fileHandle, groupHandle
real(pREAL), dimension(:,:,:,:), pointer :: F real(pREAL), dimension(:,:,:,:), pointer :: F
call DMDAVecGetArrayReadF90(da,solution_vec,F,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
if (num%update_gamma) C_minMaxAvgRestart = C_minMaxAvg if (num%update_gamma) C_minMaxAvgRestart = C_minMaxAvg
@ -427,7 +427,7 @@ subroutine grid_mechanical_spectral_basic_restartWrite()
call HDF5_closeFile(fileHandle) call HDF5_closeFile(fileHandle)
end if end if
call DMDAVecRestoreArrayReadF90(da,solution_vec,F,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,F_PETSc,F,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end subroutine grid_mechanical_spectral_basic_restartWrite end subroutine grid_mechanical_spectral_basic_restartWrite
@ -498,9 +498,9 @@ subroutine formResidual(residual_subdomain, F, &
integer(MPI_INTEGER_KIND) :: err_MPI integer(MPI_INTEGER_KIND) :: err_MPI
call SNESGetNumberFunctionEvals(SNES_mechanical,nfuncs,err_PETSc) call SNESGetNumberFunctionEvals(SNES_mech,nfuncs,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetIterationNumber(SNES_mechanical,PETScIter,err_PETSc) call SNESGetIterationNumber(SNES_mech,PETScIter,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1 ! new increment if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1 ! new increment

View File

@ -56,9 +56,9 @@ module grid_mechanical_spectral_polarisation
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! PETSc data ! PETSc data
DM :: da DM :: DM_mech
SNES :: SNES_mechanical SNES :: SNES_mech
Vec :: solution_vec Vec :: FandF_tau_PETSc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! common pointwise data ! common pointwise data
@ -183,9 +183,9 @@ subroutine grid_mechanical_spectral_polarisation_init()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! initialize solver specific parts of PETSc ! initialize solver specific parts of PETSc
call SNESCreate(PETSC_COMM_WORLD,SNES_mechanical,err_PETSc) call SNESCreate(PETSC_COMM_WORLD,SNES_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetOptionsPrefix(SNES_mechanical,'mechanical_',err_PETSc) call SNESSetOptionsPrefix(SNES_mech,'mechanical_',err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call MPI_Allgather(int(cells3,pPetscInt),1_MPI_INTEGER_KIND,MPI_INTEGER,& call MPI_Allgather(int(cells3,pPetscInt),1_MPI_INTEGER_KIND,MPI_INTEGER,&
cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI) cells3_global,1_MPI_INTEGER_KIND,MPI_INTEGER,MPI_COMM_WORLD,err_MPI)
@ -197,26 +197,26 @@ subroutine grid_mechanical_spectral_polarisation_init()
1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), & 1_pPETSCINT, 1_pPETSCINT, int(worldsize,pPETSCINT), &
18_pPETSCINT, 0_pPETSCINT, & ! #dof (2xtensor), ghost boundary width (domain overlap) 18_pPETSCINT, 0_pPETSCINT, & ! #dof (2xtensor), ghost boundary width (domain overlap)
[int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells [int(cells(1),pPETSCINT)],[int(cells(2),pPETSCINT)],int(cells3_global,pPETSCINT), & ! local cells
da,err_PETSc) ! handle, error DM_mech,err_PETSc) ! handle, error
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetFromOptions(da,err_PETSc) call DMsetFromOptions(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMsetUp(da,err_PETSc) call DMsetUp(DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMcreateGlobalVector(da,solution_vec,err_PETSc) ! global solution vector (cells x 18, i.e. every def grad tensor) call DMcreateGlobalVector(DM_mech,FandF_tau_PETSc,err_PETSc) ! global solution vector (cells x 18, i.e. every def grad tensor)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call DMDASNESsetFunctionLocal(da,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector call DMDASNESsetFunctionLocal(DM_mech,INSERT_VALUES,formResidual,PETSC_NULL_SNES,err_PETSc) ! residual vector of same shape as solution vector
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESsetConvergenceTest(SNES_mechanical,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged" call SNESsetConvergenceTest(SNES_mech,converged,PETSC_NULL_SNES,PETSC_NULL_FUNCTION,err_PETSc) ! specify custom convergence check function "converged"
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESSetDM(SNES_mechanical,da,err_PETSc) call SNESSetDM(SNES_mech,DM_mech,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESsetFromOptions(SNES_mechanical,err_PETSc) ! pull it all together with additional CLI arguments call SNESsetFromOptions(SNES_mech,err_PETSc) ! pull it all together with additional CLI arguments
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! init fields ! init fields
call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) ! places pointer on PETSc data call DMDAVecGetArrayF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc) ! places pointer on PETSc data
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
F => FandF_tau(0: 8,:,:,:) F => FandF_tau(0: 8,:,:,:)
F_tau => FandF_tau(9:17,:,:,:) F_tau => FandF_tau(9:17,:,:,:)
@ -260,7 +260,7 @@ subroutine grid_mechanical_spectral_polarisation_init()
call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 call utilities_constitutiveResponse(P,P_av,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2
reshape(F,shape(F_lastInc)), & ! target F reshape(F,shape(F_lastInc)), & ! target F
0.0_pREAL) ! time increment 0.0_pREAL) ! time increment
call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,err_PETSc) ! deassociate pointer call DMDAVecRestoreArrayF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc) ! deassociate pointer
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
restartRead2: if (CLI_restartInc > 0) then restartRead2: if (CLI_restartInc > 0) then
@ -315,9 +315,9 @@ function grid_mechanical_spectral_polarisation_solution(incInfoIn) result(soluti
S_scale = math_invSym3333(C_minMaxAvg) S_scale = math_invSym3333(C_minMaxAvg)
end if end if
call SNESSolve(SNES_mechanical,PETSC_NULL_VEC,solution_vec,err_PETSc) call SNESSolve(SNES_mech,PETSC_NULL_VEC,FandF_tau_PETSc,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetConvergedReason(SNES_mechanical,reason,err_PETSc) call SNESGetConvergedReason(SNES_mech,reason,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
solution%converged = reason > 0 solution%converged = reason > 0
@ -354,7 +354,7 @@ subroutine grid_mechanical_spectral_polarisation_forward(cutBack,guess,Delta_t,D
real(pREAL), dimension(3,3) :: F_lambda33 real(pREAL), dimension(3,3) :: F_lambda33
call DMDAVecGetArrayF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecGetArrayF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
F => FandF_tau(0: 8,:,:,:) F => FandF_tau(0: 8,:,:,:)
F_tau => FandF_tau(9:17,:,:,:) F_tau => FandF_tau(9:17,:,:,:)
@ -418,7 +418,7 @@ subroutine grid_mechanical_spectral_polarisation_forward(cutBack,guess,Delta_t,D
end do; end do; end do end do; end do; end do
end if end if
call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecRestoreArrayF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -438,10 +438,10 @@ subroutine grid_mechanical_spectral_polarisation_updateCoords()
PetscErrorCode :: err_PETSc PetscErrorCode :: err_PETSc
real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau
call DMDAVecGetArrayReadF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call utilities_updateCoords(reshape(FandF_tau(0:8,:,:,:),[3,3,size(FandF_tau,2),size(FandF_tau,3),size(FandF_tau,4)])) call utilities_updateCoords(reshape(FandF_tau(0:8,:,:,:),[3,3,size(FandF_tau,2),size(FandF_tau,3),size(FandF_tau,4)]))
call DMDAVecRestoreArrayReadF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end subroutine grid_mechanical_spectral_polarisation_updateCoords end subroutine grid_mechanical_spectral_polarisation_updateCoords
@ -456,7 +456,7 @@ subroutine grid_mechanical_spectral_polarisation_restartWrite()
integer(HID_T) :: fileHandle, groupHandle integer(HID_T) :: fileHandle, groupHandle
real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau, F, F_tau real(pREAL), dimension(:,:,:,:), pointer :: FandF_tau, F, F_tau
call DMDAVecGetArrayReadF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecGetArrayReadF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
F => FandF_tau(0: 8,:,:,:) F => FandF_tau(0: 8,:,:,:)
F_tau => FandF_tau(9:17,:,:,:) F_tau => FandF_tau(9:17,:,:,:)
@ -488,7 +488,7 @@ subroutine grid_mechanical_spectral_polarisation_restartWrite()
call HDF5_closeFile(fileHandle) call HDF5_closeFile(fileHandle)
end if end if
call DMDAVecRestoreArrayReadF90(da,solution_vec,FandF_tau,err_PETSc) call DMDAVecRestoreArrayReadF90(DM_mech,FandF_tau_PETSc,FandF_tau,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
end subroutine grid_mechanical_spectral_polarisation_restartWrite end subroutine grid_mechanical_spectral_polarisation_restartWrite
@ -576,9 +576,9 @@ subroutine formResidual(residual_subdomain, FandF_tau, &
call MPI_Allreduce(MPI_IN_PLACE,F_av,9_MPI_INTEGER_KIND,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD,err_MPI) call MPI_Allreduce(MPI_IN_PLACE,F_av,9_MPI_INTEGER_KIND,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD,err_MPI)
if (err_MPI /= 0_MPI_INTEGER_KIND) error stop 'MPI error' if (err_MPI /= 0_MPI_INTEGER_KIND) error stop 'MPI error'
call SNESGetNumberFunctionEvals(SNES_mechanical,nfuncs,err_PETSc) call SNESGetNumberFunctionEvals(SNES_mech,nfuncs,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
call SNESGetIterationNumber(SNES_mechanical,PETScIter,err_PETSc) call SNESGetIterationNumber(SNES_mech,PETScIter,err_PETSc)
CHKERRQ(err_PETSc) CHKERRQ(err_PETSc)
if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1 ! new increment if (nfuncs == 0 .and. PETScIter == 0) totalIter = -1 ! new increment