avoid code duplication

This commit is contained in:
Martin Diehl 2019-02-13 07:22:37 +01:00
parent 9574dfae2d
commit 61baa66c38
5 changed files with 57 additions and 101 deletions

View File

@ -235,6 +235,7 @@ module material
public :: & public :: &
material_init, & material_init, &
material_allocatePlasticState, & material_allocatePlasticState, &
material_allocateSourceState, &
ELASTICITY_hooke_ID ,& ELASTICITY_hooke_ID ,&
PLASTICITY_none_ID, & PLASTICITY_none_ID, &
PLASTICITY_isotropic_ID, & PLASTICITY_isotropic_ID, &
@ -966,6 +967,47 @@ subroutine material_allocatePlasticState(phase,NofMyPhase,&
end subroutine material_allocatePlasticState end subroutine material_allocatePlasticState
!--------------------------------------------------------------------------------------------------
!> @brief allocates the source state of a phase
!--------------------------------------------------------------------------------------------------
subroutine material_allocateSourceState(phase,of,NofMyPhase,sizeState)
use numerics, only: &
numerics_integrator2 => numerics_integrator ! compatibility hack
implicit none
integer(pInt), intent(in) :: &
phase, &
of, &
NofMyPhase, &
sizeState
integer(pInt) :: numerics_integrator ! compatibility hack
numerics_integrator = numerics_integrator2(1) ! compatibility hack
sourceState(phase)%p(of)%sizeState = sizeState
sourceState(phase)%p(of)%sizeDotState = sizeState
sourceState(phase)%p(of)%sizeDeltaState = 0_pInt
allocate(sourceState(phase)%p(of)%aTolState (sizeState), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%state0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%subState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%state (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%dotState (sizeState,NofMyPhase), source=0.0_pReal)
if (numerics_integrator == 1_pInt) then
allocate(sourceState(phase)%p(of)%previousDotState (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(of)%previousDotState2 (sizeState,NofMyPhase), source=0.0_pReal)
endif
if (numerics_integrator == 4_pInt) &
allocate(sourceState(phase)%p(of)%RK4dotState (sizeState,NofMyPhase), source=0.0_pReal)
if (numerics_integrator == 5_pInt) &
allocate(sourceState(phase)%p(of)%RKCK45dotState (6,sizeState,NofMyPhase), source=0.0_pReal)
allocate(plasticState(phase)%deltaState (0,NofMyPhase), source=0.0_pReal)
end subroutine material_allocateSourceState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief populates the grains !> @brief populates the grains
!> @details populates the grains by identifying active microstructure/homogenization pairs, !> @details populates the grains by identifying active microstructure/homogenization pairs,

View File

@ -102,6 +102,7 @@ subroutine source_damage_anisoBrittle_init(fileUnit)
IO_timeStamp, & IO_timeStamp, &
IO_EOF IO_EOF
use material, only: & use material, only: &
material_allocateSourceState, &
phase_source, & phase_source, &
phase_Nsources, & phase_Nsources, &
phase_Noutput, & phase_Noutput, &
@ -285,33 +286,10 @@ subroutine source_damage_anisoBrittle_init(fileUnit)
endif endif
enddo outputsLoop enddo outputsLoop
!-------------------------------------------------------------------------------------------------- call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1_pInt)
! Determine size of state array
sizeDotState = 1_pInt
sizeDeltaState = 0_pInt
sizeState = 1_pInt
sourceState(phase)%p(sourceOffset)%sizeState = sizeState
sourceState(phase)%p(sourceOffset)%sizeDotState = sizeDotState
sourceState(phase)%p(sourceOffset)%sizeDeltaState = sizeDeltaState
sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_anisoBrittle_sizePostResults(instance) sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_anisoBrittle_sizePostResults(instance)
allocate(sourceState(phase)%p(sourceOffset)%aTolState (sizeState), & sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol
source=param(instance)%aTol)
allocate(sourceState(phase)%p(sourceOffset)%state0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%subState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%state (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 1_pInt)) then
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
endif
if (any(numerics_integrator == 4_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 5_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal)
endif endif

View File

@ -105,6 +105,7 @@ subroutine source_damage_anisoDuctile_init(fileUnit)
IO_timeStamp, & IO_timeStamp, &
IO_EOF IO_EOF
use material, only: & use material, only: &
material_allocateSourceState, &
phase_source, & phase_source, &
phase_Nsources, & phase_Nsources, &
phase_Noutput, & phase_Noutput, &
@ -285,32 +286,9 @@ subroutine source_damage_anisoDuctile_init(fileUnit)
endif endif
enddo outputsLoop enddo outputsLoop
!-------------------------------------------------------------------------------------------------- call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1_pInt)
! Determine size of state array
sizeDotState = 1_pInt
sizeDeltaState = 0_pInt
sizeState = 1_pInt
sourceState(phase)%p(sourceOffset)%sizeState = sizeState
sourceState(phase)%p(sourceOffset)%sizeDotState = sizeDotState
sourceState(phase)%p(sourceOffset)%sizeDeltaState = sizeDeltaState
sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_anisoDuctile_sizePostResults(instance) sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_anisoDuctile_sizePostResults(instance)
allocate(sourceState(phase)%p(sourceOffset)%aTolState (sizeState), & sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol
source=param(instance)%aTol)
allocate(sourceState(phase)%p(sourceOffset)%state0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%subState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%state (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 1_pInt)) then
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
endif
if (any(numerics_integrator == 4_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 5_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal)
endif endif

View File

@ -86,6 +86,7 @@ subroutine source_damage_isoBrittle_init(fileUnit)
IO_timeStamp, & IO_timeStamp, &
IO_EOF IO_EOF
use material, only: & use material, only: &
material_allocateSourceState, &
phase_source, & phase_source, &
phase_Nsources, & phase_Nsources, &
phase_Noutput, & phase_Noutput, &
@ -221,32 +222,10 @@ subroutine source_damage_isoBrittle_init(fileUnit)
source_damage_isoBrittle_sizePostResults(instance) = source_damage_isoBrittle_sizePostResults(instance) + mySize source_damage_isoBrittle_sizePostResults(instance) = source_damage_isoBrittle_sizePostResults(instance) + mySize
endif endif
enddo outputsLoop enddo outputsLoop
! Determine size of state array
sizeDotState = 1_pInt call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1_pInt)
sizeDeltaState = 1_pInt
sizeState = 1_pInt
sourceState(phase)%p(sourceOffset)%sizeState = sizeState
sourceState(phase)%p(sourceOffset)%sizeDotState = sizeDotState
sourceState(phase)%p(sourceOffset)%sizeDeltaState = sizeDeltaState
sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_isoBrittle_sizePostResults(instance) sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_isoBrittle_sizePostResults(instance)
allocate(sourceState(phase)%p(sourceOffset)%aTolState (sizeState), & sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol
source=param(instance)%aTol)
allocate(sourceState(phase)%p(sourceOffset)%state0 (sizeState,NofMyPhase), source=.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%subState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%state (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 1_pInt)) then
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
endif
if (any(numerics_integrator == 4_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 5_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal)
endif endif

View File

@ -86,6 +86,7 @@ subroutine source_damage_isoDuctile_init(fileUnit)
IO_timeStamp, & IO_timeStamp, &
IO_EOF IO_EOF
use material, only: & use material, only: &
material_allocateSourceState, &
phase_source, & phase_source, &
phase_Nsources, & phase_Nsources, &
phase_Noutput, & phase_Noutput, &
@ -221,33 +222,11 @@ subroutine source_damage_isoDuctile_init(fileUnit)
source_damage_isoDuctile_sizePostResults(instance) = source_damage_isoDuctile_sizePostResults(instance) + mySize source_damage_isoDuctile_sizePostResults(instance) = source_damage_isoDuctile_sizePostResults(instance) + mySize
endif endif
enddo outputsLoop enddo outputsLoop
! Determine size of state array
sizeDotState = 1_pInt call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1_pInt)
sizeDeltaState = 0_pInt
sizeState = 1_pInt
sourceState(phase)%p(sourceOffset)%sizeState = sizeState
sourceState(phase)%p(sourceOffset)%sizeDotState = sizeDotState
sourceState(phase)%p(sourceOffset)%sizeDeltaState = sizeDeltaState
sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_isoDuctile_sizePostResults(instance) sourceState(phase)%p(sourceOffset)%sizePostResults = source_damage_isoDuctile_sizePostResults(instance)
allocate(sourceState(phase)%p(sourceOffset)%aTolState (sizeState), & sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol
source=param(instance)%aTol)
allocate(sourceState(phase)%p(sourceOffset)%state0 (sizeState,NofMyPhase), source=.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%partionedState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%subState0 (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%state (sizeState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%deltaState (sizeDeltaState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 1_pInt)) then
allocate(sourceState(phase)%p(sourceOffset)%previousDotState (sizeDotState,NofMyPhase), source=0.0_pReal)
allocate(sourceState(phase)%p(sourceOffset)%previousDotState2 (sizeDotState,NofMyPhase), source=0.0_pReal)
endif
if (any(numerics_integrator == 4_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RK4dotState (sizeDotState,NofMyPhase), source=0.0_pReal)
if (any(numerics_integrator == 5_pInt)) &
allocate(sourceState(phase)%p(sourceOffset)%RKCK45dotState (6,sizeDotState,NofMyPhase),source=0.0_pReal)
endif endif
enddo initializeInstances enddo initializeInstances