polishing

This commit is contained in:
Martin Diehl 2021-04-11 15:32:17 +02:00
parent f655a6fe5c
commit 887524bcc1
6 changed files with 52 additions and 62 deletions

View File

@ -15,8 +15,8 @@ module grid_damage_spectral
use IO
use spectral_utilities
use discretization_grid
use YAML_types
use homogenization
use YAML_types
use config
implicit none
@ -61,7 +61,7 @@ contains
!> @brief allocates all neccessary fields and fills them with data
! ToDo: Restart not implemented
!--------------------------------------------------------------------------------------------------
subroutine grid_damage_spectral_init
subroutine grid_damage_spectral_init()
PetscInt, dimension(0:worldsize-1) :: localK
DM :: damage_grid
@ -88,10 +88,10 @@ subroutine grid_damage_spectral_init
num_generic => config_numerics%get('generic',defaultVal=emptyDict)
num%residualStiffness = num_generic%get_asFloat('residualStiffness', defaultVal=1.0e-6_pReal)
if (num%residualStiffness < 0.0_pReal) call IO_error(301,ext_msg='residualStiffness')
if (num%itmax <= 1) call IO_error(301,ext_msg='itmax')
if (num%eps_damage_atol <= 0.0_pReal) call IO_error(301,ext_msg='eps_damage_atol')
if (num%eps_damage_rtol <= 0.0_pReal) call IO_error(301,ext_msg='eps_damage_rtol')
if (num%residualStiffness < 0.0_pReal) call IO_error(301,ext_msg='residualStiffness')
if (num%itmax <= 1) call IO_error(301,ext_msg='itmax')
if (num%eps_damage_atol <= 0.0_pReal) call IO_error(301,ext_msg='eps_damage_atol')
if (num%eps_damage_rtol <= 0.0_pReal) call IO_error(301,ext_msg='eps_damage_rtol')
!--------------------------------------------------------------------------------------------------
! set default and user defined options for PETSc
@ -146,6 +146,7 @@ subroutine grid_damage_spectral_init
allocate(phi_current(grid(1),grid(2),grid3), source=1.0_pReal)
allocate(phi_lastInc(grid(1),grid(2),grid3), source=1.0_pReal)
allocate(phi_stagInc(grid(1),grid(2),grid3), source=1.0_pReal)
call VecSet(solution_vec,1.0_pReal,ierr); CHKERRQ(ierr)
call updateReference
@ -216,21 +217,21 @@ end function grid_damage_spectral_solution
subroutine grid_damage_spectral_forward(cutBack)
logical, intent(in) :: cutBack
integer :: i, j, k, ce
integer :: i, j, k, ce
DM :: dm_local
PetscScalar, dimension(:,:,:), pointer :: x_scal
PetscErrorCode :: ierr
PetscScalar, dimension(:,:,:), pointer :: x_scal
PetscErrorCode :: ierr
if (cutBack) then
phi_current = phi_lastInc
phi_stagInc = phi_lastInc
!--------------------------------------------------------------------------------------------------
! reverting damage field state
ce = 0
call SNESGetDM(damage_snes,dm_local,ierr); CHKERRQ(ierr)
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) = phi_current
call DMDAVecRestoreArrayF90(dm_local,solution_vec,x_scal,ierr); CHKERRQ(ierr)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
call homogenization_set_phi(phi_current(i,j,k),ce)
@ -271,8 +272,7 @@ subroutine formResidual(in,x_scal,f_scal,dummy,ierr)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
vectorField_real(1:3,i,j,k) = matmul(homogenization_K_phi(ce) - K_ref, &
vectorField_real(1:3,i,j,k))
vectorField_real(1:3,i,j,k) = matmul(homogenization_K_phi(ce) - K_ref, vectorField_real(1:3,i,j,k))
enddo; enddo; enddo
call utilities_FFTvectorForward
call utilities_fourierVectorDivergence !< calculate damage divergence in fourier field
@ -290,6 +290,7 @@ subroutine formResidual(in,x_scal,f_scal,dummy,ierr)
call utilities_FFTscalarForward
call utilities_fourierGreenConvolution(K_ref, mu_ref, params%timeinc)
call utilities_FFTscalarBackward
where(scalarField_real(1:grid(1),1:grid(2),1:grid3) > phi_lastInc) &
scalarField_real(1:grid(1),1:grid(2),1:grid3) = phi_lastInc
where(scalarField_real(1:grid(1),1:grid(2),1:grid3) < num%residualStiffness) &
@ -305,7 +306,7 @@ end subroutine formResidual
!--------------------------------------------------------------------------------------------------
!> @brief update reference viscosity and conductivity
!--------------------------------------------------------------------------------------------------
subroutine updateReference
subroutine updateReference()
integer :: ce,ierr

View File

@ -18,7 +18,6 @@ module grid_thermal_spectral
use homogenization
use YAML_types
use config
use material
implicit none
private
@ -46,7 +45,7 @@ module grid_thermal_spectral
!--------------------------------------------------------------------------------------------------
! reference diffusion tensor, mobility etc.
integer :: totalIter = 0 !< total iteration in current increment
integer :: totalIter = 0 !< total iteration in current increment
real(pReal), dimension(3,3) :: K_ref
real(pReal) :: mu_ref
@ -68,7 +67,7 @@ subroutine grid_thermal_spectral_init(T_0)
PetscInt, dimension(0:worldsize-1) :: localK
integer :: i, j, k, ce
DM :: thermal_grid
PetscScalar, dimension(:,:,:), pointer :: x_scal
PetscScalar, dimension(:,:,:), pointer :: T_PETSc
PetscErrorCode :: ierr
class(tNode), pointer :: &
num_grid
@ -85,9 +84,9 @@ subroutine grid_thermal_spectral_init(T_0)
num%eps_thermal_atol = num_grid%get_asFloat ('eps_thermal_atol',defaultVal=1.0e-2_pReal)
num%eps_thermal_rtol = num_grid%get_asFloat ('eps_thermal_rtol',defaultVal=1.0e-6_pReal)
if (num%itmax <= 1) call IO_error(301,ext_msg='itmax')
if (num%eps_thermal_atol <= 0.0_pReal) call IO_error(301,ext_msg='eps_thermal_atol')
if (num%eps_thermal_rtol <= 0.0_pReal) call IO_error(301,ext_msg='eps_thermal_rtol')
if (num%itmax <= 1) call IO_error(301,ext_msg='itmax')
if (num%eps_thermal_atol <= 0.0_pReal) call IO_error(301,ext_msg='eps_thermal_atol')
if (num%eps_thermal_rtol <= 0.0_pReal) call IO_error(301,ext_msg='eps_thermal_rtol')
!--------------------------------------------------------------------------------------------------
! set default and user defined options for PETSc
@ -130,6 +129,7 @@ subroutine grid_thermal_spectral_init(T_0)
allocate(T_current(grid(1),grid(2),grid3), source=0.0_pReal)
allocate(T_lastInc(grid(1),grid(2),grid3), source=0.0_pReal)
allocate(T_stagInc(grid(1),grid(2),grid3), source=0.0_pReal)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
@ -138,9 +138,10 @@ subroutine grid_thermal_spectral_init(T_0)
T_stagInc(i,j,k) = T_current(i,j,k)
call homogenization_thermal_setField(T_0,0.0_pReal,ce)
enddo; enddo; enddo
call DMDAVecGetArrayF90(thermal_grid,solution_vec,x_scal,ierr); CHKERRQ(ierr) !< get the data out of PETSc to work with
x_scal(xstart:xend,ystart:yend,zstart:zend) = T_current
call DMDAVecRestoreArrayF90(thermal_grid,solution_vec,x_scal,ierr); CHKERRQ(ierr)
call DMDAVecGetArrayF90(thermal_grid,solution_vec,T_PETSc,ierr); CHKERRQ(ierr)
T_PETSc(xstart:xend,ystart:yend,zstart:zend) = T_current
call DMDAVecRestoreArrayF90(thermal_grid,solution_vec,T_PETSc,ierr); CHKERRQ(ierr)
call updateReference
@ -190,9 +191,7 @@ function grid_thermal_spectral_solution(timeinc) result(solution)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
call homogenization_thermal_setField(T_current(i,j,k), &
(T_current(i,j,k)-T_lastInc(i,j,k))/params%timeinc, &
ce)
call homogenization_thermal_setField(T_current(i,j,k),(T_current(i,j,k)-T_lastInc(i,j,k))/params%timeinc,ce)
enddo; enddo; enddo
call VecMin(solution_vec,devNull,T_min,ierr); CHKERRQ(ierr)
@ -223,16 +222,14 @@ subroutine grid_thermal_spectral_forward(cutBack)
!--------------------------------------------------------------------------------------------------
! reverting thermal field state
ce = 0
call SNESGetDM(thermal_snes,dm_local,ierr); CHKERRQ(ierr)
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) = T_current
call DMDAVecRestoreArrayF90(dm_local,solution_vec,x_scal,ierr); CHKERRQ(ierr)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
call homogenization_thermal_setField(T_current(i,j,k), &
(T_current(i,j,k)-T_lastInc(i,j,k))/params%timeinc, &
ce)
call homogenization_thermal_setField(T_current(i,j,k),(T_current(i,j,k)-T_lastInc(i,j,k))/params%timeinc,ce)
enddo; enddo; enddo
else
T_lastInc = T_current
@ -270,8 +267,7 @@ subroutine formResidual(in,x_scal,f_scal,dummy,ierr)
ce = 0
do k = 1, grid3; do j = 1, grid(2); do i = 1,grid(1)
ce = ce + 1
vectorField_real(1:3,i,j,k) = matmul(homogenization_K_T(ce) - K_ref, &
vectorField_real(1:3,i,j,k))
vectorField_real(1:3,i,j,k) = matmul(homogenization_K_T(ce) - K_ref, vectorField_real(1:3,i,j,k))
enddo; enddo; enddo
call utilities_FFTvectorForward
call utilities_fourierVectorDivergence !< calculate temperature divergence in fourier field
@ -300,7 +296,7 @@ end subroutine formResidual
!--------------------------------------------------------------------------------------------------
!> @brief update reference viscosity and conductivity
!--------------------------------------------------------------------------------------------------
subroutine updateReference
subroutine updateReference()
integer :: ce,ierr

View File

@ -38,14 +38,12 @@ module subroutine damage_init()
configHomogenizations, &
configHomogenization, &
configHomogenizationDamage, &
num_generic, &
material_homogenization
integer :: ho
integer :: Ninstances,Nmaterialpoints,h
num_generic
integer :: ho,Nmaterialpoints
print'(/,a)', ' <<<+- homogenization:damage init -+>>>'
print'(/,a)', ' <<<+- homogenization:damage:pass init -+>>>'
configHomogenizations => config_material%get('homogenization')
allocate(param(configHomogenizations%length))
@ -62,6 +60,10 @@ module subroutine damage_init()
#else
prm%output = configHomogenizationDamage%get_as1dString('output',defaultVal=emptyStringArray)
#endif
Nmaterialpoints = count(material_homogenizationAt == ho)
damageState_h(ho)%sizeState = 1
allocate(damageState_h(ho)%state0(1,Nmaterialpoints), source=1.0_pReal)
allocate(damageState_h(ho)%state (1,Nmaterialpoints), source=1.0_pReal)
else
prm%output = emptyStringArray
endif
@ -73,18 +75,7 @@ module subroutine damage_init()
num_generic => config_numerics%get('generic',defaultVal= emptyDict)
num_damage%charLength = num_generic%get_asFloat('charLength',defaultVal=1.0_pReal)
Ninstances = count(damage_type == DAMAGE_nonlocal_ID)
material_homogenization => config_material%get('homogenization')
do h = 1, material_homogenization%length
if (damage_type(h) /= DAMAGE_NONLOCAL_ID) cycle
Nmaterialpoints = count(material_homogenizationAt == h)
damageState_h(h)%sizeState = 1
allocate(damageState_h(h)%state0 (1,Nmaterialpoints), source=1.0_pReal)
allocate(damageState_h(h)%state (1,Nmaterialpoints), source=1.0_pReal)
enddo
call pass_init()
end subroutine damage_init
@ -183,13 +174,13 @@ module subroutine damage_results(ho,group)
integer :: o
associate(prm => param(ho))
outputsLoop: do o = 1,size(prm%output)
select case(prm%output(o))
case ('phi')
call results_writeDataset(group,damagestate_h(ho)%state(1,:),prm%output(o),&
'damage indicator','-')
end select
enddo outputsLoop
outputsLoop: do o = 1,size(prm%output)
select case(prm%output(o))
case ('phi')
call results_writeDataset(group,damagestate_h(ho)%state(1,:),prm%output(o),&
'damage indicator','-')
end select
enddo outputsLoop
end associate
end subroutine damage_results

View File

@ -6,8 +6,9 @@ submodule(homogenization:damage) damage_pass
contains
module subroutine pass_init
module subroutine pass_init()
print'(/,a)', ' <<<+- homogenization:damage:pass init -+>>>'
end subroutine pass_init

View File

@ -45,8 +45,6 @@ module subroutine thermal_init()
print'(/,a)', ' <<<+- homogenization:thermal init -+>>>'
print'(/,a)', ' <<<+- homogenization:thermal:pass init -+>>>'
configHomogenizations => config_material%get('homogenization')
@ -71,6 +69,8 @@ module subroutine thermal_init()
end associate
enddo
call pass_init()
end subroutine thermal_init

View File

@ -6,8 +6,9 @@ submodule(homogenization:thermal) thermal_pass
contains
module subroutine pass_init
module subroutine pass_init()
print'(/,a)', ' <<<+- homogenization:thermal:pass init -+>>>'
end subroutine pass_init