using new mappings

This commit is contained in:
Martin Diehl 2021-05-22 13:18:33 +02:00
parent d41f3a5490
commit ee80efd705
3 changed files with 27 additions and 41 deletions

View File

@ -68,8 +68,8 @@ subroutine material_init(restart)
if (.not. restart) then if (.not. restart) then
call results_openJobFile call results_openJobFile
call results_mapping_phase(material_phaseAt,material_phaseMemberAt,material_name_phase) call results_mapping_phase(material_phaseID,material_phaseMemberAt,material_name_phase)
call results_mapping_homogenization(material_homogenizationAt,material_homogenizationMemberAt,material_name_homogenization) call results_mapping_homogenization(material_homogenizationID,material_homogenizationMemberAt,material_name_homogenization)
call results_closeJobFile call results_closeJobFile
endif endif

View File

@ -92,7 +92,7 @@ module function anisobrittle_init() result(mySources)
if (any(prm%g_crit < 0.0_pReal)) extmsg = trim(extmsg)//' g_crit' if (any(prm%g_crit < 0.0_pReal)) extmsg = trim(extmsg)//' g_crit'
if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit' if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit'
Nmembers = count(material_phaseAt==p) * discretization_nIPs Nmembers = count(material_phaseID==p)
call phase_allocateState(damageState(p),Nmembers,1,1,0) call phase_allocateState(damageState(p),Nmembers,1,1,0)
damageState(p)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) damageState(p)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal)
if(any(damageState(p)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' if(any(damageState(p)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol'

View File

@ -415,16 +415,15 @@ end subroutine results_writeTensorDataset_int
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief adds the unique mapping from spatial position and constituent ID to results !> @brief adds the unique mapping from spatial position and constituent ID to results
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine results_mapping_phase(phaseAt,memberAtLocal,label) subroutine results_mapping_phase(phaseID,memberAtLocal,label)
integer, dimension(:,:), intent(in) :: phaseAt !< phase section at (constituent,element) integer, dimension(:,:), intent(in) :: phaseID !< phase ID at (constituent,ce)
integer, dimension(:,:,:), intent(in) :: memberAtLocal !< phase member at (constituent,IP,element) integer, dimension(:,:,:), intent(in) :: memberAtLocal !< phase entry at (constituent,IP,element)
character(len=*), dimension(:), intent(in) :: label !< label of each phase section character(len=*), dimension(:), intent(in) :: label !< label of each phase section
integer, dimension(size(memberAtLocal,1),size(memberAtLocal,2),size(memberAtLocal,3)) :: & integer, dimension(size(memberAtLocal,1),size(memberAtLocal,2),size(memberAtLocal,3)) :: &
phaseAtMaterialpoint, &
memberAtGlobal memberAtGlobal
integer, dimension(size(label),0:worldsize-1) :: memberOffset !< offset in member counting per process integer, dimension(size(label),0:worldsize-1) :: memberOffset !< offset in entry counting per process
integer, dimension(0:worldsize-1) :: writeSize !< amount of data written per process integer, dimension(0:worldsize-1) :: writeSize !< amount of data written per process
integer(HSIZE_T), dimension(2) :: & integer(HSIZE_T), dimension(2) :: &
myShape, & !< shape of the dataset (this process) myShape, & !< shape of the dataset (this process)
@ -451,10 +450,10 @@ subroutine results_mapping_phase(phaseAt,memberAtLocal,label)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
memberOffset = 0 memberOffset = 0
do i=1, size(label) do i=1, size(label)
memberOffset(i,worldrank) = count(phaseAt == i)*size(memberAtLocal,2) ! number of points/instance of this process memberOffset(i,worldrank) = count(phaseID == i) ! number of entries of this process
enddo enddo
writeSize = 0 writeSize = 0
writeSize(worldrank) = size(memberAtLocal(1,:,:)) ! total number of points by this process writeSize(worldrank) = size(memberAtLocal(1,:,:)) ! total number of entries of this process
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! MPI settings and communication ! MPI settings and communication
@ -469,21 +468,15 @@ subroutine results_mapping_phase(phaseAt,memberAtLocal,label)
if(ierr /= 0) error stop 'MPI error' if(ierr /= 0) error stop 'MPI error'
#endif #endif
myShape = int([size(phaseAt,1),writeSize(worldrank)], HSIZE_T) myShape = int([size(phaseID,1),writeSize(worldrank)], HSIZE_T)
myOffset = int([0,sum(writeSize(0:worldrank-1))], HSIZE_T) myOffset = int([0,sum(writeSize(0:worldrank-1))], HSIZE_T)
totalShape = int([size(phaseAt,1),sum(writeSize)], HSIZE_T) totalShape = int([size(phaseID,1),sum(writeSize)], HSIZE_T)
!---------------------------------------------------------------------------------------------------
! expand phaseAt to consider IPs (is not stored per IP)
do i = 1, size(phaseAtMaterialpoint,2)
phaseAtMaterialpoint(:,i,:) = phaseAt
enddo
!--------------------------------------------------------------------------------------------------- !---------------------------------------------------------------------------------------------------
! renumber member from my process to all processes ! renumber member from my process to all processes
do i = 1, size(label) do i = 1, size(label)
where(phaseAtMaterialpoint == i) memberAtGlobal = memberAtLocal + sum(memberOffset(i,0:worldrank-1)) -1 ! convert to 0-based where(reshape(phaseID,shape(memberAtGlobal)) == i) &
memberAtGlobal = memberAtLocal + sum(memberOffset(i,0:worldrank-1)) -1 ! 0-based
enddo enddo
!--------------------------------------------------------------------------------------------------- !---------------------------------------------------------------------------------------------------
@ -540,7 +533,7 @@ subroutine results_mapping_phase(phaseAt,memberAtLocal,label)
call h5dcreate_f(loc_id, 'phase', dtype_id, filespace_id, dset_id, hdferr) call h5dcreate_f(loc_id, 'phase', dtype_id, filespace_id, dset_id, hdferr)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
call h5dwrite_f(dset_id, label_id, reshape(label(pack(phaseAtMaterialpoint,.true.)),myShape), & call h5dwrite_f(dset_id, label_id, reshape(label(pack(phaseID,.true.)),myShape), &
myShape, hdferr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) myShape, hdferr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
call h5dwrite_f(dset_id, entry_id, reshape(pack(memberAtGlobal,.true.),myShape), & call h5dwrite_f(dset_id, entry_id, reshape(pack(memberAtGlobal,.true.),myShape), &
@ -572,16 +565,15 @@ end subroutine results_mapping_phase
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief adds the unique mapping from spatial position and constituent ID to results !> @brief adds the unique mapping from spatial position and constituent ID to results
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label) subroutine results_mapping_homogenization(homogenizationID,memberAtLocal,label)
integer, dimension(:), intent(in) :: homogenizationAt !< homogenization section at (element) integer, dimension(:), intent(in) :: homogenizationID !< homogenization ID at (cell)
integer, dimension(:,:), intent(in) :: memberAtLocal !< homogenization member at (IP,element) integer, dimension(:,:), intent(in) :: memberAtLocal !< homogenization entry at (IP,element)
character(len=*), dimension(:), intent(in) :: label !< label of each homogenization section character(len=*), dimension(:), intent(in) :: label !< label of each homogenization section
integer, dimension(size(memberAtLocal,1),size(memberAtLocal,2)) :: & integer, dimension(size(memberAtLocal,1),size(memberAtLocal,2)) :: &
homogenizationAtMaterialpoint, &
memberAtGlobal memberAtGlobal
integer, dimension(size(label),0:worldsize-1) :: memberOffset !< offset in member counting per process integer, dimension(size(label),0:worldsize-1) :: memberOffset !< offset in entry counting per process
integer, dimension(0:worldsize-1) :: writeSize !< amount of data written per process integer, dimension(0:worldsize-1) :: writeSize !< amount of data written per process
integer(HSIZE_T), dimension(1) :: & integer(HSIZE_T), dimension(1) :: &
myShape, & !< shape of the dataset (this process) myShape, & !< shape of the dataset (this process)
@ -609,10 +601,10 @@ subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
memberOffset = 0 memberOffset = 0
do i=1, size(label) do i=1, size(label)
memberOffset(i,worldrank) = count(homogenizationAt == i)*size(memberAtLocal,1) ! number of points/instance of this process memberOffset(i,worldrank) = count(homogenizationID == i) ! number of entries of this process
enddo enddo
writeSize = 0 writeSize = 0
writeSize(worldrank) = size(memberAtLocal) ! total number of points by this process writeSize(worldrank) = size(memberAtLocal) ! total number of entries of this process
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! MPI settings and communication ! MPI settings and communication
@ -631,17 +623,11 @@ subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label)
myOffset = int([sum(writeSize(0:worldrank-1))], HSIZE_T) myOffset = int([sum(writeSize(0:worldrank-1))], HSIZE_T)
totalShape = int([sum(writeSize)], HSIZE_T) totalShape = int([sum(writeSize)], HSIZE_T)
!---------------------------------------------------------------------------------------------------
! expand phaseAt to consider IPs (is not stored per IP)
do i = 1, size(homogenizationAtMaterialpoint,1)
homogenizationAtMaterialpoint(i,:) = homogenizationAt
enddo
!--------------------------------------------------------------------------------------------------- !---------------------------------------------------------------------------------------------------
! renumber member from my process to all processes ! renumber member from my process to all processes
do i = 1, size(label) do i = 1, size(label)
where(homogenizationAtMaterialpoint == i) memberAtGlobal = memberAtLocal + sum(memberOffset(i,0:worldrank-1)) - 1 ! convert to 0-based where(reshape(homogenizationID,shape(memberAtGlobal)) == i) &
memberAtGlobal = memberAtLocal + sum(memberOffset(i,0:worldrank-1)) - 1 ! 0-based
enddo enddo
!--------------------------------------------------------------------------------------------------- !---------------------------------------------------------------------------------------------------
@ -698,7 +684,7 @@ subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label)
call h5dcreate_f(loc_id, 'homogenization', dtype_id, filespace_id, dset_id, hdferr) call h5dcreate_f(loc_id, 'homogenization', dtype_id, filespace_id, dset_id, hdferr)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
call h5dwrite_f(dset_id, label_id, reshape(label(pack(homogenizationAtMaterialpoint,.true.)),myShape), & call h5dwrite_f(dset_id, label_id, reshape(label(pack(homogenizationID,.true.)),myShape), &
myShape, hdferr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) myShape, hdferr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id)
if(hdferr < 0) error stop 'HDF5 error' if(hdferr < 0) error stop 'HDF5 error'
call h5dwrite_f(dset_id, entry_id, reshape(pack(memberAtGlobal,.true.),myShape), & call h5dwrite_f(dset_id, entry_id, reshape(pack(memberAtGlobal,.true.),myShape), &