From 887524bcc1cd53bc6a842f0da5799d9082656ede Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 11 Apr 2021 15:32:17 +0200 Subject: [PATCH] polishing --- src/grid/grid_damage_spectral.f90 | 27 ++++++++++---------- src/grid/grid_thermal_spectral.f90 | 34 +++++++++++-------------- src/homogenization_damage.f90 | 39 +++++++++++------------------ src/homogenization_damage_pass.f90 | 5 ++-- src/homogenization_thermal.f90 | 4 +-- src/homogenization_thermal_pass.f90 | 5 ++-- 6 files changed, 52 insertions(+), 62 deletions(-) diff --git a/src/grid/grid_damage_spectral.f90 b/src/grid/grid_damage_spectral.f90 index 7fe9a4a25..967ddb73a 100644 --- a/src/grid/grid_damage_spectral.f90 +++ b/src/grid/grid_damage_spectral.f90 @@ -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 diff --git a/src/grid/grid_thermal_spectral.f90 b/src/grid/grid_thermal_spectral.f90 index 60d8d99db..1ccb3b901 100644 --- a/src/grid/grid_thermal_spectral.f90 +++ b/src/grid/grid_thermal_spectral.f90 @@ -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 diff --git a/src/homogenization_damage.f90 b/src/homogenization_damage.f90 index b8aecccae..0546834fd 100644 --- a/src/homogenization_damage.f90 +++ b/src/homogenization_damage.f90 @@ -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 diff --git a/src/homogenization_damage_pass.f90 b/src/homogenization_damage_pass.f90 index cbb7ec79f..bdb44b822 100644 --- a/src/homogenization_damage_pass.f90 +++ b/src/homogenization_damage_pass.f90 @@ -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 diff --git a/src/homogenization_thermal.f90 b/src/homogenization_thermal.f90 index 0effe09ed..f681036de 100644 --- a/src/homogenization_thermal.f90 +++ b/src/homogenization_thermal.f90 @@ -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 diff --git a/src/homogenization_thermal_pass.f90 b/src/homogenization_thermal_pass.f90 index 940a15024..87020d8d2 100644 --- a/src/homogenization_thermal_pass.f90 +++ b/src/homogenization_thermal_pass.f90 @@ -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