Merge branch 'development' into polishing-for-release

This commit is contained in:
Martin Diehl 2021-01-14 20:40:34 +01:00
commit 0e3462f3e0
3 changed files with 15 additions and 33 deletions

@ -1 +1 @@
Subproject commit 1d0c95c5c1c0e7e6f57bdfc94b695e47a6ad6c60 Subproject commit b1a31a79cc90d458494068a96cfd3e9497aa330c

View File

@ -1 +1 @@
v3.0.0-alpha2-258-g715504ee5 v3.0.0-alpha2-279-g8182c9c54

View File

@ -83,6 +83,7 @@ subroutine discretization_mesh_init(restart)
num_mesh num_mesh
integer :: integrationOrder !< order of quadrature rule required integer :: integrationOrder !< order of quadrature rule required
print'(/,a)', ' <<<+- discretization_mesh init -+>>>' print'(/,a)', ' <<<+- discretization_mesh init -+>>>'
!-------------------------------------------------------------------------------- !--------------------------------------------------------------------------------
@ -95,13 +96,15 @@ subroutine discretization_mesh_init(restart)
debug_element = config_debug%get_asInt('element',defaultVal=1) debug_element = config_debug%get_asInt('element',defaultVal=1)
debug_ip = config_debug%get_asInt('integrationpoint',defaultVal=1) debug_ip = config_debug%get_asInt('integrationpoint',defaultVal=1)
call DMPlexCreateFromFile(PETSC_COMM_WORLD,interface_geomFile,PETSC_TRUE,globalMesh,ierr) call DMPlexCreateFromFile(PETSC_COMM_WORLD,interface_geomFile,PETSC_TRUE,globalMesh,ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
call DMGetDimension(globalMesh,dimPlex,ierr) call DMGetDimension(globalMesh,dimPlex,ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
call DMGetStratumSize(globalMesh,'depth',dimPlex,mesh_NcpElemsGlobal,ierr) call DMGetStratumSize(globalMesh,'depth',dimPlex,mesh_NcpElemsGlobal,ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
call DMView(globalMesh, PETSC_VIEWER_STDOUT_WORLD,ierr)
CHKERRQ(ierr)
! get number of IDs in face sets (for boundary conditions?) ! get number of IDs in face sets (for boundary conditions?)
call DMGetLabelSize(globalMesh,'Face Sets',mesh_Nboundaries,ierr) call DMGetLabelSize(globalMesh,'Face Sets',mesh_Nboundaries,ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
@ -109,6 +112,13 @@ subroutine discretization_mesh_init(restart)
call MPI_Bcast(mesh_NcpElemsGlobal,1,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr) call MPI_Bcast(mesh_NcpElemsGlobal,1,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr)
call MPI_Bcast(dimPlex,1,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr) call MPI_Bcast(dimPlex,1,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr)
if (worldrank == 0) then
call DMClone(globalMesh,geomMesh,ierr)
else
call DMPlexDistribute(globalMesh,0,sf,geomMesh,ierr)
endif
CHKERRQ(ierr)
allocate(mesh_boundaries(mesh_Nboundaries), source = 0) allocate(mesh_boundaries(mesh_Nboundaries), source = 0)
call DMGetLabelSize(globalMesh,'Face Sets',nFaceSets,ierr) call DMGetLabelSize(globalMesh,'Face Sets',nFaceSets,ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
@ -123,35 +133,6 @@ subroutine discretization_mesh_init(restart)
endif endif
call MPI_Bcast(mesh_boundaries,mesh_Nboundaries,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr) call MPI_Bcast(mesh_boundaries,mesh_Nboundaries,MPI_INTEGER,0,PETSC_COMM_WORLD,ierr)
if (worldrank == 0) then
fileContent = IO_readlines(interface_geomFile)
l = 0
do
l = l + 1
if (IO_isBlank(fileContent(l))) cycle ! need also to ignore C and C++ style comments?
if (trim(fileContent(l)) == '$Elements') then
j = 0
l = l + 1
do
l = l + 1
if (trim(fileContent(l)) == '$EndElements') exit
chunkPos = IO_stringPos(fileContent(l))
if (chunkPos(1) == 3+IO_intValue(fileContent(l),chunkPos,3)+dimPlex+1) then
call DMSetLabelValue(globalMesh,'material',j,IO_intValue(fileContent(l),chunkPos,4),ierr)
CHKERRQ(ierr)
j = j + 1
endif
enddo
exit
endif
enddo
call DMClone(globalMesh,geomMesh,ierr)
CHKERRQ(ierr)
else
call DMPlexDistribute(globalMesh,0,sf,geomMesh,ierr)
CHKERRQ(ierr)
endif
call DMDestroy(globalMesh,ierr); CHKERRQ(ierr) call DMDestroy(globalMesh,ierr); CHKERRQ(ierr)
call DMGetStratumSize(geomMesh,'depth',dimPlex,mesh_NcpElems,ierr) call DMGetStratumSize(geomMesh,'depth',dimPlex,mesh_NcpElems,ierr)
@ -166,9 +147,10 @@ subroutine discretization_mesh_init(restart)
allocate(materialAt(mesh_NcpElems)) allocate(materialAt(mesh_NcpElems))
do j = 1, mesh_NcpElems do j = 1, mesh_NcpElems
call DMGetLabelValue(geomMesh,'material',j-1,materialAt(j),ierr) call DMGetLabelValue(geomMesh,'Cell Sets',j-1,materialAt(j),ierr)
CHKERRQ(ierr) CHKERRQ(ierr)
end do end do
materialAt = materialAt + 1
if (debug_element < 1 .or. debug_element > mesh_NcpElems) call IO_error(602,ext_msg='element') if (debug_element < 1 .or. debug_element > mesh_NcpElems) call IO_error(602,ext_msg='element')
if (debug_ip < 1 .or. debug_ip > mesh_maxNips) call IO_error(602,ext_msg='IP') if (debug_ip < 1 .or. debug_ip > mesh_maxNips) call IO_error(602,ext_msg='IP')