223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

less global variables

This commit is contained in:
Martin Diehl 2019-06-06 13:08:10 +02:00
parent 9e8bc7d9b1
commit 73d41ffaf7
1 changed files with 213 additions and 210 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

175
src/mesh_marc.f90
View File

@ -319,13 +319,17 @@ subroutine mesh_init(ip,el)
mesh_node = mesh_node0
if (myDebug) write(6,'(a)') ' Built nodes'; flush(6)
elemType = mesh_marc_count_cpSizes(FILEUNIT)
elemType = mesh_marc_getElemType(mesh_nElems,FILEUNIT)
if (myDebug) write(6,'(a)') ' Counted CP sizes'; flush(6)
call theMesh%init(elemType,mesh_node0)
call theMesh%setNelems(mesh_NcpElems)
call mesh_marc_build_elements(initialcondTableStyle,FILEUNIT)
allocate(mesh_element(4+theMesh%elem%nNodes,theMesh%nElems), source=0)
mesh_element(1,:) = -1 ! DEPRECATED
mesh_element(2,:) = elemType ! DEPRECATED
call mesh_marc_buildElements(initialcondTableStyle,FILEUNIT)
if (myDebug) write(6,'(a)') ' Built elements'; flush(6)
close (FILEUNIT)
@ -560,6 +564,7 @@ subroutine mesh_marc_map_elementSets(nameElemSet,mapElemSet,fileUnit)
640 end subroutine mesh_marc_map_elementSets
!--------------------------------------------------------------------------------------------------
!> @brief Maps elements from FE ID to internal (consecutive) representation.
!--------------------------------------------------------------------------------------------------
@ -665,7 +670,7 @@ subroutine mesh_marc_build_nodes(fileUnit)
integer, intent(in) :: fileUnit
integer, dimension(5), parameter :: node_ends = int([0,10,30,50,70],pInt)
integer, dimension(5), parameter :: node_ends = [0,10,30,50,70]
integer, allocatable, dimension(:) :: chunkPos
character(len=300) :: line
integer :: i,j,m
@ -695,18 +700,18 @@ end subroutine mesh_marc_build_nodes
!--------------------------------------------------------------------------------------------------
!> @brief Gets maximum count of nodes, IPs, IP neighbors, and cellnodes among cpElements.
!! Sets global values 'mesh_maxNnodes', 'mesh_maxNips', 'mesh_maxNipNeighbors',
!! and 'mesh_maxNcellnodes'
!> @brief Gets element type (and checks if the whole mesh comprises of only one type)
!--------------------------------------------------------------------------------------------------
integer function mesh_marc_count_cpSizes(fileUnit)
integer function mesh_marc_getElemType(nElem,fileUnit)
integer, intent(in) :: fileUnit
integer, intent(in) :: &
nElem, &
fileUnit
type(tElement) :: tempEl
integer, allocatable, dimension(:) :: chunkPos
character(len=300) :: line
integer :: i,t,g,e,c
integer :: i,t
t = -1
@ -716,15 +721,15 @@ integer function mesh_marc_count_cpSizes(fileUnit)
chunkPos = IO_stringPos(line)
if( IO_lc(IO_stringValue(line,chunkPos,1)) == 'connectivity' ) then
read (fileUnit,'(A300)') line ! Garbage line
do i=1,mesh_Nelems ! read all elements
do i=1,nElem ! read all elements
read (fileUnit,'(A300)') line
chunkPos = IO_stringPos(line) ! limit to id and type
chunkPos = IO_stringPos(line)
if (t == -1) then
t = FE_mapElemtype(IO_stringValue(line,chunkPos,2))
t = mapElemtype(IO_stringValue(line,chunkPos,2))
call tempEl%init(t)
mesh_marc_count_cpSizes = t
mesh_marc_getElemType = t
else
if (t /= FE_mapElemtype(IO_stringValue(line,chunkPos,2))) call IO_error(0) !ToDo: error message
if (t /= mapElemtype(IO_stringValue(line,chunkPos,2))) call IO_error(191,el=t,ip=i)
endif
call IO_skipChunks(fileUnit,tempEl%nNodes-(chunkPos(1)-2))
enddo
@ -732,16 +737,63 @@ integer function mesh_marc_count_cpSizes(fileUnit)
endif
enddo
630 end function mesh_marc_count_cpSizes
contains
!--------------------------------------------------------------------------------------------------
!> @brief mapping of Marc element types to internal representation
!--------------------------------------------------------------------------------------------------
integer function mapElemtype(what)
character(len=*), intent(in) :: what
select case (IO_lc(what))
case ( '6')
mapElemtype = 1 ! Two-dimensional Plane Strain Triangle
case ( '155', &
'125', &
'128')
mapElemtype = 2 ! Two-dimensional Plane Strain triangle (155: cubic shape function, 125/128: second order isoparametric)
case ( '11')
mapElemtype = 3 ! Arbitrary Quadrilateral Plane-strain
case ( '27')
mapElemtype = 4 ! Plane Strain, Eight-node Distorted Quadrilateral
case ( '54')
mapElemtype = 5 ! Plane Strain, Eight-node Distorted Quadrilateral with reduced integration
case ( '134')
mapElemtype = 6 ! Three-dimensional Four-node Tetrahedron
case ( '157')
mapElemtype = 7 ! Three-dimensional, Low-order, Tetrahedron, Herrmann Formulations
case ( '127')
mapElemtype = 8 ! Three-dimensional Ten-node Tetrahedron
case ( '136')
mapElemtype = 9 ! Three-dimensional Arbitrarily Distorted Pentahedral
case ( '117', &
'123')
mapElemtype = 10 ! Three-dimensional Arbitrarily Distorted linear hexahedral with reduced integration
case ( '7')
mapElemtype = 11 ! Three-dimensional Arbitrarily Distorted Brick
case ( '57')
mapElemtype = 12 ! Three-dimensional Arbitrarily Distorted quad hexahedral with reduced integration
case ( '21')
mapElemtype = 13 ! Three-dimensional Arbitrarily Distorted quadratic hexahedral
case default
call IO_error(error_ID=190,ext_msg=IO_lc(what))
end select
end function mapElemtype
630 end function mesh_marc_getElemType
!--------------------------------------------------------------------------------------------------
!> @brief Store FEid, type, mat, tex, and node list per element.
!! Allocates global array 'mesh_element'
!> @brief Stores node IDs and homogenization and microstructure ID
!--------------------------------------------------------------------------------------------------
subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
subroutine mesh_marc_buildElements(initialcondTableStyle,fileUnit)
integer, intent(in) :: initialcondTableStyle,fileUnit
integer, intent(in) :: &
initialcondTableStyle, &
fileUnit
integer, allocatable, dimension(:) :: chunkPos
character(len=300) line
@ -749,10 +801,6 @@ subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
integer, dimension(1+theMesh%nElems) :: contInts
integer :: i,j,t,sv,myVal,e,nNodesAlreadyRead
allocate(mesh_element(4+theMesh%elem%nNodes,theMesh%nElems), source=0)
mesh_elemType = -1
rewind(fileUnit)
do
read (fileUnit,'(A300)',END=620) line
@ -764,12 +812,6 @@ subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
chunkPos = IO_stringPos(line)
e = mesh_FEasCP('elem',IO_intValue(line,chunkPos,1))
if (e /= 0) then ! disregard non CP elems
mesh_element(1,e) = -1 ! DEPRECATED
t = FE_mapElemtype(IO_StringValue(line,chunkPos,2)) ! elem type
if (mesh_elemType /= t .and. mesh_elemType /= -1) &
call IO_error(191,el=t,ip=mesh_elemType)
mesh_elemType = t
mesh_element(2,e) = t
nNodesAlreadyRead = 0
do j = 1,chunkPos(1)-2
mesh_element(4+j,e) = mesh_FEasCP('node',IO_IntValue(line,chunkPos,j+2)) ! CP ids of nodes
@ -779,8 +821,7 @@ subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
read (fileUnit,'(A300)',END=620) line
chunkPos = IO_stringPos(line)
do j = 1,chunkPos(1)
mesh_element(4+nNodesAlreadyRead+j,e) &
= mesh_FEasCP('node',IO_IntValue(line,chunkPos,j)) ! CP ids of nodes
mesh_element(4+nNodesAlreadyRead+j,e) = mesh_FEasCP('node',IO_IntValue(line,chunkPos,j)) ! CP ids of nodes
enddo
nNodesAlreadyRead = nNodesAlreadyRead + chunkPos(1)
enddo
@ -790,7 +831,17 @@ subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
endif
enddo
620 rewind(fileUnit) ! just in case "initial state" appears before "connectivity"
call calcCells(theMesh,mesh_element(5:,:))
#if defined(DAMASK_HDF5)
call results_openJobFile
call HDF5_closeGroup(results_addGroup('geometry'))
call results_writeDataset('geometry',mesh_element(5:,:),'C',&
'connectivity of the elements','-')
call results_closeJobFile
#endif
call calcCells(theMesh,theMesh%elem,mesh_element(5:,:))
read (fileUnit,'(A300)',END=630) line
do
chunkPos = IO_stringPos(line)
@ -825,13 +876,14 @@ subroutine mesh_marc_build_elements(initialcondTableStyle,fileUnit)
endif
enddo
630 end subroutine mesh_marc_build_elements
630 end subroutine mesh_marc_buildElements
subroutine calcCells(thisMesh,connectivity_elem)
subroutine calcCells(thisMesh,elem,connectivity_elem)
class(tMesh) :: thisMesh
integer(pInt),dimension(:,:), intent(inout) :: connectivity_elem
type(tElement) :: elem
integer(pInt),dimension(:,:), intent(in) :: connectivity_elem
integer(pInt),dimension(:,:), allocatable :: con_elem,temp,con,parentsAndWeights,candidates_global
integer(pInt),dimension(:), allocatable :: l, nodes, candidates_local
integer(pInt),dimension(:,:,:), allocatable :: con_cell,connectivity_cell
@ -839,13 +891,6 @@ subroutine calcCells(thisMesh,connectivity_elem)
real(pReal), dimension(:,:), allocatable :: coordinates,nodes5
integer(pInt) :: e, n, c, p, s,u,i,m,j,nParentNodes,nCellNode,ierr
#if defined(DAMASK_HDF5)
call results_openJobFile
call HDF5_closeGroup(results_addGroup('geometry'))
call results_writeDataset('geometry',connectivity_elem,'connectivity_element',&
'connectivity of the elements','-')
#endif
!---------------------------------------------------------------------------------------------------
! initialize global connectivity to negative local connectivity
allocate(connectivity_cell(thisMesh%elem%NcellNodesPerCell,thisMesh%elem%nIPs,thisMesh%Nelems))
@ -971,7 +1016,8 @@ subroutine calcCells(thisMesh,connectivity_elem)
connectivity_cell_reshape = reshape(connectivity_cell,[thisMesh%elem%NcellNodesPerCell,thisMesh%elem%nIPs*thisMesh%Nelems])
#if defined(DAMASK_HDF5)
call results_writeDataset('geometry',connectivity_cell_reshape,'connectivity_cell',&
call results_openJobFile
call results_writeDataset('geometry',connectivity_cell_reshape,'c',&
'connectivity of the cells','-')
call results_closeJobFile
#endif
@ -1714,49 +1760,6 @@ end subroutine mesh_faceMatch
end subroutine mesh_build_ipNeighborhood
!--------------------------------------------------------------------------------------------------
!> @brief mapping of FE element types to internal representation
!--------------------------------------------------------------------------------------------------
integer function FE_mapElemtype(what)
character(len=*), intent(in) :: what
select case (IO_lc(what))
case ( '6')
FE_mapElemtype = 1 ! Two-dimensional Plane Strain Triangle
case ( '155', &
'125', &
'128')
FE_mapElemtype = 2 ! Two-dimensional Plane Strain triangle (155: cubic shape function, 125/128: second order isoparametric)
case ( '11')
FE_mapElemtype = 3 ! Arbitrary Quadrilateral Plane-strain
case ( '27')
FE_mapElemtype = 4 ! Plane Strain, Eight-node Distorted Quadrilateral
case ( '54')
FE_mapElemtype = 5 ! Plane Strain, Eight-node Distorted Quadrilateral with reduced integration
case ( '134')
FE_mapElemtype = 6 ! Three-dimensional Four-node Tetrahedron
case ( '157')
FE_mapElemtype = 7 ! Three-dimensional, Low-order, Tetrahedron, Herrmann Formulations
case ( '127')
FE_mapElemtype = 8 ! Three-dimensional Ten-node Tetrahedron
case ( '136')
FE_mapElemtype = 9 ! Three-dimensional Arbitrarily Distorted Pentahedral
case ( '117', &
'123')
FE_mapElemtype = 10 ! Three-dimensional Arbitrarily Distorted linear hexahedral with reduced integration
case ( '7')
FE_mapElemtype = 11 ! Three-dimensional Arbitrarily Distorted Brick
case ( '57')
FE_mapElemtype = 12 ! Three-dimensional Arbitrarily Distorted quad hexahedral with reduced integration
case ( '21')
FE_mapElemtype = 13 ! Three-dimensional Arbitrarily Distorted quadratic hexahedral
case default
call IO_error(error_ID=190,ext_msg=IO_lc(what))
end select
end function FE_mapElemtype
!--------------------------------------------------------------------------------------------------
!> @brief get properties of different types of finite elements