This commit is contained in:
Martin Diehl 2019-03-11 23:51:04 +01:00
parent 011af1f686
commit 94cf563141
2 changed files with 60 additions and 55 deletions

View File

@ -79,7 +79,7 @@ program DAMASK_spectral
FIELD_DAMAGE_ID FIELD_DAMAGE_ID
use spectral_mech_Basic use spectral_mech_Basic
use spectral_mech_Polarisation use spectral_mech_Polarisation
use spectral_damage use grid_damage_spectral
use grid_thermal_spectral use grid_thermal_spectral
use results use results
@ -368,7 +368,7 @@ program DAMASK_spectral
call grid_thermal_spectral_init call grid_thermal_spectral_init
case(FIELD_DAMAGE_ID) case(FIELD_DAMAGE_ID)
call spectral_damage_init call grid_damage_spectral_init
end select end select
enddo enddo
@ -511,7 +511,7 @@ program DAMASK_spectral
rotation_BC = loadCases(currentLoadCase)%rotation) rotation_BC = loadCases(currentLoadCase)%rotation)
case(FIELD_THERMAL_ID); call grid_thermal_spectral_forward case(FIELD_THERMAL_ID); call grid_thermal_spectral_forward
case(FIELD_DAMAGE_ID); call spectral_damage_forward case(FIELD_DAMAGE_ID); call grid_damage_spectral_forward
end select end select
enddo enddo
@ -532,7 +532,7 @@ program DAMASK_spectral
solres(field) = grid_thermal_spectral_solution(timeinc,timeIncOld,remainingLoadCaseTime) solres(field) = grid_thermal_spectral_solution(timeinc,timeIncOld,remainingLoadCaseTime)
case(FIELD_DAMAGE_ID) case(FIELD_DAMAGE_ID)
solres(field) = spectral_damage_solution(timeinc,timeIncOld,remainingLoadCaseTime) solres(field) = grid_damage_spectral_solution(timeinc,timeIncOld,remainingLoadCaseTime)
end select end select

View File

@ -3,13 +3,12 @@
!> @author Shaokang Zhang, Max-Planck-Institut für Eisenforschung GmbH !> @author Shaokang Zhang, Max-Planck-Institut für Eisenforschung GmbH
!> @brief Spectral solver for nonlocal damage !> @brief Spectral solver for nonlocal damage
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module spectral_damage module grid_damage_spectral
#include <petsc/finclude/petscsnes.h> #include <petsc/finclude/petscsnes.h>
#include <petsc/finclude/petscdmda.h> #include <petsc/finclude/petscdmda.h>
use PETScdmda use PETScdmda
use PETScsnes use PETScsnes
use prec, only: & use prec, only: &
pInt, &
pReal pReal
use spectral_utilities, only: & use spectral_utilities, only: &
tSolutionState, & tSolutionState, &
@ -18,9 +17,6 @@ module spectral_damage
implicit none implicit none
private private
character (len=*), parameter, public :: &
spectral_damage_label = 'spectraldamage'
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! derived types ! derived types
type(tSolutionParams), private :: params type(tSolutionParams), private :: params
@ -28,7 +24,7 @@ module spectral_damage
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! PETSc data ! PETSc data
SNES, private :: damage_snes SNES, private :: damage_snes
Vec, private :: solution Vec, private :: solution_vec
PetscInt, private :: xstart, xend, ystart, yend, zstart, zend PetscInt, private :: xstart, xend, ystart, yend, zstart, zend
real(pReal), private, dimension(:,:,:), allocatable :: & real(pReal), private, dimension(:,:,:), allocatable :: &
damage_current, & !< field of current damage damage_current, & !< field of current damage
@ -37,21 +33,21 @@ module spectral_damage
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! reference diffusion tensor, mobility etc. ! reference diffusion tensor, mobility etc.
integer(pInt), private :: totalIter = 0 !< total iteration in current increment integer, private :: totalIter = 0 !< total iteration in current increment
real(pReal), dimension(3,3), private :: D_ref real(pReal), dimension(3,3), private :: D_ref
real(pReal), private :: mobility_ref real(pReal), private :: mobility_ref
public :: & public :: &
spectral_damage_init, & grid_damage_spectral_init, &
spectral_damage_solution, & grid_damage_spectral_solution, &
spectral_damage_forward grid_damage_spectral_forward
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief allocates all neccessary fields and fills them with data, potentially from restart info !> @brief allocates all neccessary fields and fills them with data, potentially from restart info
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine spectral_damage_init() subroutine grid_damage_spectral_init
use IO, only: & use IO, only: &
IO_intOut IO_intOut
use spectral_utilities, only: & use spectral_utilities, only: &
@ -64,7 +60,8 @@ subroutine spectral_damage_init()
damage_nonlocal_getMobility damage_nonlocal_getMobility
use numerics, only: & use numerics, only: &
worldrank, & worldrank, &
worldsize worldsize, &
petsc_options
implicit none implicit none
PetscInt, dimension(worldsize) :: localK PetscInt, dimension(worldsize) :: localK
@ -74,11 +71,18 @@ subroutine spectral_damage_init()
PetscErrorCode :: ierr PetscErrorCode :: ierr
character(len=100) :: snes_type character(len=100) :: snes_type
write(6,'(/,a)') ' <<<+- spectral_damage init -+>>>' write(6,'(/,a)') ' <<<+- grid_spectral_damage init -+>>>'
write(6,'(/,a)') ' Shanthraj et al., Handbook of Mechanics of Materials, 2019' write(6,'(/,a)') ' Shanthraj et al., Handbook of Mechanics of Materials, 2019'
write(6,'(a)') ' https://doi.org/10.1007/978-981-10-6855-3_80' write(6,'(a)') ' https://doi.org/10.1007/978-981-10-6855-3_80'
!--------------------------------------------------------------------------------------------------
! set default and user defined options for PETSc
call PETScOptionsInsertString(PETSC_NULL_OPTIONS,'-damage_snes_type ngmres',ierr)
CHKERRQ(ierr)
call PETScOptionsInsertString(PETSC_NULL_OPTIONS,trim(petsc_options),ierr)
CHKERRQ(ierr)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! initialize solver specific parts of PETSc ! initialize solver specific parts of PETSc
call SNESCreate(PETSC_COMM_WORLD,damage_snes,ierr); CHKERRQ(ierr) call SNESCreate(PETSC_COMM_WORLD,damage_snes,ierr); CHKERRQ(ierr)
@ -98,8 +102,8 @@ subroutine spectral_damage_init()
call SNESSetDM(damage_snes,damage_grid,ierr); CHKERRQ(ierr) !< connect snes to da call SNESSetDM(damage_snes,damage_grid,ierr); CHKERRQ(ierr) !< connect snes to da
call DMsetFromOptions(damage_grid,ierr); CHKERRQ(ierr) call DMsetFromOptions(damage_grid,ierr); CHKERRQ(ierr)
call DMsetUp(damage_grid,ierr); CHKERRQ(ierr) call DMsetUp(damage_grid,ierr); CHKERRQ(ierr)
call DMCreateGlobalVector(damage_grid,solution,ierr); CHKERRQ(ierr) !< global solution vector (grid x 1, i.e. every def grad tensor) call DMCreateGlobalVector(damage_grid,solution_vec,ierr); CHKERRQ(ierr) !< global solution vector (grid x 1, i.e. every def grad tensor)
call DMDASNESSetFunctionLocal(damage_grid,INSERT_VALUES,spectral_damage_formResidual,& call DMDASNESSetFunctionLocal(damage_grid,INSERT_VALUES,grid_damage_spectral_formResidual,&
PETSC_NULL_SNES,ierr) !< residual vector of same shape as solution vector PETSC_NULL_SNES,ierr) !< residual vector of same shape as solution vector
CHKERRQ(ierr) CHKERRQ(ierr)
call SNESSetFromOptions(damage_snes,ierr); CHKERRQ(ierr) !< pull it all together with additional CLI arguments call SNESSetFromOptions(damage_snes,ierr); CHKERRQ(ierr) !< pull it all together with additional CLI arguments
@ -122,14 +126,14 @@ subroutine spectral_damage_init()
xend = xstart + xend - 1 xend = xstart + xend - 1
yend = ystart + yend - 1 yend = ystart + yend - 1
zend = zstart + zend - 1 zend = zstart + zend - 1
call VecSet(solution,1.0_pReal,ierr); CHKERRQ(ierr) call VecSet(solution_vec,1.0_pReal,ierr); CHKERRQ(ierr)
allocate(damage_current(grid(1),grid(2),grid3), source=1.0_pReal) allocate(damage_current(grid(1),grid(2),grid3), source=1.0_pReal)
allocate(damage_lastInc(grid(1),grid(2),grid3), source=1.0_pReal) allocate(damage_lastInc(grid(1),grid(2),grid3), source=1.0_pReal)
allocate(damage_stagInc(grid(1),grid(2),grid3), source=1.0_pReal) allocate(damage_stagInc(grid(1),grid(2),grid3), source=1.0_pReal)
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! damage reference diffusion update ! damage reference diffusion update
cell = 0_pInt cell = 0
D_ref = 0.0_pReal D_ref = 0.0_pReal
mobility_ref = 0.0_pReal mobility_ref = 0.0_pReal
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1) do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
@ -142,12 +146,12 @@ subroutine spectral_damage_init()
mobility_ref = mobility_ref*wgt mobility_ref = mobility_ref*wgt
call MPI_Allreduce(MPI_IN_PLACE,mobility_ref,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,mobility_ref,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
end subroutine spectral_damage_init end subroutine grid_damage_spectral_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief solution for the spectral damage scheme with internal iterations !> @brief solution for the spectral damage scheme with internal iterations
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
type(tSolutionState) function spectral_damage_solution(timeinc,timeinc_old,loadCaseTime) function grid_damage_spectral_solution(timeinc,timeinc_old,loadCaseTime) result(solution)
use numerics, only: & use numerics, only: &
itmax, & itmax, &
err_damage_tolAbs, & err_damage_tolAbs, &
@ -168,32 +172,33 @@ type(tSolutionState) function spectral_damage_solution(timeinc,timeinc_old,loadC
PetscInt ::position PetscInt ::position
PetscReal :: minDamage, maxDamage, stagNorm, solnNorm PetscReal :: minDamage, maxDamage, stagNorm, solnNorm
PetscErrorCode :: ierr PetscErrorCode :: ierr
type(tSolutionState) :: &
solution
SNESConvergedReason :: reason SNESConvergedReason :: reason
spectral_damage_solution%converged =.false. solution%converged =.false.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! set module wide availabe data ! set module wide availabe data
params%timeinc = timeinc params%timeinc = timeinc
params%timeincOld = timeinc_old params%timeincOld = timeinc_old
call SNESSolve(damage_snes,PETSC_NULL_VEC,solution,ierr); CHKERRQ(ierr) call SNESSolve(damage_snes,PETSC_NULL_VEC,solution_vec,ierr); CHKERRQ(ierr)
call SNESGetConvergedReason(damage_snes,reason,ierr); CHKERRQ(ierr) call SNESGetConvergedReason(damage_snes,reason,ierr); CHKERRQ(ierr)
if (reason < 1) then if (reason < 1) then
spectral_damage_solution%converged = .false. solution%converged = .false.
spectral_damage_solution%iterationsNeeded = itmax solution%iterationsNeeded = itmax
else else
spectral_damage_solution%converged = .true. solution%converged = .true.
spectral_damage_solution%iterationsNeeded = totalIter solution%iterationsNeeded = totalIter
endif endif
stagNorm = maxval(abs(damage_current - damage_stagInc)) stagNorm = maxval(abs(damage_current - damage_stagInc))
solnNorm = maxval(abs(damage_current)) solnNorm = maxval(abs(damage_current))
call MPI_Allreduce(MPI_IN_PLACE,stagNorm,1,MPI_DOUBLE,MPI_MAX,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,stagNorm,1,MPI_DOUBLE,MPI_MAX,PETSC_COMM_WORLD,ierr)
call MPI_Allreduce(MPI_IN_PLACE,solnNorm,1,MPI_DOUBLE,MPI_MAX,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,solnNorm,1,MPI_DOUBLE,MPI_MAX,PETSC_COMM_WORLD,ierr)
damage_stagInc = damage_current damage_stagInc = damage_current
spectral_damage_solution%stagConverged = stagNorm < err_damage_tolAbs & solution%stagConverged = stagNorm < min(err_damage_tolAbs,err_damage_tolRel*solnNorm)
.or. stagNorm < err_damage_tolRel*solnNorm
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! updating damage state ! updating damage state
@ -205,20 +210,20 @@ type(tSolutionState) function spectral_damage_solution(timeinc,timeinc_old,loadC
call VecMin(solution,position,minDamage,ierr); CHKERRQ(ierr) call VecMin(solution,position,minDamage,ierr); CHKERRQ(ierr)
call VecMax(solution,position,maxDamage,ierr); CHKERRQ(ierr) call VecMax(solution,position,maxDamage,ierr); CHKERRQ(ierr)
if (spectral_damage_solution%converged) & if (solution%converged) &
write(6,'(/,a)') ' ... nonlocal damage converged .....................................' write(6,'(/,a)') ' ... nonlocal damage converged .....................................'
write(6,'(/,a,f8.6,2x,f8.6,2x,f8.6,/)',advance='no') ' Minimum|Maximum|Delta Damage = ',& write(6,'(/,a,f8.6,2x,f8.6,2x,f8.6,/)',advance='no') ' Minimum|Maximum|Delta Damage = ',&
minDamage, maxDamage, stagNorm minDamage, maxDamage, stagNorm
write(6,'(/,a)') ' ===========================================================================' write(6,'(/,a)') ' ==========================================================================='
flush(6) flush(6)
end function spectral_damage_solution end function grid_damage_spectral_solution
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief forms the spectral damage residual vector !> @brief forms the spectral damage residual vector
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine spectral_damage_formResidual(in,x_scal,f_scal,dummy,ierr) subroutine grid_damage_spectral_formResidual(in,x_scal,f_scal,dummy,ierr)
use numerics, only: & use numerics, only: &
residualStiffness residualStiffness
use mesh, only: & use mesh, only: &
@ -252,7 +257,7 @@ subroutine spectral_damage_formResidual(in,x_scal,f_scal,dummy,ierr)
f_scal f_scal
PetscObject :: dummy PetscObject :: dummy
PetscErrorCode :: ierr PetscErrorCode :: ierr
integer(pInt) :: i, j, k, cell integer :: i, j, k, cell
real(pReal) :: phiDot, dPhiDot_dPhi, mobility real(pReal) :: phiDot, dPhiDot_dPhi, mobility
damage_current = x_scal damage_current = x_scal
@ -263,18 +268,18 @@ subroutine spectral_damage_formResidual(in,x_scal,f_scal,dummy,ierr)
call utilities_FFTscalarForward() call utilities_FFTscalarForward()
call utilities_fourierScalarGradient() !< calculate gradient of damage field call utilities_fourierScalarGradient() !< calculate gradient of damage field
call utilities_FFTvectorBackward() call utilities_FFTvectorBackward()
cell = 0_pInt cell = 0
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt,grid(1) do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
cell = cell + 1_pInt cell = cell + 1
vectorField_real(1:3,i,j,k) = math_mul33x3(damage_nonlocal_getDiffusion33(1,cell) - D_ref, & vectorField_real(1:3,i,j,k) = math_mul33x3(damage_nonlocal_getDiffusion33(1,cell) - D_ref, &
vectorField_real(1:3,i,j,k)) vectorField_real(1:3,i,j,k))
enddo; enddo; enddo enddo; enddo; enddo
call utilities_FFTvectorForward() call utilities_FFTvectorForward()
call utilities_fourierVectorDivergence() !< calculate damage divergence in fourier field call utilities_fourierVectorDivergence() !< calculate damage divergence in fourier field
call utilities_FFTscalarBackward() call utilities_FFTscalarBackward()
cell = 0_pInt cell = 0
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt,grid(1) do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
cell = cell + 1_pInt cell = cell + 1
call damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, damage_current(i,j,k), 1, cell) call damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, damage_current(i,j,k), 1, cell)
mobility = damage_nonlocal_getMobility(1,cell) mobility = damage_nonlocal_getMobility(1,cell)
scalarField_real(i,j,k) = params%timeinc*scalarField_real(i,j,k) + & scalarField_real(i,j,k) = params%timeinc*scalarField_real(i,j,k) + &
@ -298,12 +303,12 @@ subroutine spectral_damage_formResidual(in,x_scal,f_scal,dummy,ierr)
! constructing residual ! constructing residual
f_scal = scalarField_real(1:grid(1),1:grid(2),1:grid3) - damage_current f_scal = scalarField_real(1:grid(1),1:grid(2),1:grid3) - damage_current
end subroutine spectral_damage_formResidual end subroutine grid_damage_spectral_formResidual
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief spectral damage forwarding routine !> @brief spectral damage forwarding routine
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine spectral_damage_forward() subroutine grid_damage_spectral_forward()
use mesh, only: & use mesh, only: &
grid, & grid, &
grid3 grid3
@ -316,7 +321,7 @@ subroutine spectral_damage_forward()
damage_nonlocal_getMobility damage_nonlocal_getMobility
implicit none implicit none
integer(pInt) :: i, j, k, cell integer :: i, j, k, cell
DM :: dm_local DM :: dm_local
PetscScalar, dimension(:,:,:), pointer :: x_scal PetscScalar, dimension(:,:,:), pointer :: x_scal
PetscErrorCode :: ierr PetscErrorCode :: ierr
@ -326,24 +331,24 @@ subroutine spectral_damage_forward()
damage_stagInc = damage_lastInc damage_stagInc = damage_lastInc
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! reverting damage field state ! reverting damage field state
cell = 0_pInt cell = 0
call SNESGetDM(damage_snes,dm_local,ierr); CHKERRQ(ierr) call SNESGetDM(damage_snes,dm_local,ierr); CHKERRQ(ierr)
call DMDAVecGetArrayF90(dm_local,solution,x_scal,ierr); CHKERRQ(ierr) !< get the data out of PETSc to work with call DMDAVecGetArrayF90(dm_local,solution_vec,x_scal,ierr); CHKERRQ(ierr) !< get the data out of PETSc to work with
x_scal(xstart:xend,ystart:yend,zstart:zend) = damage_current x_scal(xstart:xend,ystart:yend,zstart:zend) = damage_current
call DMDAVecRestoreArrayF90(dm_local,solution,x_scal,ierr); CHKERRQ(ierr) call DMDAVecRestoreArrayF90(dm_local,solution_vec,x_scal,ierr); CHKERRQ(ierr)
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt,grid(1) do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
cell = cell + 1_pInt cell = cell + 1
call damage_nonlocal_putNonLocalDamage(damage_current(i,j,k),1,cell) call damage_nonlocal_putNonLocalDamage(damage_current(i,j,k),1,cell)
enddo; enddo; enddo enddo; enddo; enddo
else else
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update rate and forward last inc ! update rate and forward last inc
damage_lastInc = damage_current damage_lastInc = damage_current
cell = 0_pInt cell = 0
D_ref = 0.0_pReal D_ref = 0.0_pReal
mobility_ref = 0.0_pReal mobility_ref = 0.0_pReal
do k = 1_pInt, grid3; do j = 1_pInt, grid(2); do i = 1_pInt,grid(1) do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
cell = cell + 1_pInt cell = cell + 1
D_ref = D_ref + damage_nonlocal_getDiffusion33(1,cell) D_ref = D_ref + damage_nonlocal_getDiffusion33(1,cell)
mobility_ref = mobility_ref + damage_nonlocal_getMobility(1,cell) mobility_ref = mobility_ref + damage_nonlocal_getMobility(1,cell)
enddo; enddo; enddo enddo; enddo; enddo
@ -353,6 +358,6 @@ subroutine spectral_damage_forward()
call MPI_Allreduce(MPI_IN_PLACE,mobility_ref,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) call MPI_Allreduce(MPI_IN_PLACE,mobility_ref,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr)
endif endif
end subroutine spectral_damage_forward end subroutine grid_damage_spectral_forward
end module spectral_damage end module grid_damage_spectral