diff --git a/src/homogenization_mechanical_RGC.f90 b/src/homogenization_mechanical_RGC.f90 index ac7dc63ea..345592790 100644 --- a/src/homogenization_mechanical_RGC.f90 +++ b/src/homogenization_mechanical_RGC.f90 @@ -766,25 +766,21 @@ end function relaxationVector !-------------------------------------------------------------------------------------------------- !> @brief identify the normal of an interface !-------------------------------------------------------------------------------------------------- -pure function interfaceNormal(intFace,ho,en) - - real(pReal), dimension(3) :: interfaceNormal +pure function interfaceNormal(intFace,ho,en) result(n) + real(pReal), dimension(3) :: n integer, dimension(4), intent(in) :: intFace !< interface ID in 4D array (normal and position) integer, intent(in) :: & ho, & en - integer :: nPos + associate (dst => dependentState(ho)) -!-------------------------------------------------------------------------------------------------- -! get the normal of the interface, identified from the value of intFace(1) - interfaceNormal = 0.0_pReal - nPos = abs(intFace(1)) ! identify the position of the interface in global state array - interfaceNormal(nPos) = real(intFace(1)/abs(intFace(1)),pReal) ! get the normal vector w.r.t. cluster axis + n = 0.0_pReal + n(abs(intFace(1))) = real(intFace(1)/abs(intFace(1)),pReal) ! get the normal vector w.r.t. cluster axis - interfaceNormal = matmul(dst%orientation(1:3,1:3,en),interfaceNormal) ! map the normal vector into sample coordinate system (basis) + n = matmul(dst%orientation(1:3,1:3,en),n) ! map the normal vector into sample coordinate system (basis) end associate @@ -794,22 +790,18 @@ end function interfaceNormal !-------------------------------------------------------------------------------------------------- !> @brief collect six faces of a grain in 4D (normal and position) !-------------------------------------------------------------------------------------------------- -pure function getInterface(iFace,iGrain3) - - integer, dimension(4) :: getInterface +pure function getInterface(iFace,iGrain3) result(i) + integer, dimension(4) :: i integer, dimension(3), intent(in) :: iGrain3 !< grain ID in 3D array integer, intent(in) :: iFace !< face index (1..6) mapped like (-e1,-e2,-e3,+e1,+e2,+e3) or iDir = (-1,-2,-3,1,2,3) integer :: iDir !< direction of interface normal - iDir = (int(real(iFace-1,pReal)/2.0_pReal)+1)*(-1)**iFace - getInterface(1) = iDir -!-------------------------------------------------------------------------------------------------- -! identify the interface position by the direction of its normal - getInterface(2:4) = iGrain3 - if (iDir < 0) getInterface(1-iDir) = getInterface(1-iDir)-1 ! to have a correlation with coordinate/position in real space + iDir = (int(real(iFace-1,pReal)/2.0_pReal)+1)*(-1)**iFace + i = [iDir,iGrain3] + if (iDir < 0) i(1-iDir) = i(1-iDir)-1 ! to have a correlation with coordinate/position in real space end function getInterface diff --git a/src/phase.f90 b/src/phase.f90 index effe189d3..7b19254b7 100644 --- a/src/phase.f90 +++ b/src/phase.f90 @@ -67,8 +67,8 @@ module phase interface ! == cleaned:begin ================================================================================= - module subroutine mechanical_init(materials,phases) - class(tNode), pointer :: materials,phases + module subroutine mechanical_init(phases) + class(tNode), pointer :: phases end subroutine mechanical_init module subroutine damage_init @@ -386,7 +386,7 @@ subroutine phase_init phase_O(ph)%data = phase_O_0(ph)%data enddo - call mechanical_init(materials,phases) + call mechanical_init(phases) call damage_init call thermal_init(phases) @@ -482,7 +482,6 @@ end subroutine phase_results subroutine crystallite_init() integer :: & - ph, & ce, & co, & !< counter in integration point component loop ip, & !< counter in integration point loop diff --git a/src/phase_mechanical.f90 b/src/phase_mechanical.f90 index 52ce43c84..6ad96665c 100644 --- a/src/phase_mechanical.f90 +++ b/src/phase_mechanical.f90 @@ -199,10 +199,9 @@ contains !> @brief Initialize mechanical field related constitutive models !> @details Initialize elasticity, plasticity and stiffness degradation models. !-------------------------------------------------------------------------------------------------- -module subroutine mechanical_init(materials,phases) +module subroutine mechanical_init(phases) class(tNode), pointer :: & - materials, & phases integer :: & diff --git a/src/phase_mechanical_plastic_dislotungsten.f90 b/src/phase_mechanical_plastic_dislotungsten.f90 index 76dc4adeb..1e5e120e5 100644 --- a/src/phase_mechanical_plastic_dislotungsten.f90 +++ b/src/phase_mechanical_plastic_dislotungsten.f90 @@ -80,7 +80,7 @@ module function plastic_dislotungsten_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - ph, i, & + ph, & Nmembers, & sizeState, sizeDotState, & startIndex, endIndex diff --git a/src/phase_mechanical_plastic_nonlocal.f90 b/src/phase_mechanical_plastic_nonlocal.f90 index 1d22d35c2..1968d258e 100644 --- a/src/phase_mechanical_plastic_nonlocal.f90 +++ b/src/phase_mechanical_plastic_nonlocal.f90 @@ -1119,8 +1119,11 @@ end subroutine nonlocal_dotState !--------------------------------------------------------------------------------------------------- !> @brief calculates the rate of change of microstructure !--------------------------------------------------------------------------------------------------- +#if __INTEL_COMPILER >= 2020 +non_recursive function rhoDotFlux(timestep,ph,en,ip,el) +#else function rhoDotFlux(timestep,ph,en,ip,el) - +#endif real(pReal), intent(in) :: & timestep !< substepped crystallite time increment integer, intent(in) :: & diff --git a/src/results.f90 b/src/results.f90 index a5dccdaac..064829658 100644 --- a/src/results.f90 +++ b/src/results.f90 @@ -67,7 +67,6 @@ subroutine results_init(restart) character(len=pPathLen) :: commandLine integer :: hdferr - integer(HID_T) :: group_id character(len=:), allocatable :: date