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

only mechanics at the moment 

will be extended, but most likely differently
This commit is contained in:
Martin Diehl 2021-07-22 10:12:43 +02:00
parent eb834b635d
commit 5b66db8a39
2 changed files with 41 additions and 80 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

113
src/mesh/DAMASK_mesh.f90
View File

@ -21,7 +21,6 @@ program DAMASK_mesh
use mesh_mechanical_FEM use mesh_mechanical_FEM
implicit none implicit none
integer :: nActiveFields = 0
type :: tLoadCase type :: tLoadCase
real(pReal) :: time = 0.0_pReal !< length of increment real(pReal) :: time = 0.0_pReal !< length of increment
@ -78,7 +77,7 @@ program DAMASK_mesh
type(tLoadCase), allocatable, dimension(:) :: loadCases !< array of all load cases type(tLoadCase), allocatable, dimension(:) :: loadCases !< array of all load cases
type(tSolutionState), allocatable, dimension(:) :: solres type(tSolutionState), allocatable, dimension(:) :: solres
PetscInt :: faceSet, currentFaceSet, field, dimPlex PetscInt :: faceSet, currentFaceSet, dimPlex
PetscErrorCode :: ierr PetscErrorCode :: ierr
integer(kind(COMPONENT_UNDEFINED_ID)) :: ID integer(kind(COMPONENT_UNDEFINED_ID)) :: ID
external :: & external :: &
@ -101,8 +100,7 @@ program DAMASK_mesh
! reading basic information from load case file and allocate data structure containing load cases ! reading basic information from load case file and allocate data structure containing load cases
call DMGetDimension(geomMesh,dimPlex,ierr) !< dimension of mesh (2D or 3D) call DMGetDimension(geomMesh,dimPlex,ierr) !< dimension of mesh (2D or 3D)
CHKERRA(ierr) CHKERRA(ierr)
nActiveFields = 1 allocate(solres(1))
allocate(solres(nActiveFields))
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! reading basic information from load case file and allocate data structure containing load cases ! reading basic information from load case file and allocate data structure containing load cases
@ -124,32 +122,26 @@ program DAMASK_mesh
allocate(loadCases(N_def)) allocate(loadCases(N_def))
do i = 1, size(loadCases) do i = 1, size(loadCases)
allocate(loadCases(i)%fieldBC(nActiveFields)) allocate(loadCases(i)%fieldBC(1))
field = 1 loadCases(i)%fieldBC(1)%ID = FIELD_MECH_ID
loadCases(i)%fieldBC(field)%ID = FIELD_MECH_ID
enddo enddo
do i = 1, size(loadCases) do i = 1, size(loadCases)
do field = 1, nActiveFields loadCases(i)%fieldBC(1)%nComponents = dimPlex !< X, Y (, Z) displacements
select case (loadCases(i)%fieldBC(field)%ID) allocate(loadCases(i)%fieldBC(1)%componentBC(loadCases(i)%fieldBC(1)%nComponents))
case(FIELD_MECH_ID) do component = 1, loadCases(i)%fieldBC(1)%nComponents
loadCases(i)%fieldBC(field)%nComponents = dimPlex !< X, Y (, Z) displacements select case (component)
allocate(loadCases(i)%fieldBC(field)%componentBC(loadCases(i)%fieldBC(field)%nComponents)) case (1)
do component = 1, loadCases(i)%fieldBC(field)%nComponents loadCases(i)%fieldBC(1)%componentBC(component)%ID = COMPONENT_MECH_X_ID
select case (component) case (2)
case (1) loadCases(i)%fieldBC(1)%componentBC(component)%ID = COMPONENT_MECH_Y_ID
loadCases(i)%fieldBC(field)%componentBC(component)%ID = COMPONENT_MECH_X_ID case (3)
case (2) loadCases(i)%fieldBC(1)%componentBC(component)%ID = COMPONENT_MECH_Z_ID
loadCases(i)%fieldBC(field)%componentBC(component)%ID = COMPONENT_MECH_Y_ID
case (3)
loadCases(i)%fieldBC(field)%componentBC(component)%ID = COMPONENT_MECH_Z_ID
end select
enddo
end select end select
do component = 1, loadCases(i)%fieldBC(field)%nComponents enddo
allocate(loadCases(i)%fieldBC(field)%componentBC(component)%Value(mesh_Nboundaries), source = 0.0_pReal) do component = 1, loadCases(i)%fieldBC(1)%nComponents
allocate(loadCases(i)%fieldBC(field)%componentBC(component)%Mask (mesh_Nboundaries), source = .false.) allocate(loadCases(i)%fieldBC(1)%componentBC(component)%Value(mesh_Nboundaries), source = 0.0_pReal)
enddo allocate(loadCases(i)%fieldBC(1)%componentBC(component)%Mask (mesh_Nboundaries), source = .false.)
enddo enddo
enddo enddo
@ -194,16 +186,12 @@ program DAMASK_mesh
ID = COMPONENT_MECH_Z_ID ID = COMPONENT_MECH_Z_ID
end select end select
do field = 1, nActiveFields do component = 1, loadcases(currentLoadCase)%fieldBC(1)%nComponents
if (loadCases(currentLoadCase)%fieldBC(field)%ID == FIELD_MECH_ID) then if (loadCases(currentLoadCase)%fieldBC(1)%componentBC(component)%ID == ID) then
do component = 1, loadcases(currentLoadCase)%fieldBC(field)%nComponents loadCases(currentLoadCase)%fieldBC(1)%componentBC(component)%Mask (currentFaceSet) = &
if (loadCases(currentLoadCase)%fieldBC(field)%componentBC(component)%ID == ID) then .true.
loadCases(currentLoadCase)%fieldBC(field)%componentBC(component)%Mask (currentFaceSet) = & loadCases(currentLoadCase)%fieldBC(1)%componentBC(component)%Value(currentFaceSet) = &
.true. IO_floatValue(line,chunkPos,i+1)
loadCases(currentLoadCase)%fieldBC(field)%componentBC(component)%Value(currentFaceSet) = &
IO_floatValue(line,chunkPos,i+1)
endif
enddo
endif endif
enddo enddo
end select end select
@ -219,21 +207,16 @@ program DAMASK_mesh
print'(a,i0)', ' load case: ', currentLoadCase print'(a,i0)', ' load case: ', currentLoadCase
if (.not. loadCases(currentLoadCase)%followFormerTrajectory) & if (.not. loadCases(currentLoadCase)%followFormerTrajectory) &
print'(a)', ' drop guessing along trajectory' print'(a)', ' drop guessing along trajectory'
do field = 1, nActiveFields print'(a)', ' Field '//trim(FIELD_MECH_label)
select case (loadCases(currentLoadCase)%fieldBC(field)%ID)
case(FIELD_MECH_ID)
print'(a)', ' Field '//trim(FIELD_MECH_label)
end select do faceSet = 1, mesh_Nboundaries
do faceSet = 1, mesh_Nboundaries do component = 1, loadCases(currentLoadCase)%fieldBC(1)%nComponents
do component = 1, loadCases(currentLoadCase)%fieldBC(field)%nComponents if (loadCases(currentLoadCase)%fieldBC(1)%componentBC(component)%Mask(faceSet)) &
if (loadCases(currentLoadCase)%fieldBC(field)%componentBC(component)%Mask(faceSet)) & print'(a,i2,a,i2,a,f12.7)', ' Face ', mesh_boundaries(faceSet), &
print'(a,i2,a,i2,a,f12.7)', ' Face ', mesh_boundaries(faceSet), & ' Component ', component, &
' Component ', component, & ' Value ', loadCases(currentLoadCase)%fieldBC(1)% &
' Value ', loadCases(currentLoadCase)%fieldBC(field)% & componentBC(component)%Value(faceSet)
componentBC(component)%Value(faceSet) enddo
enddo
enddo
enddo enddo
print'(a,f12.6)', ' time: ', loadCases(currentLoadCase)%time print'(a,f12.6)', ' time: ', loadCases(currentLoadCase)%time
if (loadCases(currentLoadCase)%incs < 1) errorID = 835 ! non-positive incs count if (loadCases(currentLoadCase)%incs < 1) errorID = 835 ! non-positive incs count
@ -247,12 +230,7 @@ program DAMASK_mesh
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! doing initialization depending on active solvers ! doing initialization depending on active solvers
call FEM_Utilities_init call FEM_Utilities_init
do field = 1, nActiveFields call FEM_mechanical_init(loadCases(1)%fieldBC(1))
select case (loadCases(1)%fieldBC(field)%ID)
case(FIELD_MECH_ID)
call FEM_mechanical_init(loadCases(1)%fieldBC(field))
end select
enddo
if (worldrank == 0) then if (worldrank == 0) then
open(newunit=statUnit,file=trim(getSolverJobName())//'.sta',form='FORMATTED',status='REPLACE') open(newunit=statUnit,file=trim(getSolverJobName())//'.sta',form='FORMATTED',status='REPLACE')
@ -295,33 +273,16 @@ program DAMASK_mesh
'-',stepFraction, '/', subStepFactor**cutBackLevel '-',stepFraction, '/', subStepFactor**cutBackLevel
flush(IO_STDOUT) flush(IO_STDOUT)
!-------------------------------------------------------------------------------------------------- call FEM_mechanical_forward(guess,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(1))
! forward fields
do field = 1, nActiveFields
select case (loadCases(currentLoadCase)%fieldBC(field)%ID)
case(FIELD_MECH_ID)
call FEM_mechanical_forward (&
guess,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(field))
end select
enddo
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! solve fields ! solve fields
stagIter = 0 stagIter = 0
stagIterate = .true. stagIterate = .true.
do while (stagIterate) do while (stagIterate)
do field = 1, nActiveFields solres(1) = FEM_mechanical_solution(incInfo,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(1))
select case (loadCases(currentLoadCase)%fieldBC(field)%ID) if(.not. solres(1)%converged) exit
case(FIELD_MECH_ID)
solres(field) = FEM_mechanical_solution (&
incInfo,timeinc,timeIncOld,loadCases(currentLoadCase)%fieldBC(field))
end select
if(.not. solres(field)%converged) exit ! no solution found
enddo
stagIter = stagIter + 1 stagIter = stagIter + 1
stagIterate = stagIter < stagItMax & stagIterate = stagIter < stagItMax &
.and. all(solres(:)%converged) & .and. all(solres(:)%converged) &

8
src/mesh/FEM_utilities.f90
View File

@ -53,15 +53,15 @@ module FEM_utilities
end type tSolutionState end type tSolutionState
type, public :: tComponentBC type, public :: tComponentBC
integer(kind(COMPONENT_UNDEFINED_ID)) :: ID integer(kind(COMPONENT_UNDEFINED_ID)) :: ID
real(pReal), allocatable, dimension(:) :: Value real(pReal), allocatable, dimension(:) :: Value
logical, allocatable, dimension(:) :: Mask logical, allocatable, dimension(:) :: Mask
end type tComponentBC end type tComponentBC
type, public :: tFieldBC type, public :: tFieldBC
integer(kind(FIELD_UNDEFINED_ID)) :: ID integer(kind(FIELD_UNDEFINED_ID)) :: ID
integer :: nComponents = 0 integer :: nComponents = 0
type(tComponentBC), allocatable :: componentBC(:) type(tComponentBC), allocatable, dimension(:) :: componentBC
end type tFieldBC end type tFieldBC
public :: & public :: &