From 8d6c82e704ce3e1d58afd38ce39f85209a0915b5 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 14 Feb 2020 09:00:14 +0100 Subject: [PATCH] no need for enums they just complicate the code, any performance gain should be negligible --- src/constitutive_plastic_disloUCLA.f90 | 255 +++++++++------------ src/constitutive_plastic_none.f90 | 12 +- src/constitutive_plastic_phenopowerlaw.f90 | 2 +- 3 files changed, 114 insertions(+), 155 deletions(-) diff --git a/src/constitutive_plastic_disloUCLA.f90 b/src/constitutive_plastic_disloUCLA.f90 index 2b353ba26..00a1b8e21 100644 --- a/src/constitutive_plastic_disloUCLA.f90 +++ b/src/constitutive_plastic_disloUCLA.f90 @@ -6,21 +6,10 @@ !> @brief crystal plasticity model for bcc metals, especially Tungsten !-------------------------------------------------------------------------------------------------- submodule(constitutive) plastic_disloUCLA - + real(pReal), parameter :: & kB = 1.38e-23_pReal !< Boltzmann constant in J/Kelvin - - enum, bind(c) - enumerator :: & - undefined_ID, & - rho_mob_ID, & - rho_dip_ID, & - dot_gamma_sl_ID, & - gamma_sl_ID, & - Lambda_sl_ID, & - tau_pass_ID - end enum - + type :: tParameters real(pReal) :: & aTol_rho, & @@ -28,7 +17,7 @@ submodule(constitutive) plastic_disloUCLA mu, & D_0, & !< prefactor for self-diffusion coefficient Q_cl !< activation energy for dislocation climb - real(pReal), dimension(:), allocatable :: & + real(pReal), allocatable, dimension(:) :: & rho_mob_0, & !< initial dislocation density rho_dip_0, & !< initial dipole density b_sl, & !< magnitude of burgers vector [m] @@ -46,30 +35,30 @@ submodule(constitutive) plastic_disloUCLA kink_height, & !< height of the kink pair w, & !< width of the kink pair omega !< attempt frequency for kink pair nucleation - real(pReal), dimension(:,:), allocatable :: & + real(pReal), allocatable, dimension(:,:) :: & h_sl_sl, & !< slip resistance from slip activity forestProjectionEdge - real(pReal), dimension(:,:,:), allocatable :: & + real(pReal), allocatable, dimension(:,:,:) :: & Schmid, & nonSchmid_pos, & nonSchmid_neg integer :: & sum_N_sl !< total number of active slip system - integer, dimension(:), allocatable :: & + integer, allocatable, dimension(:) :: & N_sl !< number of active slip systems for each family - integer(kind(undefined_ID)), dimension(:),allocatable :: & - outputID !< ID of each post result output + character(len=pStringLen), allocatable, dimension(:) :: & + output logical :: & dipoleFormation !< flag indicating consideration of dipole formation end type !< container type for internal constitutive parameters - + type :: tDisloUCLAState real(pReal), dimension(:,:), pointer :: & rho_mob, & rho_dip, & gamma_sl end type tDisloUCLAState - + type :: tDisloUCLAdependentState real(pReal), dimension(:,:), allocatable :: & Lambda_sl, & @@ -88,41 +77,36 @@ contains !-------------------------------------------------------------------------------------------------- -!> @brief module initialization +!> @brief Perform module initialization. !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- module subroutine plastic_disloUCLA_init - + integer :: & Ninstance, & p, i, & NipcMyPhase, & sizeState, sizeDotState, & startIndex, endIndex - - integer(kind(undefined_ID)) :: & - outputID - + character(len=pStringLen) :: & extmsg = '' - character(len=pStringLen), dimension(:), allocatable :: & - outputs - - write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_DISLOUCLA_label//' init -+>>>' - - write(6,'(/,a)') ' Cereceda et al., International Journal of Plasticity 78:242–256, 2016' - write(6,'(a)') ' https://dx.doi.org/10.1016/j.ijplas.2015.09.002' - + + write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_DISLOUCLA_label//' init -+>>>'; flush(6) + + write(6,'(/,a)') ' Cereceda et al., International Journal of Plasticity 78:242–256, 2016' + write(6,'(a)') ' https://dx.doi.org/10.1016/j.ijplas.2015.09.002' + Ninstance = count(phase_plasticity == PLASTICITY_DISLOUCLA_ID) if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) & write(6,'(a16,1x,i5,/)') '# instances:',Ninstance - + allocate(param(Ninstance)) allocate(state(Ninstance)) allocate(dotState(Ninstance)) allocate(dependentState(Ninstance)) - - + + do p = 1, size(phase_plasticity) if (phase_plasticity(p) /= PLASTICITY_DISLOUCLA_ID) cycle associate(prm => param(phase_plasticityInstance(p)), & @@ -134,9 +118,9 @@ module subroutine plastic_disloUCLA_init !-------------------------------------------------------------------------------------------------- ! optional parameters that need to be defined prm%mu = lattice_mu(p) - + prm%aTol_rho = config%getFloat('atol_rho') - + ! sanity checks if (prm%aTol_rho <= 0.0_pReal) extmsg = trim(extmsg)//' atol_rho' @@ -147,7 +131,7 @@ module subroutine plastic_disloUCLA_init slipActive: if (prm%sum_N_sl > 0) then prm%Schmid = lattice_SchmidMatrix_slip(prm%N_sl,config%getString('lattice_structure'),& config%getFloat('c/a',defaultVal=0.0_pReal)) - + if(trim(config%getString('lattice_structure')) == 'bcc') then prm%nonSchmidCoeff = config%getFloats('nonschmid_coefficients',& defaultVal = emptyRealArray) @@ -158,20 +142,20 @@ module subroutine plastic_disloUCLA_init prm%nonSchmid_pos = prm%Schmid prm%nonSchmid_neg = prm%Schmid endif - + prm%h_sl_sl = lattice_interaction_SlipBySlip(prm%N_sl, & config%getFloats('interaction_slipslip'), & config%getString('lattice_structure')) prm%forestProjectionEdge = lattice_forestProjection_edge(prm%N_sl,config%getString('lattice_structure'),& config%getFloat('c/a',defaultVal=0.0_pReal)) prm%forestProjectionEdge = transpose(prm%forestProjectionEdge) - + prm%rho_mob_0 = config%getFloats('rhoedge0', requiredSize=size(prm%N_sl)) prm%rho_dip_0 = config%getFloats('rhoedgedip0', requiredSize=size(prm%N_sl)) prm%v0 = config%getFloats('v0', requiredSize=size(prm%N_sl)) prm%b_sl = config%getFloats('slipburgers', requiredSize=size(prm%N_sl)) prm%delta_F = config%getFloats('qedge', requiredSize=size(prm%N_sl)) - + prm%i_sl = config%getFloats('clambdaslip', requiredSize=size(prm%N_sl)) prm%tau_0 = config%getFloats('tau_peierls', requiredSize=size(prm%N_sl)) prm%p = config%getFloats('p_slip', requiredSize=size(prm%N_sl), & @@ -182,14 +166,14 @@ module subroutine plastic_disloUCLA_init prm%w = config%getFloats('kink_width', requiredSize=size(prm%N_sl)) prm%omega = config%getFloats('omega', requiredSize=size(prm%N_sl)) prm%B = config%getFloats('friction_coeff', requiredSize=size(prm%N_sl)) - + prm%D = config%getFloat('grainsize') prm%D_0 = config%getFloat('d0') prm%Q_cl = config%getFloat('qsd') prm%atomicVolume = config%getFloat('catomicvolume') * prm%b_sl**3.0_pReal prm%D_a = config%getFloat('cedgedipmindistance') * prm%b_sl prm%dipoleformation = config%getFloat('dipoleformationfactor') > 0.0_pReal !should be on by default, ToDo: change to /key/-type key - + ! expand: family => system prm%rho_mob_0 = math_expand(prm%rho_mob_0, prm%N_sl) prm%rho_dip_0 = math_expand(prm%rho_dip_0, prm%N_sl) @@ -206,8 +190,8 @@ module subroutine plastic_disloUCLA_init prm%i_sl = math_expand(prm%i_sl, prm%N_sl) prm%atomicVolume = math_expand(prm%atomicVolume, prm%N_sl) prm%D_a = math_expand(prm%D_a, prm%N_sl) - - + + ! sanity checks if ( prm%D_0 <= 0.0_pReal) extmsg = trim(extmsg)//' D_0' if ( prm%Q_cl <= 0.0_pReal) extmsg = trim(extmsg)//' Q_cl' @@ -219,7 +203,7 @@ module subroutine plastic_disloUCLA_init if (any(prm%tau_0 < 0.0_pReal)) extmsg = trim(extmsg)//' tau_0' if (any(prm%D_a <= 0.0_pReal)) extmsg = trim(extmsg)//' cedgedipmindistance or slipb_sl' if (any(prm%atomicVolume <= 0.0_pReal)) extmsg = trim(extmsg)//' catomicvolume or slipb_sl' - + else slipActive allocate(prm%rho_mob_0(0)) allocate(prm%rho_dip_0(0)) @@ -232,39 +216,14 @@ module subroutine plastic_disloUCLA_init !-------------------------------------------------------------------------------------------------- ! output pararameters - outputs = config%getStrings('(output)',defaultVal=emptyStringArray) - allocate(prm%outputID(0)) - do i=1, size(outputs) - outputID = undefined_ID - select case(trim(outputs(i))) - - case ('edge_density') - outputID = merge(rho_mob_ID,undefined_ID,prm%sum_N_sl>0) - case ('dipole_density') - outputID = merge(rho_dip_ID,undefined_ID,prm%sum_N_sl>0) - case ('shear_rate','shearrate','shear_rate_slip','shearrate_slip') - outputID = merge(dot_gamma_sl_ID,undefined_ID,prm%sum_N_sl>0) - case ('accumulated_shear','accumulatedshear','accumulated_shear_slip') - outputID = merge(gamma_sl_ID,undefined_ID,prm%sum_N_sl>0) - case ('mfp','mfp_slip') - outputID = merge(Lambda_sl_ID,undefined_ID,prm%sum_N_sl>0) - case ('threshold_stress','threshold_stress_slip') - outputID = merge(tau_pass_ID,undefined_ID,prm%sum_N_sl>0) - - end select - - if (outputID /= undefined_ID) then - prm%outputID = [prm%outputID, outputID] - endif - - enddo + prm%output = config%getStrings('(output)',defaultVal=emptyStringArray) !-------------------------------------------------------------------------------------------------- ! allocate state arrays NipcMyPhase = count(material_phaseAt == p) * discretization_nIP sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeState = sizeDotState - + call material_allocatePlasticState(p,NipcMyPhase,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- @@ -275,14 +234,14 @@ module subroutine plastic_disloUCLA_init stt%rho_mob= spread(prm%rho_mob_0,2,NipcMyPhase) dot%rho_mob=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%aTolState(startIndex:endIndex) = prm%aTol_rho - + startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%rho_dip=>plasticState(p)%state(startIndex:endIndex,:) stt%rho_dip= spread(prm%rho_dip_0,2,NipcMyPhase) dot%rho_dip=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%aTolState(startIndex:endIndex) = prm%aTol_rho - + startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%gamma_sl=>plasticState(p)%state(startIndex:endIndex,:) @@ -290,21 +249,21 @@ module subroutine plastic_disloUCLA_init plasticState(p)%aTolState(startIndex:endIndex) = 1.0e6_pReal ! Don't use for convergence check ! global alias plasticState(p)%slipRate => plasticState(p)%dotState(startIndex:endIndex,:) - + allocate(dst%Lambda_sl(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) allocate(dst%threshold_stress(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) - + plasticState(p)%state0 = plasticState(p)%state ! ToDo: this could be done centrally - + end associate - + enddo end subroutine plastic_disloUCLA_init !-------------------------------------------------------------------------------------------------- -!> @brief calculates plastic velocity gradient and its tangent +!> @brief Calculate plastic velocity gradient and its tangent. !-------------------------------------------------------------------------------------------------- pure module subroutine plastic_disloUCLA_LpAndItsTangent(Lp,dLp_dMp, & Mp,T,instance,of) @@ -312,7 +271,7 @@ pure module subroutine plastic_disloUCLA_LpAndItsTangent(Lp,dLp_dMp, & Lp !< plastic velocity gradient real(pReal), dimension(3,3,3,3), intent(out) :: & dLp_dMp !< derivative of Lp with respect to the Mandel stress - + real(pReal), dimension(3,3), intent(in) :: & Mp !< Mandel stress real(pReal), intent(in) :: & @@ -320,18 +279,18 @@ pure module subroutine plastic_disloUCLA_LpAndItsTangent(Lp,dLp_dMp, & integer, intent(in) :: & instance, & of - + integer :: & i,k,l,m,n real(pReal), dimension(param(instance)%sum_N_sl) :: & dot_gamma_pos,dot_gamma_neg, & ddot_gamma_dtau_pos,ddot_gamma_dtau_neg - + Lp = 0.0_pReal dLp_dMp = 0.0_pReal - + associate(prm => param(instance)) - + call kinetics(Mp,T,instance,of,dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg) do i = 1, prm%sum_N_sl Lp = Lp + (dot_gamma_pos(i)+dot_gamma_neg(i))*prm%Schmid(1:3,1:3,i) @@ -340,17 +299,17 @@ pure module subroutine plastic_disloUCLA_LpAndItsTangent(Lp,dLp_dMp, & + ddot_gamma_dtau_pos(i) * prm%Schmid(k,l,i) * prm%nonSchmid_pos(m,n,i) & + ddot_gamma_dtau_neg(i) * prm%Schmid(k,l,i) * prm%nonSchmid_neg(m,n,i) enddo - + end associate end subroutine plastic_disloUCLA_LpAndItsTangent !-------------------------------------------------------------------------------------------------- -!> @brief calculates the rate of change of microstructure +!> @brief Calculate the rate of change of microstructure. !-------------------------------------------------------------------------------------------------- module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of) - + real(pReal), dimension(3,3), intent(in) :: & Mp !< Mandel stress real(pReal), intent(in) :: & @@ -358,7 +317,7 @@ module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of) integer, intent(in) :: & instance, & of - + real(pReal) :: & VacancyDiffusion real(pReal), dimension(param(instance)%sum_N_sl) :: & @@ -369,16 +328,16 @@ module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of) dot_rho_dip_formation, & dot_rho_dip_climb, & dip_distance - + associate(prm => param(instance), stt => state(instance),dot => dotState(instance), dst => dependentState(instance)) - + call kinetics(Mp,T,instance,of,& gdot_pos,gdot_neg, & tau_pos_out = tau_pos,tau_neg_out = tau_neg) - + dot%gamma_sl(:,of) = (gdot_pos+gdot_neg) ! ToDo: needs to be abs VacancyDiffusion = prm%D_0*exp(-prm%Q_cl/(kB*T)) - + where(dEq0(tau_pos)) ! ToDo: use avg of pos and neg dot_rho_dip_formation = 0.0_pReal dot_rho_dip_climb = 0.0_pReal @@ -393,73 +352,73 @@ module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of) * (1.0_pReal/(dip_distance+prm%D_a)) dot_rho_dip_climb = (4.0_pReal*v_cl*stt%rho_dip(:,of))/(dip_distance-prm%D_a) ! ToDo: Discuss with Franz: Stress dependency? end where - + dot%rho_mob(:,of) = abs(dot%gamma_sl(:,of))/(prm%b_sl*dst%Lambda_sl(:,of)) & ! multiplication - dot_rho_dip_formation & - (2.0_pReal*prm%D_a)/prm%b_sl*stt%rho_mob(:,of)*abs(dot%gamma_sl(:,of)) ! Spontaneous annihilation of 2 single edge dislocations dot%rho_dip(:,of) = dot_rho_dip_formation & - (2.0_pReal*prm%D_a)/prm%b_sl*stt%rho_dip(:,of)*abs(dot%gamma_sl(:,of)) & ! Spontaneous annihilation of a single edge dislocation with a dipole constituent - dot_rho_dip_climb - + end associate end subroutine plastic_disloUCLA_dotState !-------------------------------------------------------------------------------------------------- -!> @brief calculates derived quantities from state +!> @brief Calculate derived quantities from state. !-------------------------------------------------------------------------------------------------- module subroutine plastic_disloUCLA_dependentState(instance,of) - + integer, intent(in) :: & instance, & of - + real(pReal), dimension(param(instance)%sum_N_sl) :: & dislocationSpacing - + associate(prm => param(instance), stt => state(instance),dst => dependentState(instance)) - + dislocationSpacing = sqrt(matmul(prm%forestProjectionEdge, & stt%rho_mob(:,of)+stt%rho_dip(:,of))) dst%threshold_stress(:,of) = prm%mu*prm%b_sl & * sqrt(matmul(prm%h_sl_sl,stt%rho_mob(:,of)+stt%rho_dip(:,of))) - + dst%Lambda_sl(:,of) = prm%D/(1.0_pReal+prm%D*dislocationSpacing/prm%i_sl) - + end associate end subroutine plastic_disloUCLA_dependentState !-------------------------------------------------------------------------------------------------- -!> @brief writes results to HDF5 output file +!> @brief Write results to HDF5 output file. !-------------------------------------------------------------------------------------------------- module subroutine plastic_disloUCLA_results(instance,group) integer, intent(in) :: instance character(len=*), intent(in) :: group - + integer :: o associate(prm => param(instance), stt => state(instance), dst => dependentState(instance)) - outputsLoop: do o = 1,size(prm%outputID) - select case(prm%outputID(o)) - case (rho_mob_ID) - call results_writeDataset(group,stt%rho_mob,'rho_mob',& - 'mobile dislocation density','1/m²') - case (rho_dip_ID) - call results_writeDataset(group,stt%rho_dip,'rho_dip',& - 'dislocation dipole density''1/m²') - case (dot_gamma_sl_ID) - call results_writeDataset(group,stt%gamma_sl,'dot_gamma_sl',& ! this is not dot!! - 'plastic shear','1') - case (Lambda_sl_ID) - call results_writeDataset(group,dst%Lambda_sl,'Lambda_sl',& - 'mean free path for slip','m') - case (tau_pass_ID) - call results_writeDataset(group,dst%threshold_stress,'tau_pass',& - 'threshold stress for slip','Pa') + outputsLoop: do o = 1,size(prm%output) + select case(trim(prm%output(o))) + case('edge_density') ! ToDo: should be rho_mob + if(prm%sum_N_sl>0) call results_writeDataset(group,stt%rho_mob,'rho_mob',& + 'mobile dislocation density','1/m²') + case('dipole_density') ! ToDo: should be rho_dip + if(prm%sum_N_sl>0) call results_writeDataset(group,stt%rho_dip,'rho_dip',& + 'dislocation dipole density''1/m²') + case('shear_rate_slip') ! should be gamma + if(prm%sum_N_sl>0) call results_writeDataset(group,stt%gamma_sl,'dot_gamma_sl',& ! this is not dot!! + 'plastic shear','1') + case('mfp_slip') !ToDo: should be Lambda + if(prm%sum_N_sl>0) call results_writeDataset(group,dst%Lambda_sl,'Lambda_sl',& + 'mean free path for slip','m') + case('threshold_stress_slip') !ToDo: should be tau_pass + if(prm%sum_N_sl>0) call results_writeDataset(group,dst%threshold_stress,'tau_pass',& + 'threshold stress for slip','Pa') end select enddo outputsLoop end associate @@ -468,15 +427,15 @@ end subroutine plastic_disloUCLA_results !-------------------------------------------------------------------------------------------------- -!> @brief Shear rates on slip systems, their derivatives with respect to resolved stress and the -! resolved stresss +!> @brief Calculate shear rates on slip systems, their derivatives with respect to resolved +! stress, and the resolved stress. !> @details Derivatives and resolved stress are calculated only optionally. ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to ! have the optional arguments at the end !-------------------------------------------------------------------------------------------------- pure subroutine kinetics(Mp,T,instance,of, & dot_gamma_pos,dot_gamma_neg,ddot_gamma_dtau_pos,ddot_gamma_dtau_neg,tau_pos_out,tau_neg_out) - + real(pReal), dimension(3,3), intent(in) :: & Mp !< Mandel stress real(pReal), intent(in) :: & @@ -484,7 +443,7 @@ pure subroutine kinetics(Mp,T,instance,of, & integer, intent(in) :: & instance, & of - + real(pReal), intent(out), dimension(param(instance)%sum_N_sl) :: & dot_gamma_pos, & dot_gamma_neg @@ -501,82 +460,82 @@ pure subroutine kinetics(Mp,T,instance,of, & t_n, t_k, dtk,dtn, & needsGoodName ! ToDo: @Karo: any idea? integer :: j - + associate(prm => param(instance), stt => state(instance), dst => dependentState(instance)) - + do j = 1, prm%sum_N_sl tau_pos(j) = math_mul33xx33(Mp,prm%nonSchmid_pos(1:3,1:3,j)) tau_neg(j) = math_mul33xx33(Mp,prm%nonSchmid_neg(1:3,1:3,j)) enddo - - + + if (present(tau_pos_out)) tau_pos_out = tau_pos if (present(tau_neg_out)) tau_neg_out = tau_neg - + associate(BoltzmannRatio => prm%delta_F/(kB*T), & dot_gamma_0 => stt%rho_mob(:,of)*prm%b_sl*prm%v0, & effectiveLength => dst%Lambda_sl(:,of) - prm%w) - + significantPositiveTau: where(abs(tau_pos)-dst%threshold_stress(:,of) > tol_math_check) StressRatio = (abs(tau_pos)-dst%threshold_stress(:,of))/prm%tau_0 StressRatio_p = StressRatio** prm%p StressRatio_pminus1 = StressRatio**(prm%p-1.0_pReal) needsGoodName = exp(-BoltzmannRatio*(1-StressRatio_p) ** prm%q) - + t_n = prm%b_sl/(needsGoodName*prm%omega*effectiveLength) t_k = effectiveLength * prm%B /(2.0_pReal*prm%b_sl*tau_pos) - + vel = prm%kink_height/(t_n + t_k) - + dot_gamma_pos = dot_gamma_0 * sign(vel,tau_pos) * 0.5_pReal else where significantPositiveTau dot_gamma_pos = 0.0_pReal end where significantPositiveTau - + if (present(ddot_gamma_dtau_pos)) then significantPositiveTau2: where(abs(tau_pos)-dst%threshold_stress(:,of) > tol_math_check) dtn = -1.0_pReal * t_n * BoltzmannRatio * prm%p * prm%q * (1.0_pReal-StressRatio_p)**(prm%q - 1.0_pReal) & * (StressRatio)**(prm%p - 1.0_pReal) / prm%tau_0 dtk = -1.0_pReal * t_k / tau_pos - + dvel = -1.0_pReal * prm%kink_height * (dtk + dtn) / (t_n + t_k)**2.0_pReal - + ddot_gamma_dtau_pos = dot_gamma_0 * dvel* 0.5_pReal else where significantPositiveTau2 ddot_gamma_dtau_pos = 0.0_pReal end where significantPositiveTau2 endif - + significantNegativeTau: where(abs(tau_neg)-dst%threshold_stress(:,of) > tol_math_check) StressRatio = (abs(tau_neg)-dst%threshold_stress(:,of))/prm%tau_0 StressRatio_p = StressRatio** prm%p StressRatio_pminus1 = StressRatio**(prm%p-1.0_pReal) needsGoodName = exp(-BoltzmannRatio*(1-StressRatio_p) ** prm%q) - + t_n = prm%b_sl/(needsGoodName*prm%omega*effectiveLength) t_k = effectiveLength * prm%B /(2.0_pReal*prm%b_sl*tau_pos) - + vel = prm%kink_height/(t_n + t_k) - + dot_gamma_neg = dot_gamma_0 * sign(vel,tau_neg) * 0.5_pReal else where significantNegativeTau dot_gamma_neg = 0.0_pReal end where significantNegativeTau - + if (present(ddot_gamma_dtau_neg)) then significantNegativeTau2: where(abs(tau_neg)-dst%threshold_stress(:,of) > tol_math_check) dtn = -1.0_pReal * t_n * BoltzmannRatio * prm%p * prm%q * (1.0_pReal-StressRatio_p)**(prm%q - 1.0_pReal) & * (StressRatio)**(prm%p - 1.0_pReal) / prm%tau_0 dtk = -1.0_pReal * t_k / tau_neg - + dvel = -1.0_pReal * prm%kink_height * (dtk + dtn) / (t_n + t_k)**2.0_pReal - + ddot_gamma_dtau_neg = dot_gamma_0 * dvel * 0.5_pReal else where significantNegativeTau2 ddot_gamma_dtau_neg = 0.0_pReal end where significantNegativeTau2 end if - + end associate end associate diff --git a/src/constitutive_plastic_none.f90 b/src/constitutive_plastic_none.f90 index 578fb4149..4840c9c09 100644 --- a/src/constitutive_plastic_none.f90 +++ b/src/constitutive_plastic_none.f90 @@ -18,19 +18,19 @@ module subroutine plastic_none_init Ninstance, & p, & NipcMyPhase - - write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_NONE_label//' init -+>>>' - + + write(6,'(/,a)') ' <<<+- plastic_'//PLASTICITY_NONE_label//' init -+>>>'; flush(6) + 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_phaseAt == p) * discretization_nIP call material_allocatePlasticState(p,NipcMyPhase,0,0,0) - + enddo end subroutine plastic_none_init diff --git a/src/constitutive_plastic_phenopowerlaw.f90 b/src/constitutive_plastic_phenopowerlaw.f90 index ade8a177f..6dfa3fb16 100644 --- a/src/constitutive_plastic_phenopowerlaw.f90 +++ b/src/constitutive_plastic_phenopowerlaw.f90 @@ -416,7 +416,7 @@ end subroutine plastic_phenopowerlaw_results !-------------------------------------------------------------------------------------------------- -!> @brief Calculate shear rates on slip systems and their derivatives with respect to resolved. +!> @brief Calculate shear rates on slip systems and their derivatives with respect to resolved ! stress. !> @details Derivatives are calculated only optionally. ! NOTE: Against the common convention, the result (i.e. intent(out)) variables are the last to