consistent data structure

This commit is contained in:
Martin Diehl 2021-05-22 18:42:18 +02:00
parent 6d31f77deb
commit d41f3a5490
4 changed files with 81 additions and 97 deletions

View File

@ -20,7 +20,7 @@ module material
type(Rotation), dimension(:), allocatable :: data type(Rotation), dimension(:), allocatable :: data
end type end type
type(tRotationContainer), dimension(:), allocatable :: material_orientation0_2 type(tRotationContainer), dimension(:), allocatable :: material_orientation0
integer, dimension(:), allocatable, public, protected :: & integer, dimension(:), allocatable, public, protected :: &
homogenization_Nconstituents !< number of grains in each homogenization homogenization_Nconstituents !< number of grains in each homogenization
@ -46,7 +46,7 @@ module material
public :: & public :: &
tRotationContainer, & tRotationContainer, &
material_orientation0_2, & material_orientation0, &
material_init material_init
contains contains
@ -161,14 +161,15 @@ subroutine parse()
enddo enddo
allocate(material_orientation0_2(materials%length)) allocate(material_orientation0(materials%length))
do ma = 1, materials%length do ma = 1, materials%length
material => materials%get(ma) material => materials%get(ma)
constituents => material%get('constituents') constituents => material%get('constituents')
allocate(material_orientation0_2(ma)%data(constituents%length)) allocate(material_orientation0(ma)%data(constituents%length))
do co = 1, constituents%length do co = 1, constituents%length
call material_orientation0_2(ma)%data(co)%fromQuaternion(constituent%get_as1dFloat('O',requiredSize=4)) constituent => constituents%get(co)
call material_orientation0(ma)%data(co)%fromQuaternion(constituent%get_as1dFloat('O',requiredSize=4))
enddo enddo
enddo enddo

View File

@ -23,9 +23,6 @@ module phase
phase_orientation0, & phase_orientation0, &
phase_orientation phase_orientation
type(rotation), dimension(:,:,:), allocatable :: &
crystallite_orientation !< current orientation
type :: tTensorContainer type :: tTensorContainer
real(pReal), dimension(:,:,:), allocatable :: data real(pReal), dimension(:,:,:), allocatable :: data
end type end type
@ -259,8 +256,7 @@ module phase
ph, & ph, &
i, & i, &
e e
type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & type(tRotationContainer), dimension(:), intent(in) :: orientation
orientation !< crystal orientation
end subroutine plastic_nonlocal_updateCompatibility end subroutine plastic_nonlocal_updateCompatibility
module subroutine plastic_dependentState(co,ip,el) module subroutine plastic_dependentState(co,ip,el)
@ -378,7 +374,7 @@ subroutine phase_init
ma = discretization_materialAt((ce-1)/discretization_nIPs+1) ma = discretization_materialAt((ce-1)/discretization_nIPs+1)
do co = 1,homogenization_Nconstituents(material_homogenizationID(ce)) do co = 1,homogenization_Nconstituents(material_homogenizationID(ce))
ph = material_phaseID(co,ce) ph = material_phaseID(co,ce)
phase_orientation0(ph)%data(material_phaseEntry(co,ce)) = material_orientation0_2(ma)%data(co) phase_orientation0(ph)%data(material_phaseEntry(co,ce)) = material_orientation0(ma)%data(co)
enddo enddo
enddo enddo
@ -523,8 +519,6 @@ subroutine crystallite_init()
iMax = discretization_nIPs iMax = discretization_nIPs
eMax = discretization_Nelems eMax = discretization_Nelems
allocate(crystallite_orientation(cMax,iMax,eMax))
num_crystallite => config_numerics%get('crystallite',defaultVal=emptyDict) num_crystallite => config_numerics%get('crystallite',defaultVal=emptyDict)
num%subStepMinCryst = num_crystallite%get_asFloat ('subStepMin', defaultVal=1.0e-3_pReal) num%subStepMinCryst = num_crystallite%get_asFloat ('subStepMin', defaultVal=1.0e-3_pReal)
@ -594,13 +588,16 @@ subroutine crystallite_orientations(co,ip,el)
ip, & !< counter in integration point loop ip, & !< counter in integration point loop
el !< counter in element loop el !< counter in element loop
integer :: ph, en
call crystallite_orientation(co,ip,el)%fromMatrix(transpose(math_rotationalPart(&
mechanical_F_e(material_phaseAt(co,el),material_phaseMemberAt(co,ip,el))))) ph = material_phaseID(co,(el-1)*discretization_nIPs + ip)
en = material_phaseEntry(co,(el-1)*discretization_nIPs + ip)
call phase_orientation(ph)%data(en)%fromMatrix(transpose(math_rotationalPart(mechanical_F_e(ph,en))))
if (plasticState(material_phaseAt(1,el))%nonlocal) & if (plasticState(material_phaseAt(1,el))%nonlocal) &
call plastic_nonlocal_updateCompatibility(crystallite_orientation, & call plastic_nonlocal_updateCompatibility(phase_orientation,material_phaseAt(1,el),ip,el)
material_phaseAt(1,el),ip,el)
end subroutine crystallite_orientations end subroutine crystallite_orientations

View File

@ -248,18 +248,14 @@ module subroutine mechanical_init(materials,phases)
phase => phases%get(ph) phase => phases%get(ph)
mech => phase%get('mechanical') mech => phase%get('mechanical')
#if defined(__GFORTRAN__) #if defined(__GFORTRAN__)
output_constituent(ph)%label = output_as1dString(mech) output_constituent(ph)%label = output_as1dString(mech)
#else #else
output_constituent(ph)%label = mech%get_as1dString('output',defaultVal=emptyStringArray) output_constituent(ph)%label = mech%get_as1dString('output',defaultVal=emptyStringArray)
#endif #endif
enddo enddo
do el = 1, size(material_phaseMemberAt,3); do ip = 1, size(material_phaseMemberAt,2) do ph = 1, phases%length
do co = 1, homogenization_Nconstituents(material_homogenizationAt(el)) do en = 1, count(material_phaseID == ph)
material => materials%get(discretization_materialAt(el))
ph = material_phaseID(co,(el-1)*discretization_nIPs + ip)
en = material_phaseEntry(co,(el-1)*discretization_nIPs + ip)
phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_orientation0(ph)%data(en)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_orientation0(ph)%data(en)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005)
phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_mechanical_Fp0(ph)%data(1:3,1:3,en) & phase_mechanical_Fp0(ph)%data(1:3,1:3,en) = phase_mechanical_Fp0(ph)%data(1:3,1:3,en) &
@ -268,13 +264,12 @@ module subroutine mechanical_init(materials,phases)
phase_mechanical_F0(ph)%data(1:3,1:3,en) = math_I3 phase_mechanical_F0(ph)%data(1:3,1:3,en) = math_I3
phase_mechanical_Fe(ph)%data(1:3,1:3,en) = math_inv33(matmul(phase_mechanical_Fi0(ph)%data(1:3,1:3,en), & phase_mechanical_Fe(ph)%data(1:3,1:3,en) = math_inv33(matmul(phase_mechanical_Fi0(ph)%data(1:3,1:3,en), &
phase_mechanical_Fp0(ph)%data(1:3,1:3,en))) ! assuming that euler angles are given in internal strain free configuration phase_mechanical_Fp0(ph)%data(1:3,1:3,en))) ! assuming that euler angles are given in internal strain free configuration
phase_mechanical_Fp(ph)%data(1:3,1:3,en) = phase_mechanical_Fp0(ph)%data(1:3,1:3,en)
phase_mechanical_Fi(ph)%data(1:3,1:3,en) = phase_mechanical_Fi0(ph)%data(1:3,1:3,en)
phase_mechanical_F(ph)%data(1:3,1:3,en) = phase_mechanical_F0(ph)%data(1:3,1:3,en)
enddo enddo
enddo; enddo phase_mechanical_Fp(ph)%data = phase_mechanical_Fp0(ph)%data
phase_mechanical_Fi(ph)%data = phase_mechanical_Fi0(ph)%data
phase_mechanical_F(ph)%data = phase_mechanical_F0(ph)%data
enddo
! initialize elasticity ! initialize elasticity
@ -903,80 +898,71 @@ subroutine crystallite_results(group,ph)
character(len=:), allocatable :: structureLabel character(len=:), allocatable :: structureLabel
call results_closeGroup(results_addGroup(group//'/mechanical')) call results_closeGroup(results_addGroup(group//'/mechanical'))
do ou = 1, size(output_constituent(ph)%label) do ou = 1, size(output_constituent(ph)%label)
select case (output_constituent(ph)%label(ou)) select case (output_constituent(ph)%label(ou))
case('F') case('F')
call results_writeDataset(group//'/mechanical/',phase_mechanical_F(ph)%data,'F',& call results_writeDataset(group//'/mechanical/',phase_mechanical_F(ph)%data,'F',&
'deformation gradient','1') 'deformation gradient','1')
case('F_e') case('F_e')
call results_writeDataset(group//'/mechanical/',phase_mechanical_Fe(ph)%data,'F_e',& call results_writeDataset(group//'/mechanical/',phase_mechanical_Fe(ph)%data,'F_e',&
'elastic deformation gradient','1') 'elastic deformation gradient','1')
case('F_p') case('F_p')
call results_writeDataset(group//'/mechanical/',phase_mechanical_Fp(ph)%data,'F_p', & call results_writeDataset(group//'/mechanical/',phase_mechanical_Fp(ph)%data,'F_p', &
'plastic deformation gradient','1') 'plastic deformation gradient','1')
case('F_i') case('F_i')
call results_writeDataset(group//'/mechanical/',phase_mechanical_Fi(ph)%data,'F_i', & call results_writeDataset(group//'/mechanical/',phase_mechanical_Fi(ph)%data,'F_i', &
'inelastic deformation gradient','1') 'inelastic deformation gradient','1')
case('L_p') case('L_p')
call results_writeDataset(group//'/mechanical/',phase_mechanical_Lp(ph)%data,'L_p', & call results_writeDataset(group//'/mechanical/',phase_mechanical_Lp(ph)%data,'L_p', &
'plastic velocity gradient','1/s') 'plastic velocity gradient','1/s')
case('L_i') case('L_i')
call results_writeDataset(group//'/mechanical/',phase_mechanical_Li(ph)%data,'L_i', & call results_writeDataset(group//'/mechanical/',phase_mechanical_Li(ph)%data,'L_i', &
'inelastic velocity gradient','1/s') 'inelastic velocity gradient','1/s')
case('P') case('P')
call results_writeDataset(group//'/mechanical/',phase_mechanical_P(ph)%data,'P', & call results_writeDataset(group//'/mechanical/',phase_mechanical_P(ph)%data,'P', &
'first Piola-Kirchhoff stress','Pa') 'first Piola-Kirchhoff stress','Pa')
case('S') case('S')
call results_writeDataset(group//'/mechanical/',phase_mechanical_S(ph)%data,'S', & call results_writeDataset(group//'/mechanical/',phase_mechanical_S(ph)%data,'S', &
'second Piola-Kirchhoff stress','Pa') 'second Piola-Kirchhoff stress','Pa')
case('O') case('O')
select case(lattice_structure(ph)) select case(lattice_structure(ph))
case(lattice_ISO_ID) case(lattice_ISO_ID)
structureLabel = 'aP' structureLabel = 'aP'
case(lattice_FCC_ID) case(lattice_FCC_ID)
structureLabel = 'cF' structureLabel = 'cF'
case(lattice_BCC_ID) case(lattice_BCC_ID)
structureLabel = 'cI' structureLabel = 'cI'
case(lattice_BCT_ID) case(lattice_BCT_ID)
structureLabel = 'tI' structureLabel = 'tI'
case(lattice_HEX_ID) case(lattice_HEX_ID)
structureLabel = 'hP' structureLabel = 'hP'
case(lattice_ORT_ID) case(lattice_ORT_ID)
structureLabel = 'oP' structureLabel = 'oP'
end select end select
selected_rotations = select_rotations(crystallite_orientation,ph) selected_rotations = select_rotations(phase_orientation(ph)%data)
call results_writeDataset(group//'/mechanical',selected_rotations,output_constituent(ph)%label(ou),& call results_writeDataset(group//'/mechanical',selected_rotations,output_constituent(ph)%label(ou),&
'crystal orientation as quaternion','q_0 (q_1 q_2 q_3)') 'crystal orientation as quaternion','q_0 (q_1 q_2 q_3)')
call results_addAttribute('lattice',structureLabel,group//'/mechanical/'//output_constituent(ph)%label(ou)) call results_addAttribute('lattice',structureLabel,group//'/mechanical/'//output_constituent(ph)%label(ou))
end select end select
enddo enddo
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief select rotations for output !> @brief Convert orientation for output: ToDo: implement in HDF5/results
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
function select_rotations(dataset,ph) function select_rotations(dataset)
integer, intent(in) :: ph type(rotation), dimension(:), intent(in) :: dataset
type(rotation), dimension(:,:,:), intent(in) :: dataset real(pReal), dimension(4,size(dataset,1)) :: select_rotations
real(pReal), dimension(4,count(material_phaseID==ph)) :: select_rotations integer :: en
integer :: el,ip,co,j
j=0 do en = 1, size(dataset,1)
do el = 1, size(material_phaseAt,2) select_rotations(:,en) = dataset(en)%asQuaternion()
do ip = 1, discretization_nIPs
do co = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains
if (material_phaseAt(co,el) == ph) then
j = j + 1
select_rotations(1:4,j) = dataset(co,ip,el)%asQuaternion()
endif
enddo
enddo
enddo enddo
end function select_rotations end function select_rotations

View File

@ -1394,7 +1394,7 @@ end function rhoDotFlux
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine plastic_nonlocal_updateCompatibility(orientation,ph,i,e) module subroutine plastic_nonlocal_updateCompatibility(orientation,ph,i,e)
type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & type(tRotationContainer), dimension(:), intent(in) :: &
orientation ! crystal orientation orientation ! crystal orientation
integer, intent(in) :: & integer, intent(in) :: &
ph, & ph, &
@ -1464,7 +1464,7 @@ module subroutine plastic_nonlocal_updateCompatibility(orientation,ph,i,e)
!* the number of compatible slip systems is minimized with the sum of the original compatibility values exceeding one. !* the number of compatible slip systems is minimized with the sum of the original compatibility values exceeding one.
!* Finally the smallest compatibility value is decreased until the sum is exactly equal to one. !* Finally the smallest compatibility value is decreased until the sum is exactly equal to one.
!* All values below the threshold are set to zero. !* All values below the threshold are set to zero.
mis = orientation(1,i,e)%misorientation(orientation(1,neighbor_i,neighbor_e)) mis = orientation(ph)%data(en)%misorientation(orientation(neighbor_phase)%data(neighbor_me))
mySlipSystems: do s1 = 1,ns mySlipSystems: do s1 = 1,ns
neighborSlipSystems: do s2 = 1,ns neighborSlipSystems: do s2 = 1,ns
my_compatibility(1,s2,s1,n) = math_inner(prm%slip_normal(1:3,s1), & my_compatibility(1,s2,s1,n) = math_inner(prm%slip_normal(1:3,s1), &