From 519cd840bde2610117ca80cf6298817028db2c41 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 15:19:23 +0200 Subject: [PATCH 01/24] cleaning --- src/meshFEM.f90 | 12 +++++------- 1 file changed, 5 insertions(+), 7 deletions(-) diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 7d79dd46d..8553fb2a6 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -17,6 +17,8 @@ use PETScis implicit none private + integer(pInt), public, parameter :: & + mesh_ElemType=1_pInt !< Element type of the mesh (only support homogeneous meshes) integer(pInt), public, protected :: & mesh_Nboundaries, & @@ -25,8 +27,7 @@ use PETScis mesh_Nnodes, & !< total number of nodes in mesh mesh_maxNnodes, & !< max number of nodes in any CP element mesh_maxNips, & !< max number of IPs in any CP element - mesh_maxNipNeighbors, & - mesh_Nelems !< total number of elements in mesh + mesh_maxNipNeighbors real(pReal), public, protected :: charLength @@ -212,11 +213,10 @@ subroutine mesh_init(ip,el) endif call DMDestroy(globalMesh,ierr); CHKERRQ(ierr) - call DMGetStratumSize(geomMesh,'depth',dimPlex,mesh_Nelems,ierr) + call DMGetStratumSize(geomMesh,'depth',dimPlex,mesh_NcpElems,ierr) CHKERRQ(ierr) call DMGetStratumSize(geomMesh,'depth',0,mesh_Nnodes,ierr) CHKERRQ(ierr) - mesh_NcpElems = mesh_Nelems FE_Nips(FE_geomtype(1_pInt)) = FEM_Zoo_nQuadrature(dimPlex,integrationOrder) mesh_maxNnodes = FE_Nnodes(1_pInt) @@ -227,14 +227,12 @@ subroutine mesh_init(ip,el) allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems)); mesh_element = 0_pInt do j = 1, mesh_NcpElems mesh_element( 1,j) = j - mesh_element( 2,j) = 1_pInt ! elem type + mesh_element( 2,j) = mesh_elemType ! elem type mesh_element( 3,j) = 1_pInt ! homogenization call DMGetLabelValue(geomMesh,'material',j-1,mesh_element(4,j),ierr) CHKERRQ(ierr) end do - if (usePingPong .and. (mesh_Nelems /= mesh_NcpElems)) & - call IO_error(600_pInt) ! ping-pong must be disabled when having non-DAMASK elements if (debug_e < 1 .or. debug_e > mesh_NcpElems) & call IO_error(602_pInt,ext_msg='element') ! selected element does not exist if (debug_i < 1 .or. debug_i > FE_Nips(FE_geomtype(mesh_element(2_pInt,debug_e)))) & From db45b7615a4f762ba0c67ee4f69f3c07e1afeacf Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 15:19:47 +0200 Subject: [PATCH 02/24] drop support for heterogeneous meshes heterogeneous meshes are neither advisable nor typically used --- src/IO.f90 | 2 ++ src/mesh.f90 | 33 +++++++++++++++++++++------------ 2 files changed, 23 insertions(+), 12 deletions(-) diff --git a/src/IO.f90 b/src/IO.f90 index c97dcfa9c..b12245ec1 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -1487,6 +1487,8 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) msg = 'no microstructure specified via State Variable 3' case (190_pInt) msg = 'unknown element type:' + case (191_pInt) + msg = 'mesh consists of more than one element type' !-------------------------------------------------------------------------------------------------- ! plasticity error messages diff --git a/src/mesh.f90 b/src/mesh.f90 index 4e72ba73e..138a9ecd3 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -3,7 +3,6 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @author Christoph Koords, Max-Planck-Institut für Eisenforschung GmbH !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH -!> @author Krishna Komerla, Max-Planck-Institut für Eisenforschung GmbH !> @brief Sets up the mesh for the solvers MSC.Marc, Abaqus and the spectral solver !-------------------------------------------------------------------------------------------------- module mesh @@ -15,6 +14,7 @@ module mesh integer(pInt), public, protected :: & mesh_NcpElems, & !< total number of CP elements in local mesh mesh_NelemSets, & + mesh_ElemType, & !< Element type of the mesh (only support homogeneous meshes) mesh_maxNelemInSet, & mesh_Nmaterials, & mesh_Nnodes, & !< total number of nodes in mesh @@ -2023,7 +2023,8 @@ subroutine mesh_marc_build_elements(fileUnit) IO_skipChunks, & IO_stringPos, & IO_intValue, & - IO_continuousIntValues + IO_continuousIntValues, & + IO_error implicit none integer(pInt), intent(in) :: fileUnit @@ -2034,7 +2035,8 @@ subroutine mesh_marc_build_elements(fileUnit) integer(pInt), dimension(1_pInt+mesh_NcpElems) :: contInts integer(pInt) :: i,j,t,sv,myVal,e,nNodesAlreadyRead - allocate (mesh_element(4_pInt+mesh_maxNnodes,mesh_NcpElems)) ; mesh_element = 0_pInt + allocate (mesh_element(4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) + mesh_elemType = -1_pInt 610 FORMAT(A300) @@ -2049,8 +2051,11 @@ subroutine mesh_marc_build_elements(fileUnit) chunkPos = IO_stringPos(line) e = mesh_FEasCP('elem',IO_intValue(line,chunkPos,1_pInt)) if (e /= 0_pInt) then ! disregard non CP elems - mesh_element(1,e) = IO_IntValue (line,chunkPos,1_pInt) ! FE id - t = FE_mapElemtype(IO_StringValue(line,chunkPos,2_pInt)) ! elem type + mesh_element(1,e) = IO_IntValue (line,chunkPos,1_pInt) ! FE id + t = FE_mapElemtype(IO_StringValue(line,chunkPos,2_pInt)) ! elem type + if (mesh_elemType /= t .and. mesh_elemType /= -1_pInt) & + call IO_error(191,el=t,ip=mesh_elemType) + mesh_elemType = t mesh_element(2,e) = t nNodesAlreadyRead = 0_pInt do j = 1_pInt,chunkPos(1)-2_pInt @@ -2688,8 +2693,8 @@ subroutine mesh_abaqus_build_elements(fileUnit) IO_intValue, & IO_extractValue, & IO_floatValue, & - IO_error, & - IO_countDataLines + IO_countDataLines, & + IO_error implicit none integer(pInt), intent(in) :: fileUnit @@ -2701,7 +2706,8 @@ subroutine mesh_abaqus_build_elements(fileUnit) character (len=64) :: materialName,elemSetName character(len=300) :: line - allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems)) ; mesh_element = 0_pInt + allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) + mesh_elemType = -1_pInt 610 FORMAT(A300) @@ -2720,17 +2726,20 @@ subroutine mesh_abaqus_build_elements(fileUnit) IO_lc(IO_stringValue(line,chunkPos,2_pInt)) /= 'matrix' .and. & IO_lc(IO_stringValue(line,chunkPos,2_pInt)) /= 'response' ) & ) then - t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,chunkPos,2_pInt)),'type')) ! remember elem type + t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,chunkPos,2_pInt)),'type')) ! remember elem type c = IO_countDataLines(fileUnit) do i = 1_pInt,c backspace(fileUnit) enddo do i = 1_pInt,c read (fileUnit,610,END=620) line - chunkPos = IO_stringPos(line) ! limit to 64 nodes max + chunkPos = IO_stringPos(line) ! limit to 64 nodes max e = mesh_FEasCP('elem',IO_intValue(line,chunkPos,1_pInt)) - if (e /= 0_pInt) then ! disregard non CP elems - mesh_element(1,e) = IO_intValue(line,chunkPos,1_pInt) ! FE id + if (e /= 0_pInt) then ! disregard non CP elems + mesh_element(1,e) = IO_intValue(line,chunkPos,1_pInt) ! FE id + if (mesh_elemType /= t .and. mesh_elemType /= -1_pInt) & + call IO_error(191,el=t,ip=mesh_elemType) + mesh_elemType = t mesh_element(2,e) = t ! elem type nNodesAlreadyRead = 0_pInt do j = 1_pInt,chunkPos(1)-1_pInt From 4862aca34016f59376ab1d87037013f227758036 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 15:27:51 +0200 Subject: [PATCH 03/24] grouping solver specific variables for better readability --- src/mesh.f90 | 106 ++++++++++++++++++--------------------------------- 1 file changed, 38 insertions(+), 68 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 138a9ecd3..727937055 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -27,20 +27,6 @@ module mesh mesh_maxNcellnodes, & !< max number of cell nodes in any CP element mesh_Nelems !< total number of elements in mesh -#ifdef Spectral - integer(pInt), dimension(3), public, protected :: & - grid !< (global) grid - integer(pInt), public, protected :: & - mesh_NcpElemsGlobal, & !< total number of CP elements in global mesh - grid3, & !< (local) grid in 3rd direction - grid3Offset !< (local) grid offset in 3rd direction - real(pReal), dimension(3), public, protected :: & - geomSize - real(pReal), public, protected :: & - size3, & !< (local) size in 3rd direction - size3offset !< (local) size offset in 3rd direction -#endif - integer(pInt), dimension(:,:), allocatable, public, protected :: & mesh_element, & !< FEid, type(internal representation), material, texture, node indices as CP IDs mesh_sharedElem, & !< entryCount and list of elements containing node @@ -71,36 +57,12 @@ module mesh logical, dimension(3), public, protected :: mesh_periodicSurface !< flag indicating periodic outer surfaces (used for fluxes) -#ifdef Marc4DAMASK - integer(pInt), private :: & - MarcVersion, & !< Version of input file format (Marc only) - hypoelasticTableStyle, & !< Table style (Marc only) - initialcondTableStyle !< Table style (Marc only) - integer(pInt), dimension(:), allocatable, private :: & - Marc_matNumber !< array of material numbers for hypoelastic material (Marc only) -#endif - integer(pInt), dimension(2), private :: & mesh_maxValStateVar = 0_pInt -#ifndef Spectral - character(len=64), dimension(:), allocatable, private :: & - mesh_nameElemSet, & !< names of elementSet - mesh_nameMaterial, & !< names of material in solid section - mesh_mapMaterial !< name of elementSet for material - - integer(pInt), dimension(:,:), allocatable, private :: & - mesh_mapElemSet !< list of elements in elementSet -#endif - integer(pInt), dimension(:,:), allocatable, private :: & +integer(pInt), dimension(:,:), allocatable, private :: & mesh_cellnodeParent !< cellnode's parent element ID, cellnode's intra-element ID -#if defined(Marc4DAMASK) || defined(Abaqus) - integer(pInt), dimension(:,:), allocatable, target, private :: & - mesh_mapFEtoCPelem, & !< [sorted FEid, corresponding CPid] - mesh_mapFEtoCPnode !< [sorted FEid, corresponding CPid] -#endif - integer(pInt),dimension(:,:,:), allocatable, private :: & mesh_cell !< cell connectivity for each element,ip/cell @@ -116,10 +78,6 @@ module mesh integer(pInt), dimension(:,:,:,:), allocatable, private :: & FE_subNodeOnIPFace -#ifdef Abaqus - logical, private :: noPart !< for cases where the ABAQUS input file does not use part/assembly information -#endif - ! These definitions should actually reside in the FE-solver specific part (different for MARC/ABAQUS) ! Hence, I suggest to prefix with "FE_" @@ -375,32 +333,44 @@ module mesh 4 & ! element 21 (3D 20node 27ip) ],pInt) +#ifdef Spectral + integer(pInt), dimension(3), public, protected :: & + grid !< (global) grid + integer(pInt), public, protected :: & + mesh_NcpElemsGlobal, & !< total number of CP elements in global mesh + grid3, & !< (local) grid in 3rd direction + grid3Offset !< (local) grid offset in 3rd direction + real(pReal), dimension(3), public, protected :: & + geomSize + real(pReal), public, protected :: & + size3, & !< (local) size in 3rd direction + size3offset !< (local) size offset in 3rd direction +#endif -! integer(pInt), dimension(FE_Nelemtypes), parameter, private :: MESH_VTKELEMTYPE = & -! int([ & -! 5, & ! element 6 (2D 3node 1ip) -! 22, & ! element 125 (2D 6node 3ip) -! 9, & ! element 11 (2D 4node 4ip) -! 23, & ! element 27 (2D 8node 9ip) -! 23, & ! element 54 (2D 8node 4ip) -! 10, & ! element 134 (3D 4node 1ip) -! 10, & ! element 157 (3D 5node 4ip) -! 24, & ! element 127 (3D 10node 4ip) -! 13, & ! element 136 (3D 6node 6ip) -! 12, & ! element 117 (3D 8node 1ip) -! 12, & ! element 7 (3D 8node 8ip) -! 25, & ! element 57 (3D 20node 8ip) -! 25 & ! element 21 (3D 20node 27ip) -! ],pInt) -! -! integer(pInt), dimension(FE_Ncelltypes), parameter, private :: MESH_VTKCELLTYPE = & -! int([ & -! 5, & ! (2D 3node) -! 9, & ! (2D 4node) -! 10, & ! (3D 4node) -! 12 & ! (3D 8node) -! ],pInt) +#if defined(Marc4DAMASK) || defined(Abaqus) + character(len=64), dimension(:), allocatable, private :: & + mesh_nameElemSet, & !< names of elementSet + mesh_nameMaterial, & !< names of material in solid section + mesh_mapMaterial !< name of elementSet for material + integer(pInt), dimension(:,:), allocatable, private :: & + mesh_mapElemSet !< list of elements in elementSet + integer(pInt), dimension(:,:), allocatable, target, private :: & + mesh_mapFEtoCPelem, & !< [sorted FEid, corresponding CPid] + mesh_mapFEtoCPnode !< [sorted FEid, corresponding CPid] +#endif +#ifdef Marc4DAMASK + integer(pInt), private :: & + MarcVersion, & !< Version of input file format (Marc only) + hypoelasticTableStyle, & !< Table style (Marc only) + initialcondTableStyle !< Table style (Marc only) + integer(pInt), dimension(:), allocatable, private :: & + Marc_matNumber !< array of material numbers for hypoelastic material (Marc only) +#endif + +#ifdef Abaqus + logical, private :: noPart !< for cases where the ABAQUS input file does not use part/assembly information +#endif public :: & mesh_init, & @@ -454,7 +424,7 @@ module mesh mesh_abaqus_count_cpSizes, & mesh_abaqus_build_elements, & #endif -#ifndef Spectral +#if defined(Marc4DAMASK) || defined(Abaqus) mesh_build_nodeTwins, & mesh_build_sharedElems, & mesh_build_ipNeighborhood, & From c1b8854132ae3f3f75db066110ef789a656e0804 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 15:31:19 +0200 Subject: [PATCH 04/24] only needed for commercial solvers --- src/mesh.f90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 727937055..c0c3e7c6f 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -13,7 +13,6 @@ module mesh private integer(pInt), public, protected :: & mesh_NcpElems, & !< total number of CP elements in local mesh - mesh_NelemSets, & mesh_ElemType, & !< Element type of the mesh (only support homogeneous meshes) mesh_maxNelemInSet, & mesh_Nmaterials, & @@ -348,6 +347,8 @@ integer(pInt), dimension(:,:), allocatable, private :: & #endif #if defined(Marc4DAMASK) || defined(Abaqus) + integer(pInt), private :: & + mesh_NelemSets character(len=64), dimension(:), allocatable, private :: & mesh_nameElemSet, & !< names of elementSet mesh_nameMaterial, & !< names of material in solid section From bd600185139e7b224963843415249929a8adc77b Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 15:57:21 +0200 Subject: [PATCH 05/24] not needed for spectral --- src/mesh.f90 | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index c0c3e7c6f..b353f759e 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -19,7 +19,6 @@ module mesh mesh_Nnodes, & !< total number of nodes in mesh mesh_Ncellnodes, & !< total number of cell nodes in mesh (including duplicates) mesh_Ncells, & !< total number of cells in mesh - mesh_maxNnodes, & !< max number of nodes in any CP element mesh_maxNips, & !< max number of IPs in any CP element mesh_maxNipNeighbors, & !< max number of IP neighbors in any CP element mesh_maxNsharedElems, & !< max number of CP elements sharing a node @@ -348,6 +347,7 @@ integer(pInt), dimension(:,:), allocatable, private :: & #if defined(Marc4DAMASK) || defined(Abaqus) integer(pInt), private :: & + mesh_maxNnodes, & !< max number of nodes in any CP element mesh_NelemSets character(len=64), dimension(:), allocatable, private :: & mesh_nameElemSet, & !< names of elementSet @@ -1166,7 +1166,7 @@ end subroutine mesh_spectral_count !-------------------------------------------------------------------------------------------------- !> @brief Gets maximum count of nodes, IPs, IP neighbors, and subNodes among cpElements. -!! Sets global values 'mesh_maxNnodes', 'mesh_maxNips', 'mesh_maxNipNeighbors', +!! Sets global values 'mesh_maxNips', 'mesh_maxNipNeighbors', !! and 'mesh_maxNcellnodes' !-------------------------------------------------------------------------------------------------- subroutine mesh_spectral_count_cpSizes @@ -1178,7 +1178,6 @@ subroutine mesh_spectral_count_cpSizes g = FE_geomtype(t) c = FE_celltype(g) - mesh_maxNnodes = FE_Nnodes(t) mesh_maxNips = FE_Nips(g) mesh_maxNipNeighbors = FE_NipNeighbors(c) mesh_maxNcellnodes = FE_Ncellnodes(g) @@ -1282,7 +1281,7 @@ subroutine mesh_spectral_build_elements(fileUnit) i = IO_countContinuousIntValues(fileUnit) maxIntCount = max(maxIntCount, i) enddo - allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems), source = 0_pInt) + allocate (mesh_element (4_pInt+8_pInt,mesh_NcpElems), source = 0_pInt) allocate (microstructures (1_pInt+maxIntCount), source = 1_pInt) allocate (mesh_microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) From e3f2ad34b23c7225c7e5977493c19a0cfb6e3fca Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 16:32:13 +0200 Subject: [PATCH 06/24] not needed --- src/meshFEM.f90 | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 8553fb2a6..22568c99d 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -25,7 +25,6 @@ use PETScis mesh_NcpElems, & !< total number of CP elements in mesh mesh_NcpElemsGlobal, & mesh_Nnodes, & !< total number of nodes in mesh - mesh_maxNnodes, & !< max number of nodes in any CP element mesh_maxNips, & !< max number of IPs in any CP element mesh_maxNipNeighbors @@ -219,12 +218,11 @@ subroutine mesh_init(ip,el) CHKERRQ(ierr) FE_Nips(FE_geomtype(1_pInt)) = FEM_Zoo_nQuadrature(dimPlex,integrationOrder) - mesh_maxNnodes = FE_Nnodes(1_pInt) mesh_maxNips = FE_Nips(1_pInt) call mesh_FEM_build_ipCoordinates(dimPlex,FEM_Zoo_QuadraturePoints(dimPlex,integrationOrder)%p) call mesh_FEM_build_ipVolumes(dimPlex) - allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems)); mesh_element = 0_pInt + allocate (mesh_element (4_pInt+FE_nodes(1_pInt),mesh_NcpElems)); mesh_element = 0_pInt do j = 1, mesh_NcpElems mesh_element( 1,j) = j mesh_element( 2,j) = mesh_elemType ! elem type From 5814e07021a04d5161898588c4258b2d6f75a450 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 16:55:03 +0200 Subject: [PATCH 07/24] simplified --- src/mesh.f90 | 4 ---- src/meshFEM.f90 | 20 +++++--------------- 2 files changed, 5 insertions(+), 19 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index b353f759e..753fe705c 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -3250,10 +3250,6 @@ subroutine mesh_tell_statistics write(6,*) mesh_Nelems, ' : total number of elements in mesh' write(6,*) mesh_NcpElems, ' : total number of CP elements in mesh' write(6,*) mesh_Nnodes, ' : total number of nodes in mesh' - write(6,*) mesh_maxNnodes, ' : max number of nodes in any CP element' - write(6,*) mesh_maxNips, ' : max number of IPs in any CP element' - write(6,*) mesh_maxNipNeighbors, ' : max number of IP neighbors in any CP element' - write(6,*) mesh_maxNsharedElems, ' : max number of CP elements sharing a node' write(6,'(/,a,/)') ' Input Parser: HOMOGENIZATION/MICROSTRUCTURE' write(6,*) mesh_maxValStateVar(1), ' : maximum homogenization index' write(6,*) mesh_maxValStateVar(2), ' : maximum microstructure index' diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 22568c99d..6abfb83d4 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -59,27 +59,17 @@ use PETScis integer(pInt), dimension(:), allocatable, public, protected :: & mesh_boundaries - - integer(pInt), parameter, public :: & - FE_Nelemtypes = 1_pInt, & - FE_Ngeomtypes = 1_pInt, & - FE_Ncelltypes = 1_pInt, & - FE_maxNnodes = 1_pInt, & - FE_maxNips = 14_pInt - integer(pInt), dimension(FE_Nelemtypes), parameter, public :: FE_geomtype = & !< geometry type of particular element type + integer(pInt), dimension(1_pInt), parameter, public :: FE_geomtype = & !< geometry type of particular element type int([1],pInt) - integer(pInt), dimension(FE_Ngeomtypes), parameter, public :: FE_celltype = & !< cell type that is used by each geometry type + integer(pInt), dimension(1_pInt), parameter, public :: FE_celltype = & !< cell type that is used by each geometry type int([1],pInt) - integer(pInt), dimension(FE_Nelemtypes), parameter, public :: FE_Nnodes = & !< number of nodes that constitute a specific type of element - int([0],pInt) - - integer(pInt), dimension(FE_Ngeomtypes), public :: FE_Nips = & !< number of IPs in a specific type of element + integer(pInt), dimension(1_pInt), public :: FE_Nips = & !< number of IPs in a specific type of element int([0],pInt) - integer(pInt), dimension(FE_Ncelltypes), parameter, public :: FE_NipNeighbors = & !< number of ip neighbors / cell faces in a specific cell type + integer(pInt), dimension(1_pInt), parameter, public :: FE_NipNeighbors = & !< number of ip neighbors / cell faces in a specific cell type int([6],pInt) @@ -222,7 +212,7 @@ subroutine mesh_init(ip,el) call mesh_FEM_build_ipCoordinates(dimPlex,FEM_Zoo_QuadraturePoints(dimPlex,integrationOrder)%p) call mesh_FEM_build_ipVolumes(dimPlex) - allocate (mesh_element (4_pInt+FE_nodes(1_pInt),mesh_NcpElems)); mesh_element = 0_pInt + allocate (mesh_element (4_pInt,mesh_NcpElems)); mesh_element = 0_pInt do j = 1, mesh_NcpElems mesh_element( 1,j) = j mesh_element( 2,j) = mesh_elemType ! elem type From 51390b1acfb027b164f33215f5b33a15b24c0422 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 17:05:01 +0200 Subject: [PATCH 08/24] Nelems /= NcpElems only in case of Abaqus/MSC.Marc --- src/mesh.f90 | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 753fe705c..8bf655d8f 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -22,8 +22,7 @@ module mesh mesh_maxNips, & !< max number of IPs in any CP element mesh_maxNipNeighbors, & !< max number of IP neighbors in any CP element mesh_maxNsharedElems, & !< max number of CP elements sharing a node - mesh_maxNcellnodes, & !< max number of cell nodes in any CP element - mesh_Nelems !< total number of elements in mesh + mesh_maxNcellnodes !< max number of cell nodes in any CP element integer(pInt), dimension(:,:), allocatable, public, protected :: & mesh_element, & !< FEid, type(internal representation), material, texture, node indices as CP IDs @@ -347,6 +346,7 @@ integer(pInt), dimension(:,:), allocatable, private :: & #if defined(Marc4DAMASK) || defined(Abaqus) integer(pInt), private :: & + mesh_Nelems, & !< total number of elements in mesh (including non-DAMASK elements) mesh_maxNnodes, & !< max number of nodes in any CP element mesh_NelemSets character(len=64), dimension(:), allocatable, private :: & @@ -618,8 +618,10 @@ subroutine mesh_init(ip,el) call mesh_tell_statistics endif +#if defined(Marc4DAMASK) || defined(Abaqus) if (usePingPong .and. (mesh_Nelems /= mesh_NcpElems)) & call IO_error(600_pInt) ! ping-pong must be disabled when having non-DAMASK elements +#endif if (debug_e < 1 .or. debug_e > mesh_NcpElems) & call IO_error(602_pInt,ext_msg='element') ! selected element does not exist if (debug_i < 1 .or. debug_i > FE_Nips(FE_geomtype(mesh_element(2_pInt,debug_e)))) & @@ -1155,8 +1157,7 @@ subroutine mesh_spectral_count() implicit none - mesh_Nelems = product(grid(1:2))*grid3 - mesh_NcpElems= mesh_Nelems + mesh_NcpElems= product(grid(1:2))*grid3 mesh_Nnodes = product(grid(1:2) + 1_pInt)*(grid3 + 1_pInt) mesh_NcpElemsGlobal = product(grid) @@ -3247,7 +3248,6 @@ subroutine mesh_tell_statistics !$OMP CRITICAL (write2out) if (iand(myDebug,debug_levelBasic) /= 0_pInt) then write(6,'(/,a,/)') ' Input Parser: STATISTICS' - write(6,*) mesh_Nelems, ' : total number of elements in mesh' write(6,*) mesh_NcpElems, ' : total number of CP elements in mesh' write(6,*) mesh_Nnodes, ' : total number of nodes in mesh' write(6,'(/,a,/)') ' Input Parser: HOMOGENIZATION/MICROSTRUCTURE' From 2fe2c4ca45167ae67d8bd309ac83f3e9bc00edcd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 17:26:13 +0200 Subject: [PATCH 09/24] leaner syntax with sourced allocation --- src/mesh.f90 | 47 +++++++++++++++++++++++------------------------ 1 file changed, 23 insertions(+), 24 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 8bf655d8f..0443f1fe8 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -628,7 +628,7 @@ subroutine mesh_init(ip,el) call IO_error(602_pInt,ext_msg='IP') ! selected element does not have requested IP FEsolving_execElem = [ 1_pInt,mesh_NcpElems ] ! parallel loop bounds set to comprise all DAMASK elements - allocate(FEsolving_execIP(2_pInt,mesh_NcpElems)); FEsolving_execIP = 1_pInt ! parallel loop bounds set to comprise from first IP... + allocate(FEsolving_execIP(2_pInt,mesh_NcpElems), source=1_pInt) ! parallel loop bounds set to comprise from first IP... forall (j = 1_pInt:mesh_NcpElems) FEsolving_execIP(2,j) = FE_Nips(FE_geomtype(mesh_element(2,j))) ! ...up to own IP count for each element allocate(calcMode(mesh_maxNips,mesh_NcpElems)) @@ -1686,8 +1686,8 @@ subroutine mesh_marc_map_elementSets(fileUnit) character(len=300) :: line integer(pInt) :: elemSet = 0_pInt - allocate (mesh_nameElemSet(mesh_NelemSets)) ; mesh_nameElemSet = '' - allocate (mesh_mapElemSet(1_pInt+mesh_maxNelemInSet,mesh_NelemSets)) ; mesh_mapElemSet = 0_pInt + allocate (mesh_nameElemSet(mesh_NelemSets)); mesh_nameElemSet = '' + allocate (mesh_mapElemSet(1_pInt+mesh_maxNelemInSet,mesh_NelemSets), source=0_pInt) 610 FORMAT(A300) @@ -1784,7 +1784,7 @@ subroutine mesh_marc_map_elements(fileUnit) integer(pInt), dimension (1_pInt+mesh_NcpElems) :: contInts integer(pInt) :: i,cpElem = 0_pInt - allocate (mesh_mapFEtoCPelem(2,mesh_NcpElems)) ; mesh_mapFEtoCPelem = 0_pInt + allocate (mesh_mapFEtoCPelem(2,mesh_NcpElems), source = 0_pInt) 610 FORMAT(A300) @@ -1854,7 +1854,7 @@ subroutine mesh_marc_map_nodes(fileUnit) integer(pInt), dimension (mesh_Nnodes) :: node_count integer(pInt) :: i - allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes)) ; mesh_mapFEtoCPnode = 0_pInt + allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes) source=0_pInt) 610 FORMAT(A300) @@ -1901,8 +1901,8 @@ subroutine mesh_marc_build_nodes(fileUnit) character(len=300) :: line integer(pInt) :: i,j,m - allocate ( mesh_node0 (3,mesh_Nnodes) ); mesh_node0 = 0.0_pReal - allocate ( mesh_node (3,mesh_Nnodes) ); mesh_node = 0.0_pReal + allocate ( mesh_node0 (3,mesh_Nnodes), source=0.0_pReal) + allocate ( mesh_node (3,mesh_Nnodes), source=0.0_pReal) 610 FORMAT(A300) @@ -2256,8 +2256,8 @@ subroutine mesh_abaqus_map_elementSets(fileUnit) integer(pInt) :: elemSet = 0_pInt,i logical :: inPart = .false. - allocate (mesh_nameElemSet(mesh_NelemSets)) ; mesh_nameElemSet = '' - allocate (mesh_mapElemSet(1_pInt+mesh_maxNelemInSet,mesh_NelemSets)) ; mesh_mapElemSet = 0_pInt + allocate (mesh_nameElemSet(mesh_NelemSets)); mesh_nameElemSet = '' + allocate (mesh_mapElemSet(1_pInt+mesh_maxNelemInSet,mesh_NelemSets),source=0_pInt) 610 FORMAT(A300) @@ -2308,8 +2308,8 @@ subroutine mesh_abaqus_map_materials(fileUnit) logical :: inPart = .false. character(len=64) :: elemSetName,materialName - allocate (mesh_nameMaterial(mesh_Nmaterials)) ; mesh_nameMaterial = '' - allocate (mesh_mapMaterial(mesh_Nmaterials)) ; mesh_mapMaterial = '' + allocate (mesh_nameMaterial(mesh_Nmaterials)); mesh_nameMaterial = '' + allocate (mesh_mapMaterial(mesh_Nmaterials)); mesh_mapMaterial = '' 610 FORMAT(A300) @@ -2426,7 +2426,7 @@ subroutine mesh_abaqus_map_elements(fileUnit) logical :: materialFound = .false. character (len=64) materialName,elemSetName ! why limited to 64? ABAQUS? - allocate (mesh_mapFEtoCPelem(2,mesh_NcpElems)) ; mesh_mapFEtoCPelem = 0_pInt + allocate (mesh_mapFEtoCPelem(2,mesh_NcpElems), source = 0_pInt) 610 FORMAT(A300) @@ -2489,7 +2489,7 @@ subroutine mesh_abaqus_map_nodes(fileUnit) integer(pInt) :: i,c,cpNode = 0_pInt logical :: inPart = .false. - allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes)) ; mesh_mapFEtoCPnode = 0_pInt + allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes), source=0_pInt) 610 FORMAT(A300) @@ -2551,8 +2551,8 @@ subroutine mesh_abaqus_build_nodes(fileUnit) integer(pInt) :: i,j,m,c logical :: inPart - allocate ( mesh_node0 (3,mesh_Nnodes) ); mesh_node0 = 0.0_pReal - allocate ( mesh_node (3,mesh_Nnodes) ); mesh_node = 0.0_pReal + allocate ( mesh_node0 (3,mesh_Nnodes), source=0.0_pReal) + allocate ( mesh_node (3,mesh_Nnodes), source=0.0_pReal) 610 FORMAT(A300) @@ -2990,7 +2990,7 @@ subroutine mesh_build_sharedElems myDim, & ! dimension index nodeTwin ! node twin in the specified dimension integer(pInt), dimension (mesh_Nnodes) :: node_count - integer(pInt), dimension (:), allocatable :: node_seen + integer(pInt), dimension(:), allocatable :: node_seen allocate(node_seen(maxval(FE_NmatchingNodes))) @@ -3015,8 +3015,7 @@ subroutine mesh_build_sharedElems mesh_maxNsharedElems = int(maxval(node_count),pInt) ! most shared node - allocate(mesh_sharedElem(1+mesh_maxNsharedElems,mesh_Nnodes)) - mesh_sharedElem = 0_pInt + allocate(mesh_sharedElem(1+mesh_maxNsharedElems,mesh_Nnodes),source=0_pInt) do e = 1_pInt,mesh_NcpElems g = FE_geomtype(mesh_element(2,e)) ! get elemGeomType @@ -3238,7 +3237,7 @@ subroutine mesh_tell_statistics if (mesh_maxValStateVar(1) < 1_pInt) call IO_error(error_ID=170_pInt) ! no homogenization specified if (mesh_maxValStateVar(2) < 1_pInt) call IO_error(error_ID=180_pInt) ! no microstructure specified - allocate (mesh_HomogMicro(mesh_maxValStateVar(1),mesh_maxValStateVar(2))); mesh_HomogMicro = 0_pInt + allocate (mesh_HomogMicro(mesh_maxValStateVar(1),mesh_maxValStateVar(2)),source = 0_pInt) do e = 1_pInt,mesh_NcpElems if (mesh_element(3,e) < 1_pInt) call IO_error(error_ID=170_pInt,el=e) ! no homogenization specified if (mesh_element(4,e) < 1_pInt) call IO_error(error_ID=180_pInt,el=e) ! no microstructure specified @@ -3502,11 +3501,11 @@ subroutine mesh_build_FEdata implicit none integer(pInt) :: me - allocate(FE_nodesAtIP(FE_maxmaxNnodesAtIP,FE_maxNips,FE_Ngeomtypes)); FE_nodesAtIP = 0_pInt - allocate(FE_ipNeighbor(FE_maxNipNeighbors,FE_maxNips,FE_Ngeomtypes)); FE_ipNeighbor = 0_pInt - allocate(FE_cell(FE_maxNcellnodesPerCell,FE_maxNips,FE_Ngeomtypes)); FE_cell = 0_pInt - allocate(FE_cellnodeParentnodeWeights(FE_maxNnodes,FE_maxNcellnodes,FE_Nelemtypes)); FE_cellnodeParentnodeWeights = 0.0_pReal - allocate(FE_cellface(FE_maxNcellnodesPerCellface,FE_maxNcellfaces,FE_Ncelltypes)); FE_cellface = 0_pInt + allocate(FE_nodesAtIP(FE_maxmaxNnodesAtIP,FE_maxNips,FE_Ngeomtypes), source=0_pInt) + allocate(FE_ipNeighbor(FE_maxNipNeighbors,FE_maxNips,FE_Ngeomtypes), source=0_pInt) + allocate(FE_cell(FE_maxNcellnodesPerCell,FE_maxNips,FE_Ngeomtypes), source=0_pInt) + allocate(FE_cellnodeParentnodeWeights(FE_maxNnodes,FE_maxNcellnodes,FE_Nelemtypes), source=0.0_pReal) + allocate(FE_cellface(FE_maxNcellnodesPerCellface,FE_maxNcellfaces,FE_Ncelltypes), source=0_pInt) !*** fill FE_nodesAtIP with data *** From 5936397ae7a92e8fdb7086db7d863f901b0adeaa Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 17:37:57 +0200 Subject: [PATCH 10/24] introducing better names allows further simplifications as we do not store max and per elem values any more for number of integration points and number of cell nodes --- src/mesh.f90 | 16 ++++++++++++++-- 1 file changed, 14 insertions(+), 2 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 0443f1fe8..ae28bc7e3 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -19,10 +19,15 @@ module mesh mesh_Nnodes, & !< total number of nodes in mesh mesh_Ncellnodes, & !< total number of cell nodes in mesh (including duplicates) mesh_Ncells, & !< total number of cells in mesh - mesh_maxNips, & !< max number of IPs in any CP element + mesh_NipsPerElem, & !< number of IPs in per element + mesh_NcellnodesPerElem, & !< number of cell nodes per element mesh_maxNipNeighbors, & !< max number of IP neighbors in any CP element - mesh_maxNsharedElems, & !< max number of CP elements sharing a node + mesh_maxNsharedElems !< max number of CP elements sharing a node +!!!! BEGIN DEPRECATED !!!!! + integer(pInt), public, protected :: & + mesh_maxNips, & !< max number of IPs in any CP element mesh_maxNcellnodes !< max number of cell nodes in any CP element +!!!! BEGIN DEPRECATED !!!!! integer(pInt), dimension(:,:), allocatable, public, protected :: & mesh_element, & !< FEid, type(internal representation), material, texture, node indices as CP IDs @@ -639,6 +644,13 @@ subroutine mesh_init(ip,el) calcMode(ip,el) = .true. ! first ip,el needs to be already pingponged to "calc" #endif +!!!! COMPATIBILITY HACK !!!! +! for a homogeneous mesh, all elements have the same number of IPs and and cell nodes. +! hence, xxPerElem instead of maxXX + mesh_NipsPerElem = mesh_maxNips + mesh_NcellnodesPerElem = mesh_maxNcellnodes +!!!!!!!!!!!!!!!!!!!!!!!! + end subroutine mesh_init From cf6d388a6bebc655d491d68322e0a5d1c1e5513b Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 17:57:48 +0200 Subject: [PATCH 11/24] consistent order of solver-specific functionality --- src/mesh.f90 | 59 +++++++++++++++++++++++----------------------------- 1 file changed, 26 insertions(+), 33 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index ae28bc7e3..d8a9eefdc 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -335,7 +335,7 @@ integer(pInt), dimension(:,:), allocatable, private :: & 4 & ! element 21 (3D 20node 27ip) ],pInt) -#ifdef Spectral +#if defined(Spectral) integer(pInt), dimension(3), public, protected :: & grid !< (global) grid integer(pInt), public, protected :: & @@ -347,9 +347,7 @@ integer(pInt), dimension(:,:), allocatable, private :: & real(pReal), public, protected :: & size3, & !< (local) size in 3rd direction size3offset !< (local) size offset in 3rd direction -#endif - -#if defined(Marc4DAMASK) || defined(Abaqus) +#elif defined(Marc4DAMASK) || defined(Abaqus) integer(pInt), private :: & mesh_Nelems, & !< total number of elements in mesh (including non-DAMASK elements) mesh_maxNnodes, & !< max number of nodes in any CP element @@ -364,47 +362,54 @@ integer(pInt), dimension(:,:), allocatable, private :: & mesh_mapFEtoCPelem, & !< [sorted FEid, corresponding CPid] mesh_mapFEtoCPnode !< [sorted FEid, corresponding CPid] #endif - -#ifdef Marc4DAMASK +#if defined(Marc4DAMASK) integer(pInt), private :: & MarcVersion, & !< Version of input file format (Marc only) hypoelasticTableStyle, & !< Table style (Marc only) initialcondTableStyle !< Table style (Marc only) integer(pInt), dimension(:), allocatable, private :: & Marc_matNumber !< array of material numbers for hypoelastic material (Marc only) -#endif - -#ifdef Abaqus +#elif defined(Abaqus) logical, private :: noPart !< for cases where the ABAQUS input file does not use part/assembly information #endif public :: & mesh_init, & -#if defined(Marc4DAMASK) || defined(Abaqus) - mesh_FEasCP, & -#endif mesh_build_cellnodes, & mesh_build_ipVolumes, & mesh_build_ipCoordinates, & mesh_cellCenterCoordinates, & mesh_get_Ncellnodes, & mesh_get_unitlength, & - mesh_get_nodeAtIP -#ifdef Spectral - public :: & + mesh_get_nodeAtIP, & +#if defined(Spectral) mesh_spectral_getGrid, & mesh_spectral_getSize +#elif defined(Marc4DAMASK) || defined(Abaqus) + mesh_FEasCP #endif private :: & -#ifdef Spectral + mesh_get_damaskOptions, & + mesh_build_cellconnectivity, & + mesh_build_ipAreas, & + mesh_tell_statistics, & + FE_mapElemtype, & + mesh_faceMatch, & + mesh_build_FEdata, & +#if defined(Spectral) mesh_spectral_getHomogenization, & mesh_spectral_count, & mesh_spectral_count_cpSizes, & mesh_spectral_build_nodes, & mesh_spectral_build_elements, & - mesh_spectral_build_ipNeighborhood, & -#elif defined Marc4DAMASK + mesh_spectral_build_ipNeighborhood +#elif defined(Marc4DAMASK) || defined(Abaqus) + mesh_build_nodeTwins, & + mesh_build_sharedElems, & + mesh_build_ipNeighborhood, & +#endif +#if defined(Marc4DAMASK) mesh_marc_get_fileFormat, & mesh_marc_get_tableStyles, & mesh_marc_get_matNumber, & @@ -416,8 +421,8 @@ integer(pInt), dimension(:,:), allocatable, private :: & mesh_marc_map_nodes, & mesh_marc_build_nodes, & mesh_marc_count_cpSizes, & - mesh_marc_build_elements, & -#elif defined Abaqus + mesh_marc_build_elements +#elif defined(Abaqus) mesh_abaqus_count_nodesAndElements, & mesh_abaqus_count_elementSets, & mesh_abaqus_count_materials, & @@ -428,20 +433,8 @@ integer(pInt), dimension(:,:), allocatable, private :: & mesh_abaqus_map_nodes, & mesh_abaqus_build_nodes, & mesh_abaqus_count_cpSizes, & - mesh_abaqus_build_elements, & + mesh_abaqus_build_elements #endif -#if defined(Marc4DAMASK) || defined(Abaqus) - mesh_build_nodeTwins, & - mesh_build_sharedElems, & - mesh_build_ipNeighborhood, & -#endif - mesh_get_damaskOptions, & - mesh_build_cellconnectivity, & - mesh_build_ipAreas, & - mesh_tell_statistics, & - FE_mapElemtype, & - mesh_faceMatch, & - mesh_build_FEdata contains From 4b14cc5560573d777e484d64b15ced19caea9106 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 18:06:18 +0200 Subject: [PATCH 12/24] calcmode only needed for Abaqus and MSC.Marc --- src/CPFEM.f90 | 2 +- src/mesh.f90 | 6 ++---- src/meshFEM.f90 | 18 +++++++++++------- 3 files changed, 14 insertions(+), 12 deletions(-) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index 6caeaf57c..dc6ba07dc 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -98,7 +98,7 @@ subroutine CPFEM_initAll(el,ip) call config_init call math_init call FE_init - call mesh_init(ip, el) ! pass on coordinates to alter calcMode of first ip + call mesh_init(ip, el) call lattice_init call material_init call constitutive_init diff --git a/src/mesh.f90 b/src/mesh.f90 index d8a9eefdc..953d29eef 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -481,9 +481,9 @@ subroutine mesh_init(ip,el) FEsolving_execElem, & #ifndef Spectral modelName, & + calcMode #endif FEsolving_execIP, & - calcMode implicit none #ifdef Spectral @@ -629,12 +629,10 @@ subroutine mesh_init(ip,el) allocate(FEsolving_execIP(2_pInt,mesh_NcpElems), source=1_pInt) ! parallel loop bounds set to comprise from first IP... forall (j = 1_pInt:mesh_NcpElems) FEsolving_execIP(2,j) = FE_Nips(FE_geomtype(mesh_element(2,j))) ! ...up to own IP count for each element +#if defined(Marc4DAMASK) || defined(Abaqus) allocate(calcMode(mesh_maxNips,mesh_NcpElems)) calcMode = .false. ! pretend to have collected what first call is asking (F = I) -#if defined(Marc4DAMASK) || defined(Abaqus) calcMode(ip,mesh_FEasCP('elem',el)) = .true. ! first ip,el needs to be already pingponged to "calc" -#else - calcMode(ip,el) = .true. ! first ip,el needs to be already pingponged to "calc" #endif !!!! COMPATIBILITY HACK !!!! diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 6abfb83d4..ca666f047 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -25,8 +25,12 @@ use PETScis mesh_NcpElems, & !< total number of CP elements in mesh mesh_NcpElemsGlobal, & mesh_Nnodes, & !< total number of nodes in mesh - mesh_maxNips, & !< max number of IPs in any CP element + mesh_maxNipNeighbors, & !< max number of IP neighbors in any CP element mesh_maxNipNeighbors +!!!! BEGIN DEPRECATED !!!!! + integer(pInt), public, protected :: & + mesh_maxNips, & !< max number of IPs in any CP element +!!!! BEGIN DEPRECATED !!!!! real(pReal), public, protected :: charLength @@ -117,8 +121,7 @@ subroutine mesh_init(ip,el) worldsize use FEsolving, only: & FEsolving_execElem, & - FEsolving_execIP, & - calcMode + FEsolving_execIP use FEM_Zoo, only: & FEM_Zoo_nQuadrature, & FEM_Zoo_QuadraturePoints @@ -231,10 +234,11 @@ subroutine mesh_init(ip,el) allocate(FEsolving_execIP(2_pInt,mesh_NcpElems)); FEsolving_execIP = 1_pInt ! parallel loop bounds set to comprise from first IP... forall (j = 1_pInt:mesh_NcpElems) FEsolving_execIP(2,j) = FE_Nips(FE_geomtype(mesh_element(2,j))) ! ...up to own IP count for each element - if (allocated(calcMode)) deallocate(calcMode) - allocate(calcMode(mesh_maxNips,mesh_NcpElems)) - calcMode = .false. ! pretend to have collected what first call is asking (F = I) - calcMode(ip,el) = .true. ! first ip,el needs to be already pingponged to "calc" +!!!! COMPATIBILITY HACK !!!! +! for a homogeneous mesh, all elements have the same number of IPs and and cell nodes. +! hence, xxPerElem instead of maxXX + mesh_NipsPerElem = mesh_maxNips +!!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init From 854d99250c9f18793c57d6e5e72667924d332211 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 18:40:14 +0200 Subject: [PATCH 13/24] mesh_element(1,:) only used for debug output set to -1 at the moment to indicate that it is not used. Re-implementation should be done for MSC.Marc and Abaqus only. --- src/mesh.f90 | 6 +++--- src/meshFEM.f90 | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 953d29eef..384804ef6 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -1310,7 +1310,7 @@ subroutine mesh_spectral_build_elements(fileUnit) e = 0_pInt do while (e < mesh_NcpElems) ! fill expected number of elements, stop at end of data (or blank line!) e = e+1_pInt ! valid element entry - mesh_element( 1,e) = e ! FE id + mesh_element( 1,e) = -1_pInt ! DEPRECATED mesh_element( 2,e) = elemType ! elem type mesh_element( 3,e) = homog ! homogenization mesh_element( 4,e) = mesh_microGlobal(e+elemOffset) ! microstructure @@ -2025,7 +2025,7 @@ subroutine mesh_marc_build_elements(fileUnit) chunkPos = IO_stringPos(line) e = mesh_FEasCP('elem',IO_intValue(line,chunkPos,1_pInt)) if (e /= 0_pInt) then ! disregard non CP elems - mesh_element(1,e) = IO_IntValue (line,chunkPos,1_pInt) ! FE id + mesh_element(1,e) = -1_pInt ! DEPRECATED t = FE_mapElemtype(IO_StringValue(line,chunkPos,2_pInt)) ! elem type if (mesh_elemType /= t .and. mesh_elemType /= -1_pInt) & call IO_error(191,el=t,ip=mesh_elemType) @@ -2710,7 +2710,7 @@ subroutine mesh_abaqus_build_elements(fileUnit) chunkPos = IO_stringPos(line) ! limit to 64 nodes max e = mesh_FEasCP('elem',IO_intValue(line,chunkPos,1_pInt)) if (e /= 0_pInt) then ! disregard non CP elems - mesh_element(1,e) = IO_intValue(line,chunkPos,1_pInt) ! FE id + mesh_element(1,e) = -1_pInt ! DEPRECATED if (mesh_elemType /= t .and. mesh_elemType /= -1_pInt) & call IO_error(191,el=t,ip=mesh_elemType) mesh_elemType = t diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index ca666f047..5202746a9 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -217,7 +217,7 @@ subroutine mesh_init(ip,el) allocate (mesh_element (4_pInt,mesh_NcpElems)); mesh_element = 0_pInt do j = 1, mesh_NcpElems - mesh_element( 1,j) = j + mesh_element( 1,j) = -1_pInt ! DEPRECATED mesh_element( 2,j) = mesh_elemType ! elem type mesh_element( 3,j) = 1_pInt ! homogenization call DMGetLabelValue(geomMesh,'material',j-1,mesh_element(4,j),ierr) From 67483487eac41e4784360cfc824d30796ccc0176 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 18:50:54 +0200 Subject: [PATCH 14/24] more descriptive names --- src/mesh.f90 | 15 +++++++++++---- src/meshFEM.f90 | 7 +++++++ 2 files changed, 18 insertions(+), 4 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 384804ef6..c32ddbce6 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -29,6 +29,10 @@ module mesh mesh_maxNcellnodes !< max number of cell nodes in any CP element !!!! BEGIN DEPRECATED !!!!! + integer(pInt), dimension(:), allocatable, public, protected :: & + mesh_homogenization, & !< homogenization ID of each element + mesh_microstructure !< homogenization ID of each element + integer(pInt), dimension(:,:), allocatable, public, protected :: & mesh_element, & !< FEid, type(internal representation), material, texture, node indices as CP IDs mesh_sharedElem, & !< entryCount and list of elements containing node @@ -640,6 +644,9 @@ subroutine mesh_init(ip,el) ! hence, xxPerElem instead of maxXX mesh_NipsPerElem = mesh_maxNips mesh_NcellnodesPerElem = mesh_maxNcellnodes +! better name + mesh_homogenization = mesh_element(3,:) + mesh_microstructure = mesh_element(4,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init @@ -1248,7 +1255,7 @@ subroutine mesh_spectral_build_elements(fileUnit) elemOffset integer(pInt), dimension(:), allocatable :: & microstructures, & - mesh_microGlobal + microGlobal integer(pInt), dimension(1,1) :: & dummySet = 0_pInt character(len=65536) :: & @@ -1287,7 +1294,7 @@ subroutine mesh_spectral_build_elements(fileUnit) enddo allocate (mesh_element (4_pInt+8_pInt,mesh_NcpElems), source = 0_pInt) allocate (microstructures (1_pInt+maxIntCount), source = 1_pInt) - allocate (mesh_microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) + allocate (microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) !-------------------------------------------------------------------------------------------------- ! read in microstructures @@ -1301,7 +1308,7 @@ subroutine mesh_spectral_build_elements(fileUnit) microstructures = IO_continuousIntValues(fileUnit,maxIntCount,dummyName,dummySet,0_pInt) ! get affected elements do i = 1_pInt,microstructures(1_pInt) e = e+1_pInt ! valid element entry - mesh_microGlobal(e) = microstructures(1_pInt+i) + microGlobal(e) = microstructures(1_pInt+i) enddo enddo @@ -1313,7 +1320,7 @@ subroutine mesh_spectral_build_elements(fileUnit) mesh_element( 1,e) = -1_pInt ! DEPRECATED mesh_element( 2,e) = elemType ! elem type mesh_element( 3,e) = homog ! homogenization - mesh_element( 4,e) = mesh_microGlobal(e+elemOffset) ! microstructure + mesh_element( 4,e) = microGlobal(e+elemOffset) ! microstructure mesh_element( 5,e) = e + (e-1_pInt)/grid(1) + & ((e-1_pInt)/(grid(1)*grid(2)))*(grid(1)+1_pInt) ! base node mesh_element( 6,e) = mesh_element(5,e) + 1_pInt diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 5202746a9..afead2b07 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -32,6 +32,10 @@ use PETScis mesh_maxNips, & !< max number of IPs in any CP element !!!! BEGIN DEPRECATED !!!!! + integer(pInt), dimension(:), allocatable, public, protected :: & + mesh_homogenization, & !< homogenization ID of each element + mesh_microstructure !< homogenization ID of each element + real(pReal), public, protected :: charLength integer(pInt), dimension(:,:), allocatable, public, protected :: & @@ -238,6 +242,9 @@ subroutine mesh_init(ip,el) ! for a homogeneous mesh, all elements have the same number of IPs and and cell nodes. ! hence, xxPerElem instead of maxXX mesh_NipsPerElem = mesh_maxNips +! better name + mesh_homogenization = mesh_element(3,:) + mesh_microstructure = mesh_element(4,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init From ebef12e44600679f85a52a05d9f56c2dc8f4163a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 19:04:17 +0200 Subject: [PATCH 15/24] syntax errors in declaration --- src/mesh.f90 | 6 +++--- src/meshFEM.f90 | 4 ++-- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index c32ddbce6..8736f10ac 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -482,12 +482,12 @@ subroutine mesh_init(ip,el) numerics_unitlength, & worldrank use FEsolving, only: & - FEsolving_execElem, & #ifndef Spectral modelName, & - calcMode + calcMode, & #endif - FEsolving_execIP, & + FEsolving_execElem, & + FEsolving_execIP implicit none #ifdef Spectral diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index afead2b07..68113eaf2 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -25,11 +25,11 @@ use PETScis mesh_NcpElems, & !< total number of CP elements in mesh mesh_NcpElemsGlobal, & mesh_Nnodes, & !< total number of nodes in mesh - mesh_maxNipNeighbors, & !< max number of IP neighbors in any CP element + mesh_NipsPerElem, & !< number of IPs in per element mesh_maxNipNeighbors !!!! BEGIN DEPRECATED !!!!! integer(pInt), public, protected :: & - mesh_maxNips, & !< max number of IPs in any CP element + mesh_maxNips !< max number of IPs in any CP element !!!! BEGIN DEPRECATED !!!!! integer(pInt), dimension(:), allocatable, public, protected :: & From 72b87b0a9bcbfffe7d8fe761e246566d37205171 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 19:58:43 +0200 Subject: [PATCH 16/24] better name --- src/mesh.f90 | 14 ++++++++------ src/meshFEM.f90 | 4 +++- 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 8736f10ac..87d5737dd 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -34,7 +34,8 @@ module mesh mesh_microstructure !< homogenization ID of each element integer(pInt), dimension(:,:), allocatable, public, protected :: & - mesh_element, & !< FEid, type(internal representation), material, texture, node indices as CP IDs + mesh_CPnodeID, & + mesh_element, & !DEPRECATED mesh_sharedElem, & !< entryCount and list of elements containing node mesh_nodeTwins !< node twins are surface nodes that lie exactly on opposite sides of the mesh (surfaces nodes with equal coordinate values in two dimensions) @@ -647,6 +648,7 @@ subroutine mesh_init(ip,el) ! better name mesh_homogenization = mesh_element(3,:) mesh_microstructure = mesh_element(4,:) + mesh_CPnodeID = mesh_element(5:4+mesh_NipsPerElem,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init @@ -1292,9 +1294,9 @@ subroutine mesh_spectral_build_elements(fileUnit) i = IO_countContinuousIntValues(fileUnit) maxIntCount = max(maxIntCount, i) enddo - allocate (mesh_element (4_pInt+8_pInt,mesh_NcpElems), source = 0_pInt) - allocate (microstructures (1_pInt+maxIntCount), source = 1_pInt) - allocate (microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) + allocate(mesh_element (4_pInt+8_pInt,mesh_NcpElems), source = 0_pInt) + allocate(microstructures (1_pInt+maxIntCount), source = 1_pInt) + allocate(microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) !-------------------------------------------------------------------------------------------------- ! read in microstructures @@ -2016,7 +2018,7 @@ subroutine mesh_marc_build_elements(fileUnit) integer(pInt), dimension(1_pInt+mesh_NcpElems) :: contInts integer(pInt) :: i,j,t,sv,myVal,e,nNodesAlreadyRead - allocate (mesh_element(4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) + allocate(mesh_element(4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) mesh_elemType = -1_pInt 610 FORMAT(A300) @@ -2687,7 +2689,7 @@ subroutine mesh_abaqus_build_elements(fileUnit) character (len=64) :: materialName,elemSetName character(len=300) :: line - allocate (mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) + allocate(mesh_element (4_pInt+mesh_maxNnodes,mesh_NcpElems), source=0_pInt) mesh_elemType = -1_pInt 610 FORMAT(A300) diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 68113eaf2..31c1932c6 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -39,7 +39,8 @@ use PETScis real(pReal), public, protected :: charLength integer(pInt), dimension(:,:), allocatable, public, protected :: & - mesh_element !< FEid, type(internal representation), material, texture, node indices as CP IDs + mesh_CPnodeID, & + mesh_element !DEPRECATED real(pReal), dimension(:,:), allocatable, public :: & mesh_node !< node x,y,z coordinates (after deformation! ONLY FOR MARC!!!) @@ -245,6 +246,7 @@ subroutine mesh_init(ip,el) ! better name mesh_homogenization = mesh_element(3,:) mesh_microstructure = mesh_element(4,:) + mesh_CPnodeID = mesh_element(5:4+mesh_NipsPerElem,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init From c42eb87a3386d9546f9839a87805ca6abfb76016 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 20:53:35 +0200 Subject: [PATCH 17/24] using arrays with new names --- src/material.f90 | 63 ++++++++++++++++++++++++++---------------------- src/mesh.f90 | 2 +- 2 files changed, 35 insertions(+), 30 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index 812b0c55d..d4600948a 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -199,6 +199,7 @@ module material integer(pInt), dimension(:,:,:), allocatable, public :: & material_phase !< phase (index) of each grain,IP,element +! DEPRECATED. DID WE EVER ALLOWED DIFFERENT HOMOGENIZATION SCHEMES WITHIN ONE ELEMENT? integer(pInt), dimension(:,:), allocatable, public :: & material_homog !< homogenization (index) of each IP,element type(tPlasticState), allocatable, dimension(:), public :: & @@ -362,10 +363,10 @@ subroutine material_init() phase_name, & texture_name use mesh, only: & + mesh_homogenization, & + mesh_NipsPerElem, & mesh_maxNips, & mesh_NcpElems, & - mesh_element, & - FE_Nips, & FE_geomtype implicit none @@ -480,11 +481,11 @@ subroutine material_init() allocate(CrystallitePosition (size(config_phase)), source=0_pInt) ElemLoop:do e = 1_pInt,mesh_NcpElems - myHomog = mesh_element(3,e) - IPloop:do i = 1_pInt,FE_Nips(FE_geomtype(mesh_element(2,e))) + myHomog = mesh_homogenization(e) + IPloop:do i = 1_pInt, mesh_NipsPerElem HomogenizationPosition(myHomog) = HomogenizationPosition(myHomog) + 1_pInt mappingHomogenization(1:2,i,e) = [HomogenizationPosition(myHomog),myHomog] - GrainLoop:do g = 1_pInt,homogenization_Ngrains(mesh_element(3,e)) + GrainLoop:do g = 1_pInt,homogenization_Ngrains(myHomog) phase = material_phase(g,i,e) ConstitutivePosition(phase) = ConstitutivePosition(phase)+1_pInt ! not distinguishing between instances of same phase phaseAt(g,i,e) = phase @@ -519,10 +520,10 @@ end subroutine material_init subroutine material_parseHomogenization use config, only : & config_homogenization + use mesh, only: & + mesh_homogenization use IO, only: & IO_error - use mesh, only: & - mesh_element implicit none integer(pInt) :: h @@ -549,7 +550,8 @@ subroutine material_parseHomogenization allocate(porosity_initialPhi(size(config_homogenization)), source=1.0_pReal) allocate(hydrogenflux_initialCh(size(config_homogenization)), source=0.0_pReal) - forall (h = 1_pInt:size(config_homogenization)) homogenization_active(h) = any(mesh_element(3,:) == h) + forall (h = 1_pInt:size(config_homogenization)) & + homogenization_active(h) = any(mesh_homogenization == h) do h=1_pInt, size(config_homogenization) @@ -685,7 +687,7 @@ subroutine material_parseMicrostructure config_microstructure, & microstructure_name use mesh, only: & - mesh_element, & + mesh_microstructure, & mesh_NcpElems implicit none @@ -701,10 +703,11 @@ subroutine material_parseMicrostructure allocate(microstructure_active(size(config_microstructure)), source=.false.) allocate(microstructure_elemhomo(size(config_microstructure)), source=.false.) - if(any(mesh_element(4,1:mesh_NcpElems) > size(config_microstructure))) & + if(any(mesh_microstructure > size(config_microstructure))) & call IO_error(155_pInt,ext_msg='More microstructures in geometry than sections in material.config') - forall (e = 1_pInt:mesh_NcpElems) microstructure_active(mesh_element(4,e)) = .true. ! current microstructure used in model? Elementwise view, maximum N operations for N elements + forall (e = 1_pInt:mesh_NcpElems) & + microstructure_active(mesh_microstructure(e)) = .true. ! current microstructure used in model? Elementwise view, maximum N operations for N elements do m=1_pInt, size(config_microstructure) microstructure_Nconstituents(m) = config_microstructure(m)%countKeys('(constituent)') @@ -1082,11 +1085,13 @@ subroutine material_populateGrains math_sampleFiberOri, & math_symmetricEulers use mesh, only: & - mesh_element, & + mesh_NipsPerElem, & + mesh_elemType, & + mesh_homogenization, & + mesh_microstructure, & mesh_maxNips, & mesh_NcpElems, & mesh_ipVolume, & - FE_Nips, & FE_geomtype use config, only: & config_homogenization, & @@ -1136,14 +1141,14 @@ subroutine material_populateGrains ! populating homogenization schemes in each !-------------------------------------------------------------------------------------------------- do e = 1_pInt, mesh_NcpElems - material_homog(1_pInt:FE_Nips(FE_geomtype(mesh_element(2,e))),e) = mesh_element(3,e) + material_homog(1_pInt:mesh_NipsPerElem,e) = mesh_homogenization(e) enddo !-------------------------------------------------------------------------------------------------- ! precounting of elements for each homog/micro pair do e = 1_pInt, mesh_NcpElems - homog = mesh_element(3,e) - micro = mesh_element(4,e) + homog = mesh_homogenization(e) + micro = mesh_microstructure(e) Nelems(homog,micro) = Nelems(homog,micro) + 1_pInt enddo allocate(elemsOfHomogMicro(size(config_homogenization),size(config_microstructure))) @@ -1160,9 +1165,9 @@ subroutine material_populateGrains ! identify maximum grain count per IP (from element) and find grains per homog/micro pair Nelems = 0_pInt ! reuse as counter elementLooping: do e = 1_pInt,mesh_NcpElems - t = FE_geomtype(mesh_element(2,e)) - homog = mesh_element(3,e) - micro = mesh_element(4,e) + t = mesh_elemType + homog = mesh_homogenization(e) + micro = mesh_microstructure(e) if (homog < 1_pInt .or. homog > size(config_homogenization)) & ! out of bounds call IO_error(154_pInt,e,0_pInt,0_pInt) if (micro < 1_pInt .or. micro > size(config_microstructure)) & ! out of bounds @@ -1170,7 +1175,7 @@ subroutine material_populateGrains if (microstructure_elemhomo(micro)) then ! how many grains are needed at this element? dGrains = homogenization_Ngrains(homog) ! only one set of Ngrains (other IPs are plain copies) else - dGrains = homogenization_Ngrains(homog) * FE_Nips(t) ! each IP has Ngrains + dGrains = homogenization_Ngrains(homog) * mesh_NipsPerElem ! each IP has Ngrains endif Ngrains(homog,micro) = Ngrains(homog,micro) + dGrains ! total grain count Nelems(homog,micro) = Nelems(homog,micro) + 1_pInt ! total element count @@ -1183,8 +1188,8 @@ subroutine material_populateGrains allocate(orientationOfGrain(3,maxval(Ngrains)),source=0.0_pReal) ! reserve memory for maximum case if (iand(myDebug,debug_levelBasic) /= 0_pInt) then - write(6,'(/,a/)') ' MATERIAL grain population' - write(6,'(a32,1x,a32,1x,a6)') 'homogenization_name','microstructure_name','grain#' + write(6,'(/,a/)') ' MATERIAL grain population' + write(6,'(a32,1x,a32,1x,a6)') 'homogenization_name','microstructure_name','grain#' endif homogenizationLoop: do homog = 1_pInt,size(config_homogenization) dGrains = homogenization_Ngrains(homog) ! grain number per material point @@ -1204,16 +1209,16 @@ subroutine material_populateGrains do hme = 1_pInt, Nelems(homog,micro) e = elemsOfHomogMicro(homog,micro)%p(hme) ! my combination of homog and micro, only perform calculations for elements with homog, micro combinations which is indexed in cpElemsindex - t = FE_geomtype(mesh_element(2,e)) + t = mesh_elemType if (microstructure_elemhomo(micro)) then ! homogeneous distribution of grains over each element's IPs - volumeOfGrain(grain+1_pInt:grain+dGrains) = sum(mesh_ipVolume(1:FE_Nips(t),e))/& + volumeOfGrain(grain+1_pInt:grain+dGrains) = sum(mesh_ipVolume(1:mesh_NipsPerElem,e))/& real(dGrains,pReal) ! each grain combines size of all IPs in that element grain = grain + dGrains ! wind forward by Ngrains@IP else - forall (i = 1_pInt:FE_Nips(t)) & ! loop over IPs + forall (i = 1_pInt:mesh_NipsPerElem) & ! loop over IPs volumeOfGrain(grain+(i-1)*dGrains+1_pInt:grain+i*dGrains) = & mesh_ipVolume(i,e)/real(dGrains,pReal) ! assign IPvolume/Ngrains@IP to all grains of IP - grain = grain + FE_Nips(t) * dGrains ! wind forward by Nips*Ngrains@IP + grain = grain + mesh_NipsPerElem * dGrains ! wind forward by Nips*Ngrains@IP endif enddo @@ -1367,11 +1372,11 @@ subroutine material_populateGrains do hme = 1_pInt, Nelems(homog,micro) e = elemsOfHomogMicro(homog,micro)%p(hme) ! only perform calculations for elements with homog, micro combinations which is indexed in cpElemsindex - t = FE_geomtype(mesh_element(2,e)) + t = mesh_elemType if (microstructure_elemhomo(micro)) then ! homogeneous distribution of grains over each element's IPs m = 1_pInt ! process only first IP else - m = FE_Nips(t) ! process all IPs + m = mesh_NipsPerElem endif do i = 1_pInt, m ! loop over necessary IPs @@ -1409,7 +1414,7 @@ subroutine material_populateGrains enddo - do i = i, FE_Nips(t) ! loop over IPs to (possibly) distribute copies from first IP + do i = i, mesh_NipsPerElem ! loop over IPs to (possibly) distribute copies from first IP material_volume (1_pInt:dGrains,i,e) = material_volume (1_pInt:dGrains,1,e) material_phase (1_pInt:dGrains,i,e) = material_phase (1_pInt:dGrains,1,e) material_texture(1_pInt:dGrains,i,e) = material_texture(1_pInt:dGrains,1,e) diff --git a/src/mesh.f90 b/src/mesh.f90 index 87d5737dd..19035b067 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -13,7 +13,7 @@ module mesh private integer(pInt), public, protected :: & mesh_NcpElems, & !< total number of CP elements in local mesh - mesh_ElemType, & !< Element type of the mesh (only support homogeneous meshes) + mesh_elemType, & !< Element type of the mesh (only support homogeneous meshes) mesh_maxNelemInSet, & mesh_Nmaterials, & mesh_Nnodes, & !< total number of nodes in mesh From ceb385ef39a09828eca0cc7b01327d19883e8030 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 23 Sep 2018 21:31:30 +0200 Subject: [PATCH 18/24] calcMode not needed for spectral and FEM --- src/CPFEM2.f90 | 8 +++----- src/DAMASK_FEM.f90 | 2 +- src/DAMASK_spectral.f90 | 2 +- src/mesh.f90 | 2 +- 4 files changed, 6 insertions(+), 8 deletions(-) diff --git a/src/CPFEM2.f90 b/src/CPFEM2.f90 index 89e65f5fd..1d1fb8342 100644 --- a/src/CPFEM2.f90 +++ b/src/CPFEM2.f90 @@ -18,7 +18,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief call (thread safe) all module initializations !-------------------------------------------------------------------------------------------------- -subroutine CPFEM_initAll(el,ip) +subroutine CPFEM_initAll() use prec, only: & pInt use prec, only: & @@ -55,10 +55,8 @@ subroutine CPFEM_initAll(el,ip) #endif implicit none - integer(pInt), intent(in) :: el, & !< FE el number - ip !< FE integration point number - call DAMASK_interface_init ! Spectral and FEM interface to commandline + call DAMASK_interface_init ! Spectral and FEM interface to commandline call prec_init call IO_init #ifdef FEM @@ -69,7 +67,7 @@ subroutine CPFEM_initAll(el,ip) call config_init call math_init call FE_init - call mesh_init(ip, el) ! pass on coordinates to alter calcMode of first ip + call mesh_init call lattice_init call material_init call constitutive_init diff --git a/src/DAMASK_FEM.f90 b/src/DAMASK_FEM.f90 index ee425585c..f561e38bd 100644 --- a/src/DAMASK_FEM.f90 +++ b/src/DAMASK_FEM.f90 @@ -142,7 +142,7 @@ program DAMASK_FEM quit !-------------------------------------------------------------------------------------------------- ! init DAMASK (all modules) - call CPFEM_initAll(el = 1_pInt, ip = 1_pInt) + call CPFEM_initAll write(6,'(/,a)') ' <<<+- DAMASK_FEM init -+>>>' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" diff --git a/src/DAMASK_spectral.f90 b/src/DAMASK_spectral.f90 index ee5af421c..bab5cfd5c 100644 --- a/src/DAMASK_spectral.f90 +++ b/src/DAMASK_spectral.f90 @@ -154,7 +154,7 @@ program DAMASK_spectral !-------------------------------------------------------------------------------------------------- ! init DAMASK (all modules) - call CPFEM_initAll(el = 1_pInt, ip = 1_pInt) + call CPFEM_initAll write(6,'(/,a)') ' <<<+- DAMASK_spectral init -+>>>' write(6,'(/,a,/)') ' Roters et al., Computational Materials Science, 2018' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() diff --git a/src/mesh.f90 b/src/mesh.f90 index 19035b067..356594e53 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -497,7 +497,7 @@ subroutine mesh_init(ip,el) integer :: ierr, worldsize #endif integer(pInt), parameter :: FILEUNIT = 222_pInt - integer(pInt), intent(in) :: el, ip + integer(pInt), intent(in), optional :: el, ip integer(pInt) :: j logical :: myDebug From e8f687a99cac63de26599c9b930377b253e8b03d Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 27 Sep 2018 20:18:37 +0200 Subject: [PATCH 19/24] typo, Abaqus/Marc did not compile --- src/mesh.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index 356594e53..18f986b33 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -1866,7 +1866,7 @@ subroutine mesh_marc_map_nodes(fileUnit) integer(pInt), dimension (mesh_Nnodes) :: node_count integer(pInt) :: i - allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes) source=0_pInt) + allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes),source=0_pInt) 610 FORMAT(A300) From 49a4202e2651df72d7b0b2aff6f28eaebdf2d49a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 29 Sep 2018 10:31:45 +0200 Subject: [PATCH 20/24] not needed here --- src/meshFEM.f90 | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 31c1932c6..5d81fa0b3 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -104,7 +104,7 @@ contains !> @brief initializes the mesh by calling all necessary private routines the mesh module !! Order and routines strongly depend on type of solver !-------------------------------------------------------------------------------------------------- -subroutine mesh_init(ip,el) +subroutine mesh_init() use DAMASK_interface use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) use IO, only: & @@ -133,7 +133,6 @@ subroutine mesh_init(ip,el) implicit none integer(pInt), parameter :: FILEUNIT = 222_pInt - integer(pInt), intent(in) :: el, ip integer(pInt) :: j integer(pInt), allocatable, dimension(:) :: chunkPos integer :: dimPlex From 97d2c2b3531acb5f72d293a1e0f24a496b848b8f Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 2 Oct 2018 22:54:44 +0200 Subject: [PATCH 21/24] mesh_element(1) has a size of 4 fixed sigsegv --- src/meshFEM.f90 | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index 0aab11779..c58cbf2e3 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -41,7 +41,7 @@ use PETScis real(pReal), public, protected :: charLength integer(pInt), dimension(:,:), allocatable, public, protected :: & - mesh_CPnodeID, & + !mesh_CPnodeID, & mesh_element !DEPRECATED real(pReal), dimension(:,:), allocatable, public :: & @@ -245,7 +245,6 @@ subroutine mesh_init() ! better name mesh_homogenization = mesh_element(3,:) mesh_microstructure = mesh_element(4,:) - mesh_CPnodeID = mesh_element(5:4+mesh_NipsPerElem,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init From df473302f4a624429e124a21ff4ed3864aa8e363 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 4 Oct 2018 06:03:48 +0200 Subject: [PATCH 22/24] consistent (and better understandable) names --- src/material.f90 | 28 ++++++++++++++-------------- src/mesh.f90 | 10 +++++----- src/meshFEM.f90 | 13 +++++-------- 3 files changed, 24 insertions(+), 27 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index d4600948a..98fa1220f 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -363,7 +363,7 @@ subroutine material_init() phase_name, & texture_name use mesh, only: & - mesh_homogenization, & + mesh_homogenizationAt, & mesh_NipsPerElem, & mesh_maxNips, & mesh_NcpElems, & @@ -481,7 +481,7 @@ subroutine material_init() allocate(CrystallitePosition (size(config_phase)), source=0_pInt) ElemLoop:do e = 1_pInt,mesh_NcpElems - myHomog = mesh_homogenization(e) + myHomog = mesh_homogenizationAt(e) IPloop:do i = 1_pInt, mesh_NipsPerElem HomogenizationPosition(myHomog) = HomogenizationPosition(myHomog) + 1_pInt mappingHomogenization(1:2,i,e) = [HomogenizationPosition(myHomog),myHomog] @@ -521,7 +521,7 @@ subroutine material_parseHomogenization use config, only : & config_homogenization use mesh, only: & - mesh_homogenization + mesh_homogenizationAt use IO, only: & IO_error @@ -551,7 +551,7 @@ subroutine material_parseHomogenization allocate(hydrogenflux_initialCh(size(config_homogenization)), source=0.0_pReal) forall (h = 1_pInt:size(config_homogenization)) & - homogenization_active(h) = any(mesh_homogenization == h) + homogenization_active(h) = any(mesh_homogenizationAt == h) do h=1_pInt, size(config_homogenization) @@ -687,7 +687,7 @@ subroutine material_parseMicrostructure config_microstructure, & microstructure_name use mesh, only: & - mesh_microstructure, & + mesh_microstructureAt, & mesh_NcpElems implicit none @@ -703,11 +703,11 @@ subroutine material_parseMicrostructure allocate(microstructure_active(size(config_microstructure)), source=.false.) allocate(microstructure_elemhomo(size(config_microstructure)), source=.false.) - if(any(mesh_microstructure > size(config_microstructure))) & + if(any(mesh_microstructureAt > size(config_microstructure))) & call IO_error(155_pInt,ext_msg='More microstructures in geometry than sections in material.config') forall (e = 1_pInt:mesh_NcpElems) & - microstructure_active(mesh_microstructure(e)) = .true. ! current microstructure used in model? Elementwise view, maximum N operations for N elements + microstructure_active(mesh_microstructureAt(e)) = .true. ! current microstructure used in model? Elementwise view, maximum N operations for N elements do m=1_pInt, size(config_microstructure) microstructure_Nconstituents(m) = config_microstructure(m)%countKeys('(constituent)') @@ -1087,8 +1087,8 @@ subroutine material_populateGrains use mesh, only: & mesh_NipsPerElem, & mesh_elemType, & - mesh_homogenization, & - mesh_microstructure, & + mesh_homogenizationAt, & + mesh_microstructureAt, & mesh_maxNips, & mesh_NcpElems, & mesh_ipVolume, & @@ -1141,14 +1141,14 @@ subroutine material_populateGrains ! populating homogenization schemes in each !-------------------------------------------------------------------------------------------------- do e = 1_pInt, mesh_NcpElems - material_homog(1_pInt:mesh_NipsPerElem,e) = mesh_homogenization(e) + material_homog(1_pInt:mesh_NipsPerElem,e) = mesh_homogenizationAt(e) enddo !-------------------------------------------------------------------------------------------------- ! precounting of elements for each homog/micro pair do e = 1_pInt, mesh_NcpElems - homog = mesh_homogenization(e) - micro = mesh_microstructure(e) + homog = mesh_homogenizationAt(e) + micro = mesh_microstructureAt(e) Nelems(homog,micro) = Nelems(homog,micro) + 1_pInt enddo allocate(elemsOfHomogMicro(size(config_homogenization),size(config_microstructure))) @@ -1166,8 +1166,8 @@ subroutine material_populateGrains Nelems = 0_pInt ! reuse as counter elementLooping: do e = 1_pInt,mesh_NcpElems t = mesh_elemType - homog = mesh_homogenization(e) - micro = mesh_microstructure(e) + homog = mesh_homogenizationAt(e) + micro = mesh_microstructureAt(e) if (homog < 1_pInt .or. homog > size(config_homogenization)) & ! out of bounds call IO_error(154_pInt,e,0_pInt,0_pInt) if (micro < 1_pInt .or. micro > size(config_microstructure)) & ! out of bounds diff --git a/src/mesh.f90 b/src/mesh.f90 index 18f986b33..b2c63004b 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -30,11 +30,11 @@ module mesh !!!! BEGIN DEPRECATED !!!!! integer(pInt), dimension(:), allocatable, public, protected :: & - mesh_homogenization, & !< homogenization ID of each element - mesh_microstructure !< homogenization ID of each element + mesh_homogenizationAt, & !< homogenization ID of each element + mesh_microstructureAt !< microstructure ID of each element integer(pInt), dimension(:,:), allocatable, public, protected :: & - mesh_CPnodeID, & + mesh_CPnodeID, & !< nodes forming an element mesh_element, & !DEPRECATED mesh_sharedElem, & !< entryCount and list of elements containing node mesh_nodeTwins !< node twins are surface nodes that lie exactly on opposite sides of the mesh (surfaces nodes with equal coordinate values in two dimensions) @@ -646,8 +646,8 @@ subroutine mesh_init(ip,el) mesh_NipsPerElem = mesh_maxNips mesh_NcellnodesPerElem = mesh_maxNcellnodes ! better name - mesh_homogenization = mesh_element(3,:) - mesh_microstructure = mesh_element(4,:) + mesh_homogenizationAt = mesh_element(3,:) + mesh_microstructureAt = mesh_element(4,:) mesh_CPnodeID = mesh_element(5:4+mesh_NipsPerElem,:) !!!!!!!!!!!!!!!!!!!!!!!! diff --git a/src/meshFEM.f90 b/src/meshFEM.f90 index c58cbf2e3..1362063f8 100644 --- a/src/meshFEM.f90 +++ b/src/meshFEM.f90 @@ -35,13 +35,10 @@ use PETScis !!!! BEGIN DEPRECATED !!!!! integer(pInt), dimension(:), allocatable, public, protected :: & - mesh_homogenization, & !< homogenization ID of each element - mesh_microstructure !< homogenization ID of each element - - real(pReal), public, protected :: charLength + mesh_homogenizationAt, & !< homogenization ID of each element + mesh_microstructureAt !< microstructure ID of each element integer(pInt), dimension(:,:), allocatable, public, protected :: & - !mesh_CPnodeID, & mesh_element !DEPRECATED real(pReal), dimension(:,:), allocatable, public :: & @@ -241,10 +238,10 @@ subroutine mesh_init() !!!! COMPATIBILITY HACK !!!! ! for a homogeneous mesh, all elements have the same number of IPs and and cell nodes. ! hence, xxPerElem instead of maxXX - mesh_NipsPerElem = mesh_maxNips + mesh_NipsPerElem = mesh_maxNips ! better name - mesh_homogenization = mesh_element(3,:) - mesh_microstructure = mesh_element(4,:) + mesh_homogenizationAt = mesh_element(3,:) + mesh_microstructureAt = mesh_element(4,:) !!!!!!!!!!!!!!!!!!!!!!!! end subroutine mesh_init From dc289a278b1d6f0617750094fbc1d00bcf9de850 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 4 Oct 2018 06:39:03 +0200 Subject: [PATCH 23/24] clearer code homogenization is defined per element, not per IP hence, use material_homogenizationAt instead of deprecated material_homog which pretends a dependency on the integration point --- src/CPFEM.f90 | 6 +++--- src/constitutive.f90 | 21 ++++++++++----------- src/material.f90 | 5 ++++- 3 files changed, 17 insertions(+), 15 deletions(-) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index dc6ba07dc..674a557b5 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -314,7 +314,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt thermal_type, & THERMAL_conduction_ID, & phase_Nsources, & - material_homog + material_homogenizationAt use config, only: & material_Nhomogenization use crystallite, only: & @@ -503,7 +503,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt if (.not. parallelExecution) then chosenThermal1: select case (thermal_type(mesh_element(3,elCP))) case (THERMAL_conduction_ID) chosenThermal1 - temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = & + temperature(material_homogenizationAt(elCP))%p(thermalMapping(material_homogenizationAt(elCP))%p(ip,elCP)) = & temperature_inp end select chosenThermal1 materialpoint_F0(1:3,1:3,ip,elCP) = ffn @@ -516,7 +516,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt CPFEM_dcsde(1:6,1:6,ip,elCP) = CPFEM_odd_jacobian * math_identity2nd(6) chosenThermal2: select case (thermal_type(mesh_element(3,elCP))) case (THERMAL_conduction_ID) chosenThermal2 - temperature(material_homog(ip,elCP))%p(thermalMapping(material_homog(ip,elCP))%p(ip,elCP)) = & + temperature(material_homogenizationAt(elCP))%p(thermalMapping(material_homogenizationAt(elCP))%p(ip,elCP)) = & temperature_inp end select chosenThermal2 materialpoint_F0(1:3,1:3,ip,elCP) = ffn diff --git a/src/constitutive.f90 b/src/constitutive.f90 index eedc3e509..85c5c2da2 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -386,7 +386,7 @@ subroutine constitutive_microstructure(orientations, Fe, Fp, ipc, ip, el) use material, only: & phase_plasticity, & material_phase, & - material_homog, & + material_homogenizationAt, & temperature, & thermalMapping, & PLASTICITY_dislotwin_ID, & @@ -413,7 +413,7 @@ subroutine constitutive_microstructure(orientations, Fe, Fp, ipc, ip, el) real(pReal), intent(in), dimension(:,:,:,:) :: & orientations !< crystal orientations as quaternions - ho = material_homog(ip,el) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) @@ -444,7 +444,7 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, S6, Fi, ipc, ip, e phase_plasticity, & phase_plasticityInstance, & material_phase, & - material_homog, & + material_homogenizationAt, & temperature, & thermalMapping, & PLASTICITY_NONE_ID, & @@ -494,7 +494,7 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, S6, Fi, ipc, ip, e integer(pInt) :: & i, j, instance, of - ho = material_homog(ip,el) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) S = math_Mandel6to33(S6) @@ -752,7 +752,7 @@ subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, Fe, Fi, ipc, ip math_I3 use material, only: & material_phase, & - material_homog, & + material_homogenizationAt, & phase_NstiffnessDegradations, & phase_stiffnessDegradation, & damage, & @@ -783,8 +783,7 @@ subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, Fe, Fi, ipc, ip integer(pInt) :: & i, j - ho = material_homog(ip,el) - + ho = material_homogenizationAt(el) C = math_Mandel66to3333(constitutive_homogenizedC(ipc,ip,el)) DegradationLoop: do d = 1_pInt, phase_NstiffnessDegradations(material_phase(ipc,ip,el)) @@ -835,7 +834,7 @@ subroutine constitutive_collectDotState(S6, FeArray, Fi, FpArray, subdt, subfrac phase_source, & phase_Nsources, & material_phase, & - material_homog, & + material_homogenizationAt, & temperature, & thermalMapping, & homogenization_maxNgrains, & @@ -895,7 +894,7 @@ subroutine constitutive_collectDotState(S6, FeArray, Fi, FpArray, subdt, subfrac s, & !< counter in source loop instance, of - ho = material_homog( ip,el) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) Mp = math_mul33x33(math_mul33x33(transpose(Fi),Fi),math_Mandel6to33(S6)) @@ -1054,7 +1053,7 @@ function constitutive_postResults(S6, Fi, FeArray, ipc, ip, el) phase_source, & phase_Nsources, & material_phase, & - material_homog, & + material_homogenizationAt, & temperature, & thermalMapping, & homogenization_maxNgrains, & @@ -1117,7 +1116,7 @@ function constitutive_postResults(S6, Fi, FeArray, ipc, ip, el) Mp = math_mul33x33(math_mul33x33(transpose(Fi),Fi),math_Mandel6to33(S6)) - ho = material_homog( ip,el) + ho = material_homogenizationAt(el) tme = thermalMapping(ho)%p(ip,el) startPos = 1_pInt diff --git a/src/material.f90 b/src/material.f90 index 98fa1220f..4b4b9e31a 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -169,6 +169,7 @@ module material homogenization_maxNgrains !< max number of grains in any USED homogenization integer(pInt), dimension(:), allocatable, public, protected :: & + material_homogenizationAt, & !< homogenization ID of each element (copy of mesh_homogenizationAt) phase_Nsources, & !< number of source mechanisms active in each phase phase_Nkinematics, & !< number of kinematic mechanisms active in each phase phase_NstiffnessDegradations, & !< number of stiffness degradation mechanisms active in each phase @@ -199,9 +200,10 @@ module material integer(pInt), dimension(:,:,:), allocatable, public :: & material_phase !< phase (index) of each grain,IP,element -! DEPRECATED. DID WE EVER ALLOWED DIFFERENT HOMOGENIZATION SCHEMES WITHIN ONE ELEMENT? +! BEGIN DEPRECATED: use material_homogenizationAt integer(pInt), dimension(:,:), allocatable, public :: & material_homog !< homogenization (index) of each IP,element +! END DEPRECATED type(tPlasticState), allocatable, dimension(:), public :: & plasticState type(tSourceState), allocatable, dimension(:), public :: & @@ -1132,6 +1134,7 @@ subroutine material_populateGrains allocate(material_volume(homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source=0.0_pReal) allocate(material_phase(homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source=0_pInt) allocate(material_homog(mesh_maxNips,mesh_NcpElems), source=0_pInt) + allocate(material_homogenizationAt,source=mesh_homogenizationAt) allocate(material_texture(homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), source=0_pInt) allocate(material_EulerAngles(3,homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems),source=0.0_pReal) From f0b9c0caf70e5904ec0625d410dc24f57e70bc4c Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Tue, 9 Oct 2018 17:52:54 -0400 Subject: [PATCH 24/24] polishing and encapsulating of Abaqus-specific local variables --- src/mesh.f90 | 28 ++++++++++++++++++---------- 1 file changed, 18 insertions(+), 10 deletions(-) diff --git a/src/mesh.f90 b/src/mesh.f90 index b2c63004b..e55165d51 100644 --- a/src/mesh.f90 +++ b/src/mesh.f90 @@ -14,8 +14,6 @@ module mesh integer(pInt), public, protected :: & mesh_NcpElems, & !< total number of CP elements in local mesh mesh_elemType, & !< Element type of the mesh (only support homogeneous meshes) - mesh_maxNelemInSet, & - mesh_Nmaterials, & mesh_Nnodes, & !< total number of nodes in mesh mesh_Ncellnodes, & !< total number of cell nodes in mesh (including duplicates) mesh_Ncells, & !< total number of cells in mesh @@ -64,6 +62,12 @@ module mesh logical, dimension(3), public, protected :: mesh_periodicSurface !< flag indicating periodic outer surfaces (used for fluxes) +#if defined(Marc4DAMASK) || defined(Abaqus) + integer(pInt), private :: & + mesh_maxNelemInSet, & + mesh_Nmaterials +#endif + integer(pInt), dimension(2), private :: & mesh_maxValStateVar = 0_pInt @@ -520,8 +524,12 @@ subroutine mesh_init(ip,el) if(worldsize>grid(3)) call IO_error(894_pInt, ext_msg='number of processes exceeds grid(3)') geomSize = mesh_spectral_getSize(fileUnit) - devNull = fftw_mpi_local_size_3d(int(grid(3),C_INTPTR_T),int(grid(2),C_INTPTR_T),& - int(grid(1),C_INTPTR_T)/2+1,PETSC_COMM_WORLD,local_K,local_K_offset) + devNull = fftw_mpi_local_size_3d(int(grid(3),C_INTPTR_T), & + int(grid(2),C_INTPTR_T), & + int(grid(1),C_INTPTR_T)/2+1, & + PETSC_COMM_WORLD, & + local_K, & ! domain grid size along z + local_K_offset) ! domain grid offset along z grid3 = int(local_K,pInt) grid3Offset = int(local_K_offset,pInt) size3 = geomSize(3)*real(grid3,pReal) /real(grid(3),pReal) @@ -1251,7 +1259,7 @@ subroutine mesh_spectral_build_elements(fileUnit) integer(pInt) :: & e, i, & headerLength = 0_pInt, & - maxIntCount, & + maxDataPerLine, & homog, & elemType, & elemOffset @@ -1287,16 +1295,16 @@ subroutine mesh_spectral_build_elements(fileUnit) read(fileUnit,'(a65536)') line enddo - maxIntCount = 0_pInt + maxDataPerLine = 0_pInt i = 1_pInt do while (i > 0_pInt) i = IO_countContinuousIntValues(fileUnit) - maxIntCount = max(maxIntCount, i) + maxDataPerLine = max(maxDataPerLine, i) ! found a longer line? enddo allocate(mesh_element (4_pInt+8_pInt,mesh_NcpElems), source = 0_pInt) - allocate(microstructures (1_pInt+maxIntCount), source = 1_pInt) - allocate(microGlobal(mesh_NcpElemsGlobal), source = 1_pInt) + allocate(microstructures (1_pInt+maxDataPerLine), source = 1_pInt) ! prepare to receive counter and max data size + allocate(microGlobal (mesh_NcpElemsGlobal), source = 1_pInt) !-------------------------------------------------------------------------------------------------- ! read in microstructures @@ -1307,7 +1315,7 @@ subroutine mesh_spectral_build_elements(fileUnit) e = 0_pInt do while (e < mesh_NcpElemsGlobal .and. microstructures(1) > 0_pInt) ! fill expected number of elements, stop at end of data (or blank line!) - microstructures = IO_continuousIntValues(fileUnit,maxIntCount,dummyName,dummySet,0_pInt) ! get affected elements + microstructures = IO_continuousIntValues(fileUnit,maxDataPerLine,dummyName,dummySet,0_pInt) ! get affected elements do i = 1_pInt,microstructures(1_pInt) e = e+1_pInt ! valid element entry microGlobal(e) = microstructures(1_pInt+i)