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

simplified

This commit is contained in:
Martin Diehl 2020-12-27 11:48:02 +01:00
parent cee04c9b5f
commit e8ea815d92
16 changed files with 49 additions and 94 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
src/CPFEM.f90
View File

@ -5,7 +5,6 @@
!--------------------------------------------------------------------------------------------------
module CPFEM
use prec
use FEsolving
use math
use rotations
use YAML_types
@ -197,11 +196,9 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS
CPFEM_dcsde(1:6,1:6,ip,elCP) = ODD_JACOBIAN * math_eye(6)
else validCalculation
FEsolving_execElem = elCP
FEsolving_execIP = ip
if (debugCPFEM%extensive) &
print'(a,i8,1x,i2)', '<< CPFEM >> calculation for elFE ip ',elFE,ip
call materialpoint_stressAndItsTangent(dt)
call materialpoint_stressAndItsTangent(dt,[ip,ip],[elCP,elCP])
terminalIllness: if (terminallyIll) then

1
src/CPFEM2.f90
View File

@ -6,7 +6,6 @@
module CPFEM2
use prec
use config
use FEsolving
use math
use rotations
use YAML_types

1
src/DAMASK_marc.f90
View File

@ -176,7 +176,6 @@ subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,nshear,disp, &
use DAMASK_interface
use config
use YAML_types
use FEsolving
use discretization_marc
use homogenization
use CPFEM

15
src/FEsolving.f90
View File

@ -1,15 +0,0 @@
!--------------------------------------------------------------------------------------------------
!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
!> @brief global variables for flow control
!--------------------------------------------------------------------------------------------------
module FEsolving
implicit none
public
integer, dimension(2) :: &
FEsolving_execElem, & !< for ping-pong scheme always whole range, otherwise one specific element
FEsolving_execIP !< for ping-pong scheme always range to max IP, otherwise one specific IP
end module FEsolving

1
src/commercialFEM_fileList.f90
View File

@ -13,7 +13,6 @@
#include "math.f90"
#include "quaternions.f90"
#include "rotations.f90"
#include "FEsolving.f90"
#include "element.f90"
#include "HDF5_utilities.f90"
#include "results.f90"

39
src/constitutive.f90
View File

@ -16,7 +16,6 @@ module constitutive
use parallelization
use HDF5_utilities
use DAMASK_interface
use FEsolving
use results
implicit none
@ -940,8 +939,8 @@ subroutine crystallite_init
flush(IO_STDOUT)
!$OMP PARALLEL DO PRIVATE(ph,me)
do el = FEsolving_execElem(1),FEsolving_execElem(2)
do ip = FEsolving_execIP(1), FEsolving_execIP(2); do co = 1, homogenization_Nconstituents(material_homogenizationAt(el))
do el = 1, size(material_phaseMemberAt,3)
do ip = 1, size(material_phaseMemberAt,2); do co = 1, homogenization_Nconstituents(material_homogenizationAt(el))
ph = material_phaseAt(co,el)
me = material_phaseMemberAt(co,ip,el)
@ -967,14 +966,14 @@ subroutine crystallite_init
crystallite_partitionedF0 = crystallite_F0
crystallite_partitionedF = crystallite_F0
call crystallite_orientations()
!$OMP PARALLEL DO PRIVATE(ph,me)
do el = FEsolving_execElem(1),FEsolving_execElem(2)
do ip = FEsolving_execIP(1),FEsolving_execIP(2)
do el = 1, size(material_phaseMemberAt,3)
do ip = 1, size(material_phaseMemberAt,2)
do co = 1,homogenization_Nconstituents(material_homogenizationAt(el))
ph = material_phaseAt(co,el)
me = material_phaseMemberAt(co,ip,el)
call crystallite_orientations(co,ip,el)
call constitutive_plastic_dependentState(crystallite_partitionedF0(1:3,1:3,co,ip,el),co,ip,el) ! update dependent state variables to be consistent with basic states
enddo
enddo
@ -1210,34 +1209,20 @@ end function crystallite_stressTangent
!--------------------------------------------------------------------------------------------------
!> @brief calculates orientations
!--------------------------------------------------------------------------------------------------
subroutine crystallite_orientations
subroutine crystallite_orientations(co,ip,el)
integer &
integer, intent(in) :: &
co, & !< counter in integration point component loop
ip, & !< counter in integration point loop
el !< counter in element loop
!$OMP PARALLEL DO
do el = FEsolving_execElem(1),FEsolving_execElem(2)
do ip = FEsolving_execIP(1),FEsolving_execIP(2)
do co = 1,homogenization_Nconstituents(material_homogenizationAt(el))
call crystallite_orientation(co,ip,el)%fromMatrix(transpose(math_rotationalPart(crystallite_Fe(1:3,1:3,co,ip,el))))
enddo; enddo; enddo
!$OMP END PARALLEL DO
call crystallite_orientation(co,ip,el)%fromMatrix(transpose(math_rotationalPart(crystallite_Fe(1:3,1:3,co,ip,el))))
if (plasticState(material_phaseAt(1,el))%nonlocal) &
call plastic_nonlocal_updateCompatibility(crystallite_orientation, &
phase_plasticityInstance(material_phaseAt(1,el)),ip,el)
nonlocalPresent: if (any(plasticState%nonlocal)) then
!$OMP PARALLEL DO
do el = FEsolving_execElem(1),FEsolving_execElem(2)
if (plasticState(material_phaseAt(1,el))%nonlocal) then
do ip = FEsolving_execIP(1),FEsolving_execIP(2)
call plastic_nonlocal_updateCompatibility(crystallite_orientation, &
phase_plasticityInstance(material_phaseAt(1,el)),ip,el)
enddo
endif
enddo
!$OMP END PARALLEL DO
endif nonlocalPresent
end subroutine crystallite_orientations

4
src/grid/discretization_grid.f90
View File

@ -19,7 +19,6 @@ module discretization_grid
use results
use discretization
use geometry_plastic_nonlocal
use FEsolving
implicit none
private
@ -117,9 +116,6 @@ subroutine discretization_grid_init(restart)
(grid(1)+1) * (grid(2)+1) * grid3,& ! ...unless not last process
worldrank+1==worldsize))
FEsolving_execElem = [1,product(myGrid)] ! parallel loop bounds set to comprise all elements
FEsolving_execIP = [1,1] ! parallel loop bounds set to comprise the only IP
!--------------------------------------------------------------------------------------------------
! store geometry information for post processing
if(.not. restart) then

1
src/grid/grid_mech_FEM.f90
View File

@ -18,7 +18,6 @@ module grid_mech_FEM
use math
use rotations
use spectral_utilities
use FEsolving
use config
use homogenization
use discretization

1
src/grid/grid_mech_spectral_basic.f90
View File

@ -18,7 +18,6 @@ module grid_mech_spectral_basic
use math
use rotations
use spectral_utilities
use FEsolving
use config
use homogenization
use discretization_grid

1
src/grid/grid_mech_spectral_polarisation.f90
View File

@ -18,7 +18,6 @@ module grid_mech_spectral_polarisation
use math
use rotations
use spectral_utilities
use FEsolving
use config
use homogenization
use discretization_grid

4
src/grid/spectral_utilities.f90
View File

@ -810,9 +810,9 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,&
print'(/,a)', ' ... evaluating constitutive response ......................................'
flush(IO_STDOUT)
homogenization_F = reshape(F,[3,3,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field
homogenization_F = reshape(F,[3,3,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field
call materialpoint_stressAndItsTangent(timeinc) ! calculate P field
call materialpoint_stressAndItsTangent(timeinc,[1,1],[1,product(grid(1:2))*grid3]) ! calculate P field
P = reshape(homogenization_P, [3,3,grid(1),grid(2),grid3])
P_av = sum(sum(sum(P,dim=5),dim=4),dim=3) * wgt ! average of P

57
src/homogenization.f90
View File

@ -11,7 +11,6 @@ module homogenization
use math
use material
use constitutive
use FEsolving
use discretization
use thermal_isothermal
use thermal_conduction
@ -144,15 +143,16 @@ end subroutine homogenization_init
!--------------------------------------------------------------------------------------------------
!> @brief parallelized calculation of stress and corresponding tangent at material points
!--------------------------------------------------------------------------------------------------
subroutine materialpoint_stressAndItsTangent(dt)
subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execElem)
real(pReal), intent(in) :: dt !< time increment
integer, dimension(2), intent(in) :: FEsolving_execElem, FEsolving_execIP
integer :: &
NiterationHomog, &
NiterationMPstate, &
ip, & !< integration point number
el, & !< element number
myNgrains, co, ce
myNgrains, co, ce, ho
real(pReal) :: &
subFrac, &
subStep
@ -162,8 +162,10 @@ subroutine materialpoint_stressAndItsTangent(dt)
doneAndHappy
!$OMP PARALLEL DO PRIVATE(ce,myNgrains,NiterationMPstate,NiterationHomog,subFrac,converged,subStep,doneAndHappy)
!$OMP PARALLEL DO PRIVATE(ce,ho,myNgrains,NiterationMPstate,NiterationHomog,subFrac,converged,subStep,doneAndHappy)
do el = FEsolving_execElem(1),FEsolving_execElem(2)
ho = material_homogenizationAt(el)
myNgrains = homogenization_Nconstituents(ho)
do ip = FEsolving_execIP(1),FEsolving_execIP(2)
!--------------------------------------------------------------------------------------------------
@ -174,18 +176,19 @@ subroutine materialpoint_stressAndItsTangent(dt)
converged = .false. ! pretend failed step ...
subStep = 1.0_pReal/num%subStepSizeHomog ! ... larger then the requested calculation
if (homogState(material_homogenizationAt(el))%sizeState > 0) &
homogState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el)) = &
homogState(material_homogenizationAt(el))%State0( :,material_homogenizationMemberAt(ip,el))
if (homogState(ho)%sizeState > 0) &
homogState(ho)%subState0(:,material_homogenizationMemberAt(ip,el)) = &
homogState(ho)%State0( :,material_homogenizationMemberAt(ip,el))
if (damageState(ho)%sizeState > 0) &
damageState(ho)%subState0(:,material_homogenizationMemberAt(ip,el)) = &
damageState(ho)%State0( :,material_homogenizationMemberAt(ip,el))
if (damageState(material_homogenizationAt(el))%sizeState > 0) &
damageState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el)) = &
damageState(material_homogenizationAt(el))%State0( :,material_homogenizationMemberAt(ip,el))
NiterationHomog = 0
cutBackLooping: do while (.not. terminallyIll .and. subStep > num%subStepMinHomog)
myNgrains = homogenization_Nconstituents(material_homogenizationAt(el))
if (converged) then
subFrac = subFrac + subStep
@ -196,12 +199,12 @@ subroutine materialpoint_stressAndItsTangent(dt)
! wind forward grain starting point
call constitutive_windForward(ip,el)
if(homogState(material_homogenizationAt(el))%sizeState > 0) &
homogState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el)) = &
homogState(material_homogenizationAt(el))%State (:,material_homogenizationMemberAt(ip,el))
if(damageState(material_homogenizationAt(el))%sizeState > 0) &
damageState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el)) = &
damageState(material_homogenizationAt(el))%State (:,material_homogenizationMemberAt(ip,el))
if(homogState(ho)%sizeState > 0) &
homogState(ho)%subState0(:,material_homogenizationMemberAt(ip,el)) = &
homogState(ho)%State (:,material_homogenizationMemberAt(ip,el))
if(damageState(ho)%sizeState > 0) &
damageState(ho)%subState0(:,material_homogenizationMemberAt(ip,el)) = &
damageState(ho)%State (:,material_homogenizationMemberAt(ip,el))
endif steppingNeeded
@ -219,12 +222,12 @@ subroutine materialpoint_stressAndItsTangent(dt)
call crystallite_restore(ip,el,subStep < 1.0_pReal)
call constitutive_restore(ip,el)
if(homogState(material_homogenizationAt(el))%sizeState > 0) &
homogState(material_homogenizationAt(el))%State( :,material_homogenizationMemberAt(ip,el)) = &
homogState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el))
if(damageState(material_homogenizationAt(el))%sizeState > 0) &
damageState(material_homogenizationAt(el))%State( :,material_homogenizationMemberAt(ip,el)) = &
damageState(material_homogenizationAt(el))%subState0(:,material_homogenizationMemberAt(ip,el))
if(homogState(ho)%sizeState > 0) &
homogState(ho)%State( :,material_homogenizationMemberAt(ip,el)) = &
homogState(ho)%subState0(:,material_homogenizationMemberAt(ip,el))
if(damageState(ho)%sizeState > 0) &
damageState(ho)%State( :,material_homogenizationMemberAt(ip,el)) = &
damageState(ho)%subState0(:,material_homogenizationMemberAt(ip,el))
endif
endif
@ -275,10 +278,14 @@ subroutine materialpoint_stressAndItsTangent(dt)
!$OMP END PARALLEL DO
if (.not. terminallyIll ) then
call crystallite_orientations() ! calculate crystal orientations
!$OMP PARALLEL DO
!$OMP PARALLEL DO PRIVATE(ho,myNgrains)
elementLooping3: do el = FEsolving_execElem(1),FEsolving_execElem(2)
ho = material_homogenizationAt(el)
myNgrains = homogenization_Nconstituents(ho)
IpLooping3: do ip = FEsolving_execIP(1),FEsolving_execIP(2)
do co = 1, myNgrains
call crystallite_orientations(co,ip,el)
enddo
call mech_homogenize(ip,el)
enddo IpLooping3
enddo elementLooping3

4
src/marc/discretization_marc.f90
View File

@ -12,7 +12,6 @@ module discretization_marc
use DAMASK_interface
use IO
use config
use FEsolving
use element
use discretization
use geometry_plastic_nonlocal
@ -89,9 +88,6 @@ subroutine discretization_marc_init
if (debug_e < 1 .or. debug_e > nElems) call IO_error(602,ext_msg='element')
if (debug_i < 1 .or. debug_i > elem%nIPs) call IO_error(602,ext_msg='IP')
FEsolving_execElem = [1,nElems]
FEsolving_execIP = [1,elem%nIPs]
allocate(cellNodeDefinition(elem%nNodes-1))
allocate(connectivity_cell(elem%NcellNodesPerCell,elem%nIPs,nElems))
call buildCells(connectivity_cell,cellNodeDefinition,&

1
src/mesh/DAMASK_mesh.f90
View File

@ -15,7 +15,6 @@ program DAMASK_mesh
use IO
use math
use CPFEM2
use FEsolving
use config
use discretization_mesh
use FEM_Utilities

2
src/mesh/FEM_utilities.f90
View File

@ -160,7 +160,7 @@ subroutine utilities_constitutiveResponse(timeinc,P_av,forwardData)
print'(/,a)', ' ... evaluating constitutive response ......................................'
call materialpoint_stressAndItsTangent(timeinc) ! calculate P field
call materialpoint_stressAndItsTangent(timeinc,[1,mesh_maxNips],[1,mesh_NcpElems]) ! calculate P field
cutBack = .false. ! reset cutBack status

6
src/mesh/discretization_mesh.f90
View File

@ -18,7 +18,6 @@ module discretization_mesh
use config
use discretization
use results
use FEsolving
use FEM_quadrature
use YAML_types
use prec
@ -30,7 +29,7 @@ module discretization_mesh
mesh_Nboundaries, &
mesh_NcpElemsGlobal
integer :: &
integer, public, protected :: &
mesh_NcpElems !< total number of CP elements in mesh
!!!! BEGIN DEPRECATED !!!!!
@ -174,9 +173,6 @@ subroutine discretization_mesh_init(restart)
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')
FEsolving_execElem = [1,mesh_NcpElems] ! parallel loop bounds set to comprise all DAMASK elements
FEsolving_execIP = [1,mesh_maxNips]
allocate(mesh_node0(3,mesh_Nnodes),source=0.0_pReal)
call discretization_init(materialAt,&