From e18f304c5a4b348b39d1e00ec2031feb91c973d4 Mon Sep 17 00:00:00 2001 From: Luc Hantcherli Date: Fri, 7 Dec 2007 12:36:17 +0000 Subject: [PATCH] Crystal.f90 contains now the twin systems (at least for fcc, not implemented for bcc and hcp). Schmid and transformations matrices for twin systemss are defined and constructed. Twin systems are specified using the type of variables as slip systems: "MaxNTwin", "MaxMaxNtwin" ... --- trunk/crystal.f90 | 71 ++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 70 insertions(+), 1 deletion(-) diff --git a/trunk/crystal.f90 b/trunk/crystal.f90 index 495517863..197300132 100644 --- a/trunk/crystal.f90 +++ b/trunk/crystal.f90 @@ -23,14 +23,28 @@ integer(pInt), parameter :: crystal_MaxCrystalStructure = 3 !* (has to be changed according the definition of slip systems) integer(pInt), dimension(crystal_MaxCrystalStructure), parameter :: crystal_MaxNslipOfStructure = & reshape((/12,48,12/),(/crystal_MaxCrystalStructure/)) +!* Total number of twin systems per crystal structure +!* (has to be changed according the definition of twin systems) +integer(pInt), dimension(crystal_MaxCrystalStructure), parameter :: crystal_MaxNtwinOfStructure = & +reshape((/12,12,6/),(/crystal_MaxCrystalStructure/)) !* Maximum number of slip systems over crystal structures integer(pInt), parameter :: crystal_MaxMaxNslipOfStructure = 48 +!* Maximum number of twin systems over crystal structures +integer(pInt), parameter :: crystal_MaxMaxNtwinOfStructure = 12 !* Slip direction, slip normales and Schmid matrices real(pReal), dimension(3,3,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: crystal_Sslip real(pReal), dimension(6,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: crystal_Sslip_v real(pReal), dimension(3,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: crystal_sn real(pReal), dimension(3,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: crystal_sd real(pReal), dimension(3,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: crystal_st +!* twin direction, twin normales, Schmid matrices and transformation matrices +real(pReal), dimension(3,3,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_Stwin +real(pReal), dimension(6,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_Stwin_v +real(pReal), dimension(3,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_tn +real(pReal), dimension(3,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_td +real(pReal), dimension(3,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_tt +real(pReal), dimension(3,3,crystal_MaxMaxNtwinOfStructure,crystal_MaxCrystalStructure) :: crystal_Qtwin + !* Slip_slip interaction matrices integer(pInt), dimension(crystal_MaxMaxNslipOfStructure,crystal_MaxMaxNslipOfStructure,crystal_MaxCrystalStructure) :: & crystal_SlipIntType @@ -50,6 +64,21 @@ data crystal_sd(:,10,1)/ 0, 1, 1/ ; data crystal_sn(:,10,1)/-1, 1,-1/ data crystal_sd(:,11,1)/ 1, 0,-1/ ; data crystal_sn(:,11,1)/-1, 1,-1/ data crystal_sd(:,12,1)/-1,-1, 0/ ; data crystal_sn(:,12,1)/-1, 1,-1/ +!*** Twin systems for FCC structures (1) *** +!* System {111}<112> Sort according Eisenlohr&Hantcherli +data crystal_td(:, 1,1)/-2, 1, 1/ ; data crystal_tn(:, 1,1)/ 1, 1, 1/ +data crystal_td(:, 2,1)/ 1,-2, 1/ ; data crystal_tn(:, 2,1)/ 1, 1, 1/ +data crystal_td(:, 3,1)/ 1, 1,-2/ ; data crystal_tn(:, 3,1)/ 1, 1, 1/ +data crystal_td(:, 4,1)/ 2,-1, 1/ ; data crystal_tn(:, 4,1)/-1,-1, 1/ +data crystal_td(:, 5,1)/-1, 2, 1/ ; data crystal_tn(:, 5,1)/-1,-1, 1/ +data crystal_td(:, 6,1)/-1,-1,-2/ ; data crystal_tn(:, 6,1)/-1,-1, 1/ +data crystal_td(:, 7,1)/-2,-1,-1/ ; data crystal_tn(:, 7,1)/ 1,-1,-1/ +data crystal_td(:, 8,1)/ 1, 2,-1/ ; data crystal_tn(:, 8,1)/ 1,-1,-1/ +data crystal_td(:, 9,1)/ 1,-1, 2/ ; data crystal_tn(:, 9,1)/ 1,-1,-1/ +data crystal_td(:,10,1)/ 2, 1,-1/ ; data crystal_tn(:,10,1)/-1, 1,-1/ +data crystal_td(:,11,1)/-1,-2,-1/ ; data crystal_tn(:,11,1)/-1, 1,-1/ +data crystal_td(:,12,1)/-1, 1, 2/ ; data crystal_tn(:,12,1)/-1, 1,-1/ + !*** Slip-Slip interactions for FCC structures (1) *** data crystal_SlipIntType( 1,1:crystal_MaxNslipOfStructure(1),1)/1,2,2,4,6,5,3,5,5,4,5,6/ data crystal_SlipIntType( 2,1:crystal_MaxNslipOfStructure(1),1)/2,1,2,6,4,5,5,4,6,5,3,5/ @@ -120,6 +149,11 @@ data crystal_sd(:,46,2)/ 1, 1, 1/ ; data crystal_sn(:,46,2)/-3, 2, 1/ data crystal_sd(:,47,2)/ 1, 1,-1/ ; data crystal_sn(:,47,2)/ 3,-2, 1/ data crystal_sd(:,48,2)/ 1,-1, 1/ ; data crystal_sn(:,48,2)/ 3, 2,-1/ +!*** Twin systems for BCC structures (2) *** +!* System {112}<111> +!* Sort? +!* Not implemented yet + !*** Slip-Slip interactions for BCC structures (2) *** data crystal_SlipIntType( 1,:,2)/1,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2/ data crystal_SlipIntType( 2,:,2)/2,1,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2/ @@ -207,6 +241,11 @@ data crystal_sd(:,10,3)/-1, 0, 0/ ; data crystal_sn(:,10,3)/ 1, 0, 1/ data crystal_sd(:,11,3)/ 0,-1, 0/ ; data crystal_sn(:,11,3)/-1, 1, 1/ data crystal_sd(:,12,3)/ 1, 1, 0/ ; data crystal_sn(:,12,3)/ 1,-1, 1/ +!*** Twin systems for HCP structures (2) *** +!* System {1012}<1011> +!* Sort? +!* Not implemented yet + !*** Slip-Slip interactions for HCP structures (3) *** data crystal_SlipIntType( 1,1:crystal_MaxNslipOfStructure(3),3)/1,2,2,2,2,2,2,2,2,2,2,2/ data crystal_SlipIntType( 2,1:crystal_MaxNslipOfStructure(3),3)/2,1,2,2,2,2,2,2,2,2,2,2/ @@ -242,6 +281,7 @@ subroutine crystal_SchmidMatrices() !* Calculation of Schmid matrices * !************************************** use prec, only: pReal,pInt +use math, only: math_identity2nd implicit none !* Definition of variables @@ -250,7 +290,7 @@ real(pReal) invNorm !* Iteration over the crystal structures do l=1,crystal_MaxCrystalStructure -!* Iteration over the systems +!* Iteration over the slip systems do k=1,crystal_MaxNslipOfStructure(l) !* Definition of transverse direction st for the frame (sd,st,sn) crystal_st(1,k,l)=crystal_sn(2,k,l)*crystal_sd(3,k,l)-crystal_sn(3,k,l)*crystal_sd(2,k,l) @@ -273,8 +313,37 @@ do l=1,crystal_MaxCrystalStructure crystal_Sslip_v(5,k,l)=(crystal_Sslip(2,3,k,l)+crystal_Sslip(3,2,k,l))/dsqrt(2.0_pReal) crystal_Sslip_v(6,k,l)=(crystal_Sslip(1,3,k,l)+crystal_Sslip(3,1,k,l))/dsqrt(2.0_pReal) enddo + +!* Iteration over the twin systems + do k=1,crystal_MaxNslipOfStructure(l) +!* Definition of transverse direction tt for the frame (td,tt,tn) + crystal_tt(1,k,l)=crystal_tn(2,k,l)*crystal_td(3,k,l)-crystal_tn(3,k,l)*crystal_td(2,k,l) + crystal_tt(2,k,l)=crystal_tn(3,k,l)*crystal_td(1,k,l)-crystal_tn(1,k,l)*crystal_td(3,k,l) + crystal_tt(3,k,l)=crystal_tn(1,k,l)*crystal_td(2,k,l)-crystal_tn(2,k,l)*crystal_td(1,k,l) +!* Defintion of Schmid matrix and transformation matrices + crystal_Qtwin(:,:,k,l)=-math_identity2nd(3) + forall (i=1:3,j=1:3) + crystal_Stwin(i,j,k,l)=crystal_td(i,k,l)*crystal_tn(j,k,l) + crystal_Qtwin(i,j,k,l)=crystal_Qtwin(i,j,k,l)+2*crystal_tn(i,k,l)*crystal_tn(j,k,l) + endforall +!* Normalization of Schmid matrix + invNorm=dsqrt(1.0_pReal/((crystal_tn(1,k,l)**2+crystal_tn(2,k,l)**2+crystal_tn(3,k,l)**2)*& + (crystal_td(1,k,l)**2+crystal_td(2,k,l)**2+crystal_td(3,k,l)**2))) + crystal_Stwin(:,:,k,l)=crystal_Stwin(:,:,k,l)*invNorm +!* Vectorization of normalized Schmid matrix + crystal_Stwin_v(1,k,l)=crystal_Stwin(1,1,k,l) + crystal_Stwin_v(2,k,l)=crystal_Stwin(2,2,k,l) + crystal_Stwin_v(3,k,l)=crystal_Stwin(3,3,k,l) + !* be compatible with Mandel notation of Tstar + crystal_Stwin_v(4,k,l)=(crystal_Stwin(1,2,k,l)+crystal_Stwin(2,1,k,l))/dsqrt(2.0_pReal) + crystal_Stwin_v(5,k,l)=(crystal_Stwin(2,3,k,l)+crystal_Stwin(3,2,k,l))/dsqrt(2.0_pReal) + crystal_Stwin_v(6,k,l)=(crystal_Stwin(1,3,k,l)+crystal_Stwin(3,1,k,l))/dsqrt(2.0_pReal) + enddo enddo end subroutine END MODULE + + +