diff --git a/src/crystallite.f90 b/src/crystallite.f90 index df3782756..b7ea191fd 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -25,6 +25,7 @@ module crystallite use future use plastic_nonlocal use geometry_plastic_nonlocal, only: & + nIPneighbors => geometry_plastic_nonlocal_nIPneighbors, & IPneighborhood => geometry_plastic_nonlocal_IPneighborhood #if defined(PETSc) || defined(DAMASK_HDF5) use HDF5_utilities @@ -345,7 +346,7 @@ subroutine crystallite_init case(elasmatrix_ID) mySize = 36 case(neighboringip_ID,neighboringelement_ID) - mySize = theMesh%elem%nIPneighbors + mySize = nIPneighbors case default mySize = 0 end select @@ -960,12 +961,12 @@ function crystallite_postResults(ipc, ip, el) mySize = 36 crystallite_postResults(c+1:c+mySize) = reshape(constitutive_homogenizedC(ipc,ip,el),[mySize]) case(neighboringelement_ID) - mySize = theMesh%elem%nIPneighbors + mySize = nIPneighbors crystallite_postResults(c+1:c+mySize) = 0.0_pReal forall (n = 1:mySize) & crystallite_postResults(c+n) = real(IPneighborhood(1,n,ip,el),pReal) case(neighboringip_ID) - mySize = theMesh%elem%nIPneighbors + mySize = nIPneighbors crystallite_postResults(c+1:c+mySize) = 0.0_pReal forall (n = 1:mySize) & crystallite_postResults(c+n) = real(IPneighborhood(2,n,ip,el),pReal) diff --git a/src/geometry_plastic_nonlocal.f90 b/src/geometry_plastic_nonlocal.f90 index 370deb014..37909a4a5 100644 --- a/src/geometry_plastic_nonlocal.f90 +++ b/src/geometry_plastic_nonlocal.f90 @@ -28,10 +28,11 @@ module geometry_plastic_nonlocal public :: & - geometry_plastic_nonlocal_set_IPneighborhood, & - geometry_plastic_nonlocal_set_IPvolume, & - geometry_plastic_nonlocal_set_IParea, & - geometry_plastic_nonlocal_set_IPareaNormal + geometry_plastic_nonlocal_setIPneighborhood, & + geometry_plastic_nonlocal_setIPvolume, & + geometry_plastic_nonlocal_setIParea, & + geometry_plastic_nonlocal_setIPareaNormal, & + geometry_plastic_nonlocal_disable contains @@ -43,7 +44,7 @@ contains ! A triangle (2D) has 3 faces, a quadrilateral (2D) had 4 faces, a tetrahedron (3D) has ! 4 faces, and a hexahedron (3D) has 6 faces. !--------------------------------------------------------------------------------------------------- -subroutine geometry_plastic_nonlocal_set_IPneighborhood(IPneighborhood) +subroutine geometry_plastic_nonlocal_setIPneighborhood(IPneighborhood) integer, dimension(:,:,:,:), intent(in) :: IPneighborhood @@ -51,45 +52,64 @@ subroutine geometry_plastic_nonlocal_set_IPneighborhood(IPneighborhood) geometry_plastic_nonlocal_nIPneighbors = size(IPneighborhood,2) -end subroutine geometry_plastic_nonlocal_set_IPneighborhood +end subroutine geometry_plastic_nonlocal_setIPneighborhood !--------------------------------------------------------------------------------------------------- !> @brief Set the initial volume associated with an integration point !--------------------------------------------------------------------------------------------------- -subroutine geometry_plastic_nonlocal_set_IPvolume(IPvolume) +subroutine geometry_plastic_nonlocal_setIPvolume(IPvolume) real(pReal), dimension(:,:), intent(in) :: IPvolume geometry_plastic_nonlocal_IPvolume0 = IPvolume -end subroutine geometry_plastic_nonlocal_set_IPvolume +end subroutine geometry_plastic_nonlocal_setIPvolume !--------------------------------------------------------------------------------------------------- !> @brief Set the initial areas of the unit triangle/unit quadrilateral/tetrahedron/hexahedron ! encompassing an integration point !--------------------------------------------------------------------------------------------------- -subroutine geometry_plastic_nonlocal_set_IParea(IParea) +subroutine geometry_plastic_nonlocal_setIParea(IParea) real(pReal), dimension(:,:,:), intent(in) :: IParea geometry_plastic_nonlocal_IParea0 = IParea -end subroutine geometry_plastic_nonlocal_set_IParea +end subroutine geometry_plastic_nonlocal_setIParea !--------------------------------------------------------------------------------------------------- !> @brief Set the direction normal of the areas of the triangle/quadrilateral/tetrahedron/hexahedron ! encompassing an integration point !--------------------------------------------------------------------------------------------------- -subroutine geometry_plastic_nonlocal_set_IPareaNormal(IPareaNormal) +subroutine geometry_plastic_nonlocal_setIPareaNormal(IPareaNormal) real(pReal), dimension(:,:,:,:), intent(in) :: IPareaNormal geometry_plastic_nonlocal_IPareaNormal0 = IPareaNormal -end subroutine geometry_plastic_nonlocal_set_IPareaNormal +end subroutine geometry_plastic_nonlocal_setIPareaNormal +!--------------------------------------------------------------------------------------------------- +!> @brief Frees memory used by variables only needed by plastic_nonlocal +!--------------------------------------------------------------------------------------------------- +subroutine geometry_plastic_nonlocal_disable + + if(allocated(geometry_plastic_nonlocal_IPneighborhood)) & + deallocate(geometry_plastic_nonlocal_IPneighborhood) + + if(allocated(geometry_plastic_nonlocal_IPvolume0)) & + deallocate(geometry_plastic_nonlocal_IPvolume0) + + if(allocated(geometry_plastic_nonlocal_IParea0)) & + deallocate(geometry_plastic_nonlocal_IParea0) + + if(allocated(geometry_plastic_nonlocal_IPareaNormal0)) & + deallocate(geometry_plastic_nonlocal_IPareaNormal0) + +end subroutine geometry_plastic_nonlocal_disable + end module geometry_plastic_nonlocal diff --git a/src/mesh_abaqus.f90 b/src/mesh_abaqus.f90 index 6e11a522c..2b6a8ee14 100644 --- a/src/mesh_abaqus.f90 +++ b/src/mesh_abaqus.f90 @@ -524,8 +524,10 @@ subroutine mesh_init(ip,el) call discretization_init(mesh_element(3,:),mesh_element(4,:),& reshape(mesh_ipCoordinates,[3,theMesh%elem%nIPs*theMesh%nElems]),& mesh_node0) - call geometry_plastic_nonlocal_set_IPvolume(mesh_ipVolume) - call geometry_plastic_nonlocal_set_IPneighborhood(mesh_ipNeighborhood) + call geometry_plastic_nonlocal_setIPvolume(mesh_ipVolume) + call geometry_plastic_nonlocal_setIPneighborhood(mesh_ipNeighborhood) + call geometry_plastic_nonlocal_setIParea(mesh_IParea) + call geometry_plastic_nonlocal_setIPareaNormal(mesh_IPareaNormal) contains diff --git a/src/mesh_grid.f90 b/src/mesh_grid.f90 index d41517daa..e9a9b3d33 100644 --- a/src/mesh_grid.f90 +++ b/src/mesh_grid.f90 @@ -77,7 +77,7 @@ module mesh mesh_init - type(tMesh), public, protected :: theMesh + type(tMesh) :: theMesh contains @@ -143,16 +143,18 @@ subroutine mesh_init(ip,el) if (myDebug) write(6,'(a)') ' Built IP coordinates'; flush(6) allocate(IPvolume(1,theMesh%nElems),source=product([geomSize(1:2),size3]/real([grid(1:2),grid3]))) - call geometry_plastic_nonlocal_set_IPvolume(IPvolume) + call geometry_plastic_nonlocal_setIPvolume(IPvolume) if (myDebug) write(6,'(a)') ' Built IP volumes'; flush(6) mesh_ipArea = mesh_build_ipAreas([geomSize(1:2),size3],[grid(1:2),grid3]) + call geometry_plastic_nonlocal_setIParea(mesh_IParea) mesh_ipAreaNormal = mesh_build_ipNormals(grid(1)*grid(2)*grid3) + call geometry_plastic_nonlocal_setIPareaNormal(mesh_ipAreaNormal) if (myDebug) write(6,'(a)') ' Built IP areas'; flush(6) - call geometry_plastic_nonlocal_set_IPneighborhood(mesh_spectral_build_ipNeighborhood([grid(1:2),grid3])) + call geometry_plastic_nonlocal_setIPneighborhood(mesh_spectral_build_ipNeighborhood([grid(1:2),grid3])) if (myDebug) write(6,'(a)') ' Built IP neighborhood'; flush(6) if (debug_e < 1 .or. debug_e > theMesh%nElems) & @@ -164,13 +166,6 @@ subroutine mesh_init(ip,el) allocate(FEsolving_execIP(2_pInt,theMesh%nElems), source=1_pInt) ! parallel loop bounds set to comprise from first IP... forall (j = 1_pInt:theMesh%nElems) FEsolving_execIP(2,j) = theMesh%elem%nIPs ! ...up to own IP count for each element - - - -!!!! COMPATIBILITY HACK !!!! - theMesh%homogenizationAt = mesh_element(3,:) - theMesh%microstructureAt = mesh_element(4,:) -!!!!!!!!!!!!!!!!!!!!!!!! call discretization_init(mesh_element(3,:),mesh_element(4,:),& reshape(mesh_ipCoordinates,[3,grid(1)*grid(2)*grid3]),& mesh_node0) diff --git a/src/mesh_marc.f90 b/src/mesh_marc.f90 index ba4e4e070..44561fa7e 100644 --- a/src/mesh_marc.f90 +++ b/src/mesh_marc.f90 @@ -375,11 +375,13 @@ subroutine mesh_init(ip,el) theMesh%homogenizationAt = mesh_element(3,:) theMesh%microstructureAt = mesh_element(4,:) - call discretization_init(mesh_element(3,:),mesh_element(4,:),& - reshape(mesh_ipCoordinates,[3,theMesh%elem%nIPs*theMesh%nElems]),& - mesh_node0) - call geometry_plastic_nonlocal_set_IPvolume(mesh_ipVolume) - call geometry_plastic_nonlocal_set_IPneighborhood(mesh_ipNeighborhood) + call discretization_init(mesh_element(3,:),mesh_element(4,:),& + reshape(mesh_ipCoordinates,[3,theMesh%elem%nIPs*theMesh%nElems]),& + mesh_node0) + call geometry_plastic_nonlocal_setIPvolume(mesh_ipVolume) + call geometry_plastic_nonlocal_setIPneighborhood(mesh_ipNeighborhood) + call geometry_plastic_nonlocal_setIParea(mesh_IParea) + call geometry_plastic_nonlocal_setIPareaNormal(mesh_IPareaNormal) end subroutine mesh_init @@ -1261,7 +1263,7 @@ subroutine IP_neighborhood2 endif f = f +1 enddo - call geometry_plastic_nonlocal_set_IPneighborhood(mesh_ipNeighborhood2) + call geometry_plastic_nonlocal_setIPneighborhood(mesh_ipNeighborhood2) do e = 1,theMesh%nElems do i = 1,theMesh%elem%nIPs diff --git a/src/plastic_nonlocal.f90 b/src/plastic_nonlocal.f90 index 7bb82e853..9083635fe 100644 --- a/src/plastic_nonlocal.f90 +++ b/src/plastic_nonlocal.f90 @@ -18,6 +18,8 @@ module plastic_nonlocal use lattice use discretization use geometry_plastic_nonlocal, only: & + geometry_plastic_nonlocal_disable, & + nIPneighbors => geometry_plastic_nonlocal_nIPneighbors, & IPneighborhood => geometry_plastic_nonlocal_IPneighborhood, & IPvolume => geometry_plastic_nonlocal_IPvolume0, & IParea => geometry_plastic_nonlocal_IParea0, & @@ -282,9 +284,11 @@ subroutine plastic_nonlocal_init write(6,'(a)') ' http://publications.rwth-aachen.de/record/229993' maxNinstances = count(phase_plasticity == PLASTICITY_NONLOCAL_ID) - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) then write(6,'(a16,1x,i5,/)') '# instances:',maxNinstances - + else + call geometry_plastic_nonlocal_disable + endif allocate(param(maxNinstances)) allocate(state(maxNinstances)) @@ -666,7 +670,7 @@ subroutine plastic_nonlocal_init enddo - allocate(compatibility(2,maxval(totalNslip),maxval(totalNslip),theMesh%elem%nIPneighbors,& + allocate(compatibility(2,maxval(totalNslip),maxval(totalNslip),nIPneighbors,& discretization_nIP,discretization_nElem), source=0.0_pReal) ! BEGIN DEPRECATED---------------------------------------------------------------------------------- @@ -852,7 +856,7 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el) invFp, & !< inverse of plastic deformation gradient connections, & invConnections - real(pReal), dimension(3,theMesh%elem%nIPneighbors) :: & + real(pReal), dimension(3,nIPneighbors) :: & connection_latticeConf real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & rhoExcess @@ -866,10 +870,10 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el) totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & myInteractionMatrix ! corrected slip interaction matrix - real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),theMesh%elem%nIPneighbors) :: & + real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),nIPneighbors) :: & rho_edg_delta_neighbor, & rho_scr_delta_neighbor - real(pReal), dimension(2,maxval(totalNslip),theMesh%elem%nIPneighbors) :: & + real(pReal), dimension(2,maxval(totalNslip),nIPneighbors) :: & neighbor_rhoExcess, & ! excess density at neighboring material point neighbor_rhoTotal ! total density at neighboring material point real(pReal), dimension(3,totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),2) :: & @@ -932,7 +936,7 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el) nRealNeighbors = 0.0_pReal neighbor_rhoTotal = 0.0_pReal - do n = 1,theMesh%elem%nIPneighbors + do n = 1,nIPneighbors neighbor_el = IPneighborhood(1,n,ip,el) neighbor_ip = IPneighborhood(2,n,ip,el) no = phasememberAt(1,neighbor_ip,neighbor_el) @@ -1625,7 +1629,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, & my_Fe = Fe(1:3,1:3,1,ip,el) my_F = matmul(my_Fe, Fp(1:3,1:3,1,ip,el)) - neighbors: do n = 1,theMesh%elem%nIPneighbors + neighbors: do n = 1,nIPneighbors neighbor_el = IPneighborhood(1,n,ip,el) neighbor_ip = IPneighborhood(2,n,ip,el) @@ -1907,7 +1911,7 @@ subroutine plastic_nonlocal_updateCompatibility(orientation,i,e) absoluteMisorientation ! absolute misorientation (without symmetry) between me and my neighbor real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phase(1,i,e))),& totalNslip(phase_plasticityInstance(material_phase(1,i,e))),& - theMesh%elem%nIPneighbors) :: & + nIPneighbors) :: & my_compatibility ! my_compatibility for current element and ip real(pReal) :: & my_compatibilitySum, & @@ -1917,7 +1921,7 @@ subroutine plastic_nonlocal_updateCompatibility(orientation,i,e) belowThreshold type(rotation) :: rot - Nneighbors = theMesh%elem%nIPneighbors + Nneighbors = nIPneighbors ph = material_phase(1,i,e) textureID = material_texture(1,i,e) instance = phase_plasticityInstance(ph)