From fd94c786f0e38bcea0c7cdd527ec2f0ac6d35658 Mon Sep 17 00:00:00 2001 From: Pratheek Shanthraj Date: Mon, 21 Jan 2013 21:57:26 +0000 Subject: [PATCH] moved stiffness tensor calculation to lattice introduced 'isotropic' and 'orthorhombic' lattice types to use corresponding symmetries in stiffness tensor. intended to be used with non-crystal plasticity models (j2, constitutive_none with isotropic, cubic or orthotropic elasticity). --- code/config/material.config | 9 +++- code/constitutive_dislotwin.f90 | 61 ++++++--------------- code/constitutive_j2.f90 | 46 ++++++++++------ code/constitutive_none.f90 | 81 +++++++--------------------- code/constitutive_nonlocal.f90 | 83 ++++++++++------------------- code/constitutive_phenopowerlaw.f90 | 63 +++++++--------------- code/constitutive_titanmod.f90 | 75 +++++++++----------------- code/crystallite.f90 | 22 ++++---- code/lattice.f90 | 77 +++++++++++++++++++++++--- 9 files changed, 226 insertions(+), 291 deletions(-) diff --git a/code/config/material.config b/code/config/material.config index c5e7af71e..98c633ff3 100644 --- a/code/config/material.config +++ b/code/config/material.config @@ -99,10 +99,16 @@ crystallite 2 elasticity hooke plasticity none -lattice_structure fcc +lattice_structure orthorhombic c11 106.75e9 +c22 106.75e9 +c33 106.75e9 c12 60.41e9 +c13 60.41e9 +c23 60.41e9 c44 28.34e9 +c55 28.34e9 +c66 28.34e9 [Aluminum_J2isotropic] @@ -112,6 +118,7 @@ plasticity j2 (output) flowstress (output) strainrate +lattice_structure isotropic c11 110.9e9 c12 58.34e9 taylorfactor 3 diff --git a/code/constitutive_dislotwin.f90 b/code/constitutive_dislotwin.f90 index 62bb872de..598a77625 100644 --- a/code/constitutive_dislotwin.f90 +++ b/code/constitutive_dislotwin.f90 @@ -57,11 +57,6 @@ integer(pInt), dimension(:), allocatable :: constitutive_dislotwin integer(pInt), dimension(:,:), allocatable :: constitutive_dislotwin_Nslip, & ! number of active slip systems for each family and instance constitutive_dislotwin_Ntwin ! number of active twin systems for each family and instance real(pReal), dimension(:), allocatable :: constitutive_dislotwin_CoverA, & ! c/a ratio for hex type lattice - constitutive_dislotwin_C11, & ! C11 element in elasticity matrix - constitutive_dislotwin_C12, & ! C12 element in elasticity matrix - constitutive_dislotwin_C13, & ! C13 element in elasticity matrix - constitutive_dislotwin_C33, & ! C33 element in elasticity matrix - constitutive_dislotwin_C44, & ! C44 element in elasticity matrix constitutive_dislotwin_Gmod, & ! shear modulus constitutive_dislotwin_CAtomicVolume, & ! atomic volume in Bugers vector unit constitutive_dislotwin_D0, & ! prefactor for self-diffusion coefficient @@ -186,16 +181,6 @@ allocate(constitutive_dislotwin_totalNtwin(maxNinstance)) constitutive_dislotwin_totalNtwin = 0_pInt allocate(constitutive_dislotwin_CoverA(maxNinstance)) constitutive_dislotwin_CoverA = 0.0_pReal -allocate(constitutive_dislotwin_C11(maxNinstance)) - constitutive_dislotwin_C11 = 0.0_pReal -allocate(constitutive_dislotwin_C12(maxNinstance)) - constitutive_dislotwin_C12 = 0.0_pReal -allocate(constitutive_dislotwin_C13(maxNinstance)) - constitutive_dislotwin_C13 = 0.0_pReal -allocate(constitutive_dislotwin_C33(maxNinstance)) - constitutive_dislotwin_C33 = 0.0_pReal -allocate(constitutive_dislotwin_C44(maxNinstance)) - constitutive_dislotwin_C44 = 0.0_pReal allocate(constitutive_dislotwin_Gmod(maxNinstance)) constitutive_dislotwin_Gmod = 0.0_pReal allocate(constitutive_dislotwin_CAtomicVolume(maxNinstance)) @@ -301,15 +286,23 @@ do ! read thru sections of case ('covera_ratio') constitutive_dislotwin_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_dislotwin_C11(i) = IO_floatValue(line,positions,2_pInt) + constitutive_dislotwin_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_dislotwin_C12(i) = IO_floatValue(line,positions,2_pInt) + constitutive_dislotwin_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_dislotwin_C13(i) = IO_floatValue(line,positions,2_pInt) + constitutive_dislotwin_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c22') + constitutive_dislotwin_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c23') + constitutive_dislotwin_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_dislotwin_C33(i) = IO_floatValue(line,positions,2_pInt) + constitutive_dislotwin_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_dislotwin_C44(i) = IO_floatValue(line,positions,2_pInt) + constitutive_dislotwin_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + case ('c55') + constitutive_dislotwin_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + case ('c66') + constitutive_dislotwin_Cslip_66(6,6,i) = IO_floatValue(line,positions,2_pInt) case ('nslip') forall (j = 1_pInt:lattice_maxNslipFamily) & constitutive_dislotwin_Nslip(j,i) = IO_intValue(line,positions,1_pInt+j) @@ -549,32 +542,12 @@ do i = 1_pInt,maxNinstance !* Elasticity matrix and shear modulus according to material.config - select case (myStructure) - case(1_pInt:2_pInt) ! cubic(s) - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) & - constitutive_dislotwin_Cslip_66(k,j,i) = constitutive_dislotwin_C12(i) - constitutive_dislotwin_Cslip_66(k,k,i) = constitutive_dislotwin_C11(i) - constitutive_dislotwin_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_dislotwin_C44(i) - end forall - case(3_pInt:) ! all hex - constitutive_dislotwin_Cslip_66(1,1,i) = constitutive_dislotwin_C11(i) - constitutive_dislotwin_Cslip_66(2,2,i) = constitutive_dislotwin_C11(i) - constitutive_dislotwin_Cslip_66(3,3,i) = constitutive_dislotwin_C33(i) - constitutive_dislotwin_Cslip_66(1,2,i) = constitutive_dislotwin_C12(i) - constitutive_dislotwin_Cslip_66(2,1,i) = constitutive_dislotwin_C12(i) - constitutive_dislotwin_Cslip_66(1,3,i) = constitutive_dislotwin_C13(i) - constitutive_dislotwin_Cslip_66(3,1,i) = constitutive_dislotwin_C13(i) - constitutive_dislotwin_Cslip_66(2,3,i) = constitutive_dislotwin_C13(i) - constitutive_dislotwin_Cslip_66(3,2,i) = constitutive_dislotwin_C13(i) - constitutive_dislotwin_Cslip_66(4,4,i) = constitutive_dislotwin_C44(i) - constitutive_dislotwin_Cslip_66(5,5,i) = constitutive_dislotwin_C44(i) - constitutive_dislotwin_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_dislotwin_C11(i)-constitutive_dislotwin_C12(i)) - end select + constitutive_dislotwin_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_dislotwin_structureName(i),& + constitutive_dislotwin_Cslip_66) + constitutive_dislotwin_Gmod(i) = & + 0.2_pReal*(constitutive_dislotwin_Cslip_66(1,1,i)-constitutive_dislotwin_Cslip_66(1,2,i))+0.3_pReal*constitutive_dislotwin_Cslip_66(4,4,i) constitutive_dislotwin_Cslip_66(:,:,i) = math_Mandel3333to66(math_Voigt66to3333(constitutive_dislotwin_Cslip_66(:,:,i))) constitutive_dislotwin_Cslip_3333(:,:,:,:,i) = math_Voigt66to3333(constitutive_dislotwin_Cslip_66(:,:,i)) - constitutive_dislotwin_Gmod(i) = & - 0.2_pReal*(constitutive_dislotwin_C11(i)-constitutive_dislotwin_C12(i))+0.3_pReal*constitutive_dislotwin_C44(i) !* Process slip related parameters ------------------------------------------------ diff --git a/code/constitutive_j2.f90 b/code/constitutive_j2.f90 index 1384415a3..938a17199 100644 --- a/code/constitutive_j2.f90 +++ b/code/constitutive_j2.f90 @@ -60,11 +60,12 @@ module constitutive_j2 constitutive_j2_output ! name of each post result output integer(pInt), dimension(:), allocatable, private :: & - constitutive_j2_Noutput + constitutive_j2_Noutput !< name of each post result output + + character(len=32), dimension(:), allocatable, private :: & + constitutive_j2_structureName real(pReal), dimension(:), allocatable, private :: & - constitutive_j2_C11, & - constitutive_j2_C12, & !* Visco-plastic constitutive_j2 parameters constitutive_j2_fTaylor, & constitutive_j2_tau0, & @@ -120,6 +121,7 @@ subroutine constitutive_j2_init(myFile) debug_level, & debug_constitutive, & debug_levelBasic + use lattice, only: lattice_symmetrizeC66 implicit none integer(pInt), intent(in) :: myFile @@ -156,10 +158,8 @@ subroutine constitutive_j2_init(myFile) constitutive_j2_output = '' allocate(constitutive_j2_Noutput(maxNinstance)) constitutive_j2_Noutput = 0_pInt - allocate(constitutive_j2_C11(maxNinstance)) - constitutive_j2_C11 = 0.0_pReal - allocate(constitutive_j2_C12(maxNinstance)) - constitutive_j2_C12 = 0.0_pReal + allocate(constitutive_j2_structureName(maxNinstance)) + constitutive_j2_structureName = '' allocate(constitutive_j2_Cslip_66(6,6,maxNinstance)) constitutive_j2_Cslip_66 = 0.0_pReal allocate(constitutive_j2_fTaylor(maxNinstance)) @@ -213,10 +213,26 @@ subroutine constitutive_j2_init(myFile) case ('(output)') constitutive_j2_Noutput(i) = constitutive_j2_Noutput(i) + 1_pInt constitutive_j2_output(constitutive_j2_Noutput(i),i) = IO_lc(IO_stringValue(line,positions,2_pInt)) + case ('lattice_structure') + constitutive_j2_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('c11') - constitutive_j2_C11(i) = IO_floatValue(line,positions,2_pInt) + constitutive_j2_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_j2_C12(i) = IO_floatValue(line,positions,2_pInt) + constitutive_j2_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c13') + constitutive_j2_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c22') + constitutive_j2_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c23') + constitutive_j2_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c33') + constitutive_j2_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c44') + constitutive_j2_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + case ('c55') + constitutive_j2_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + case ('c66') + constitutive_j2_Cslip_66(6,6,i) = IO_floatValue(line,positions,2_pInt) case ('tau0') constitutive_j2_tau0(i) = IO_floatValue(line,positions,2_pInt) case ('gdot0') @@ -250,6 +266,7 @@ subroutine constitutive_j2_init(myFile) enddo 100 do i = 1_pInt,maxNinstance ! sanity checks + if (constitutive_j2_structureName(i) == '') call IO_error(205_pInt,e=i) if (constitutive_j2_tau0(i) < 0.0_pReal) call IO_error(211_pInt,ext_msg='tau0 (' & //constitutive_j2_label//')') if (constitutive_j2_gdot0(i) <= 0.0_pReal) call IO_error(211_pInt,ext_msg='gdot0 (' & @@ -286,14 +303,9 @@ subroutine constitutive_j2_init(myFile) constitutive_j2_sizeDotState(i) = 1_pInt constitutive_j2_sizeState(i) = 1_pInt - - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) - constitutive_j2_Cslip_66(k,j,i) = constitutive_j2_C12(i) - end forall - constitutive_j2_Cslip_66(k,k,i) = constitutive_j2_C11(i) - constitutive_j2_Cslip_66(k+3,k+3,i) = 0.5_pReal*(constitutive_j2_C11(i)-constitutive_j2_C12(i)) - end forall + + constitutive_j2_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_j2_structureName(i),& + constitutive_j2_Cslip_66) constitutive_j2_Cslip_66(1:6,1:6,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_j2_Cslip_66(1:6,1:6,i))) diff --git a/code/constitutive_none.f90 b/code/constitutive_none.f90 index 9a00ee80a..c1eb16bef 100644 --- a/code/constitutive_none.f90 +++ b/code/constitutive_none.f90 @@ -46,22 +46,13 @@ module constitutive_none integer(pInt), dimension(:), allocatable, public :: & constitutive_none_sizeDotState, & constitutive_none_sizeState, & - constitutive_none_sizePostResults, & - constitutive_none_structure + constitutive_none_sizePostResults character(len=32), dimension(:), allocatable, private :: & constitutive_none_structureName integer(pInt), dimension(:,:), allocatable, target, public :: & constitutive_none_sizePostResult ! size of each post result output - - real(pReal), dimension(:), allocatable, private :: & - constitutive_none_CoverA, & - constitutive_none_C11, & - constitutive_none_C12, & - constitutive_none_C13, & - constitutive_none_C33, & - constitutive_none_C44 real(pReal), dimension(:,:,:), allocatable, private :: & constitutive_none_Cslip_66 @@ -100,7 +91,7 @@ subroutine constitutive_none_init(myFile) debug_level, & debug_constitutive, & debug_levelBasic - use lattice, only: lattice_initializeStructure, lattice_symmetryType + use lattice, only: lattice_symmetrizeC66 implicit none integer(pInt), intent(in) :: myFile @@ -132,20 +123,6 @@ subroutine constitutive_none_init(myFile) constitutive_none_sizePostResults = 0_pInt allocate(constitutive_none_structureName(maxNinstance)) constitutive_none_structureName = '' - allocate(constitutive_none_structure(maxNinstance)) - constitutive_none_structure = 0_pInt - allocate(constitutive_none_CoverA(maxNinstance)) - constitutive_none_CoverA = 0.0_pReal - allocate(constitutive_none_C11(maxNinstance)) - constitutive_none_C11 = 0.0_pReal - allocate(constitutive_none_C12(maxNinstance)) - constitutive_none_C12 = 0.0_pReal - allocate(constitutive_none_C13(maxNinstance)) - constitutive_none_C13 = 0.0_pReal - allocate(constitutive_none_C33(maxNinstance)) - constitutive_none_C33 = 0.0_pReal - allocate(constitutive_none_C44(maxNinstance)) - constitutive_none_C44 = 0.0_pReal allocate(constitutive_none_Cslip_66(6,6,maxNinstance)) constitutive_none_Cslip_66 = 0.0_pReal @@ -172,60 +149,40 @@ subroutine constitutive_none_init(myFile) cycle case ('lattice_structure') constitutive_none_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) - case ('covera_ratio') - constitutive_none_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_none_C11(i) = IO_floatValue(line,positions,2_pInt) + constitutive_none_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_none_C12(i) = IO_floatValue(line,positions,2_pInt) + constitutive_none_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_none_C13(i) = IO_floatValue(line,positions,2_pInt) + constitutive_none_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c22') + constitutive_none_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c23') + constitutive_none_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_none_C33(i) = IO_floatValue(line,positions,2_pInt) + constitutive_none_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_none_C44(i) = IO_floatValue(line,positions,2_pInt) + constitutive_none_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + case ('c55') + constitutive_none_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + case ('c66') + constitutive_none_Cslip_66(6,6,i) = IO_floatValue(line,positions,2_pInt) case default call IO_error(210_pInt,ext_msg=tag//' ('//constitutive_none_label//')') end select endif enddo -100 do i = 1_pInt,maxNinstance ! sanity checks - constitutive_none_structure(i) = lattice_initializeStructure(constitutive_none_structureName(i), & ! get structure - constitutive_none_CoverA(i)) - myStructure = constitutive_none_structure(i) - if (myStructure < 1_pInt) call IO_error(205_pInt,e=i) +100 do i = 1_pInt,maxNinstance + if (constitutive_none_structureName(i) == '') call IO_error(205_pInt,e=i) enddo do i = 1_pInt,maxNinstance constitutive_none_sizeDotState(i) = 1_pInt constitutive_none_sizeState(i) = 1_pInt - myStructure = constitutive_none_structure(i) - - select case (lattice_symmetryType(myStructure)) ! assign elasticity tensor - case(1_pInt) ! cubic(s) - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) & - constitutive_none_Cslip_66(k,j,i) = constitutive_none_C12(i) - constitutive_none_Cslip_66(k,k,i) = constitutive_none_C11(i) - constitutive_none_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_none_C44(i) - end forall - case(2_pInt) ! hex - constitutive_none_Cslip_66(1,1,i) = constitutive_none_C11(i) - constitutive_none_Cslip_66(2,2,i) = constitutive_none_C11(i) - constitutive_none_Cslip_66(3,3,i) = constitutive_none_C33(i) - constitutive_none_Cslip_66(1,2,i) = constitutive_none_C12(i) - constitutive_none_Cslip_66(2,1,i) = constitutive_none_C12(i) - constitutive_none_Cslip_66(1,3,i) = constitutive_none_C13(i) - constitutive_none_Cslip_66(3,1,i) = constitutive_none_C13(i) - constitutive_none_Cslip_66(2,3,i) = constitutive_none_C13(i) - constitutive_none_Cslip_66(3,2,i) = constitutive_none_C13(i) - constitutive_none_Cslip_66(4,4,i) = constitutive_none_C44(i) - constitutive_none_Cslip_66(5,5,i) = constitutive_none_C44(i) - constitutive_none_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_none_C11(i)- & - constitutive_none_C12(i)) - end select + constitutive_none_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_none_structureName(i),& + constitutive_none_Cslip_66) constitutive_none_Cslip_66(:,:,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_none_Cslip_66(:,:,i))) diff --git a/code/constitutive_nonlocal.f90 b/code/constitutive_nonlocal.f90 index 679900624..06fc31c3e 100644 --- a/code/constitutive_nonlocal.f90 +++ b/code/constitutive_nonlocal.f90 @@ -87,7 +87,7 @@ constitutive_nonlocal_output ! name of e integer(pInt), dimension(:), allocatable, private :: & constitutive_nonlocal_Noutput ! number of outputs per instance of this plasticity -character(len=32), dimension(:), allocatable, private :: & +character(len=32), dimension(:), allocatable, public :: & constitutive_nonlocal_structureName ! name of the lattice structure integer(pInt), dimension(:), allocatable, public :: & @@ -104,11 +104,6 @@ constitutive_nonlocal_colinearSystem ! colinear real(pReal), dimension(:), allocatable, private :: & constitutive_nonlocal_CoverA, & ! c/a ratio for hex type lattice -constitutive_nonlocal_C11, & ! C11 element in elasticity matrix -constitutive_nonlocal_C12, & ! C12 element in elasticity matrix -constitutive_nonlocal_C13, & ! C13 element in elasticity matrix -constitutive_nonlocal_C33, & ! C33 element in elasticity matrix -constitutive_nonlocal_C44, & ! C44 element in elasticity matrix constitutive_nonlocal_Gmod, & ! shear modulus constitutive_nonlocal_nu, & ! poisson's ratio constitutive_nonlocal_atomicVolume, & ! atomic volume @@ -235,6 +230,7 @@ use lattice, only: lattice_maxNslipFamily, & lattice_maxNinteraction, & lattice_NslipSystem, & lattice_initializeStructure, & + lattice_symmetrizeC66, & lattice_sd, & lattice_sn, & lattice_st, & @@ -312,11 +308,6 @@ constitutive_nonlocal_slipSystemLattice = 0_pInt constitutive_nonlocal_totalNslip = 0_pInt allocate(constitutive_nonlocal_CoverA(maxNinstance)) -allocate(constitutive_nonlocal_C11(maxNinstance)) -allocate(constitutive_nonlocal_C12(maxNinstance)) -allocate(constitutive_nonlocal_C13(maxNinstance)) -allocate(constitutive_nonlocal_C33(maxNinstance)) -allocate(constitutive_nonlocal_C44(maxNinstance)) allocate(constitutive_nonlocal_Gmod(maxNinstance)) allocate(constitutive_nonlocal_nu(maxNinstance)) allocate(constitutive_nonlocal_atomicVolume(maxNinstance)) @@ -350,11 +341,6 @@ allocate(constitutive_nonlocal_fEdgeMultiplication(maxNinstance)) allocate(constitutive_nonlocal_linetensionEffect(maxNinstance)) allocate(constitutive_nonlocal_edgeJogFactor(maxNinstance)) constitutive_nonlocal_CoverA = 0.0_pReal -constitutive_nonlocal_C11 = 0.0_pReal -constitutive_nonlocal_C12 = 0.0_pReal -constitutive_nonlocal_C13 = 0.0_pReal -constitutive_nonlocal_C33 = 0.0_pReal -constitutive_nonlocal_C44 = 0.0_pReal constitutive_nonlocal_Gmod = 0.0_pReal constitutive_nonlocal_atomicVolume = 0.0_pReal constitutive_nonlocal_Dsd0 = -1.0_pReal @@ -444,16 +430,24 @@ do constitutive_nonlocal_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('c/a_ratio','covera_ratio') constitutive_nonlocal_CoverA(i) = IO_floatValue(line,positions,2_pInt) - case ('c11') - constitutive_nonlocal_C11(i) = IO_floatValue(line,positions,2_pInt) - case ('c12') - constitutive_nonlocal_C12(i) = IO_floatValue(line,positions,2_pInt) - case ('c13') - constitutive_nonlocal_C13(i) = IO_floatValue(line,positions,2_pInt) - case ('c33') - constitutive_nonlocal_C33(i) = IO_floatValue(line,positions,2_pInt) - case ('c44') - constitutive_nonlocal_C44(i) = IO_floatValue(line,positions,2_pInt) + case ('c11') + constitutive_nonlocal_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) + case ('c12') + constitutive_nonlocal_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c13') + constitutive_nonlocal_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c22') + constitutive_nonlocal_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c23') + constitutive_nonlocal_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c33') + constitutive_nonlocal_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c44') + constitutive_nonlocal_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + case ('c55') + constitutive_nonlocal_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + case ('c66') + constitutive_nonlocal_Cslip_66(6,6,i) = IO_floatValue(line,positions,2_pInt) case ('nslip') forall (f = 1_pInt:lattice_maxNslipFamily) & constitutive_nonlocal_Nslip(f,i) = IO_intValue(line,positions,1_pInt+f) @@ -850,36 +844,17 @@ do i = 1,maxNinstance !*** elasticity matrix and shear modulus according to material.config - select case (myStructure) - case(1_pInt:2_pInt) ! cubic(s) - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) constitutive_nonlocal_Cslip_66(k,j,i) = constitutive_nonlocal_C12(i) - constitutive_nonlocal_Cslip_66(k,k,i) = constitutive_nonlocal_C11(i) - constitutive_nonlocal_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_nonlocal_C44(i) - end forall - case(3_pInt:) ! all hex - constitutive_nonlocal_Cslip_66(1,1,i) = constitutive_nonlocal_C11(i) - constitutive_nonlocal_Cslip_66(2,2,i) = constitutive_nonlocal_C11(i) - constitutive_nonlocal_Cslip_66(3,3,i) = constitutive_nonlocal_C33(i) - constitutive_nonlocal_Cslip_66(1,2,i) = constitutive_nonlocal_C12(i) - constitutive_nonlocal_Cslip_66(2,1,i) = constitutive_nonlocal_C12(i) - constitutive_nonlocal_Cslip_66(1,3,i) = constitutive_nonlocal_C13(i) - constitutive_nonlocal_Cslip_66(3,1,i) = constitutive_nonlocal_C13(i) - constitutive_nonlocal_Cslip_66(2,3,i) = constitutive_nonlocal_C13(i) - constitutive_nonlocal_Cslip_66(3,2,i) = constitutive_nonlocal_C13(i) - constitutive_nonlocal_Cslip_66(4,4,i) = constitutive_nonlocal_C44(i) - constitutive_nonlocal_Cslip_66(5,5,i) = constitutive_nonlocal_C44(i) - constitutive_nonlocal_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_nonlocal_C11(i)- constitutive_nonlocal_C12(i)) - end select + constitutive_nonlocal_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_nonlocal_structureName(i),& + constitutive_nonlocal_Cslip_66) + constitutive_nonlocal_Gmod(i) = 0.2_pReal * ( constitutive_nonlocal_Cslip_66(1,1,i) - constitutive_nonlocal_Cslip_66(1,2,i) & + + 3.0_pReal*constitutive_nonlocal_Cslip_66(4,4,i) ) ! (C11iso-C12iso)/2 with C11iso=(3*C11+2*C12+4*C44)/5 and C12iso=(C11+4*C12-2*C44)/5 + constitutive_nonlocal_nu(i) = ( constitutive_nonlocal_Cslip_66(1,1,i) + 4.0_pReal*constitutive_nonlocal_Cslip_66(1,2,i) & + - 2.0_pReal*constitutive_nonlocal_Cslip_66(1,2,i) ) & + / ( 4.0_pReal*constitutive_nonlocal_Cslip_66(1,1,i) + 6.0_pReal*constitutive_nonlocal_Cslip_66(1,2,i) & + + 2.0_pReal*constitutive_nonlocal_Cslip_66(4,4,i) ) constitutive_nonlocal_Cslip_66(1:6,1:6,i) = math_Mandel3333to66(math_Voigt66to3333(constitutive_nonlocal_Cslip_66(1:6,1:6,i))) constitutive_nonlocal_Cslip_3333(1:3,1:3,1:3,1:3,i) = math_Voigt66to3333(constitutive_nonlocal_Cslip_66(1:6,1:6,i)) - - constitutive_nonlocal_Gmod(i) = 0.2_pReal * ( constitutive_nonlocal_C11(i) - constitutive_nonlocal_C12(i) & - + 3.0_pReal*constitutive_nonlocal_C44(i) ) ! (C11iso-C12iso)/2 with C11iso=(3*C11+2*C12+4*C44)/5 and C12iso=(C11+4*C12-2*C44)/5 - constitutive_nonlocal_nu(i) = ( constitutive_nonlocal_C11(i) + 4.0_pReal*constitutive_nonlocal_C12(i) & - - 2.0_pReal*constitutive_nonlocal_C44(i) ) & - / ( 4.0_pReal*constitutive_nonlocal_C11(i) + 6.0_pReal*constitutive_nonlocal_C12(i) & - + 2.0_pReal*constitutive_nonlocal_C44(i) ) ! C12iso/(C11iso+C12iso) with C11iso=(3*C11+2*C12+4*C44)/5 and C12iso=(C11+4*C12-2*C44)/5 + ! C12iso/(C11iso+C12iso) with C11iso=(3*C11+2*C12+4*C44)/5 and C12iso=(C11+4*C12-2*C44)/5 do s1 = 1_pInt,ns f = constitutive_nonlocal_slipFamily(s1,i) diff --git a/code/constitutive_phenopowerlaw.f90 b/code/constitutive_phenopowerlaw.f90 index 3599a479a..6b3cadc07 100644 --- a/code/constitutive_phenopowerlaw.f90 +++ b/code/constitutive_phenopowerlaw.f90 @@ -53,16 +53,11 @@ module constitutive_phenopowerlaw character(len=64), dimension(:,:), allocatable, target, public :: & constitutive_phenopowerlaw_output !< name of each post result output - character(len=32), dimension(:), allocatable, private :: & + character(len=32), dimension(:), allocatable, public :: & constitutive_phenopowerlaw_structureName real(pReal), dimension(:), allocatable, private :: & constitutive_phenopowerlaw_CoverA, & !< c/a of the crystal (input parameter) - constitutive_phenopowerlaw_C11, & !< component 11 of the stiffness matrix (input parameter) - constitutive_phenopowerlaw_C12, & !< component 12 of the stiffness matrix (input parameter) - constitutive_phenopowerlaw_C13, & !< component 13 of the stiffness matrix (input parameter) - constitutive_phenopowerlaw_C33, & !< component 33 of the stiffness matrix (input parameter) - constitutive_phenopowerlaw_C44, & !< component 44 of the stiffness matrix (input parameter) constitutive_phenopowerlaw_gdot0_slip, & !< reference shear strain rate for slip (input parameter) constitutive_phenopowerlaw_gdot0_twin, & !< reference shear strain rate for twin (input parameter) constitutive_phenopowerlaw_n_slip, & !< stress exponent for slip (input parameter) @@ -127,7 +122,7 @@ subroutine constitutive_phenopowerlaw_init(myFile) use debug, only: debug_level,& debug_constitutive,& debug_levelBasic - use lattice, only: lattice_initializeStructure, lattice_symmetryType, & + use lattice, only: lattice_initializeStructure, lattice_symmetrizeC66, & lattice_maxNslipFamily, lattice_maxNtwinFamily, & lattice_maxNinteraction, lattice_NslipSystem, lattice_NtwinSystem, & lattice_interactionSlipSlip, & @@ -182,16 +177,6 @@ subroutine constitutive_phenopowerlaw_init(myFile) constitutive_phenopowerlaw_totalNtwin = 0_pInt allocate(constitutive_phenopowerlaw_CoverA(maxNinstance)) constitutive_phenopowerlaw_CoverA = 0.0_pReal - allocate(constitutive_phenopowerlaw_C11(maxNinstance)) - constitutive_phenopowerlaw_C11 = 0.0_pReal - allocate(constitutive_phenopowerlaw_C12(maxNinstance)) - constitutive_phenopowerlaw_C12 = 0.0_pReal - allocate(constitutive_phenopowerlaw_C13(maxNinstance)) - constitutive_phenopowerlaw_C13 = 0.0_pReal - allocate(constitutive_phenopowerlaw_C33(maxNinstance)) - constitutive_phenopowerlaw_C33 = 0.0_pReal - allocate(constitutive_phenopowerlaw_C44(maxNinstance)) - constitutive_phenopowerlaw_C44 = 0.0_pReal allocate(constitutive_phenopowerlaw_Cslip_66(6,6,maxNinstance)) constitutive_phenopowerlaw_Cslip_66 = 0.0_pReal allocate(constitutive_phenopowerlaw_gdot0_slip(maxNinstance)) @@ -274,15 +259,23 @@ subroutine constitutive_phenopowerlaw_init(myFile) case ('covera_ratio') constitutive_phenopowerlaw_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_phenopowerlaw_C11(i) = IO_floatValue(line,positions,2_pInt) + constitutive_phenopowerlaw_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_phenopowerlaw_C12(i) = IO_floatValue(line,positions,2_pInt) + constitutive_phenopowerlaw_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_phenopowerlaw_C13(i) = IO_floatValue(line,positions,2_pInt) + constitutive_phenopowerlaw_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + case ('c22') + constitutive_phenopowerlaw_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + case ('c23') + constitutive_phenopowerlaw_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_phenopowerlaw_C33(i) = IO_floatValue(line,positions,2_pInt) + constitutive_phenopowerlaw_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_phenopowerlaw_C44(i) = IO_floatValue(line,positions,2_pInt) + constitutive_phenopowerlaw_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + case ('c55') + constitutive_phenopowerlaw_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + case ('c66') + constitutive_phenopowerlaw_Cslip_66(6,6,i) = IO_floatValue(line,positions,2_pInt) case ('nslip') forall (j = 1_pInt:lattice_maxNslipFamily)& constitutive_phenopowerlaw_Nslip(j,i) = IO_intValue(line,positions,1_pInt+j) @@ -448,29 +441,9 @@ subroutine constitutive_phenopowerlaw_init(myFile) myStructure = constitutive_phenopowerlaw_structure(i) - select case (lattice_symmetryType(myStructure)) ! assign elasticity tensor - case(1_pInt) ! cubic(s) - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) & - constitutive_phenopowerlaw_Cslip_66(k,j,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(k,k,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_phenopowerlaw_C44(i) - end forall - case(2_pInt) ! hex - constitutive_phenopowerlaw_Cslip_66(1,1,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(2,2,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(3,3,i) = constitutive_phenopowerlaw_C33(i) - constitutive_phenopowerlaw_Cslip_66(1,2,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(2,1,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(1,3,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(3,1,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(2,3,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(3,2,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(4,4,i) = constitutive_phenopowerlaw_C44(i) - constitutive_phenopowerlaw_Cslip_66(5,5,i) = constitutive_phenopowerlaw_C44(i) - constitutive_phenopowerlaw_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_phenopowerlaw_C11(i)- & - constitutive_phenopowerlaw_C12(i)) - end select + constitutive_phenopowerlaw_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_phenopowerlaw_structureName(i),& + constitutive_phenopowerlaw_Cslip_66) + ! assign elasticity tensor constitutive_phenopowerlaw_Cslip_66(:,:,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_phenopowerlaw_Cslip_66(:,:,i))) diff --git a/code/constitutive_titanmod.f90 b/code/constitutive_titanmod.f90 index a87dfb032..83906f4fc 100644 --- a/code/constitutive_titanmod.f90 +++ b/code/constitutive_titanmod.f90 @@ -114,11 +114,6 @@ integer(pInt), dimension(:,:), allocatable :: & real(pReal), dimension(:), allocatable :: & constitutive_titanmod_CoverA, & ! c/a ratio for hex type lattice - constitutive_titanmod_C11, & ! C11 element in elasticity matrix - constitutive_titanmod_C12, & ! C12 element in elasticity matrix - constitutive_titanmod_C13, & ! C13 element in elasticity matrix - constitutive_titanmod_C33, & ! C33 element in elasticity matrix - constitutive_titanmod_C44, & ! C44 element in elasticity matrix constitutive_titanmod_debyefrequency, & !Debye frequency constitutive_titanmod_kinkf0, & !Debye frequency constitutive_titanmod_Gmod, & ! shear modulus @@ -298,16 +293,6 @@ allocate(constitutive_titanmod_totalNtwin(maxNinstance)) constitutive_titanmod_totalNtwin = 0_pInt allocate(constitutive_titanmod_CoverA(maxNinstance)) constitutive_titanmod_CoverA = 0.0_pReal -allocate(constitutive_titanmod_C11(maxNinstance)) - constitutive_titanmod_C11 = 0.0_pReal -allocate(constitutive_titanmod_C12(maxNinstance)) - constitutive_titanmod_C12 = 0.0_pReal -allocate(constitutive_titanmod_C13(maxNinstance)) - constitutive_titanmod_C13 = 0.0_pReal -allocate(constitutive_titanmod_C33(maxNinstance)) - constitutive_titanmod_C33 = 0.0_pReal -allocate(constitutive_titanmod_C44(maxNinstance)) - constitutive_titanmod_C44 = 0.0_pReal allocate(constitutive_titanmod_debyefrequency(maxNinstance)) constitutive_titanmod_debyefrequency = 0.0_pReal allocate(constitutive_titanmod_kinkf0(maxNinstance)) @@ -446,20 +431,32 @@ enddo constitutive_titanmod_CoverA(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_CoverA(i) case ('c11') - constitutive_titanmod_C11(i) = IO_floatValue(line,positions,2_pInt) - write(6,*) tag,constitutive_titanmod_C11(i) + constitutive_titanmod_Cslip_66(1,1,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(1,1,i) case ('c12') - constitutive_titanmod_C12(i) = IO_floatValue(line,positions,2_pInt) - write(6,*) tag,constitutive_titanmod_C12(i) + constitutive_titanmod_Cslip_66(1,2,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(1,1,i) case ('c13') - constitutive_titanmod_C13(i) = IO_floatValue(line,positions,2_pInt) - write(6,*) tag,constitutive_titanmod_C13(i) + constitutive_titanmod_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(1,1,i) + case ('c22') + constitutive_titanmod_Cslip_66(2,2,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(2,2,i) + case ('c23') + constitutive_titanmod_Cslip_66(2,3,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(2,3,i) case ('c33') - constitutive_titanmod_C33(i) = IO_floatValue(line,positions,2_pInt) - write(6,*) tag,constitutive_titanmod_C33(i) + constitutive_titanmod_Cslip_66(3,3,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(3,3,i) case ('c44') - constitutive_titanmod_C44(i) = IO_floatValue(line,positions,2_pInt) - write(6,*) tag,constitutive_titanmod_C44(i) + constitutive_titanmod_Cslip_66(4,4,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(4,4,i) + case ('c55') + constitutive_titanmod_Cslip_66(5,5,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(5,5,i) + case ('c66') + constitutive_titanmod_Cslip_66(1,3,i) = IO_floatValue(line,positions,2_pInt) + write(6,*) tag,constitutive_titanmod_Cslip_66(6,6,i) case ('debyefrequency') constitutive_titanmod_debyefrequency(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_debyefrequency(i) @@ -870,32 +867,12 @@ do i = 1_pInt,maxNinstance write(6,*) 'Determining elasticity matrix' !* Elasticity matrix and shear modulus according to material.config - select case (myStructure) - case(1_pInt:2_pInt) ! cubic(s) - forall(k=1_pInt:3_pInt) - forall(j=1_pInt:3_pInt) & - constitutive_titanmod_Cslip_66(k,j,i) = constitutive_titanmod_C12(i) - constitutive_titanmod_Cslip_66(k,k,i) = constitutive_titanmod_C11(i) - constitutive_titanmod_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_titanmod_C44(i) - end forall - case(3_pInt:) ! all hex - constitutive_titanmod_Cslip_66(1,1,i) = constitutive_titanmod_C11(i) - constitutive_titanmod_Cslip_66(2,2,i) = constitutive_titanmod_C11(i) - constitutive_titanmod_Cslip_66(3,3,i) = constitutive_titanmod_C33(i) - constitutive_titanmod_Cslip_66(1,2,i) = constitutive_titanmod_C12(i) - constitutive_titanmod_Cslip_66(2,1,i) = constitutive_titanmod_C12(i) - constitutive_titanmod_Cslip_66(1,3,i) = constitutive_titanmod_C13(i) - constitutive_titanmod_Cslip_66(3,1,i) = constitutive_titanmod_C13(i) - constitutive_titanmod_Cslip_66(2,3,i) = constitutive_titanmod_C13(i) - constitutive_titanmod_Cslip_66(3,2,i) = constitutive_titanmod_C13(i) - constitutive_titanmod_Cslip_66(4,4,i) = constitutive_titanmod_C44(i) - constitutive_titanmod_Cslip_66(5,5,i) = constitutive_titanmod_C44(i) - constitutive_titanmod_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_titanmod_C11(i)-constitutive_titanmod_C12(i)) - end select + constitutive_titanmod_Cslip_66(:,:,i) = lattice_symmetrizeC66(constitutive_titanmod_structureName(i),& + constitutive_titanmod_Cslip_66) + constitutive_titanmod_Gmod(i) = & + 0.2_pReal*(constitutive_titanmod_Cslip_66(1,1,i)-constitutive_titanmod_Cslip_66(1,2,i))+0.3_pReal*constitutive_titanmod_Cslip_66(4,4,i) constitutive_titanmod_Cslip_66(:,:,i) = math_Mandel3333to66(math_Voigt66to3333(constitutive_titanmod_Cslip_66(:,:,i))) constitutive_titanmod_Cslip_3333(:,:,:,:,i) = math_Voigt66to3333(constitutive_titanmod_Cslip_66(:,:,i)) - constitutive_titanmod_Gmod(i) = & - 0.2_pReal*(constitutive_titanmod_C11(i)-constitutive_titanmod_C12(i))+0.3_pReal*constitutive_titanmod_C44(i) !* Construction of the twin elasticity matrices do j=1_pInt,lattice_maxNtwinFamily diff --git a/code/crystallite.f90 b/code/crystallite.f90 index 930e52c27..0598a27f7 100644 --- a/code/crystallite.f90 +++ b/code/crystallite.f90 @@ -139,13 +139,13 @@ subroutine crystallite_init(Temperature) use constitutive, only: constitutive_microstructure use constitutive_phenopowerlaw, only: constitutive_phenopowerlaw_label, & - constitutive_phenopowerlaw_structure + constitutive_phenopowerlaw_structureName use constitutive_titanmod, only: constitutive_titanmod_label, & - constitutive_titanmod_structure + constitutive_titanmod_structureName use constitutive_dislotwin, only: constitutive_dislotwin_label, & - constitutive_dislotwin_structure + constitutive_dislotwin_structureName use constitutive_nonlocal, only: constitutive_nonlocal_label, & - constitutive_nonlocal_structure + constitutive_nonlocal_structureName implicit none integer(pInt), parameter :: myFile = 200_pInt, & @@ -169,7 +169,6 @@ subroutine crystallite_init(Temperature) p, & output, & mySize, & - myStructure, & ! lattice structure myPhase, & myMat character(len=64) tag @@ -351,19 +350,16 @@ do e = FEsolving_execElem(1),FEsolving_execElem(2) myMat = phase_plasticityInstance(myPhase) select case (phase_plasticity(myPhase)) case (constitutive_phenopowerlaw_label) - myStructure = constitutive_phenopowerlaw_structure(myMat) + crystallite_symmetryID(g,i,e) = lattice_symmetryType(constitutive_phenopowerlaw_structureName(myMat)) case (constitutive_titanmod_label) - myStructure = constitutive_titanmod_structure(myMat) + crystallite_symmetryID(g,i,e) = lattice_symmetryType(constitutive_titanmod_structureName(myMat)) case (constitutive_dislotwin_label) - myStructure = constitutive_dislotwin_structure(myMat) + crystallite_symmetryID(g,i,e) = lattice_symmetryType(constitutive_dislotwin_structureName(myMat)) case (constitutive_nonlocal_label) - myStructure = constitutive_nonlocal_structure(myMat) + crystallite_symmetryID(g,i,e) = lattice_symmetryType(constitutive_nonlocal_structureName(myMat)) case default - myStructure = -1_pInt ! does this happen for j2 material? + crystallite_symmetryID(g,i,e) = 0_pInt ! does this happen for j2 material? end select - if (myStructure > 0_pInt) then - crystallite_symmetryID(g,i,e) = lattice_symmetryType(myStructure) ! structure = 1(fcc) or 2(bcc) => 1; 3(hex)=>2 - endif enddo enddo enddo diff --git a/code/lattice.f90 b/code/lattice.f90 index 9fde22baa..6d9bb1bf9 100644 --- a/code/lattice.f90 +++ b/code/lattice.f90 @@ -681,7 +681,8 @@ module lattice public :: & lattice_init, & lattice_initializeStructure, & - lattice_symmetryType + lattice_symmetryType, & + lattice_symmetrizeC66 contains @@ -690,15 +691,15 @@ contains !> @brief Maps structure to symmetry type !> @details fcc(1) and bcc(2) are cubic(1) hex(3+) is hexagonal(2) !-------------------------------------------------------------------------------------------------- -integer(pInt) pure function lattice_symmetryType(structID) +integer(pInt) pure function lattice_symmetryType(structName) implicit none - integer(pInt), intent(in) :: structID + character(len=32), intent(in) :: structName - select case(structID) - case (1_pInt,2_pInt) + select case(structName(1:3)) + case ('fcc','bcc') lattice_symmetryType = 1_pInt - case (3_pInt:) + case ('hex') lattice_symmetryType = 2_pInt case default lattice_symmetryType = 0_pInt @@ -709,6 +710,70 @@ integer(pInt) pure function lattice_symmetryType(structID) end function lattice_symmetryType +!-------------------------------------------------------------------------------------------------- +!> @brief Symmetrizes stiffness matrix according to lattice type +!-------------------------------------------------------------------------------------------------- +pure function lattice_symmetrizeC66(structName,C66) + + implicit none + + character(len=32), intent(in) :: structName + real(pReal), dimension(6,6), intent(in) :: C66 + real(pReal), dimension(6,6) :: lattice_symmetrizeC66 + integer(pInt) :: j,k + + lattice_symmetrizeC66 = 0.0_pReal + + select case(structName(1:3)) + case ('iso') + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) + lattice_symmetrizeC66(k,j) = C66(1,2) + end forall + lattice_symmetrizeC66(k,k) = C66(1,1) + lattice_symmetrizeC66(k+3,k+3) = 0.5_pReal*(C66(1,1)-C66(1,2)) + end forall + case ('fcc','bcc') + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) + lattice_symmetrizeC66(k,j) = C66(1,2) + lattice_symmetrizeC66(k,k) = C66(1,1) + lattice_symmetrizeC66(k+3_pInt,k+3_pInt) = C66(4,4) + end forall + end forall + case ('hex') + lattice_symmetrizeC66(1,1) = C66(1,1) + lattice_symmetrizeC66(2,2) = C66(1,1) + lattice_symmetrizeC66(3,3) = C66(3,3) + lattice_symmetrizeC66(1,2) = C66(1,2) + lattice_symmetrizeC66(2,1) = C66(1,2) + lattice_symmetrizeC66(1,3) = C66(1,3) + lattice_symmetrizeC66(3,1) = C66(1,3) + lattice_symmetrizeC66(2,3) = C66(1,3) + lattice_symmetrizeC66(3,2) = C66(1,3) + lattice_symmetrizeC66(4,4) = C66(4,4) + lattice_symmetrizeC66(5,5) = C66(4,4) + lattice_symmetrizeC66(6,6) = 0.5_pReal*(C66(1,1)-C66(1,2)) + case ('ort') + lattice_symmetrizeC66(1,1) = C66(1,1) + lattice_symmetrizeC66(2,2) = C66(2,2) + lattice_symmetrizeC66(3,3) = C66(3,3) + lattice_symmetrizeC66(1,2) = C66(1,2) + lattice_symmetrizeC66(2,1) = C66(1,2) + lattice_symmetrizeC66(1,3) = C66(1,3) + lattice_symmetrizeC66(3,1) = C66(1,3) + lattice_symmetrizeC66(2,3) = C66(2,3) + lattice_symmetrizeC66(3,2) = C66(2,3) + lattice_symmetrizeC66(4,4) = C66(4,4) + lattice_symmetrizeC66(5,5) = C66(5,5) + lattice_symmetrizeC66(6,6) = C66(6,6) + end select + + return + + end function lattice_symmetrizeC66 + + !-------------------------------------------------------------------------------------------------- !> @brief Module initialization !--------------------------------------------------------------------------------------------------