From 93562d5142532cf5c410cad183a6ff2447866fb7 Mon Sep 17 00:00:00 2001
From: Martin Diehl <m.diehl@mpie.de>
Date: Mon, 20 Aug 2018 16:42:47 +0200
Subject: [PATCH] mapping of elements etc not needed for PETSc-based FEM and
 spectral solvers

---
 src/mesh.f90    | 33 ++++++++++-----------------------
 src/meshFEM.f90 | 14 +++-----------
 2 files changed, 13 insertions(+), 34 deletions(-)

diff --git a/src/mesh.f90 b/src/mesh.f90
index 5606b656b..4e72ba73e 100644
--- a/src/mesh.f90
+++ b/src/mesh.f90
@@ -95,9 +95,11 @@ module mesh
  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
@@ -402,7 +404,9 @@ module mesh
 
  public :: &
    mesh_init, &
+#if defined(Marc4DAMASK) || defined(Abaqus)
    mesh_FEasCP, &
+#endif
    mesh_build_cellnodes, &
    mesh_build_ipVolumes, &
    mesh_build_ipCoordinates, &
@@ -420,7 +424,6 @@ module mesh
 #ifdef Spectral
    mesh_spectral_getHomogenization, &
    mesh_spectral_count, &
-   mesh_spectral_mapNodesAndElems, &
    mesh_spectral_count_cpSizes, &
    mesh_spectral_build_nodes, &
    mesh_spectral_build_elements, &
@@ -552,8 +555,6 @@ subroutine mesh_init(ip,el)
  if (myDebug) write(6,'(a)') ' Grid partitioned'; flush(6)
  call mesh_spectral_count()
  if (myDebug) write(6,'(a)') ' Counted nodes/elements'; flush(6)
- call mesh_spectral_mapNodesAndElems
- if (myDebug) write(6,'(a)') ' Mapped nodes and elements'; flush(6)
  call mesh_spectral_count_cpSizes
  if (myDebug) write(6,'(a)') ' Built CP statistics'; flush(6)
  call mesh_spectral_build_nodes()
@@ -659,12 +660,16 @@ subroutine mesh_init(ip,el)
 
  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
 
 end subroutine mesh_init
 
 
+#if defined(Marc4DAMASK) || defined(Abaqus)
 !--------------------------------------------------------------------------------------------------
 !> @brief Gives the FE to CP ID mapping by binary search through lookup array
 !! valid questions (what) are 'elem', 'node'
@@ -713,7 +718,7 @@ integer(pInt) function mesh_FEasCP(what,myID)
  enddo binarySearch
 
 end function mesh_FEasCP
-
+#endif
 
 !--------------------------------------------------------------------------------------------------
 !> @brief Split CP elements into cells.
@@ -1188,24 +1193,6 @@ subroutine mesh_spectral_count()
 end subroutine mesh_spectral_count
 
 
-!--------------------------------------------------------------------------------------------------
-!> @brief fake map node from FE ID to internal (consecutive) representation for node and element
-!! Allocates global array 'mesh_mapFEtoCPnode' and 'mesh_mapFEtoCPelem'
-!--------------------------------------------------------------------------------------------------
-subroutine mesh_spectral_mapNodesAndElems
- use math, only: &
-   math_range
-
- implicit none
- allocate (mesh_mapFEtoCPnode(2_pInt,mesh_Nnodes), source = 0_pInt)
- allocate (mesh_mapFEtoCPelem(2_pInt,mesh_NcpElems), source = 0_pInt)
-
- mesh_mapFEtoCPnode = spread(math_range(mesh_Nnodes),1,2)
- mesh_mapFEtoCPelem = spread(math_range(mesh_NcpElems),1,2)
-
-end subroutine mesh_spectral_mapNodesAndElems
-
-
 !--------------------------------------------------------------------------------------------------
 !> @brief Gets maximum count of nodes, IPs, IP neighbors, and subNodes among cpElements.
 !! Sets global values 'mesh_maxNnodes', 'mesh_maxNips', 'mesh_maxNipNeighbors',
diff --git a/src/meshFEM.f90 b/src/meshFEM.f90
index 141b1b0a9..7d79dd46d 100644
--- a/src/meshFEM.f90
+++ b/src/meshFEM.f90
@@ -154,12 +154,6 @@ subroutine mesh_init(ip,el)
  write(6,'(a15,a)') ' Current time: ',IO_timeStamp()
 #include "compilation_info.f90"
 
- if (allocated(mesh_node0))                   deallocate(mesh_node0)
- if (allocated(mesh_node))                    deallocate(mesh_node)
- if (allocated(mesh_element))                 deallocate(mesh_element)
- if (allocated(mesh_ipCoordinates))           deallocate(mesh_ipCoordinates)
- if (allocated(mesh_ipVolume))                deallocate(mesh_ipVolume)
-
  call DMPlexCreateFromFile(PETSC_COMM_WORLD,geometryFile,PETSC_TRUE,globalMesh,ierr)
  CHKERRQ(ierr)
  call DMGetDimension(globalMesh,dimPlex,ierr)
@@ -334,11 +328,9 @@ subroutine mesh_FEM_build_ipCoordinates(dimPlex,qPoints)
  PetscInt                  :: cellStart, cellEnd, cell, qPt, dirI, dirJ, qOffset
  PetscErrorCode            :: ierr
 
- if (.not. allocated(mesh_ipCoordinates)) then
-   allocate(mesh_ipCoordinates(3,mesh_maxNips,mesh_NcpElems))
-   mesh_ipCoordinates = 0.0_pReal
- endif
- 
+
+ allocate(mesh_ipCoordinates(3,mesh_maxNips,mesh_NcpElems),source=0.0_pReal)
+
  pV0 => v0
  pCellJ => cellJ
  pInvcellJ => invcellJ