1 /*! 2 * \file src/NUMODIS/AnalyseJunction.cxx 3 * \brief 4 * \author Laurent Dupuy 5 * \date 9/06/2017 6 * \copyright Copyright (C) 2006-2018 CEA/DEN, EDF R&D. All rights 7 * reserved. 8 * This project is publicly released under either the GNU GPL Licence 9 * or the CECILL-A licence. A copy of thoses licences are delivered 10 * with the sources of TFEL. CEA or EDF may also distribute this 11 * project under specific licensing conditions. 12 */ 13 14 #include<algorithm> 15 #include<iostream> 16 #include "NUMODIS/AnalyseJunction.hxx" 17 #include "NUMODIS/GSystem.hxx" 18 #include "NUMODIS/Crystallo.hxx" 19 #include "NUMODIS/TripleNode.hxx" 20 #include "NUMODIS/IsotropicLineTensionModel.hxx" 21 #include "NUMODIS/FrankRead.hxx" 22 23 namespace numodis 24 { 25 26 namespace Dupuy2017 27 { 28 //====================================================================================================================== 29 // precomputeJunctionParameters 30 //---------------------------------------------------------------------------------------------------------------------- 31 //! analyse the interaction between two glide systems and return the parameters to be used for a line tension analysis 32 //---------------------------------------------------------------------------------------------------------------------- 33 /*! 34 Considering two interacting glide systems (gsystem1 = studied system, gsystem2 = forest), pre-compute the parameters 35 that will be later used within the line tension model to estimate the strenght of the junction and the corresponding 36 hardening coefficient. 37 38 \param pcrystallo pointer to the relevant crystallography of the material 39 \param gsystem1 glide system under consideration 40 \param gsystem2 glide system of the forest 41 \param Burgers vectors of the gsystem1, gsystem2 and the junction 42 \param phi angles between the Burgers vectors and the junction direction 43 44 */ 45 //====================================================================================================================== precomputeJunctionParameters(const Crystallo * pcrystallo,const GSystem & gsystem1,const GSystem & gsystem2,std::vector<double> & burgers,std::vector<double> & phi)46 void precomputeJunctionParameters(const Crystallo* pcrystallo, 47 const GSystem& gsystem1, 48 const GSystem& gsystem2, 49 std::vector<double>& burgers, 50 std::vector<double>& phi) 51 { 52 53 // same glide plane? 54 if(gsystem1.getIPlane()==gsystem2.getIPlane()) 55 { 56 // same Burgers vector? 57 if(Coincide(gsystem1.getIBurgers(),gsystem2.getIBurgers())!=0) 58 { 59 std::cout << "self"; 60 } 61 else 62 { 63 std::cout << "coplanar"; 64 } 65 throw(-1); 66 } 67 else 68 { 69 // compute junction direction 70 IDirection ijunction(pcrystallo->getNindices()); 71 pcrystallo->CrossProduct(gsystem1.getIPlane(), 72 gsystem2.getIPlane(), 73 ijunction); 74 Vect3 xjunction(pcrystallo->direction(ijunction)); 75 76 // first Burgers vector 77 Vect3 xburgers1(pcrystallo->burgers_vector(gsystem1.getIBurgers())); 78 burgers[0] = xburgers1.Length(); 79 double cosphi1 = xburgers1.Dot(xjunction)/burgers[0]; 80 cosphi1 = std::max(std::min(cosphi1,1.0),-1.0); 81 phi[0] = acos(cosphi1); 82 83 // second Burgers vector 84 Vect3 xburgers2(pcrystallo->burgers_vector(gsystem2.getIBurgers())); 85 burgers[1] = xburgers2.Length(); 86 double cosphi2 = xburgers2.Dot(xjunction)/burgers[1]; 87 cosphi2 = std::max(std::min(cosphi2,1.0),-1.0); 88 phi[1] = acos(cosphi2); 89 90 // Burgers vector of the junction 91 burgers[2] = 0.0; 92 double cosphi3 = 0.0; 93 phi[2] = 0.0; 94 if(Coincide(gsystem1.getIBurgers(),gsystem2.getIBurgers())==0) 95 { 96 GSystem gsystem3(pcrystallo->ComputeJunctionGSystem(gsystem1,gsystem2)); 97 Vect3 xburgers3(pcrystallo->burgers_vector(gsystem3.getIBurgers())); 98 burgers[2] = xburgers3.Length(); 99 cosphi3 = xburgers3.Dot(xjunction)/burgers[2]; 100 cosphi3 = std::max(std::min(cosphi3,1.0),-1.0); 101 phi[2] = acos(cosphi3); 102 } 103 104 } 105 106 } 107 108 //==================================================================================== 109 // computeNormalizedStress 110 //------------------------------------------------------------------------------------ 111 //! compute the critical/normalized stress of a junction based on line tension 112 //------------------------------------------------------------------------------------ 113 /*! 114 115 In this method, we assume that the normalized stress is given by the stress 116 required to unstabilize the weakest triple node of a junction. A very simple, 117 yet, straightforward average is made considering all possible symmetric geometries 118 of these junctions (see Dupuy et al. 2017). 119 120 \param nu Poisson coefficient 121 \param burgers Burgers vectors of the junction (b1 = glide, b2 = forest and b3 = junction) 122 \param phi angles between the Burgers vectors and the junction direction (rad) 123 \return normalized junction strength. 124 125 */ 126 //==================================================================================== computeNormalizedStress(double nu,const std::vector<double> & burgers,const std::vector<double> & phi)127 double computeNormalizedStress(double nu, 128 const std::vector<double>& burgers, 129 const std::vector<double>& phi) 130 { 131 132 // line tension model for all three segments 133 double mu = 1.0; 134 IsotropicLineTensionModel lt1(mu,burgers[0],nu); 135 FrankRead fr1(lt1); 136 IsotropicLineTensionModel lt2(mu,burgers[1],nu); 137 FrankRead fr2(lt2); 138 IsotropicLineTensionModel lt3(mu,burgers[2],nu); 139 FrankRead fr3(lt3); 140 141 // consider all the 4 possibles symmetries 142 double averageStress = 0; 143 for(int i=0; i!=4; i++) 144 { 145 int sgn1 = 1 - 2*(i%2); 146 int sgn2 = 1 - 2*(i/2); 147 148 TripleNode tn(fr1,sgn1*phi[0],fr2,sgn2*phi[1],fr3,phi[2]); 149 double anglerad; 150 try 151 { 152 anglerad = tn.computeInitialAngle(); 153 } 154 catch(int error) 155 { 156 if( error == -1 ) 157 std::cerr << "=> no critical angle <=> no junction " << std::endl; 158 continue; 159 } 160 161 double stress = 0; 162 try 163 { 164 stress = tn.computeCriticalStress(anglerad,0.001,0.00001,20); 165 } 166 catch(int) 167 { 168 std::cerr << "=>could not compute critical stress of the triple node" << std::endl; 169 } 170 averageStress += stress; 171 } 172 173 // normalized average stress for the 4 symmetries 174 averageStress /= 4.0 * burgers[0] *burgers[0]; 175 176 return averageStress; 177 178 } 179 180 //========================================================================= 181 // computeHeuristicNonColinearHardeningCoefficient 182 //------------------------------------------------------------------------- 183 //! Compute an heuristic hardening coefficient for non colinear junctions 184 //------------------------------------------------------------------------- 185 /*! 186 187 see Dupuy et al. 2017 188 189 \param nu Poisson coefficient 190 \param sij normalized junction strength computed using linetension model 191 \return hardening coefficient 192 193 */ 194 //========================================================================= computeHeuristicNonColinearHardeningCoefficient(double,double sij)195 double computeHeuristicNonColinearHardeningCoefficient(double, 196 double sij) 197 { 198 return pow( 0.21 * sij + 0.37 , 2.0 ); 199 } 200 201 //========================================================================= 202 // computeHeuristicColinearHardeningCoefficient 203 //------------------------------------------------------------------------- 204 //! Compute an heuristic hardening coefficient for colinear junctions 205 //------------------------------------------------------------------------- 206 /*! 207 208 see Dupuy et al. 2017 209 210 \param pcrystallo pointer to a Crystallo object 211 \param nu Poisson coefficient 212 \param iplane1 glide plane of the first dislocation 213 \param iplane2 glide plane of the second dislocation 214 \return hardening coefficient 215 216 */ 217 //========================================================================= computeHeuristicColinearHardeningCoefficient(const Crystallo * pcrystallo,double nu,const IPlane & iplane1,const IPlane & iplane2)218 double computeHeuristicColinearHardeningCoefficient(const Crystallo* pcrystallo, 219 double nu, 220 const IPlane& iplane1, 221 const IPlane& iplane2) 222 { 223 Vect3 n1(pcrystallo->normal(iplane1)); 224 Vect3 n2(pcrystallo->normal(iplane2)); 225 return 0.4 * nu + 0.8 - 0.3* std::abs(n1.Dot(n2)); 226 } 227 228 //========================================================================= 229 // computeHeuristicHardeningCoefficient 230 //------------------------------------------------------------------------- 231 //! Compute an heuristic hardening coefficient 232 //------------------------------------------------------------------------- 233 /*! 234 235 see Dupuy et al. 2017 236 237 \param pcrystallo pointer to a Crystallo object 238 \param nu Poisson coefficient 239 \param iplane1 glide plane of the first dislocation 240 \param iplane2 glide plane of the second dislocation 241 \param bj norm of the Burgers vector of the junction 242 \param sij normalized junction strength computed using linetension model 243 \return hardening coefficient 244 245 */ 246 //========================================================================= computeHeuristicHardeningCoefficient(const Crystallo * pcrystallo,double nu,const IPlane & iplane1,const IPlane & iplane2,double bj,double sij)247 double computeHeuristicHardeningCoefficient(const Crystallo* pcrystallo, 248 double nu, 249 const IPlane& iplane1, 250 const IPlane& iplane2, 251 double bj, 252 double sij) 253 { 254 if(bj<0.001) 255 return computeHeuristicColinearHardeningCoefficient(pcrystallo,nu,iplane1,iplane2); 256 else 257 return computeHeuristicNonColinearHardeningCoefficient(nu,sij); 258 } 259 260 } // end of namespace Dupuy2017 261 262 } // end of namespace numodis 263