diff --git a/src/element.f90 b/src/element.f90 index 208f6e718..dca5cbd28 100644 --- a/src/element.f90 +++ b/src/element.f90 @@ -386,15 +386,15 @@ module element - real(pReal), dimension(nNode(1),NcellNode(geomType(1))), parameter :: cellNodeParentNodeWeights1 = & - reshape(real([& + integer, dimension(nNode(1),NcellNode(geomType(1))), parameter :: cellNodeParentNodeWeights1 = & + reshape([& 1, 0, 0, & 0, 1, 0, & 0, 0, 1 & - ],pReal),[nNode(1),NcellNode(geomType(1))]) ! 2D 3node 1ip + ],[nNode(1),NcellNode(geomType(1))]) !< 2D 3node 1ip - real(pReal), dimension(nNode(2),NcellNode(geomType(2))), parameter :: cellNodeParentNodeWeights2 = & - reshape(real([& + integer, dimension(nNode(2),NcellNode(geomType(2))), parameter :: cellNodeParentNodeWeights2 = & + reshape([& 1, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, & @@ -402,10 +402,10 @@ module element 0, 0, 0, 0, 1, 0, & 0, 0, 0, 0, 0, 1, & 1, 1, 1, 2, 2, 2 & - ],pReal),[nNode(2),NcellNode(geomType(2))]) ! 2D 6node 3ip + ],[nNode(2),NcellNode(geomType(2))]) !< 2D 6node 3ip - real(pReal), dimension(nNode(3),NcellNode(geomType(3))), parameter :: cellNodeParentNodeWeights3 = & - reshape(real([& + integer, dimension(nNode(3),NcellNode(geomType(3))), parameter :: cellNodeParentNodeWeights3 = & + reshape([& 1, 0, 0, 0, & 0, 1, 0, 0, & 0, 0, 1, 0, & @@ -415,10 +415,10 @@ module element 0, 0, 1, 1, & 1, 0, 0, 1, & 1, 1, 1, 1 & - ],pReal),[nNode(3),NcellNode(geomType(3))]) ! 2D 6node 3ip + ],[nNode(3),NcellNode(geomType(3))]) !< 2D 6node 3ip - real(pReal), dimension(nNode(4),NcellNode(geomType(4))), parameter :: cellNodeParentNodeWeights4 = & - reshape(real([& + integer, dimension(nNode(4),NcellNode(geomType(4))), parameter :: cellNodeParentNodeWeights4 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, 0, 0, & @@ -435,10 +435,10 @@ module element 1, 4, 1, 1, 8, 8, 2, 2, & 1, 1, 4, 1, 2, 8, 8, 2, & 1, 1, 1, 4, 2, 2, 8, 8 & - ],pReal),[nNode(4),NcellNode(geomType(4))]) ! 2D 8node 9ip + ],[nNode(4),NcellNode(geomType(4))]) !< 2D 8node 9ip - real(pReal), dimension(nNode(5),NcellNode(geomType(5))), parameter :: cellNodeParentNodeWeights5 = & - reshape(real([& + integer, dimension(nNode(5),NcellNode(geomType(5))), parameter :: cellNodeParentNodeWeights5 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, 0, 0, & @@ -448,18 +448,18 @@ module element 0, 0, 0, 0, 0, 0, 1, 0, & 0, 0, 0, 0, 0, 0, 0, 1, & 1, 1, 1, 1, 2, 2, 2, 2 & - ],pReal),[nNode(5),NcellNode(geomType(5))]) ! 2D 8node 4ip + ],[nNode(5),NcellNode(geomType(5))]) !< 2D 8node 4ip - real(pReal), dimension(nNode(6),NcellNode(geomType(6))), parameter :: cellNodeParentNodeWeights6 = & - reshape(real([& + integer, dimension(nNode(6),NcellNode(geomType(6))), parameter :: cellNodeParentNodeWeights6 = & + reshape([& 1, 0, 0, 0, & 0, 1, 0, 0, & 0, 0, 1, 0, & 0, 0, 0, 1 & - ],pReal),[nNode(6),NcellNode(geomType(6))]) ! 3D 4node 1ip + ],[nNode(6),NcellNode(geomType(6))]) !< 3D 4node 1ip - real(pReal), dimension(nNode(7),NcellNode(geomType(7))), parameter :: cellNodeParentNodeWeights7 = & - reshape(real([& + integer, dimension(nNode(7),NcellNode(geomType(7))), parameter :: cellNodeParentNodeWeights7 = & + reshape([& 1, 0, 0, 0, 0, & 0, 1, 0, 0, 0, & 0, 0, 1, 0, 0, & @@ -475,10 +475,10 @@ module element 0, 1, 1, 1, 0, & 1, 0, 1, 1, 0, & 0, 0, 0, 0, 1 & - ],pReal),[nNode(7),NcellNode(geomType(7))]) ! 3D 5node 4ip + ],[nNode(7),NcellNode(geomType(7))]) !< 3D 5node 4ip - real(pReal), dimension(nNode(8),NcellNode(geomType(8))), parameter :: cellNodeParentNodeWeights8 = & - reshape(real([& + integer, dimension(nNode(8),NcellNode(geomType(8))), parameter :: cellNodeParentNodeWeights8 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, & @@ -494,10 +494,10 @@ module element 0, 1, 1, 1, 0, 2, 0, 0, 2, 2, & 1, 0, 1, 1, 0, 0, 2, 2, 0, 2, & 3, 3, 3, 3, 4, 4, 4, 4, 4, 4 & - ],pReal),[nNode(8),NcellNode(geomType(8))]) ! 3D 10node 4ip + ],[nNode(8),NcellNode(geomType(8))]) !< 3D 10node 4ip - real(pReal), dimension(nNode(9),NcellNode(geomType(9))), parameter :: cellNodeParentNodeWeights9 = & - reshape(real([& + integer, dimension(nNode(9),NcellNode(geomType(9))), parameter :: cellNodeParentNodeWeights9 = & + reshape([& 1, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, & @@ -519,10 +519,10 @@ module element 1, 0, 1, 1, 0, 1, & 0, 0, 0, 1, 1, 1, & 1, 1, 1, 1, 1, 1 & - ],pReal),[nNode(9),NcellNode(geomType(9))]) ! 3D 6node 6ip + ],[nNode(9),NcellNode(geomType(9))]) !< 3D 6node 6ip - real(pReal), dimension(nNode(10),NcellNode(geomType(10))), parameter :: cellNodeParentNodeWeights10 = & - reshape(real([& + integer, dimension(nNode(10),NcellNode(geomType(10))), parameter :: cellNodeParentNodeWeights10 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, & 0, 1, 0, 0, 0, 0, 0, 0, & 0, 0, 1, 0, 0, 0, 0, 0, & @@ -531,10 +531,10 @@ module element 0, 0, 0, 0, 0, 1, 0, 0, & 0, 0, 0, 0, 0, 0, 1, 0, & 0, 0, 0, 0, 0, 0, 0, 1 & - ],pReal),[nNode(10),NcellNode(geomType(10))]) ! 3D 8node 1ip + ],[nNode(10),NcellNode(geomType(10))]) !< 3D 8node 1ip - real(pReal), dimension(nNode(11),NcellNode(geomType(11))), parameter :: cellNodeParentNodeWeights11 = & - reshape(real([& + integer, dimension(nNode(11),NcellNode(geomType(11))), parameter :: cellNodeParentNodeWeights11 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, & ! 0, 1, 0, 0, 0, 0, 0, 0, & ! 0, 0, 1, 0, 0, 0, 0, 0, & ! @@ -562,10 +562,10 @@ module element 1, 0, 0, 1, 1, 0, 0, 1, & ! 25 0, 0, 0, 0, 1, 1, 1, 1, & ! 1, 1, 1, 1, 1, 1, 1, 1 & ! - ],pReal),[nNode(11),NcellNode(geomType(11))]) ! 3D 8node 8ip + ],[nNode(11),NcellNode(geomType(11))]) !< 3D 8node 8ip - real(pReal), dimension(nNode(12),NcellNode(geomType(12))), parameter :: cellNodeParentNodeWeights12 = & - reshape(real([& + integer, dimension(nNode(12),NcellNode(geomType(12))), parameter :: cellNodeParentNodeWeights12 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! @@ -593,10 +593,10 @@ module element 1, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 2, 0, 0, 0, 2, 2, 0, 0, 2, & ! 25 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 2, 2, 2, 2, 0, 0, 0, 0, & ! 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4 & ! - ],pReal),[nNode(12),NcellNode(geomType(12))]) ! 3D 20node 8ip + ],[nNode(12),NcellNode(geomType(12))]) !< 3D 20node 8ip - real(pReal), dimension(nNode(13),NcellNode(geomType(13))), parameter :: cellNodeParentNodeWeights13 = & - reshape(real([& + integer, dimension(nNode(13),NcellNode(geomType(13))), parameter :: cellNodeParentNodeWeights13 = & + reshape([& 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, & ! @@ -661,7 +661,7 @@ module element 4, 8, 4, 3, 8,24, 8, 4, 12,12, 4, 4, 32,32,12,12, 12,32,12, 4, & ! 3, 4, 8, 4, 4, 8,24, 8, 4,12,12, 4, 12,32,32,12, 4,12,32,12, & ! 4, 3, 4, 8, 8, 4, 8,24, 4, 4,12,12, 12,12,32,32, 12, 4,12,32 & ! - ],pReal),[nNode(13),NcellNode(geomType(13))]) ! 3D 20node 27ip + ],[nNode(13),NcellNode(geomType(13))]) !< 3D 20node 27ip integer, dimension(NCELLNODEPERCELL(CELLTYPE(1)),NIP(1)), parameter :: CELL1 = & diff --git a/src/results.f90 b/src/results.f90 index 0b9bec9f1..0304fd30b 100644 --- a/src/results.f90 +++ b/src/results.f90 @@ -17,8 +17,7 @@ module results private #if defined(PETSc) || defined(DAMASK_HDF5) - integer(HID_T), public, protected :: tempCoordinates, tempResults - integer(HID_T), private :: resultsFile, currentIncID, plist_id + integer(HID_T) :: resultsFile interface results_writeDataset diff --git a/src/source_damage_anisoDuctile.f90 b/src/source_damage_anisoDuctile.f90 index a4b4561e1..aa8c9b8fa 100644 --- a/src/source_damage_anisoDuctile.f90 +++ b/src/source_damage_anisoDuctile.f90 @@ -5,55 +5,55 @@ !> @details to be done !-------------------------------------------------------------------------------------------------- module source_damage_anisoDuctile - use prec - use debug - use IO - use math - use material - use config - - implicit none - private - - integer, dimension(:), allocatable, public, protected :: & - source_damage_anisoDuctile_offset, & !< which source is my current damage mechanism? - source_damage_anisoDuctile_instance !< instance of damage source mechanism - - integer, dimension(:,:), allocatable, target, public :: & - source_damage_anisoDuctile_sizePostResult !< size of each post result output - - character(len=64), dimension(:,:), allocatable, target, public :: & - source_damage_anisoDuctile_output !< name of each post result output + use prec + use debug + use IO + use math + use material + use config + + implicit none + private - - enum, bind(c) - enumerator :: undefined_ID, & - damage_drivingforce_ID - end enum - - - type, private :: tParameters !< container type for internal constitutive parameters - real(pReal) :: & - aTol, & - N - real(pReal), dimension(:), allocatable :: & - critPlasticStrain - integer :: & - totalNslip - integer, dimension(:), allocatable :: & - Nslip - integer(kind(undefined_ID)), allocatable, dimension(:) :: & - outputID - end type tParameters - - type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) - - - public :: & - source_damage_anisoDuctile_init, & - source_damage_anisoDuctile_dotState, & - source_damage_anisoDuctile_getRateAndItsTangent, & - source_damage_anisoDuctile_postResults + integer, dimension(:), allocatable, public, protected :: & + source_damage_anisoDuctile_offset, & !< which source is my current damage mechanism? + source_damage_anisoDuctile_instance !< instance of damage source mechanism + + integer, dimension(:,:), allocatable, target, public :: & + source_damage_anisoDuctile_sizePostResult !< size of each post result output + + character(len=64), dimension(:,:), allocatable, target, public :: & + source_damage_anisoDuctile_output !< name of each post result output + + + enum, bind(c) + enumerator :: undefined_ID, & + damage_drivingforce_ID + end enum + + + type, private :: tParameters !< container type for internal constitutive parameters + real(pReal) :: & + aTol, & + N + real(pReal), dimension(:), allocatable :: & + critPlasticStrain + integer :: & + totalNslip + integer, dimension(:), allocatable :: & + Nslip + integer(kind(undefined_ID)), allocatable, dimension(:) :: & + outputID + end type tParameters + + type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) + + + public :: & + source_damage_anisoDuctile_init, & + source_damage_anisoDuctile_dotState, & + source_damage_anisoDuctile_getRateAndItsTangent, & + source_damage_anisoDuctile_postResults contains @@ -64,193 +64,196 @@ contains !-------------------------------------------------------------------------------------------------- subroutine source_damage_anisoDuctile_init - integer :: Ninstance,phase,instance,source,sourceOffset - integer :: NofMyPhase,p ,i - - integer, dimension(0), parameter :: emptyIntArray = [integer::] - character(len=65536), dimension(0), parameter :: emptyStringArray = [character(len=65536)::] - integer(kind(undefined_ID)) :: & - outputID - - character(len=pStringLen) :: & - extmsg = '' - character(len=65536), dimension(:), allocatable :: & - outputs - - write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ANISODUCTILE_LABEL//' init -+>>>' - - Ninstance = count(phase_source == SOURCE_damage_anisoDuctile_ID) - if (Ninstance == 0) return + integer :: Ninstance,phase,instance,source,sourceOffset + integer :: NofMyPhase,p ,i - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & - write(6,'(a16,1x,i5,/)') '# instances:',Ninstance + integer, dimension(0), parameter :: emptyIntArray = [integer::] + character(len=65536), dimension(0), parameter :: emptyStringArray = [character(len=65536)::] + integer(kind(undefined_ID)) :: & + outputID - allocate(source_damage_anisoDuctile_offset(size(config_phase)), source=0) - allocate(source_damage_anisoDuctile_instance(size(config_phase)), source=0) - do phase = 1, size(config_phase) - source_damage_anisoDuctile_instance(phase) = count(phase_source(:,1:phase) == source_damage_anisoDuctile_ID) - do source = 1, phase_Nsources(phase) - if (phase_source(source,phase) == source_damage_anisoDuctile_ID) & - source_damage_anisoDuctile_offset(phase) = source - enddo - enddo - - allocate(source_damage_anisoDuctile_sizePostResult(maxval(phase_Noutput),Ninstance),source=0) - allocate(source_damage_anisoDuctile_output(maxval(phase_Noutput),Ninstance)) - source_damage_anisoDuctile_output = '' - - - allocate(param(Ninstance)) + character(len=pStringLen) :: & + extmsg = '' + character(len=65536), dimension(:), allocatable :: & + outputs - do p=1, size(config_phase) - if (all(phase_source(:,p) /= SOURCE_DAMAGE_ANISODUCTILE_ID)) cycle - associate(prm => param(source_damage_anisoDuctile_instance(p)), & - config => config_phase(p)) - - prm%aTol = config%getFloat('anisoductile_atol',defaultVal = 1.0e-3_pReal) - - prm%N = config%getFloat('anisoductile_ratesensitivity') - prm%totalNslip = sum(prm%Nslip) - ! sanity checks - if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_atol' - - if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_ratesensitivity' - - prm%Nslip = config%getInts('nslip',defaultVal=emptyIntArray) - - prm%critPlasticStrain = config%getFloats('anisoductile_criticalplasticstrain',requiredSize=size(prm%Nslip)) - - ! expand: family => system - prm%critPlasticStrain = math_expand(prm%critPlasticStrain, prm%Nslip) - - if (any(prm%critPlasticStrain < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_criticalplasticstrain' + write(6,'(/,a)') ' <<<+- source_'//SOURCE_DAMAGE_ANISODUCTILE_LABEL//' init -+>>>' + + Ninstance = count(phase_source == SOURCE_damage_anisoDuctile_ID) + if (Ninstance == 0) return + + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + write(6,'(a16,1x,i5,/)') '# instances:',Ninstance + + allocate(source_damage_anisoDuctile_offset(size(config_phase)), source=0) + allocate(source_damage_anisoDuctile_instance(size(config_phase)), source=0) + do phase = 1, size(config_phase) + source_damage_anisoDuctile_instance(phase) = count(phase_source(:,1:phase) == source_damage_anisoDuctile_ID) + do source = 1, phase_Nsources(phase) + if (phase_source(source,phase) == source_damage_anisoDuctile_ID) & + source_damage_anisoDuctile_offset(phase) = source + enddo + enddo + + allocate(source_damage_anisoDuctile_sizePostResult(maxval(phase_Noutput),Ninstance),source=0) + allocate(source_damage_anisoDuctile_output(maxval(phase_Noutput),Ninstance)) + source_damage_anisoDuctile_output = '' + + + allocate(param(Ninstance)) + + do p=1, size(config_phase) + if (all(phase_source(:,p) /= SOURCE_DAMAGE_ANISODUCTILE_ID)) cycle + associate(prm => param(source_damage_anisoDuctile_instance(p)), & + config => config_phase(p)) + + prm%aTol = config%getFloat('anisoductile_atol',defaultVal = 1.0e-3_pReal) + + prm%N = config%getFloat('anisoductile_ratesensitivity') + prm%totalNslip = sum(prm%Nslip) + ! sanity checks + if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_atol' + + if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_ratesensitivity' + + prm%Nslip = config%getInts('nslip',defaultVal=emptyIntArray) + + prm%critPlasticStrain = config%getFloats('anisoductile_criticalplasticstrain',requiredSize=size(prm%Nslip)) + + ! expand: family => system + prm%critPlasticStrain = math_expand(prm%critPlasticStrain, prm%Nslip) + + if (any(prm%critPlasticStrain < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_criticalplasticstrain' !-------------------------------------------------------------------------------------------------- ! exit if any parameter is out of range - if (extmsg /= '') & - call IO_error(211,ext_msg=trim(extmsg)//'('//SOURCE_DAMAGE_ANISODUCTILE_LABEL//')') + if (extmsg /= '') call IO_error(211,ext_msg=trim(extmsg)//'('//SOURCE_DAMAGE_ANISODUCTILE_LABEL//')') !-------------------------------------------------------------------------------------------------- ! output pararameters - outputs = config%getStrings('(output)',defaultVal=emptyStringArray) - allocate(prm%outputID(0)) - do i=1, size(outputs) - outputID = undefined_ID - select case(outputs(i)) - - case ('anisoductile_drivingforce') - source_damage_anisoDuctile_sizePostResult(i,source_damage_anisoDuctile_instance(p)) = 1 - source_damage_anisoDuctile_output(i,source_damage_anisoDuctile_instance(p)) = outputs(i) - prm%outputID = [prm%outputID, damage_drivingforce_ID] - - end select - - enddo - - end associate - - phase = p - - NofMyPhase=count(material_phase==phase) - instance = source_damage_anisoDuctile_instance(phase) - sourceOffset = source_damage_anisoDuctile_offset(phase) - - call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1,1,0) - sourceState(phase)%p(sourceOffset)%sizePostResults = sum(source_damage_anisoDuctile_sizePostResult(:,instance)) - sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol - - enddo + outputs = config%getStrings('(output)',defaultVal=emptyStringArray) + allocate(prm%outputID(0)) + do i=1, size(outputs) + outputID = undefined_ID + select case(outputs(i)) + + case ('anisoductile_drivingforce') + source_damage_anisoDuctile_sizePostResult(i,source_damage_anisoDuctile_instance(p)) = 1 + source_damage_anisoDuctile_output(i,source_damage_anisoDuctile_instance(p)) = outputs(i) + prm%outputID = [prm%outputID, damage_drivingforce_ID] + + end select + + enddo + + end associate + + phase = p + + NofMyPhase=count(material_phase==phase) + instance = source_damage_anisoDuctile_instance(phase) + sourceOffset = source_damage_anisoDuctile_offset(phase) + + call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1,1,0) + sourceState(phase)%p(sourceOffset)%sizePostResults = sum(source_damage_anisoDuctile_sizePostResult(:,instance)) + sourceState(phase)%p(sourceOffset)%aTolState=param(instance)%aTol + + enddo end subroutine source_damage_anisoDuctile_init + !-------------------------------------------------------------------------------------------------- !> @brief calculates derived quantities from state !-------------------------------------------------------------------------------------------------- subroutine source_damage_anisoDuctile_dotState(ipc, ip, el) - integer, intent(in) :: & - ipc, & !< component-ID of integration point - ip, & !< integration point - el !< element - integer :: & - phase, & - constituent, & - sourceOffset, & - homog, damageOffset, & - instance, & - f, i - - phase = phaseAt(ipc,ip,el) - constituent = phasememberAt(ipc,ip,el) - instance = source_damage_anisoDuctile_instance(phase) - sourceOffset = source_damage_anisoDuctile_offset(phase) - homog = material_homogenizationAt(el) - damageOffset = damageMapping(homog)%p(ip,el) - - - do i = 1, param(instance)%totalNslip - sourceState(phase)%p(sourceOffset)%dotState(1,constituent) = & - sourceState(phase)%p(sourceOffset)%dotState(1,constituent) + & - plasticState(phase)%slipRate(i,constituent)/ & - ((damage(homog)%p(damageOffset))**param(instance)%N)/param(instance)%critPlasticStrain(i) - enddo + integer, intent(in) :: & + ipc, & !< component-ID of integration point + ip, & !< integration point + el !< element + integer :: & + phase, & + constituent, & + sourceOffset, & + homog, damageOffset, & + instance, & + i + + phase = phaseAt(ipc,ip,el) + constituent = phasememberAt(ipc,ip,el) + instance = source_damage_anisoDuctile_instance(phase) + sourceOffset = source_damage_anisoDuctile_offset(phase) + homog = material_homogenizationAt(el) + damageOffset = damageMapping(homog)%p(ip,el) + + + do i = 1, param(instance)%totalNslip + sourceState(phase)%p(sourceOffset)%dotState(1,constituent) = & + sourceState(phase)%p(sourceOffset)%dotState(1,constituent) + & + plasticState(phase)%slipRate(i,constituent)/ & + ((damage(homog)%p(damageOffset))**param(instance)%N)/param(instance)%critPlasticStrain(i) + enddo end subroutine source_damage_anisoDuctile_dotState + !-------------------------------------------------------------------------------------------------- !> @brief returns local part of nonlocal damage driving force !-------------------------------------------------------------------------------------------------- subroutine source_damage_anisoDuctile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) - integer, intent(in) :: & - phase, & - constituent - real(pReal), intent(in) :: & - phi - real(pReal), intent(out) :: & - localphiDot, & - dLocalphiDot_dPhi - integer :: & - sourceOffset - - sourceOffset = source_damage_anisoDuctile_offset(phase) + integer, intent(in) :: & + phase, & + constituent + real(pReal), intent(in) :: & + phi + real(pReal), intent(out) :: & + localphiDot, & + dLocalphiDot_dPhi + integer :: & + sourceOffset - localphiDot = 1.0_pReal & - - sourceState(phase)%p(sourceOffset)%state(1,constituent) * phi - - dLocalphiDot_dPhi = -sourceState(phase)%p(sourceOffset)%state(1,constituent) + sourceOffset = source_damage_anisoDuctile_offset(phase) + + localphiDot = 1.0_pReal & + - sourceState(phase)%p(sourceOffset)%state(1,constituent) * phi + + dLocalphiDot_dPhi = -sourceState(phase)%p(sourceOffset)%state(1,constituent) end subroutine source_damage_anisoDuctile_getRateAndItsTangent - + + !-------------------------------------------------------------------------------------------------- !> @brief return array of local damage results !-------------------------------------------------------------------------------------------------- function source_damage_anisoDuctile_postResults(phase, constituent) - integer, intent(in) :: & - phase, & - constituent - real(pReal), dimension(sum(source_damage_anisoDuctile_sizePostResult(:, & - source_damage_anisoDuctile_instance(phase)))) :: & - source_damage_anisoDuctile_postResults + integer, intent(in) :: & + phase, & + constituent + real(pReal), dimension(sum(source_damage_anisoDuctile_sizePostResult(:, & + source_damage_anisoDuctile_instance(phase)))) :: & + source_damage_anisoDuctile_postResults + + integer :: & + instance, sourceOffset, o, c + + instance = source_damage_anisoDuctile_instance(phase) + sourceOffset = source_damage_anisoDuctile_offset(phase) + + c = 0 + + do o = 1,size(param(instance)%outputID) + select case(param(instance)%outputID(o)) + case (damage_drivingforce_ID) + source_damage_anisoDuctile_postResults(c+1) = & + sourceState(phase)%p(sourceOffset)%state(1,constituent) + c = c + 1 + + end select + enddo - integer :: & - instance, sourceOffset, o, c - - instance = source_damage_anisoDuctile_instance(phase) - sourceOffset = source_damage_anisoDuctile_offset(phase) - - c = 0 - - do o = 1,size(param(instance)%outputID) - select case(param(instance)%outputID(o)) - case (damage_drivingforce_ID) - source_damage_anisoDuctile_postResults(c+1) = & - sourceState(phase)%p(sourceOffset)%state(1,constituent) - c = c + 1 - - end select - enddo end function source_damage_anisoDuctile_postResults end module source_damage_anisoDuctile