to transpose already when constructing the interaction matrix

This commit is contained in:
Martin Diehl 2019-04-13 00:46:27 +02:00
parent ffdd3955a3
commit 1036c76ae0
8 changed files with 2334 additions and 2403 deletions

File diff suppressed because it is too large Load Diff

View File

@ -134,7 +134,6 @@ subroutine plastic_disloUCLA_init()
config_phase
use lattice
implicit none
integer :: &
Ninstance, &
p, i, &
@ -208,9 +207,9 @@ subroutine plastic_disloUCLA_init()
prm%nonSchmid_neg = prm%Schmid
endif
prm%h_sl_sl = transpose(lattice_interaction_SlipBySlip(prm%N_sl, &
prm%h_sl_sl = lattice_interaction_SlipBySlip(prm%N_sl, &
config%getFloats('interaction_slipslip'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%forestProjectionEdge = lattice_forestProjection(prm%N_sl,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
@ -361,7 +360,6 @@ end subroutine plastic_disloUCLA_init
!--------------------------------------------------------------------------------------------------
pure subroutine plastic_disloUCLA_LpAndItsTangent(Lp,dLp_dMp, &
Mp,T,instance,of)
implicit none
real(pReal), dimension(3,3), intent(out) :: &
Lp !< plastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
@ -411,7 +409,6 @@ subroutine plastic_disloUCLA_dotState(Mp,T,instance,of)
PI, &
math_clip
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
@ -472,7 +469,6 @@ end subroutine plastic_disloUCLA_dotState
!--------------------------------------------------------------------------------------------------
subroutine plastic_disloUCLA_dependentState(instance,of)
implicit none
integer, intent(in) :: &
instance, &
of
@ -507,7 +503,6 @@ function plastic_disloUCLA_postResults(Mp,T,instance,of) result(postResults)
PI, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
@ -616,7 +611,6 @@ pure subroutine kinetics(Mp,T,instance,of, &
PI, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &

View File

@ -198,7 +198,6 @@ subroutine plastic_dislotwin_init
config_phase
use lattice
implicit none
integer :: &
Ninstance, &
p, i, &
@ -268,9 +267,9 @@ subroutine plastic_dislotwin_init
slipActive: if (prm%sum_N_sl > 0) then
prm%P_sl = lattice_SchmidMatrix_slip(prm%N_sl,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
prm%h_sl_sl = transpose(lattice_interaction_SlipBySlip(prm%N_sl, &
prm%h_sl_sl = lattice_interaction_SlipBySlip(prm%N_sl, &
config%getFloats('interaction_slipslip'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%forestProjection = lattice_forestProjection (prm%N_sl,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
@ -332,9 +331,9 @@ subroutine plastic_dislotwin_init
if (prm%sum_N_tw > 0) then
prm%P_tw = lattice_SchmidMatrix_twin(prm%N_tw,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
prm%h_tw_tw = transpose(lattice_interaction_TwinByTwin(prm%N_tw,&
prm%h_tw_tw = lattice_interaction_TwinByTwin(prm%N_tw,&
config%getFloats('interaction_twintwin'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%b_tw = config%getFloats('twinburgers', requiredSize=size(prm%N_tw))
prm%t_tw = config%getFloats('twinsize', requiredSize=size(prm%N_tw))
@ -374,15 +373,15 @@ subroutine plastic_dislotwin_init
prm%b_tr = config%getFloats('transburgers')
prm%b_tr = math_expand(prm%b_tr,prm%N_tr)
prm%h = config%getFloat('transstackheight', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%i_tr = config%getFloat('cmfptrans', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%gamma_fcc_hex = config%getFloat('deltag')
prm%xc_trans = config%getFloat('xc_trans', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%L_tr = config%getFloat('l0_trans')
prm%h = config%getFloat('transstackheight', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%i_tr = config%getFloat('cmfptrans', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%gamma_fcc_hex = config%getFloat('deltag')
prm%xc_trans = config%getFloat('xc_trans', defaultVal=0.0_pReal) ! ToDo: How to handle that???
prm%L_tr = config%getFloat('l0_trans')
prm%h_tr_tr = transpose(lattice_interaction_TransByTrans(prm%N_tr,&
config%getFloats('interaction_transtrans'), &
config%getString('lattice_structure')))
prm%h_tr_tr = lattice_interaction_TransByTrans(prm%N_tr,&
config%getFloats('interaction_transtrans'), &
config%getString('lattice_structure'))
prm%C66_tr = lattice_C66_trans(prm%N_tr,prm%C66, &
config%getString('trans_lattice_structure'), &
@ -390,7 +389,7 @@ subroutine plastic_dislotwin_init
config%getFloat('a_bcc', defaultVal=0.0_pReal), &
config%getFloat('a_fcc', defaultVal=0.0_pReal))
prm%P_tr = lattice_SchmidMatrix_trans(prm%N_tr, &
prm%P_tr = lattice_SchmidMatrix_trans(prm%N_tr, &
config%getString('trans_lattice_structure'), &
0.0_pReal, &
config%getFloat('a_bcc', defaultVal=0.0_pReal), &
@ -416,16 +415,16 @@ subroutine plastic_dislotwin_init
endif
if (prm%sum_N_sl > 0 .and. prm%sum_N_tw > 0) then
prm%h_sl_tw = transpose(lattice_interaction_SlipByTwin(prm%N_sl,prm%N_tw,&
prm%h_sl_tw = lattice_interaction_SlipByTwin(prm%N_sl,prm%N_tw,&
config%getFloats('interaction_sliptwin'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
if (prm%fccTwinTransNucleation .and. prm%sum_N_tw > 12) write(6,*) 'mist' ! ToDo: implement better test. The model will fail also if N_tw is [6,6]
endif
if (prm%sum_N_sl > 0 .and. prm%sum_N_tr > 0) then
prm%h_sl_tr = transpose(lattice_interaction_SlipByTrans(prm%N_sl,prm%N_tr,&
prm%h_sl_tr = lattice_interaction_SlipByTrans(prm%N_sl,prm%N_tr,&
config%getFloats('interaction_sliptrans'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
if (prm%fccTwinTransNucleation .and. prm%sum_N_tr > 12) write(6,*) 'mist' ! ToDo: implement better test. The model will fail also if N_tr is [6,6]
endif
@ -605,7 +604,6 @@ function plastic_dislotwin_homogenizedC(ipc,ip,el) result(homogenizedC)
phase_plasticityInstance, &
phasememberAt
implicit none
real(pReal), dimension(6,6) :: &
homogenizedC
integer, intent(in) :: &
@ -653,7 +651,6 @@ subroutine plastic_dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,instance,of)
math_symmetric33, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(out) :: Lp
real(pReal), dimension(3,3,3,3), intent(out) :: dLp_dMp
real(pReal), dimension(3,3), intent(in) :: Mp
@ -776,7 +773,6 @@ subroutine plastic_dislotwin_dotState(Mp,T,instance,of)
math_mul33xx33, &
PI
implicit none
real(pReal), dimension(3,3), intent(in):: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
@ -869,7 +865,6 @@ subroutine plastic_dislotwin_dependentState(T,instance,of)
use math, only: &
PI
implicit none
integer, intent(in) :: &
instance, &
of
@ -987,7 +982,6 @@ function plastic_dislotwin_postResults(Mp,T,instance,of) result(postResults)
PI, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3),intent(in) :: &
Mp !< 2nd Piola Kirchhoff stress tensor in Mandel notation
real(pReal), intent(in) :: &
@ -1133,7 +1127,6 @@ pure subroutine kinetics_slip(Mp,T,instance,of, &
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
@ -1212,7 +1205,6 @@ pure subroutine kinetics_twin(Mp,T,dot_gamma_sl,instance,of,&
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
@ -1284,7 +1276,6 @@ pure subroutine kinetics_trans(Mp,T,dot_gamma_sl,instance,of,&
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &

View File

@ -108,7 +108,6 @@ subroutine plastic_isotropic_init
config_phase
use lattice
implicit none
integer :: &
Ninstance, &
p, i, &
@ -259,7 +258,6 @@ subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of)
math_deviatoric33, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(out) :: &
Lp !< plastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
@ -326,7 +324,6 @@ subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dTstar,Tstar,instance,of)
math_spherical33, &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(out) :: &
Li !< inleastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
@ -383,7 +380,6 @@ subroutine plastic_isotropic_dotState(Mp,instance,of)
math_mul33xx33, &
math_deviatoric33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -436,7 +432,6 @@ function plastic_isotropic_postResults(Mp,instance,of) result(postResults)
math_mul33xx33, &
math_deviatoric33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &

View File

@ -129,7 +129,6 @@ subroutine plastic_kinehardening_init
config_phase
use lattice
implicit none
integer :: &
Ninstance, &
p, i, o, &
@ -204,9 +203,9 @@ subroutine plastic_kinehardening_init
prm%nonSchmid_pos = prm%Schmid
prm%nonSchmid_neg = prm%Schmid
endif
prm%interaction_SlipSlip = transpose(lattice_interaction_SlipBySlip(prm%Nslip, &
prm%interaction_SlipSlip = lattice_interaction_SlipBySlip(prm%Nslip, &
config%getFloats('interaction_slipslip'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%crss0 = config%getFloats('crss0', requiredSize=size(prm%Nslip))
prm%tau1 = config%getFloats('tau1', requiredSize=size(prm%Nslip))
@ -347,7 +346,6 @@ end subroutine plastic_kinehardening_init
!--------------------------------------------------------------------------------------------------
pure subroutine plastic_kinehardening_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of)
implicit none
real(pReal), dimension(3,3), intent(out) :: &
Lp !< plastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
@ -390,7 +388,6 @@ end subroutine plastic_kinehardening_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
subroutine plastic_kinehardening_dotState(Mp,instance,of)
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -443,7 +440,6 @@ subroutine plastic_kinehardening_deltaState(Mp,instance,of)
debug_levelSelective
#endif
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -494,7 +490,6 @@ function plastic_kinehardening_postResults(Mp,instance,of) result(postResults)
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -608,7 +603,6 @@ pure subroutine kinetics(Mp,instance,of, &
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &

View File

@ -5,13 +5,13 @@
!> @brief Dummy plasticity for purely elastic material
!--------------------------------------------------------------------------------------------------
module plastic_none
implicit none
private
public :: &
plastic_none_init
implicit none
private
public :: &
plastic_none_init
contains
!--------------------------------------------------------------------------------------------------
@ -19,39 +19,39 @@ contains
!> @details reads in material parameters, allocates arrays, and does sanity checks
!--------------------------------------------------------------------------------------------------
subroutine plastic_none_init
use debug, only: &
debug_level, &
debug_constitutive, &
debug_levelBasic
use material, only: &
phase_plasticity, &
material_allocatePlasticState, &
PLASTICITY_NONE_label, &
PLASTICITY_NONE_ID, &
material_phase, &
plasticState
implicit none
integer :: &
Ninstance, &
p, &
NipcMyPhase
write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_NONE_label//' init -+>>>'
Ninstance = count(phase_plasticity == PLASTICITY_NONE_ID)
if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
do p = 1, size(phase_plasticity)
if (phase_plasticity(p) /= PLASTICITY_NONE_ID) cycle
NipcMyPhase = count(material_phase == p)
call material_allocatePlasticState(p,NipcMyPhase,0,0,0, &
0,0,0)
plasticState(p)%sizePostResults = 0
enddo
use debug, only: &
debug_level, &
debug_constitutive, &
debug_levelBasic
use material, only: &
phase_plasticity, &
material_allocatePlasticState, &
PLASTICITY_NONE_label, &
PLASTICITY_NONE_ID, &
material_phase, &
plasticState
implicit none
integer :: &
Ninstance, &
p, &
NipcMyPhase
write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_NONE_label//' init -+>>>'
Ninstance = count(phase_plasticity == PLASTICITY_NONE_ID)
if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
do p = 1, size(phase_plasticity)
if (phase_plasticity(p) /= PLASTICITY_NONE_ID) cycle
NipcMyPhase = count(material_phase == p)
call material_allocatePlasticState(p,NipcMyPhase,0,0,0, &
0,0,0)
plasticState(p)%sizePostResults = 0
enddo
end subroutine plastic_none_init

View File

@ -260,7 +260,6 @@ subroutine plastic_nonlocal_init
use config
use lattice
implicit none
character(len=65536), dimension(0), parameter :: emptyStringArray = [character(len=65536)::]
integer, dimension(0), parameter :: emptyIntArray = [integer::]
real(pReal), dimension(0), parameter :: emptyRealArray = [real(pReal)::]
@ -751,7 +750,6 @@ subroutine plastic_nonlocal_init
material_phase, &
phase_plasticityInstance, &
phasememberAt
implicit none
integer,intent(in) ::&
phase, &
@ -867,7 +865,6 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el)
LATTICE_fcc_ID, &
lattice_structure
implicit none
integer, intent(in) :: &
ip, &
el
@ -1090,7 +1087,6 @@ end subroutine plastic_nonlocal_dependentState
!--------------------------------------------------------------------------------------------------
subroutine plastic_nonlocal_kinetics(v, dv_dtau, dv_dtauNS, tau, tauNS, &
tauThreshold, c, Temperature, instance, of)
implicit none
integer, intent(in) :: &
c, & !< dislocation character (1:edge, 2:screw)
instance, of
@ -1239,7 +1235,6 @@ subroutine plastic_nonlocal_LpAndItsTangent(Lp, dLp_dMp, &
phaseAt, phasememberAt, &
phase_plasticityInstance
implicit none
integer, intent(in) :: &
ip, & !< current integration point
el !< current element number
@ -1392,7 +1387,6 @@ subroutine plastic_nonlocal_deltaState(Mp,ip,el)
phaseAt, phasememberAt, &
phase_plasticityInstance
implicit none
integer, intent(in) :: &
ip, &
el
@ -1553,7 +1547,6 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
LATTICE_bcc_ID, &
LATTICE_fcc_ID
implicit none
integer, intent(in) :: &
ip, & !< current integration point
el !< current element number
@ -2027,7 +2020,6 @@ subroutine plastic_nonlocal_updateCompatibility(orientation,i,e)
use lattice, only: &
lattice_qDisorientation
implicit none
integer, intent(in) :: &
i, &
e
@ -2175,7 +2167,6 @@ function plastic_nonlocal_postResults(ph,instance,of) result(postResults)
use material, only: &
plasticState
implicit none
integer, intent(in) :: &
ph, &
instance, &
@ -2378,7 +2369,6 @@ end function plastic_nonlocal_postResults
function getRho(instance,of,ip,el)
use mesh
implicit none
integer, intent(in) :: instance, of,ip,el
real(pReal), dimension(param(instance)%totalNslip,10) :: getRho

View File

@ -129,7 +129,6 @@ subroutine plastic_phenopowerlaw_init
config_phase
use lattice
implicit none
integer :: &
Ninstance, &
p, i, &
@ -203,9 +202,9 @@ subroutine plastic_phenopowerlaw_init
prm%nonSchmid_pos = prm%Schmid_slip
prm%nonSchmid_neg = prm%Schmid_slip
endif
prm%interaction_SlipSlip = transpose(lattice_interaction_SlipBySlip(prm%Nslip, &
prm%interaction_SlipSlip = lattice_interaction_SlipBySlip(prm%Nslip, &
config%getFloats('interaction_slipslip'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%xi_slip_0 = config%getFloats('tau0_slip', requiredSize=size(prm%Nslip))
prm%xi_slip_sat = config%getFloats('tausat_slip', requiredSize=size(prm%Nslip))
@ -240,9 +239,9 @@ subroutine plastic_phenopowerlaw_init
twinActive: if (prm%totalNtwin > 0) then
prm%Schmid_twin = lattice_SchmidMatrix_twin(prm%Ntwin,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
prm%interaction_TwinTwin = transpose(lattice_interaction_TwinByTwin(prm%Ntwin,&
prm%interaction_TwinTwin = lattice_interaction_TwinByTwin(prm%Ntwin,&
config%getFloats('interaction_twintwin'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
prm%gamma_twin_char = lattice_characteristicShear_twin(prm%Ntwin,config%getString('lattice_structure'),&
config%getFloat('c/a'))
@ -268,12 +267,12 @@ subroutine plastic_phenopowerlaw_init
!--------------------------------------------------------------------------------------------------
! slip-twin related parameters
slipAndTwinActive: if (prm%totalNslip > 0 .and. prm%totalNtwin > 0) then
prm%interaction_SlipTwin = transpose(lattice_interaction_SlipByTwin(prm%Nslip,prm%Ntwin,&
prm%interaction_SlipTwin = lattice_interaction_SlipByTwin(prm%Nslip,prm%Ntwin,&
config%getFloats('interaction_sliptwin'), &
config%getString('lattice_structure')))
prm%interaction_TwinSlip = transpose(lattice_interaction_TwinBySlip(prm%Ntwin,prm%Nslip,&
config%getString('lattice_structure'))
prm%interaction_TwinSlip = lattice_interaction_TwinBySlip(prm%Ntwin,prm%Nslip,&
config%getFloats('interaction_twinslip'), &
config%getString('lattice_structure')))
config%getString('lattice_structure'))
else slipAndTwinActive
allocate(prm%interaction_SlipTwin(prm%TotalNslip,prm%TotalNtwin)) ! at least one dimension is 0
allocate(prm%interaction_TwinSlip(prm%TotalNtwin,prm%TotalNslip)) ! at least one dimension is 0
@ -387,7 +386,6 @@ end subroutine plastic_phenopowerlaw_init
!--------------------------------------------------------------------------------------------------
pure subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of)
implicit none
real(pReal), dimension(3,3), intent(out) :: &
Lp !< plastic velocity gradient
real(pReal), dimension(3,3,3,3), intent(out) :: &
@ -439,7 +437,6 @@ end subroutine plastic_phenopowerlaw_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
subroutine plastic_phenopowerlaw_dotState(Mp,instance,of)
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -498,7 +495,6 @@ function plastic_phenopowerlaw_postResults(Mp,instance,of) result(postResults)
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -616,7 +612,6 @@ pure subroutine kinetics_slip(Mp,instance,of, &
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
@ -693,7 +688,6 @@ pure subroutine kinetics_twin(Mp,instance,of,&
use math, only: &
math_mul33xx33
implicit none
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &