DAMASK_EICMD/src/material.f90

228 lines
9.4 KiB
Fortran
Raw Normal View History

!--------------------------------------------------------------------------------------------------
!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
2018-06-26 22:39:46 +05:30
!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
2020-09-20 23:16:33 +05:30
!> @brief Defines phase and homogenization
!--------------------------------------------------------------------------------------------------
module material
2019-06-05 13:32:55 +05:30
use prec
use config
use results
2021-07-24 18:33:26 +05:30
use math
2019-06-05 13:32:55 +05:30
use IO
2019-06-07 18:01:42 +05:30
use rotations
2019-06-07 09:48:42 +05:30
use discretization
2021-02-12 20:01:43 +05:30
use YAML_types
2019-09-20 00:10:59 +05:30
implicit none
private
2020-03-16 22:59:15 +05:30
2021-05-22 20:51:07 +05:30
type :: tRotationContainer
type(Rotation), dimension(:), allocatable :: data
end type
2021-07-24 18:33:26 +05:30
type :: tTensorContainer
real(pReal), dimension(:,:,:), allocatable :: data
end type
2021-05-22 20:51:07 +05:30
2021-07-24 18:33:26 +05:30
type(tRotationContainer), dimension(:), allocatable :: material_O_0
2021-05-22 20:51:07 +05:30
2021-02-12 20:01:43 +05:30
integer, dimension(:), allocatable, public, protected :: &
2021-03-01 10:46:16 +05:30
homogenization_Nconstituents !< number of grains in each homogenization
2021-05-22 13:03:58 +05:30
integer, public, protected :: &
homogenization_maxNconstituents !< max number of grains in any homogenization
2020-03-16 22:59:15 +05:30
character(len=:), public, protected, allocatable, dimension(:) :: &
2020-08-15 19:32:10 +05:30
material_name_phase, & !< name of each phase
material_name_homogenization !< name of each homogenization
2019-09-20 00:10:59 +05:30
integer, dimension(:), allocatable, public, protected :: & ! (elem)
2021-05-23 03:40:46 +05:30
material_homogenizationID, & !< per cell TODO: material_ID_homogenization
material_homogenizationEntry !< per cell TODO: material_entry_homogenization
2020-10-24 18:26:03 +05:30
integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem)
2021-05-23 03:40:46 +05:30
material_phaseAt, & !< phase ID of each element TODO: remove
material_phaseID, & !< per (constituent,cell) TODO: material_ID_phase
material_phaseEntry !< per (constituent,cell) TODO: material_entry_phase
integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,IP,elem)
2021-05-23 13:40:25 +05:30
material_phaseMemberAt !TODO: remove
2019-09-20 00:10:59 +05:30
public :: &
2021-07-24 18:33:26 +05:30
tTensorContainer, &
2021-05-22 20:51:07 +05:30
tRotationContainer, &
2021-07-24 18:33:26 +05:30
material_O_0, &
2021-02-12 20:01:43 +05:30
material_init
contains
!--------------------------------------------------------------------------------------------------
!> @brief Parse material configuration file (material.yaml).
!--------------------------------------------------------------------------------------------------
subroutine material_init(restart)
logical, intent(in) :: restart
2021-02-12 20:01:43 +05:30
2020-09-22 16:39:12 +05:30
print'(/,a)', ' <<<+- material init -+>>>'; flush(IO_STDOUT)
2020-03-16 22:59:15 +05:30
2020-08-15 19:32:10 +05:30
2021-05-22 13:03:58 +05:30
call parse
print*, 'parsed material.yaml'
2020-03-16 22:59:15 +05:30
if (.not. restart) then
call results_openJobFile
2021-05-22 17:06:23 +05:30
call results_mapping_phase(material_phaseID,material_phaseEntry,material_name_phase)
call results_mapping_homogenization(material_homogenizationID,material_homogenizationEntry,material_name_homogenization)
call results_closeJobFile
endif
2019-04-04 16:55:29 +05:30
end subroutine material_init
!--------------------------------------------------------------------------------------------------
2021-05-22 13:03:58 +05:30
!> @brief Parse material.yaml to get the global structure
!--------------------------------------------------------------------------------------------------
2021-05-22 13:03:58 +05:30
subroutine parse()
2015-10-14 00:22:01 +05:30
2020-10-01 16:13:05 +05:30
class(tNode), pointer :: materials, & !> list of materials
material, & !> material definition
constituents, & !> list of constituents
constituent, & !> constituent definition
2020-08-15 19:32:10 +05:30
phases, &
2020-10-24 20:00:48 +05:30
homogenizations, &
homogenization
2020-03-16 22:59:15 +05:30
2020-08-15 19:32:10 +05:30
integer, dimension(:), allocatable :: &
2020-09-19 13:30:49 +05:30
counterPhase, &
counterHomogenization
2020-03-16 22:59:15 +05:30
2020-09-19 13:30:49 +05:30
real(pReal) :: &
frac
2020-08-15 19:32:10 +05:30
integer :: &
2021-05-22 20:51:07 +05:30
el, ip, co, ma, &
2021-01-17 19:08:12 +05:30
h, ce
2020-10-24 20:00:48 +05:30
materials => config_material%get('material')
phases => config_material%get('phase')
homogenizations => config_material%get('homogenization')
2020-08-15 19:32:10 +05:30
2020-10-24 21:30:17 +05:30
call sanityCheck(materials, homogenizations)
material_name_phase = getKeys(phases)
material_name_homogenization = getKeys(homogenizations)
2020-10-24 18:18:07 +05:30
allocate(homogenization_Nconstituents(homogenizations%length))
2020-10-24 21:30:17 +05:30
do h=1, homogenizations%length
2020-10-24 20:00:48 +05:30
homogenization => homogenizations%get(h)
homogenization_Nconstituents(h) = homogenization%get_asInt('N_constituents')
2020-10-24 20:00:48 +05:30
enddo
homogenization_maxNconstituents = maxval(homogenization_Nconstituents)
2020-10-24 20:00:48 +05:30
allocate(counterPhase(phases%length),source=0)
allocate(counterHomogenization(homogenizations%length),source=0)
allocate(material_phaseAt(homogenization_maxNconstituents,discretization_Nelems),source=0)
allocate(material_phaseMemberAt(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems),source=0)
2021-01-17 19:08:12 +05:30
2021-04-06 15:08:44 +05:30
allocate(material_homogenizationID(discretization_nIPs*discretization_Nelems),source=0)
allocate(material_homogenizationEntry(discretization_nIPs*discretization_Nelems),source=0)
allocate(material_phaseID(homogenization_maxNconstituents,discretization_nIPs*discretization_Nelems),source=0)
allocate(material_phaseEntry(homogenization_maxNconstituents,discretization_nIPs*discretization_Nelems),source=0)
2021-01-17 19:08:12 +05:30
do el = 1, discretization_Nelems
material => materials%get(discretization_materialAt(el))
constituents => material%get('constituents')
2021-01-17 19:08:12 +05:30
do ip = 1, discretization_nIPs
ce = (el-1)*discretization_nIPs + ip
2021-05-23 13:40:25 +05:30
material_homogenizationID(ce) = homogenizations%getIndex(material%get_asString('homogenization'))
counterHomogenization(material_homogenizationID(ce)) = counterHomogenization(material_homogenizationID(ce)) + 1
material_homogenizationEntry(ce) = counterHomogenization(material_homogenizationID(ce))
2020-08-15 19:32:10 +05:30
enddo
2020-09-19 13:30:49 +05:30
frac = 0.0_pReal
2021-01-17 19:08:12 +05:30
do co = 1, constituents%length
constituent => constituents%get(co)
2021-02-04 18:16:01 +05:30
frac = frac + constituent%get_asFloat('v')
2021-01-17 19:08:12 +05:30
material_phaseAt(co,el) = phases%getIndex(constituent%get_asString('phase'))
do ip = 1, discretization_nIPs
ce = (el-1)*discretization_nIPs + ip
counterPhase(material_phaseAt(co,el)) = counterPhase(material_phaseAt(co,el)) + 1
2021-05-23 03:40:46 +05:30
material_phaseMemberAt(co,ip,el) = counterPhase(material_phaseAt(co,el))
material_phaseEntry(co,ce) = counterPhase(material_phaseAt(co,el))
material_phaseID(co,ce) = material_phaseAt(co,el)
enddo
2020-08-15 19:32:10 +05:30
enddo
if (dNeq(frac,1.0_pReal,1.e-12_pReal)) call IO_error(153,ext_msg='constituent')
2020-03-16 22:59:15 +05:30
enddo
2021-07-24 18:33:26 +05:30
allocate(material_O_0(materials%length))
2021-05-22 20:51:07 +05:30
do ma = 1, materials%length
material => materials%get(ma)
constituents => material%get('constituents')
2021-07-24 18:33:26 +05:30
allocate(material_O_0(ma)%data(constituents%length))
2021-05-22 20:51:07 +05:30
do co = 1, constituents%length
2021-05-22 22:12:18 +05:30
constituent => constituents%get(co)
2021-07-24 18:33:26 +05:30
call material_O_0(ma)%data(co)%fromQuaternion(constituent%get_as1dFloat('O',requiredSize=4))
2021-05-22 20:51:07 +05:30
enddo
enddo
2021-05-22 13:03:58 +05:30
end subroutine parse
2018-10-14 23:46:30 +05:30
2020-10-24 21:30:17 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief Check if material.yaml is consistent and contains sufficient # of materials
!--------------------------------------------------------------------------------------------------
subroutine sanityCheck(materials,homogenizations)
class(tNode), intent(in) :: materials, &
homogenizations
class(tNode), pointer :: material, &
homogenization, &
constituents
integer :: m
if(maxval(discretization_materialAt) > materials%length) &
call IO_error(155,ext_msg='More materials requested than found in material.yaml')
do m = 1, materials%length
material => materials%get(m)
constituents => material%get('constituents')
homogenization => homogenizations%get(material%get_asString('homogenization'))
if(constituents%length /= homogenization%get_asInt('N_constituents')) call IO_error(148)
enddo
end subroutine sanityCheck
!--------------------------------------------------------------------------------------------------
2020-11-19 02:21:37 +05:30
!> @brief Get all keys from a dictionary
2020-10-24 21:30:17 +05:30
!--------------------------------------------------------------------------------------------------
function getKeys(dict)
class(tNode), intent(in) :: dict
character(len=:), dimension(:), allocatable :: getKeys
character(len=pStringLen), dimension(:), allocatable :: temp
2020-10-24 21:30:17 +05:30
integer :: i,l
2020-10-24 21:30:17 +05:30
allocate(temp(dict%length))
l = 0
2020-10-24 21:30:17 +05:30
do i=1, dict%length
temp(i) = dict%getKey(i)
l = max(len_trim(temp(i)),l)
enddo
2020-12-10 05:01:58 +05:30
allocate(character(l)::getKeys(dict%length))
do i=1, dict%length
getKeys(i) = trim(temp(i))
2020-10-24 21:30:17 +05:30
enddo
end function getKeys
end module material