223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

capitalize constants 

not sure whether we should make exceptions for k_B and T_room
This commit is contained in:
Martin Diehl 2021-11-25 11:52:52 +01:00
parent e941867290
commit fb8c515a98
5 changed files with 24 additions and 25 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
python/damask/_colormap.py
View File

@ -5,7 +5,6 @@ import colorsys
from pathlib import Path from pathlib import Path
from typing import Sequence, Union, TextIO from typing import Sequence, Union, TextIO
import numpy as np import numpy as np
import matplotlib as mpl import matplotlib as mpl
if os.name == 'posix' and 'DISPLAY' not in os.environ: if os.name == 'posix' and 'DISPLAY' not in os.environ:
@ -17,9 +16,9 @@ from PIL import Image
from . import util from . import util
from . import Table from . import Table
_eps = 216./24389. _EPS = 216./24389.
_kappa = 24389./27. _KAPPA = 24389./27.
_ref_white = np.array([.95047, 1.00000, 1.08883]) # Observer = 2, Illuminant = D65 _REF_WHITE = np.array([.95047, 1.00000, 1.08883]) # Observer = 2, Illuminant = D65
# ToDo (if needed) # ToDo (if needed)
# - support alpha channel (paraview/ASCII/input) # - support alpha channel (paraview/ASCII/input)
@ -522,10 +521,10 @@ class Colormap(mpl.colors.ListedColormap):
f_z = (lab[0]+16.)/116. - lab[2]/200. f_z = (lab[0]+16.)/116. - lab[2]/200.
return np.array([ return np.array([
f_x**3. if f_x**3. > _eps else (116.*f_x-16.)/_kappa, f_x**3. if f_x**3. > _EPS else (116.*f_x-16.)/_KAPPA,
((lab[0]+16.)/116.)**3 if lab[0]>_kappa*_eps else lab[0]/_kappa, ((lab[0]+16.)/116.)**3 if lab[0]>_KAPPA*_EPS else lab[0]/_KAPPA,
f_z**3. if f_z**3. > _eps else (116.*f_z-16.)/_kappa f_z**3. if f_z**3. > _EPS else (116.*f_z-16.)/_KAPPA
])*(ref_white if ref_white is not None else _ref_white) ])*(ref_white if ref_white is not None else _REF_WHITE)
@staticmethod @staticmethod
def _xyz2lab(xyz: np.ndarray, ref_white: np.ndarray = None) -> np.ndarray: def _xyz2lab(xyz: np.ndarray, ref_white: np.ndarray = None) -> np.ndarray:
@ -537,8 +536,8 @@ class Colormap(mpl.colors.ListedColormap):
http://www.brucelindbloom.com/index.html?Eqn_Lab_to_XYZ.html http://www.brucelindbloom.com/index.html?Eqn_Lab_to_XYZ.html
""" """
ref_white = ref_white if ref_white is not None else _ref_white ref_white = ref_white if ref_white is not None else _REF_WHITE
f = np.where(xyz/ref_white > _eps,(xyz/ref_white)**(1./3.),(_kappa*xyz/ref_white+16.)/116.) f = np.where(xyz/ref_white > _EPS,(xyz/ref_white)**(1./3.),(_KAPPA*xyz/ref_white+16.)/116.)
return np.array([ return np.array([
116.0 * f[1] - 16.0, 116.0 * f[1] - 16.0,

2
src/constants.f90
View File

@ -10,6 +10,6 @@ module constants
real(pReal), parameter :: & real(pReal), parameter :: &
T_ROOM = 300.0_pReal, & !< Room temperature in K T_ROOM = 300.0_pReal, & !< Room temperature in K
kB = 1.38e-23_pReal !< Boltzmann constant in J/Kelvin K_B = 1.38e-23_pReal !< Boltzmann constant in J/Kelvin
end module constants end module constants

4
src/phase_mechanical_plastic_dislotungsten.f90
View File

@ -341,7 +341,7 @@ module subroutine dislotungsten_dotState(Mp,T,ph,en)
dot_rho_dip_formation = merge(2.0_pReal*(d_hat-prm%d_caron)*stt%rho_mob(:,en)*dot%gamma_sl(:,en)/prm%b_sl, & dot_rho_dip_formation = merge(2.0_pReal*(d_hat-prm%d_caron)*stt%rho_mob(:,en)*dot%gamma_sl(:,en)/prm%b_sl, &
0.0_pReal, & 0.0_pReal, &
prm%dipoleformation) prm%dipoleformation)
v_cl = (3.0_pReal*mu*prm%D_0*exp(-prm%Q_cl/(kB*T))*prm%f_at/(2.0_pReal*PI*kB*T)) & v_cl = (3.0_pReal*mu*prm%D_0*exp(-prm%Q_cl/(K_B*T))*prm%f_at/(2.0_pReal*PI*K_B*T)) &
* (1.0_pReal/(d_hat+prm%d_caron)) * (1.0_pReal/(d_hat+prm%d_caron))
dot_rho_dip_climb = (4.0_pReal*v_cl*stt%rho_dip(:,en))/(d_hat-prm%d_caron) ! ToDo: Discuss with Franz: Stress dependency? dot_rho_dip_climb = (4.0_pReal*v_cl*stt%rho_dip(:,en))/(d_hat-prm%d_caron) ! ToDo: Discuss with Franz: Stress dependency?
end where end where
@ -472,7 +472,7 @@ pure subroutine kinetics(Mp,T,ph,en, &
if (present(tau_pos_out)) tau_pos_out = tau_pos if (present(tau_pos_out)) tau_pos_out = tau_pos
if (present(tau_neg_out)) tau_neg_out = tau_neg if (present(tau_neg_out)) tau_neg_out = tau_neg
associate(BoltzmannRatio => prm%Q_s/(kB*T), & associate(BoltzmannRatio => prm%Q_s/(K_B*T), &
b_rho_half => stt%rho_mob(:,en) * prm%b_sl * 0.5_pReal, & b_rho_half => stt%rho_mob(:,en) * prm%b_sl * 0.5_pReal, &
effectiveLength => dst%Lambda_sl(:,en) - prm%w) effectiveLength => dst%Lambda_sl(:,en) - prm%w)

12
src/phase_mechanical_plastic_dislotwin.f90
View File

@ -594,7 +594,7 @@ module subroutine dislotwin_LpAndItsTangent(Lp,dLp_dMp,Mp,T,ph,en)
shearBandingContribution: if (dNeq0(prm%v_sb)) then shearBandingContribution: if (dNeq0(prm%v_sb)) then
E_kB_T = prm%E_sb/(kB*T) E_kB_T = prm%E_sb/(K_B*T)
call math_eigh33(eigValues,eigVectors,Mp) ! is Mp symmetric by design? call math_eigh33(eigValues,eigVectors,Mp) ! is Mp symmetric by design?
do i = 1,6 do i = 1,6
@ -691,8 +691,8 @@ module subroutine dislotwin_dotState(Mp,T,ph,en)
* (prm%Gamma_sf(1) + prm%Gamma_sf(2) * T) / (mu*prm%b_sl(i)), & * (prm%Gamma_sf(1) + prm%Gamma_sf(2) * T) / (mu*prm%b_sl(i)), &
1.0_pReal, & 1.0_pReal, &
prm%ExtendedDislocations) prm%ExtendedDislocations)
v_cl = 2.0_pReal*prm%omega*b_d**2.0_pReal*exp(-prm%Q_cl/(kB*T)) & v_cl = 2.0_pReal*prm%omega*b_d**2.0_pReal*exp(-prm%Q_cl/(K_B*T)) &
* (exp(abs(sigma_cl)*prm%b_sl(i)**3.0_pReal/(kB*T)) - 1.0_pReal) * (exp(abs(sigma_cl)*prm%b_sl(i)**3.0_pReal/(K_B*T)) - 1.0_pReal)
dot_rho_dip_climb(i) = 4.0_pReal*v_cl*stt%rho_dip(i,en) & dot_rho_dip_climb(i) = 4.0_pReal*v_cl*stt%rho_dip(i,en) &
/ (d_hat-prm%d_caron(i)) / (d_hat-prm%d_caron(i))
@ -904,7 +904,7 @@ pure subroutine kinetics_sl(Mp,T,ph,en, &
significantStress: where(tau_eff > tol_math_check) significantStress: where(tau_eff > tol_math_check)
stressRatio = tau_eff/prm%tau_0 stressRatio = tau_eff/prm%tau_0
StressRatio_p = stressRatio** prm%p StressRatio_p = stressRatio** prm%p
Q_kB_T = prm%Q_sl/(kB*T) Q_kB_T = prm%Q_sl/(K_B*T)
v_wait_inverse = exp(Q_kB_T*(1.0_pReal-StressRatio_p)** prm%q) & v_wait_inverse = exp(Q_kB_T*(1.0_pReal-StressRatio_p)** prm%q) &
/ prm%v_0 / prm%v_0
v_run_inverse = prm%B/(tau_eff*prm%b_sl) v_run_inverse = prm%B/(tau_eff*prm%b_sl)
@ -977,7 +977,7 @@ pure subroutine kinetics_tw(Mp,T,dot_gamma_sl,ph,en,&
Ndot0=(abs(dot_gamma_sl(s1))*(stt%rho_mob(s2,en)+stt%rho_dip(s2,en))+& Ndot0=(abs(dot_gamma_sl(s1))*(stt%rho_mob(s2,en)+stt%rho_dip(s2,en))+&
abs(dot_gamma_sl(s2))*(stt%rho_mob(s1,en)+stt%rho_dip(s1,en)))/& abs(dot_gamma_sl(s2))*(stt%rho_mob(s1,en)+stt%rho_dip(s1,en)))/&
(prm%L_tw*prm%b_sl(i))*& (prm%L_tw*prm%b_sl(i))*&
(1.0_pReal-exp(-prm%V_cs/(kB*T)*(dst%tau_r_tw(i,en)-tau(i)))) (1.0_pReal-exp(-prm%V_cs/(K_B*T)*(dst%tau_r_tw(i,en)-tau(i))))
else else
Ndot0=0.0_pReal Ndot0=0.0_pReal
end if end if
@ -1046,7 +1046,7 @@ pure subroutine kinetics_tr(Mp,T,dot_gamma_sl,ph,en,&
Ndot0=(abs(dot_gamma_sl(s1))*(stt%rho_mob(s2,en)+stt%rho_dip(s2,en))+& Ndot0=(abs(dot_gamma_sl(s1))*(stt%rho_mob(s2,en)+stt%rho_dip(s2,en))+&
abs(dot_gamma_sl(s2))*(stt%rho_mob(s1,en)+stt%rho_dip(s1,en)))/& abs(dot_gamma_sl(s2))*(stt%rho_mob(s1,en)+stt%rho_dip(s1,en)))/&
(prm%L_tr*prm%b_sl(i))*& (prm%L_tr*prm%b_sl(i))*&
(1.0_pReal-exp(-prm%V_cs/(kB*T)*(dst%tau_r_tr(i,en)-tau(i)))) (1.0_pReal-exp(-prm%V_cs/(K_B*T)*(dst%tau_r_tr(i,en)-tau(i))))
else else
Ndot0=0.0_pReal Ndot0=0.0_pReal
end if end if

12
src/phase_mechanical_plastic_nonlocal.f90
View File

@ -1091,9 +1091,9 @@ module subroutine nonlocal_dotState(Mp, Temperature,timestep, &
! thermally activated annihilation of edge dipoles by climb ! thermally activated annihilation of edge dipoles by climb
rhoDotThermalAnnihilation = 0.0_pReal rhoDotThermalAnnihilation = 0.0_pReal
D_SD = prm%D_0 * exp(-prm%Q_cl / (kB * Temperature)) ! eq. 3.53 D_SD = prm%D_0 * exp(-prm%Q_cl / (K_B * Temperature)) ! eq. 3.53
v_climb = D_SD * mu * prm%V_at & v_climb = D_SD * mu * prm%V_at &
/ (PI * (1.0_pReal-nu) * (dUpper(:,1) + dLower(:,1)) * kB * Temperature) ! eq. 3.54 / (PI * (1.0_pReal-nu) * (dUpper(:,1) + dLower(:,1)) * K_B * Temperature) ! eq. 3.54
forall (s = 1:prm%sum_N_sl, dUpper(s,1) > dLower(s,1)) & forall (s = 1:prm%sum_N_sl, dUpper(s,1) > dLower(s,1)) &
rhoDotThermalAnnihilation(s,9) = max(- 4.0_pReal * rhoDip(s,1) * v_climb(s) / (dUpper(s,1) - dLower(s,1)), & rhoDotThermalAnnihilation(s,9) = max(- 4.0_pReal * rhoDip(s,1) * v_climb(s) / (dUpper(s,1) - dLower(s,1)), &
- rhoDip(s,1) / timestep - rhoDotAthermalAnnihilation(s,9) & - rhoDip(s,1) / timestep - rhoDotAthermalAnnihilation(s,9) &
@ -1668,9 +1668,9 @@ pure subroutine kinetics(v, dv_dtau, dv_dtauNS, tau, tauNS, tauThreshold, c, T,
activationEnergy_P = criticalStress_P * activationVolume_P activationEnergy_P = criticalStress_P * activationVolume_P
tauRel_P = min(1.0_pReal, tauEff / criticalStress_P) tauRel_P = min(1.0_pReal, tauEff / criticalStress_P)
tPeierls = 1.0_pReal / prm%nu_a & tPeierls = 1.0_pReal / prm%nu_a &
* exp(activationEnergy_P / (kB * T) & * exp(activationEnergy_P / (K_B * T) &
* (1.0_pReal - tauRel_P**prm%p)**prm%q) * (1.0_pReal - tauRel_P**prm%p)**prm%q)
dtPeierls_dtau = merge(tPeierls * prm%p * prm%q * activationVolume_P / (kB * T) & dtPeierls_dtau = merge(tPeierls * prm%p * prm%q * activationVolume_P / (K_B * T) &
* (1.0_pReal - tauRel_P**prm%p)**(prm%q-1.0_pReal) * tauRel_P**(prm%p-1.0_pReal), & * (1.0_pReal - tauRel_P**prm%p)**(prm%q-1.0_pReal) * tauRel_P**(prm%p-1.0_pReal), &
0.0_pReal, & 0.0_pReal, &
tauEff < criticalStress_P) tauEff < criticalStress_P)
@ -1682,8 +1682,8 @@ pure subroutine kinetics(v, dv_dtau, dv_dtauNS, tau, tauNS, tauThreshold, c, T,
criticalStress_S = prm%Q_sol / activationVolume_S criticalStress_S = prm%Q_sol / activationVolume_S
tauRel_S = min(1.0_pReal, tauEff / criticalStress_S) tauRel_S = min(1.0_pReal, tauEff / criticalStress_S)
tSolidSolution = 1.0_pReal / prm%nu_a & tSolidSolution = 1.0_pReal / prm%nu_a &
* exp(prm%Q_sol / (kB * T)* (1.0_pReal - tauRel_S**prm%p)**prm%q) * exp(prm%Q_sol / (K_B * T)* (1.0_pReal - tauRel_S**prm%p)**prm%q)
dtSolidSolution_dtau = merge(tSolidSolution * prm%p * prm%q * activationVolume_S / (kB * T) & dtSolidSolution_dtau = merge(tSolidSolution * prm%p * prm%q * activationVolume_S / (K_B * T) &
* (1.0_pReal - tauRel_S**prm%p)**(prm%q-1.0_pReal)* tauRel_S**(prm%p-1.0_pReal), & * (1.0_pReal - tauRel_S**prm%p)**(prm%q-1.0_pReal)* tauRel_S**(prm%p-1.0_pReal), &
0.0_pReal, & 0.0_pReal, &
tauEff < criticalStress_S) tauEff < criticalStress_S)