From 88df7f2957f793f5f384360d0d9e5f30d48a7d64 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 17 Oct 2019 07:48:57 +0200 Subject: [PATCH] store cell normal directions --- src/geometry_plastic_nonlocal.f90 | 23 ++++++++++++++++------- src/mesh_marc.f90 | 11 ++++++----- 2 files changed, 22 insertions(+), 12 deletions(-) diff --git a/src/geometry_plastic_nonlocal.f90 b/src/geometry_plastic_nonlocal.f90 index b042943d0..056c7cf57 100644 --- a/src/geometry_plastic_nonlocal.f90 +++ b/src/geometry_plastic_nonlocal.f90 @@ -120,26 +120,35 @@ end subroutine geometry_plastic_nonlocal_disable !--------------------------------------------------------------------------------------------------- subroutine geometry_plastic_nonlocal_results - integer, dimension(:), allocatable :: s + integer, dimension(:), allocatable :: shp #if defined(DAMASK_HDF5) call results_openJobFile writeVolume: block real(pReal), dimension(:), allocatable :: temp - s = shape(geometry_plastic_nonlocal_IPvolume0) - temp = reshape(geometry_plastic_nonlocal_IPvolume0,[s(1)*s(2)]) + shp = shape(geometry_plastic_nonlocal_IPvolume0) + temp = reshape(geometry_plastic_nonlocal_IPvolume0,[shp(1)*shp(2)]) call results_writeDataset('geometry',temp,'v_0',& 'initial cell volume','m³') end block writeVolume - writeArea: block + writeAreas: block real(pReal), dimension(:,:), allocatable :: temp - s = shape(geometry_plastic_nonlocal_IParea0) - temp = reshape(geometry_plastic_nonlocal_IParea0,[s(1),s(2)*s(3)]) + shp = shape(geometry_plastic_nonlocal_IParea0) + temp = reshape(geometry_plastic_nonlocal_IParea0,[shp(1),shp(2)*shp(3)]) call results_writeDataset('geometry',temp,'a_0',& 'initial cell face area','m²') - end block writeArea + end block writeAreas + + writeNormals: block + real(pReal), dimension(:,:,:), allocatable :: temp + shp = shape(geometry_plastic_nonlocal_IPareaNormal0) + temp = reshape(geometry_plastic_nonlocal_IPareaNormal0,[shp(1),shp(2),shp(3)*shp(4)]) + call results_writeDataset('geometry',temp,'n_0',& + 'initial cell face normals','-',transposed=.false.) + end block writeNormals + call results_closeJobFile #endif diff --git a/src/mesh_marc.f90 b/src/mesh_marc.f90 index 4154e75d2..66906207a 100644 --- a/src/mesh_marc.f90 +++ b/src/mesh_marc.f90 @@ -77,8 +77,8 @@ subroutine mesh_init(ip,el) connectivity_cell !< cell connectivity for each element,ip/cell integer, dimension(:,:), allocatable :: & connectivity_elem - - real(pReal), dimension(:,:,:,:),allocatable :: x + real(pReal), dimension(:,:,:,:),allocatable :: & + unscaledNormals write(6,'(/,a)') ' <<<+- mesh init -+>>>' @@ -120,8 +120,9 @@ subroutine mesh_init(ip,el) node0_cell,ip_reshaped) call geometry_plastic_nonlocal_setIPvolume(IPvolume(elem,node0_cell,connectivity_cell)) - x = IPareaNormal(elem,nElems,connectivity_cell,node0_cell) - call geometry_plastic_nonlocal_setIParea(norm2(x,1)) + unscaledNormals = IPareaNormal(elem,nElems,connectivity_cell,node0_cell) + call geometry_plastic_nonlocal_setIParea(norm2(unscaledNormals,1)) + call geometry_plastic_nonlocal_setIPareaNormal(unscaledNormals/spread(norm2(unscaledNormals,1),1,3)) call geometry_plastic_nonlocal_results @@ -1034,7 +1035,7 @@ function IPareaNormal(elem,nElem,connectivity,node) do e = 1,nElem do i = 1,elem%nIPs - do f = 1,size(elem%cellFace,2) + do f = 1,elem%nIPneighbors nodePos = node(1:3,connectivity(elem%cellface(1:m,f),i,e)) select case (elem%cellType)