#@title
# HIDDEN
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from ipywidgets import interactive, widgets
from matplotlib.ticker import FormatStrFormatter
from scipy import stats
from IPython.core.display import HTML, clear_output
import pandas as pd
pd.set_option('display.max_rows', None)
%matplotlib inline
HTML("""
<style>
.output_png {
display: table-cell;
text-align: center;
vertical-align: middle;
}
.output{align-items: center;}
</style>
""")
# HIDDEN
mpl.rcParams.update({'font.size': 22})
#@title
def run_test(epsilon_i, E, sigma_f):
if epsilon_i * E < sigma_f:
epsilons = np.linspace(0, epsilon_i, 1000)
sigmas = E*epsilons
else:
epsilons = np.linspace(0, epsilon_i, 1000)
sigmas = E*epsilons
sigmas[sigmas>sigma_f] = np.zeros(sigmas[sigmas>sigma_f].shape[0])
return epsilons, sigmas
def hide_code_in_slideshow():
from IPython import display
import binascii
import os
uid = binascii.hexlify(os.urandom(8)).decode()
html = """<div id="%s"></div>
<script type="text/javascript">
$(function(){
var p = $("#%s");
if (p.length==0) return;
while (!p.hasClass("cell")) {
p=p.parent();
if (p.prop("tagName") =="body") return;
}
var cell = p;
cell.find(".input").addClass("hide-in-slideshow")
});
</script>""" % (uid, uid)
display.display_html(html, raw=True)
Estatística de Weibull - Uma Ferramenta de Conversão e Comparação de Comportamentos Mecânicos
Referências:¶
Wachtman, Mechanical Properties of Ceramics - Páginas 119-136
Kendall et al., Influence of toughness on Weibull modulus of ceramic bending strength
Zwaag et al., The concept of filament strength and the Weibull modulus
Glandus et al., Uncertainty on the mean strength and Weibull's modulus of an alumina batch as a function of the number of samples
Introdução¶
Alguns exemplos de aplicação da estatística de Weibull:
- Usada para a definição de apólice de seguros e planos de saúde
- Na previsão climática (na distrbiuição de velocidade dos ventos)
- Na Teoria de Valores Extremos (os "Cisnes Negros")
- Estimativa do tempo de vida de componentes de elevada importância (Oleodutos, dutos e brocas para extração de petróleo e gás e até mesmo componentes da Estação Espacial Internacional)
- Em situações onde elevada reprodutibilidade do comportamento é fundamental (próteses, componentes aeronáuticos, transportes de alta velocidade)
O comportamento mecânico dos materiais é fundamental para sua aplicação segura, seja em componentes de máquinas, ou quando usado como material estrutural.
Sendo assim, caracterizá-lo é de grande importância:
Em um ensaio, aplica-se a carga (ou deslocamento) no material, até atingir sua falha. Esta é caracterizada (para materiais frágeis fraturando de forma catastrófica) como a queda repentina da tensão.
#@title
# HIDDEN
m=3
E=400e9
sigma_0=220e6
sigma_f_1 = 220e6 # np.random.weibull(m, 1)*sigma_0
mpl.rcParams.update({'font.size': 22})
def f(epsilon_i):
epsilons_1, sigmas_1 = run_test(epsilon_i, E, sigma_f_1)
epsilon_max = sigma_f_1/E
bar = np.array([0.09, 0.1, 0.11])
plt.figure(figsize=(6.5, 6.5))
# plt.text(0.0925, 215, 'Amostra')
# plt.text(bar[0]*(1-epsilon_i*100)*0.9, 205, 'F')
# plt.text(bar[2]*(1+epsilon_i*100)*1.06, 205, 'F')
# plt.arrow(bar[0]*(1-epsilon_i*100), 200, -0.01, 0, head_width=10, head_length=0.005,
# fc='k', ec='k', zorder=4)
# plt.arrow(bar[2]*(1+epsilon_i*100), 200, 0.01, 0, head_width=10, head_length=0.005,
# fc='k', ec='k', zorder=4)
# if epsilon_i < epsilon_max:
# plt.plot([bar[0]*(1-epsilon_i*100), bar[1],
# bar[2]*(1+epsilon_i*100)], [200, 200, 200], lw=10, c='k', label='Amostra')
# if epsilon_i >= epsilon_max:
# plt.plot([bar[0]*(1-epsilon_i*100),
# bar[1]*(1-(epsilon_i-epsilon_max)*100)-0.0021], [200, 200], lw=10, c='k')
# plt.plot([bar[1]*(1+(epsilon_i-epsilon_max)*100)+0.0021,
# bar[2]*(1+epsilon_i*100)], [200, 200], lw=10, c='k', label='Amostra')
plt.plot(epsilons_1*100, sigmas_1*1e-6, c='Navy', lw=5)
plt.ylim(0, 300)
plt.xlim(0, 0.14e-2*100)
plt.ylabel('Tensão, $\sigma$ [MPa]')
plt.xlabel('Deformação, $\epsilon$ [%]')
plt.grid(lw=0.25)
# plt.legend(borderpad=1.2)
plt.show()
int_plot = interactive(f,
# epsilon_i=(0, 1e-4, 0.1e-5),
epsilon_i=widgets.FloatSlider(value=0,
min=0,
max=1e-3,
step=2e-5,
description='$\epsilon$',
readout_format='.2e',
orientation='horizontal')
)
output = int_plot.children[-1]
output.layout.height = '630px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Mas onde o material falhou?¶
Uma maneira de responder essa pergunta é através da análise das tensões em um corpo sujeito a flexão em 3 pontos.
Isto pode ser feito de maneira analítica (como na disciplina Mecânica dos Sólidos) ou de maneira numérica como através da simulação por elementos finitos abaixo:
Este resultado pode ser avaliado também pela deformação no corpo de prova, como também é possível observar na malha deformada abaixo:
Isso vale para um material perfeito (sem defeitos), porém, devemos lembrar que os defeitos são intrínsecos ao material! Sendo assim, a pergunta contínua:
Onde o material falhou?
Lembrando que defeitos concentram tensões, não é só o local da falha, mas também a resistência mecânica que dependerá das suas posições e características. Como exemplo, consideremos 3 possíveis defeitos:
O material pode romper em alguns possíveis cenários como:
- Romper na superfície inferior, onde a tração é máxima, assim que alcançar o limite
- Romper no defeito 1, caso a combinação de tensão trativa local, geometria e tamanho do defeito resultarem no σth
- Romper no defeito 2
- Romper no defeito 3
Em síntese:¶
O valor de resistência mecânica vai depender dos seguintes fatores:
- Tamanho (e geometria) do defeito
- Posição do defeito no corpo
- Orientação do defeito relativa ao corpo
Todas estas características serão dependentes de fatores determinísticos (parâmetros de processamento, matérias-primas, aditivos, ...), e estocásticos e assim, dentro de uma amostragem de n corpos de um mesmo material teremos uma distribuição de valores de resistência mecânica.
Assim, temos três consequências principais:
- O uso do σf médio não é apropriado para projetos de estruturas ou componentes
- A probabilidade de falha de uma peça maior é maior do que o de uma peça menor (comparação entre amostras de diferentes áreas)
- O valor de resistência mecânica obtido depende da técnica utilizada (já que diferentes técnicas solicitam diferentes volumes do material em tração)
Dessa forma, apresentaremos ferramentas estatísticas para tratar estes dados e poder quantificar o comportamento geral de um dado material.
Começaremos com uma breve revisão estatística e em seguida introduziremos análises passíveis de serem feitas usando estas ferramentas e os conceitos já vistos no curso.
Distribuições Estatísticas¶
Começemos considerando uma distribuição de uma variável contínua, no caso uma amostragem de 10000 pessoas nas quais se perguntou a sua altura. Nós podemos visualizar o resultado através de um histograma:
# HIDDEN
fig = plt.figure(figsize=(8, 8))
mu_3 = 173.1
sigma_3 = 2.74
plt.hist(stats.norm.rvs(loc=mu_3, scale=sigma_3, size=10000), 100, color='green', zorder=3, label=f'$\mu={mu_3}$cm, $\sigma={sigma_3}$cm')
plt.ylabel('Frequência')
plt.xlabel('Altura [cm]')
plt.legend(fontsize=14)
plt.grid()
plt.show()
Esta distribuição pode ser descrita por dois valores, a média, μ e o desvio padrão, σ. Seus efeitos no perfil do histograma de distribuição pode ser visto abaixo:
#@title
# HIDDEN
mu_3 = 173.1
sigma_3 = 7.5
dist_3 = stats.norm.rvs(loc=mu_3, scale=sigma_3, size=10000)
def f(sigma_1, mu_1):
plt.figure(figsize=(6.5, 6.5))
plt.hist(dist_3, color='green', zorder=3,
bins=np.arange(110, 250 + 2, 2),
label=f'$\mu={mu_3}$cm, $\sigma={round(sigma_3**0.5, 2)}$cm')
plt.hist(stats.norm.rvs(loc=mu_1, scale=sigma_1**2, size=10000),
color='navy', zorder=3,
bins=np.arange(110, 250 + 2, 2),
label=f'$\mu={mu_1}$cm, $\sigma={round(sigma_1, 2)}$cm', alpha=0.5)
plt.ylabel('Frequência')
plt.xlabel('Altura [cm]')
plt.legend(fontsize=14)
plt.grid()
plt.show()
int_plot = interactive(f,
# epsilon_i=(0, 1e-4, 0.1e-5),
mu_1=widgets.FloatSlider(value=173.1,
min=120,
max=210,
step=5,
description='$\mu$',
orientation='horizontal'),
sigma_1=widgets.FloatSlider(value=2.74,
min=0.25,
max=5.5,
step=0.25,
description='$\sigma$',
orientation='horizontal')
)
output = int_plot.children[-1]
output.layout.height = '630px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Nós podemos aproximar os dados coletados por uma curva de distribuição (ou função) de densidade de probabilidade, a FDP, levando a espessura da barra para zero.
Outra forma de visualizar estes resultados é através da função de distribuição acumulada, a FDA, que apresenta qual a fração da amostra que tem um valor igual ou menor do que determinado valor.
# HIDDEN
xs = np.linspace(150, 205, 1000)
fig, axs = plt.subplots(1, 2, figsize=(16, 6))
mu_3 = 173.1
sigma_3 = 7.5
normal_3 = stats.norm(loc=mu_3, scale=sigma_3)
axs[1].plot(xs, normal_3.cdf(xs), c='green', lw=4,
label=f'$\mu={mu_3}$cm, $\sigma={round(sigma_3**0.5, 2)}$cm')
axs[0].hist(normal_3.rvs(1000), density=True, bins=80, zorder=0)
axs[0].plot(xs, normal_3.pdf(xs), c='green', lw=4,
label=f'$\mu={mu_3}$cm, $\sigma={round(sigma_3**0.5, 2)}$cm')
axs[1].set_ylabel('Função de Distribuição \n Acumulada (FDA)')
axs[0].set_ylabel('Função Densidade de \n Probabilidade (FDP)')
for ax in axs:
ax.set_xlabel('Altura [cm]')
ax.grid()
ax.legend(fontsize=12)
plt.tight_layout(w_pad=2)
plt.show()
Efeito dos Parâmetros em uma Distribuição Normal¶
#@title
# HIDDEN
mu_3 = 173.1
sigma_3 = 7.5
xs = np.linspace(100, 250, 1000)
cdf_3 = stats.norm.cdf(xs, loc=mu_3, scale=sigma_3)
pdf_3 = stats.norm.pdf(xs, loc=mu_3, scale=sigma_3)
def f(sigma_1, mu_1):
fig, axs = plt.subplots(1, 2, figsize=(18, 8))
axs[1].plot(xs, cdf_3, c='Navy', lw=2,
label=f'$\mu={mu_3}$cm, $\sigma={round(sigma_3**0.5, 2)}$cm')
axs[0].plot(xs, pdf_3, c='Navy', lw=2,
label=f'$\mu={mu_3}$cm, $\sigma={round(sigma_3**0.5, 2)}$cm')
axs[1].plot(xs, stats.norm.cdf(xs, loc=mu_1, scale=sigma_1**2), c='red', lw=2,
label=f'$\mu={mu_1}$cm, $\sigma={round(sigma_1, 2)}$cm')
axs[0].plot(xs, stats.norm.pdf(xs, loc=mu_1, scale=sigma_1**2), c='red', lw=2,
label=f'$\mu={mu_1}$cm, $\sigma={round(sigma_1, 2)}$cm')
axs[1].set_ylabel('Função de Distribuição \n Acumulada (FDA)')
axs[0].set_ylabel('Função Densidade de \n Probabilidade (FDP)')
for ax in axs:
ax.set_xlabel('Altura [cm]')
ax.legend(loc='upper left', fontsize=14)
ax.grid()
plt.tight_layout()
plt.show()
int_plot = interactive(f,
# epsilon_i=(0, 1e-4, 0.1e-5),
sigma_1=widgets.FloatSlider(value=2.74,
min=0.25,
max=5.5,
step=0.25,
description='$\sigma$',
orientation='horizontal', layout={'fontsize': 20}),
mu_1=widgets.FloatSlider(value=173.1,
min=120,
max=210,
step=5,
description='$\mu$',
orientation='horizontal')
)
output = int_plot.children[-1]
output.layout.height = '450px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Estas funções podem ser interpoladas dos dados experimentais
Devido o fato de que um grande número de fenômenos serem representados por um pequeno número de funções analíticas, estas ditribuições possuem nomes próprios:
- Distribuição normal (Gaussiana)
- Distribuição de Pareto
- Distribuição Logística
- Distribuição de Cauchy
- Distribuição de Boltzmann
- Distribuição de Weilbull
A Distribuição de Weibull¶
A função de distribuição de Weibull é definida pela seguinte função de distribuição acumulada:
Pf(σ)={1−e−(σ−σuσ0)mσ>σu0σ≤σu
E sua função densidade de probabilidade é dada por:
p(σ)={−mσ0(σ−σuσ0)m−1e−(σ−σuσ0)mσ>σu0σ≤σu
Os gráficos a seguir demonstram a sua FDP e FDA:
#@title
# HIDDEN
sigma = np.linspace(0,50, 1000)
m = 15
sigma_0 = 20
weibull_pdf = lambda sigma, m, sigma_0: (m/sigma_0)*(sigma/sigma_0)**(m-1)*np.exp(-(sigma/sigma_0)**m)
weibull_cdf = lambda sigma, m, sigma_0: 1-np.exp(-(sigma/sigma_0)**m)
fig, axs = plt.subplots(1, 2, figsize=(18, 8))
ax1, ax2 = axs
ax1.plot(sigma, weibull_pdf(sigma, 50, 20), c='Navy', label='$\sigma_0$ =20, m = 50')
ax2.plot(sigma, weibull_cdf(sigma, 50, 20), c='Navy', label='$\sigma_0$ =20, m = 50', lw=1, ls='--')
ax1.plot(sigma, weibull_pdf(sigma, 5, 20), c='red', label='$\sigma_0$ =20, m = 5')
ax2.plot(sigma, weibull_cdf(sigma, 5, 20), c='red', label='$\sigma_0$ =20, m = 5', lw=1, ls='--')
# ax1.set_ylabel('Probability Density Function')
# ax2.set_ylabel('Cumulative Distribution Function')
ax2.set_ylabel('Função de Distribuição\n Acumulada (FDA)')
ax1.set_ylabel('Função Densidade \n de Probabilidade (FDP)')
for ax in [ax1, ax2]:
ax.set_xlim(0,50)
ax.set_ylim(0,1.05)
ax.legend(fontsize=14)
ax.grid()
plt.tight_layout()
plt.show()
Por que usar a distribuição de Weibull?¶
A distribuição de Weibull é qualitativamente próxima de uma distribuição normal:
#@title
# HIDDEN
sigma = np.linspace(0,50, 1000)
m = 15
sigma_0 = 20
weibull_pdf = lambda sigma, m, sigma_0: (m/sigma_0)*(sigma/sigma_0)**(m-1)*np.exp(-(sigma/sigma_0)**m)
weibull_cdf = lambda sigma, m, sigma_0: 1-np.exp(-(sigma/sigma_0)**m)
fig, axs = plt.subplots(1, 2, figsize=(18, 6))
ax1, ax2 = axs
ax1.plot(sigma, stats.weibull_min.pdf(sigma, 20, 0, 20), c='Navy', label='Weibull - $\sigma_0$ =20 MPa, m = 20')
ax2.plot(sigma, stats.weibull_min.cdf(sigma, 20, 0, 20), c='Navy', label='Weibull - $\sigma_0$ =20 MPa, m = 20', lw=1, ls='--')
ax1.plot(sigma, stats.norm.pdf(sigma, 20, 1.206**0.5), c='red', label='Normal - $\mu$ =20 MPa, $\sigma$ = 1.099')
ax2.plot(sigma, stats.norm.cdf(sigma, 20, 1.206**0.5), c='red', label='Normal - $\mu$ =20 MPa, $\sigma$ = 1.099', lw=1, ls='--')
# ax1.set_ylabel('Probability Density Function')
# ax2.set_ylabel('Cumulative Distribution Function')
ax2.set_ylabel('Função de Distribuição\n Acumulada (FDA)')
ax1.set_ylabel('Função Densidade \n de Probabilidade (FDP)')
for ax in [ax1, ax2]:
ax.set_xlim(12.5, 25)
ax.set_ylim(0,1.05)
ax.legend(loc='upper left', fontsize=16)
ax.grid()
ax1.set_ylim(0,0.450)
plt.tight_layout()
plt.show()
#
Então, por que usar a distribuição de Weibull?¶
- Em primeiro lugar, a distribuição normal é consideravelmente mais complexa:
| Distribuição | FDP | FDA |
|---|---|---|
| Weibull | p(σ)={mσ0(σ−σuσ0)m−1e−(σ−σuσ0)mσ>σu0σ≤σu |
Pf(σ)={1−e−(σ−σuσ0)mσ>σu0σ≤σu |
| Normal | p(x)=12[1+erf(x−μς√2)] |
P(x)=1ς√2πe−12(x−μς)2 |
Isto seria problemático para as etapas posteriores da Análise de Weibull!
- Em segundo lugar, observe que a distribuição normal admite valores negativos para a resistência mecânica!
Efeito dos Parâmetros na Distribuição de Weibull¶
#@title
# HIDDEN
m = 15
sigma_0 = 20
def weibull_view_1(m, sigma_0):
weibull_pdf = lambda sigma, m, sigma_0: (m/sigma_0)*(sigma/sigma_0)**(m-1)*np.exp(-(sigma/sigma_0)**m)
weibull_cdf = lambda sigma, m, sigma_0: 1-np.exp(-(sigma/sigma_0)**m)
sigmas = np.linspace(0, 150, 1000)
fig, axs = plt.subplots(1, 2, figsize=(18, 8))
axs[0].plot(sigmas, weibull_pdf(sigmas, m, sigma_0), c='Navy', label=f'$\sigma_0$ = {sigma_0}, m = {m}')
axs[0].set_xlim(0, 50)
# axs[0].set_ylim(0, 5)
axs[0].set_xlabel('Tensão, $\sigma$ [MPa]')
axs[0].set_ylabel('Função Densidade \n de Probabilidade (FDP)')
axs[1].plot(sigmas, weibull_cdf(sigmas, m, sigma_0), c='Navy', label=f'$\sigma_0$ = {sigma_0}, m = {m}')
axs[1].set_xlim(0, 50)
# axs[1].set_ylim(0, 1.25)
axs[1].set_xlabel('Tensão, $\sigma$ [MPa]')
axs[1].set_ylabel('Função de Distribuição \n Acumulada (FDA)')
plt.tight_layout()
for ax in axs:
ax.yaxis.set_major_formatter(FormatStrFormatter('%.2f'))
ax.grid()
ax.legend()
plt.show()
int_plot = interactive(weibull_view_1,
sigma_0=widgets.FloatSlider(value=5,
min=5,
max=45,
step=1,
description='$\sigma_0$',
readout_format='.2f'),
m=widgets.FloatSlider(value=2,
min=2,
max=200,
step=5,
description='$m$',
readout_format='.2f')
)
output = int_plot.children[-1]
output.layout.height = '450px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Qual o sentido físico dos parâmetros de Weibull?
Análise de Weibull Aplicada na Resistência Mecânica do Material¶
A análise de Weibull é inspirada na Teoria do Elo Mais Fraco (proposta por Aristóteles e revisitada por Leonardo Da Vinci):
"Uma corrente é tão forte quanto o for o seu elo mais fraco."
"Quanto maior a corrente, maior a é a probabilidade de que ela falhe quando uma carga é aplicada."
Esta teoria considera fraturas frágeis em uma escala macroscópica.
Todo material solicitado por uma carga está sujeito a uma probabilidade de falha, Pf
Esta, pode ser representada matematicamente como:
Pf=1−PSOnde PS é a probabilidade de sobrevivência do material.
Nós podemos assumir estas quantias como funções da tensão ao qual o material está sujeito, σ e devido a teoria do elo mais fraco, do volume V:
Pf(σ,V)=1−Ps(σ,V)Para conduzir a análise, consideremos o que ocorre em um determinado corpo de volume V+δV. Este corpo apenas sobreviverá se ambas as partes, V e δV sobreviverem, isto é, precisamos considerar a probabilidade combinada de ambos os cenários. Matematicamente:
Ps(σ,V+δV)=Ps(σ,V) Ps(σ,δV)Se considerarmos que a distribuição de probabilidades de uma pequena parte do corpo δV falhar é descrita por uma densidade de porbabilidade φ(σ) (isto é, uma função que retorna o valor de probabilidade de um volume infinitesimal falhar para uma dada aplicação de uma tensão local σ), temos:
Ps(σ,δV)=1−Pf(σ,δV)=1−φ(σ) δVSubstituindo na Equação 5, e rearranjando temos: Ps(σ,V+δV)−Ps(σ,V)δV=−φ(σ) Ps(σ,V)
Se reduzirmos o valor de δV podemos aproximar a derivada da probabilidade de sobrevivência em função do volume:
∂Ps(σ,V)∂V=−φ(σ) Ps(σ,V)Lembrando que a tensão depende da posição e integrando é possível obter Ps:
Ps(σ,V)=exp{−∫Vφ[σ(r)]dV}Usando a Equação 3, podemos achar a probabilidade de falha como sendo: Pf(σ,V)=1−Ps(σ,V)=1−exp{−∫Vφ[σ(r)]dV}
Para um caso mais simples no qual a tensão independe da posição (por exemplo na caso de uma amostra em tração uniaxial), a Equação 10 se torna: Pf(σ,V)=1−Ps(σ,V)=1−exp{−V φ(σ)}
Weibull assumiu que φ(σ) tinha a seguinte forma:
φ(σ)=(σΣ0)mFinalmente chegando que Pf é dado por:
Pf=1−exp[−V(σΣ0)m]Que é a função de distribuição acumulada de Weibull com parâmetro de escala σ0 dado por:
σ0=Σ0V1mEfeito do Volume¶
Nós iremos avaliar o efeito do volume de uma amostra, usando a probabilidade de falha. Da Equação 13, para uma amostra a de um determinado material com volume Va e tensão de ruptura σa, temos:
Pf,a=1−exp[−Va(σaΣ0)m]De maneira análoga para uma amostra b do mesmo material temos: Pf,b=1−exp[−Vb(σbΣ0)m]
Podemos considerar um mesmo valor de probabilidade de falha igualando as Equações 15 e 16:
exp[−Va(σaΣ0)m]=exp[−Vb(σbΣ0)m]Se aplicarmos −log em ambos os lados temos:
Va(σaΣ0)m=Vb(σbΣ0)mComo Σ0 independe do volume, podemos cortar o (1Σ0)m de ambos os lados, chegando assim em:
Vaσma=VbσmbE rearranjando temos: σaσb=(VbVa)1m
É possível comparar tais quantias em um gráfico, podemos ver o efeito de reduzir o tamanho da amostra:
# HIDDEN
VbVa = np.logspace(1, -5, 1000)
m_1 = 1
m_2 = 10
m_3 = 100
m_4 = 500
fig = plt.figure(figsize=(12, 9))
plt.plot(VbVa, VbVa**(1/m_1), c='navy', label=f'm={m_1}', lw=3)
plt.plot(VbVa, VbVa**(1/m_2), c='orangered', label=f'm={m_2}', lw=3)
plt.plot(VbVa, VbVa**(1/m_3), c='y', label=f'm={m_3}', lw=3)
plt.plot(VbVa, VbVa**(1/m_4), c='k', label=f'm={m_4}', lw=3)
plt.ylim(-0.2, 1.05)
plt.xlim(1, 1e-5)
plt.gca().set_xscale('log')
plt.text(0.002, 0.15, 'Aumentando o valor de $m$ $\longrightarrow$', rotation=50)
plt.ylabel(r'$\frac{\sigma_a}{\sigma_b}$')
plt.xlabel(r'$\frac{V_b}{V_a}$')
plt.text(0.02, -0.15, 'Diminuindo $V_b$ →')
plt.text(0.8, -0.175, '← Aumentando $\sigma_b$', rotation=90)
plt.legend(loc='lower right', fontsize=14)
plt.grid()
plt.show()
Efeito do Ensaio¶
Nós podemos considerar o efeito do Ensaio na resistência mecânica considerando a forma mais geral da Probabilidade de falha, dada na Equação 10, e considerando a que a tensão pode ser descrita como:
σ(r)=σmaxg(r)Substituindo na Equação 10, temos: Pf=1−exp[−(σmaxΣ0)m∫Vgm(r)dV]
E assim podemos definir o fator do ensaio k≤1:
k=1V∫Vgm(r)dV ou k=1V∫V(σ(r)σmax)mdVDessa maneira temos que a Equação 22 pode ser reescrita como: Pf=1−exp[−kV(σmaxΣ0)m]
Calculando o fator do ensaio de tração para um paralelepípedo com volume V=a b c temos que σ(r)=σmax e portanto kt é dado por:
kt=1V∫V1dV=1V∭dxdydz=a b cV=1Para uma amostra em Flexão em três pontos, um resultado da Mecânica Aplicada indica que: σ(y)σmax=4 x yL h
Assim, considerando que o volume da amostra solicitado (isto é, o volume entre os apoios) é V=b L h:
k3pt=1V∫L0(2xL)mdx∫h/20(2yh)mdy∫b0dz=L h bV 2 (m+1)2=12(m+1)2Adotando uma estratégia similar àquela utilizada para a análise do efeito do volume, igualando as probabilidade de falha para um ensaio em tração uniaxial e em flexão três pontos chegamos em:
σf,3ptσf,t=(ktk3pt)1/m=[2(m+1)2]1/mUsando esta metodologia analítica para o ensaio de quatro pontos resulta em: σf,4ptσf,t=(ktk4pt)1/m=[4(m+1)2m+2]1/m
Esses comportamentos podem ser visualizados a seguir:
# HIDDEN
ms = np.linspace(1, 500, 10000)
fig = plt.figure(figsize=(8, 6))
three_t = lambda m: (2*(m+1)**2)**(1/m)
four_t = lambda m: (4*(m+1)**2/(m+2))**(1/m)
plt.plot(ms, three_t(ms), lw=3, label=r'$\frac{\sigma_{f, 3pt}}{\sigma_{f,t}}$', c='navy')
plt.plot(ms, four_t(ms), lw=3, label=r'$\frac{\sigma_{f, 4pt}}{\sigma_{f,t}}$', c='orangered')
plt.ylabel(r'$\frac{\sigma_{f, b}}{\sigma_{f, t}}$')
plt.xlabel(r'$m$')
plt.gca().set_xscale('log')
plt.grid()
plt.legend()
plt.show()
Como Obter os Parâmetros de Weibull?¶
Uma forma é através do Gráfico de Weibull. Ele é uma linearização da probabilidade de falha (como a teoria de cristalização de Avrami ou a lei de potências de um fluído não-newtoniano).
Para tanto, utiliza-se a aplicação da função log duplamente (uma vez para reduzir a exponenciação e uma segunda para a potência m):
#@title
# HIDDEN
weibull_pdf = lambda sigma, m, sigma_0: np.exp(-(sigma/sigma_0)**m)
sigma_0 = 20e6
sigmas_f_m5 = np.sort(np.random.weibull(5, 60)*sigma_0)
sigmas_f_m15 = np.sort(np.random.weibull(15, 60)*sigma_0)
S_m5 = weibull_pdf(sigmas_f_m5, 5, sigma_0)
S_m15 = weibull_pdf(sigmas_f_m15, 15, sigma_0)
fig, axs = plt.subplots(1, 2, figsize=(16, 6))
axs[0].plot(sigmas_f_m5*1e-6, S_m5, '-o', c='navy', label='m = 5')
axs[0].plot(sigmas_f_m15*1e-6, S_m15, '-o', c='coral', label='m = 15')
axs[0].axhline(0, c='k', zorder=0, lw=1)
axs[0].axvline(0, c='k', zorder=0, lw=1)
axs[0].set_ylabel(r'$P_S$', fontsize=15)
axs[0].set_xlabel(r'$\frac{\sigma}{\sigma_0}$', fontsize=15)
axs[1].scatter(np.log(sigmas_f_m5 / sigma_0), np.log(np.log(1 / S_m5)), c='navy', label='m = 5')
axs[1].scatter(np.log(sigmas_f_m15 / sigma_0), np.log(np.log(1 / S_m15)), c='coral', label='m = 15')
axs[1].axhline(0, c='k', zorder=0, lw=1)
axs[1].axvline(0, c='k', zorder=0, lw=1)
axs[1].set_ylabel(r'$\ln{\ln{\ \frac{1}{P_S}}}$', fontsize=15)
axs[1].set_xlabel(r'$\ln{\ \frac{\sigma}{\sigma_0}}$', fontsize=15)
for ax in axs:
for side in ['top','right','bottom','left']:
ax.spines[side].set_visible(False)
ax.tick_params(axis='both',which='both',labelbottom=True,bottom=False,left=False)
ax.legend(fontsize=15)
axs[0].set_xlim(0, 30)
plt.tight_layout(w_pad=2)
plt.show()
Como queremos obter o módulo de Weibull a partir de uma interpolação linear, podemos tirar proveito da propriedade de divisão do argumento da função log para separar o coeficiente angular (a) e o coeficiente linear (b) da reta a ser interpolada:
E assim podemos encontrar o módulo de Weibull, m bem como o parâmetro σ0 da seguinte maneira:
#@title
# HIDDEN
def weibull_view_2(m_1, sigma_1):
weibull_cdf = lambda sigma, m, sigma_0: np.exp(-(sigma/sigma_0)**m)
sigma_0 = 20e6
sigmas_f_m5 = np.sort(np.random.weibull(5, 60)*sigma_0)
sigmas_f_m15 = np.sort(np.random.weibull(m_1, 60)*sigma_1*1e6)
S_m5 = weibull_cdf(sigmas_f_m5, 5, sigma_0)
S_m15 = weibull_cdf(sigmas_f_m15, m_1, sigma_1*1e6)
fig, axs = plt.subplots(1, 2, figsize=(8, 4.5))
axs[0].plot(sigmas_f_m5*1e-6, S_m5, '-o', c='navy', label='$\sigma_0$ = 20MPa, m = 5')
axs[0].plot(sigmas_f_m15*1e-6, S_m15, '-o', c='coral', label=f'$\sigma_0$ = {round(sigma_1, 2)}MPa, m = {m_1}')
axs[0].axhline(0, c='k', zorder=0, lw=1)
axs[0].axvline(0, c='k', zorder=0, lw=1)
axs[0].set_ylabel(r'$P_S$', fontsize=15)
axs[0].set_xlabel(r'$\frac{\sigma}{\sigma_0}$', fontsize=15)
axs[1].scatter(np.log(sigmas_f_m5 / sigma_0), np.log(np.log(1 / S_m5)), c='navy', label='$\sigma_0$ = 20MPa, m = 5')
axs[1].scatter(np.log(sigmas_f_m15 / (sigma_1*1e6)), np.log(np.log(1 / S_m15)), c='coral', label=f'$\sigma_0$ = {round(sigma_1, 2)}MPa, m = {m_1}')
axs[1].axhline(0, c='k', zorder=0, lw=1)
axs[1].axvline(0, c='k', zorder=0, lw=1)
axs[1].set_ylabel(r'$\ln{\ln{\ \frac{1}{P_S}}}$', fontsize=15)
axs[1].set_xlabel(r'$\ln{\ \frac{\sigma}{\sigma_0}}$', fontsize=15)
for ax in axs:
for side in ['top','right','bottom','left']:
ax.spines[side].set_visible(False)
ax.tick_params(axis='both',which='both',labelbottom=True,bottom=False,left=False)
ax.legend(loc='upper left', fontsize=10)
# axs[0].set_xlim(0, 30)
plt.tight_layout(w_pad=2)
plt.show()
int_plot = interactive(weibull_view_2,
sigma_1=widgets.FloatSlider(value=5,
min=5,
max=45,
step=1,
description='$\sigma_0$ [MPa]',
readout_format='.2f'),
m_1=widgets.FloatSlider(value=2,
min=1,
max=80,
step=4,
description='$m$ [-]',
readout_format='.2f')
)
output = int_plot.children[-1]
output.layout.height = '580px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
#@title
# HIDDEN
def weibull_view_3(m_1, sigma_1):
weibull_cdf = lambda sigma, m, sigma_0: np.exp(-(sigma/sigma_0)**m)
sigma_0 = 20e6
sigmas_f_m5 = np.sort(np.random.weibull(5, 60)*sigma_0)
sigmas_f_m15 = np.sort(np.random.weibull(m_1, 60)*sigma_1*1e6)
S_m5 = weibull_cdf(sigmas_f_m5, 5, sigma_0)
S_m15 = weibull_cdf(sigmas_f_m15, m_1, sigma_1*1e6)
fig, axs = plt.subplots(1, 2, figsize=(8, 4.5))
axs[0].plot(sigmas_f_m5*1e-6, S_m5, '-o', c='navy', label='$\sigma_0$ = 20MPa, m = 5')
axs[0].plot(sigmas_f_m15*1e-6, S_m15, '-o', c='coral', label=f'$\sigma_0$ = {round(sigma_1, 2)}MPa, m = {m_1}')
axs[0].axhline(0, c='k', zorder=0, lw=1)
axs[0].axvline(0, c='k', zorder=0, lw=1)
axs[0].set_ylabel(r'$P_S$', fontsize=15)
axs[0].set_xlabel(r'$\frac{\sigma}{\sigma_0}$', fontsize=15)
axs[1].scatter(np.log(sigmas_f_m5), np.log(np.log(1 / S_m5)), c='navy', label='$\sigma_0$ = 20MPa, m = 5')
axs[1].scatter(np.log(sigmas_f_m15), np.log(np.log(1 / S_m15)), c='coral', label=f'$\sigma_0$ = {round(sigma_1, 2)}MPa, m = {m_1}')
axs[1].axhline(0, c='k', zorder=0, lw=1)
axs[1].axvline(0, c='k', zorder=0, lw=1)
axs[1].set_ylabel(r'$\ln{\ln{\ \frac{1}{P_S}}}$', fontsize=15)
axs[1].set_xlabel(r'$\ln{\sigma}$', fontsize=15)
for ax in axs:
for side in ['top','right','bottom','left']:
ax.spines[side].set_visible(False)
ax.tick_params(axis='both',which='both',labelbottom=True,bottom=False,left=False)
ax.legend(loc='upper left', fontsize=10)
axs[1].set_xlim(10, 20)
plt.tight_layout(w_pad=2)
plt.show()
int_plot = interactive(weibull_view_3,
sigma_1=widgets.FloatSlider(value=5,
min=5,
max=45,
step=1,
description='$\sigma_0$ [MPa]',
readout_format='.2f'),
m_1=widgets.FloatSlider(value=5,
min=1,
max=80,
step=4,
description='$m$ [-]',
readout_format='.2f')
)
output = int_plot.children[-1]
output.layout.height = '580px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Entretanto, para este plot é necessário obter a probabilidade de sobrevivência cada corpo. Vejamos como obter tal quantia. Considere um ensaio com um número n de amostras a partir do qual será obtido o módulo de Weibull. Iniciamos listando os dados brutos em uma tabela e em seguida ordenamos do menor para o maior, ranqueando os dados:
# HIDDEN
# nbi:left
pd.set_option('display.max_rows', None)
sigma_0 = 20e6
weibull_m5 = stats.weibull_min(5, 0, sigma_0)
N=60
sigmas_f_m5 = weibull_m5.rvs(size=N)
# sigmas_f_m5 = pd.read_clipboard(header=None).values.flatten()
df_1 = pd.DataFrame([sigmas_f_m5*1e-6]).T
df_1.columns = ['$\sigma_f$ [MPa]']
df_1.index.name = 'Amostra'
print('Dados Brutos:')
df_1.round(decimals=2)
# HIDDEN
# nbi:right
print('Dados Ordenados:')
df_1 = df_1.sort_values('$\sigma_f$ [MPa]')
df_1['Rank, i'] = np.arange(1, N+1)
df_1.round(decimals=2)
Assim é possível estimar Pf a partir do rank i seguindo a Equação abaixo:
Pf(i)=i−0.3n+0.4Note que é a partir desta etapa em que podemos categoricamente afirmar qual é a probabilidade do material falhar para um dado estado e valor de tensão, conforme ilustrado no Widget a seguir:
# HIDDEN
m_1 = 15
sigma_1 = 20
weibull_cdf = lambda sigma, m, sigma_0: np.exp(-(sigma/sigma_0)**m)
sigma_0 = 20e6
sigmas_f_m5 = np.sort(np.random.weibull(5, 100)*sigma_0)
sigmas_f_m15 = np.sort(np.random.weibull(15, 100)*sigma_1*1e6)
def weibull_view_2(percent):
S_m5 = 1-weibull_cdf(sigmas_f_m5, 5, sigma_0)
S_m15 = 1-weibull_cdf(sigmas_f_m15, 15, sigma_1*1e6)
fig, axs = plt.subplots(1, 1, figsize=(14*0.85, 8*0.85))
axs.plot(sigmas_f_m5*1e-6, S_m5 * 100, '-o', c='navy',
label='$\sigma_0$ = 20MPa, m = 5')
axs.plot(sigmas_f_m15*1e-6, S_m15 * 100, '-o', c='coral',
label=f'$\sigma_0$ = {round(sigma_1, 2)}MPa, m = {15}')
sig_percent_m5 = np.interp(percent, S_m5[::1] * 100, sigmas_f_m5[::1])
sig_percent_m15 = np.interp(percent, S_m15[::1] * 100, sigmas_f_m15[::1])
axs.plot(sig_percent_m5 * 1e-6, percent, 'o',
markerfacecolor='w', markersize=12, markeredgewidth=4, c='gold', zorder=150)
axs.plot(sig_percent_m15 * 1e-6, percent, 'o',
markerfacecolor='w', markersize=12, markeredgewidth=4, c='gold', zorder=150)
axs.axvline(sig_percent_m5 * 1e-6, lw=2, color='navy')
axs.axvline(sig_percent_m15 * 1e-6, lw=2, color='coral')
axs.axhline(0, c='k', zorder=0, lw=1)
axs.axvline(0, c='k', zorder=0, lw=1)
axs.set_ylabel(r'$P_f$ [%]', fontsize=20)
axs.set_xlabel(r'$\sigma$ [MPa]', fontsize=20)
axs.axhline(percent, 0.03, color='k', lw=2, ls='-')
axs.set_xlim(-1.5, 30)
axs.set_ylim(-5, 102.5)
bb = mpl.transforms.Bbox([[2, 13.5], [10, 50]])
fancy = mpl.patches.FancyBboxPatch(bb.p0, bb.width, bb.height,
fc='w', ec='0.8', alpha=0.8, zorder=100)
axs.add_patch(fancy)
axs.text(6, 32, f'Probabilidade de\nFalha:\n$P_f$ = {round(percent, 2)}%',
horizontalalignment='center', zorder=100, fontsize=20)
for i in range(3):
axs.text(6.1, 24, f'$\sigma_f$ = {round(sig_percent_m5 * 1e-6, 2)} MPa',
horizontalalignment='center', zorder=100, color='navy', fontsize=20)
axs.text(6.1, 16, f'$\sigma_f$ = {round(sig_percent_m15 * 1e-6, 2)} MPa',
horizontalalignment='center', zorder=100, color='coral', fontsize=20)
for ax in [axs]:
for side in ['top','right','bottom','left']:
ax.spines[side].set_visible(False)
ax.tick_params(axis='both', which='both', labelbottom=True, bottom=False, left=False)
axs.legend(loc='upper left', bbox_to_anchor=(0.0907, 0.7), fontsize=16)
plt.tight_layout(w_pad=2)
plt.show()
int_plot = interactive(weibull_view_2,
percent=widgets.FloatSlider(value=50,
min=0,
max=100,
step=5,
description='$P_f$',
readout_format='.2f')
)
output = int_plot.children[-1]
output.layout.height = '550px'
int_plot.layout = widgets.Layout(display='flex',
flex_flow='column-reverse',
align_items='center',
align_content='center',
justify_content='center',
width='100%')
int_plot
Observe também como o módulo de Weibull afeta o valor da tensão máxima (isto é, a tensão para Pf = 100%) suportada para materiais de mesmo valor de σ0: a maior dispersão do material com menor módulo resulta em maiores valores de tensão máxima. Esse "melhor desempenho", porém, é acompanhado de uma maior dispersão dos resultados, o que também afeta negativamente o valor de tensão que garanta máxima sobrevivência do material (por exemplo quando PS = 100% - Pf = 99%)
Uma vez obtida a probabilidade de falha é possível realizar a linearização do gráfico conforme ilustrado a seguir:
# HIDDEN
df_1['P_s'] = 1 - (df_1['Rank, i']-0.3)/(N+0.4)
df_1[r'lnln($\frac{1}{P_s}$)'] = np.log(np.log(1/df_1['P_s']))
df_1[r'ln($\sigma_f$)'] = np.log(df_1['$\sigma_f$ [MPa]'])
df_1_r = df_1.round(decimals=2)
df_1_r['P_s'] = df_1['P_s'].round(decimals=3)
df_1_r
# HIDDEN
def r_squared(y, y_hat):
y_bar = y.mean()
ss_tot = ((y-y_bar)**2).sum()
ss_res = ((y-y_hat)**2).sum()
return 1 - (ss_res/ss_tot)
plt.figure(figsize=(7, 7))
plt.plot(df_1[r'ln($\sigma_f$)'], df_1[r'lnln($\frac{1}{P_s}$)'], 'o', label='Dados Experimentais', c='navy', lw=3)
m_interp, sig_0_interp = np.polyfit(df_1[r'ln($\sigma_f$)'], df_1[r'lnln($\frac{1}{P_s}$)'], 1)
plt.plot(df_1[r'ln($\sigma_f$)'], df_1[r'ln($\sigma_f$)']*m_interp+sig_0_interp, '-', label='Interpolação', c='orangered', lw=3)
sigs = np.linspace(1e6, 40e6, 1000)
# plt.plot(np.log(sigs*1e-6), np.log(np.log(1/(1-weibull_m5.cdf(sigs)))), '-k', label='Resultado Analítico', lw=3)
plt.text(2.8, -1.95, '$m_{exp}$ = ' + str(round(m_interp, 2)))
plt.text(2.8, -2.3, '$\sigma_{0, exp}$ = ' + str(round(np.exp(-sig_0_interp/m_interp), 2)) + ' MPa')
plt.text(2.8, -1.6, '$r^2$ = ' + str(round(r_squared(df_1[r'lnln($\frac{1}{P_s}$)'], df_1[r'ln($\sigma_f$)']*m_interp+sig_0_interp), 2)))
plt.ylabel(r'$\log{\log{\ \frac{1}{P_S}}}$', fontsize=15)
plt.xlabel(r'$\log{\ \sigma}$', fontsize=15)
plt.xlim(2.25, 3.5)
plt.ylim(-5, 4)
plt.legend(fontsize=16)
plt.grid()
plt.show()
Alguns Comentários Sobre o Módulo de Weibull:¶
- Qual o número de corpos de prova que compõe uma amostra para se estimar o módulo de Weibull?
O gráfico a seguir ilustra como o erro na obtenção do módulo de Weibull diminui com o aumento do número de amostras consideradas:
# HIDDEN
plt.figure(figsize=(10, 5.5))
Ns = np.linspace(1, 100, 1000)
m_1 = 1
m_5 = 5
plt.plot(Ns, 1 / (2*Ns)**0.5*100, c='navy', lw=2)
plt.text(20, 50, r'Módulo de Weibull = $m \pm \dfrac{m}{\sqrt{2 \ N}}$')
plt.text(99, 10, str(int((1 / (2*Ns)**0.5*100)[-1]))+'%')
plt.grid()
plt.ylabel('$e$ [%]')
plt.xlabel('Número de corpos de prova')
plt.show()
Note que a partir de 100 amostras o erro cai para 7%, enquanto o uso de 50 amostras leva a um erro de 10%. Estas estimativas teóricas também são confirmadas experimentalmente conforme o trabalho de Glandus et al.
- O módulo de Weibull não está relacionado apenas ao material, mas é influenciado também pela repetibilidade dos CPs
- É necessário se atentar que o comportamento da amostra (de seguir uma distribuição de probabilidade de falha igual ao de Weibull) foi ASSUMIDO como sendo o comportamento real. Este caráter fenomenológico é uma das principais críticas à analise de Weibull
Qual Material Possui Maior m?¶
# HIDDEN
fig, axs = plt.subplots(3, 3, figsize=(20, 4))
for ax in axs.flatten():
ax.tick_params(left=False,
bottom=False,
labelleft=False,
labelbottom=False)
axs[0, 0].set_title('Material A', fontsize=16)
axs[0, 1].set_title('Material B', fontsize=16)
axs[0, 2].set_title('Material C', fontsize=16)
for i in range(7):
axs[0,0].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
axs[1,0].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
axs[2,0].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
for i in range(10):
axs[0,1].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
for i in range(5):
axs[1,1].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
for i in range(30):
axs[2,1].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
for i in range(5):
axs[0,2].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
for i in range(2):
axs[2,2].scatter(np.random.normal(), np.random.normal(), marker='$/$', s=200, c='k')
plt.show()
O Material A é aquele cuja dispersão de defeitos é a mais homogênea tanto em quantidade quanto em distribuição espacial.
Logo, é aquele de maior reprodutibilidade, e assim, maior módulo de Weibull.
Observações Experimentais¶
Valores comuns de módulo de Weibull¶
| Materiais Cerâmicos | Materiais Metálicos* | Materiais Poliméricos* |
|---|---|---|
| 5<m<15 | 100<m<500 | 50<m<150 |
Exemplo da influência do processamento no módulo de Weibull¶
| Concreto com Cimento Portland | Concreto com Cimento Portland + Dispersante |
|---|---|
| m=3,5 | m=20 |
O concreto (seja ele de construção civil ou refratário) é um material composto por diferentes fases destribuídas em sua microestrutura em diversas escalas de tamanho, como revelado pela tomografia de raios-X apresentada abaixo:
É possível visualizar poros, agregados e a matriz. Adicionalmente, em escalas menores, dentro da própria matriz seria possível encontrar microporos, agregados finos e fases hidratadas do cimento, ou ainda particulados que não reagiram com a água adicionada. Nos poros, também identifica-se filmes líquidos de água que não foi consumida durante sua cura. Assim, dada sua complexidade e heterogeneidade, medidas que promovam uma distribuição mais homogenêa destas fases em todo seu volume, como o uso de dispersantes, garantem uma maior reprodutibilidade e, assim, um maior módulo de Weibull.
O que acontece com o m quando uma amostra a verde é queimada?¶
O trabalho de Kendall et al. demonstra que a queima de corpos a verde não altera significativamente o módulo de Weibull, embora o KIC do material tenha sido ampliado mais de sete vezes. Esse resultado é de extrema importância pois além de demonstrar que processos que atuem de forma homogênea no material (como uma queima realizada em um forno com baixo gradiente térmico) não são capazes de alterar a distribuição dos defeitos nos materiais, ilustram que não há nenhuma relação entre o módulo de Weibull, m, e a tenacidade a fratura dos materiais, KIC.
Materiais de maior tenacidade à fratura possuem a fama de ser mais confiáveis. Por quê?¶
Efeito dos danos ao longo do tempo de vida do material alteram o módulo de Weibull, especialmente para aqueles materiais onde tais danos ocorrem de maneira mais fácil (menor tenacidade).
Considerações Finais¶
- O comportamento mecânico dos materiais é sucetível a diversas fontes de variabilidade (até mesmo erros de medida!)
- Representações matemáticas podem descrever este comportamento (médias, desvios, distribuições estatísticas)
- A probabilidade de falha de um material pode ser descrita por tais distribuições
- Em sua derivação, conceitos clássicos como a teoria do elo mais fraco geraram insights de grande importância como o efeito do volume e do ensaio
- O módulo de Weibull descreve a variabilidade de um determinado conjunto de amostras e serve para:
- Poder comparar materiais distintos e obtidos por processamentos diferentes
- Converter resultados entre corpos de prova de diferentes dimensões e de ensaios distintos
Visualização Prática do Módulo de Weibull¶
O Widget abaixo apresenta todas as etapas do processo de obtenção do Módulo de Weibull. Tente explorar os parâmetros (mova os sliders para os resultados aparecerem)!
Esta análise pode ser realizada em ensaios de flexão 3 ou 4 pontos com amostras prismáticas considerando uma aplicação final para uma barra sob tração.
Lembre-se das equações de conversão:
σf,t=[2(m+1)2]−1m σf,3pt# HIDDEN
import warnings
warnings.filterwarnings("ignore", category=RuntimeWarning)
pd.set_option('display.max_rows', 10)
# Wayne Nelson (2004) Applied Life Data Analysis
# ASTM Pf = (i-0.5) / N
np.random.seed(seed=233423)
estimador = 'Wayne'
def weibull_experiment(N, m, sigma_0):
global weibull_m
weibull_m = stats.weibull_min(m, 0, sigma_0)
global sigmas_f_m
sigmas_f_m = weibull_m.rvs(size=int(N), random_state=233423)
global df_1
df_1 = pd.DataFrame([sigmas_f_m]).T
df_1.columns = ['σ_f [MPa]']
df_1.index.name = 'Amostra'
df_1 = df_1.sort_values('σ_f [MPa]')
df_1['Rank, i'] = np.arange(1, N+1)
if estimador=='Wayne':
df_1['P_s'] = 1 - (df_1['Rank, i']-0.3)/(N+0.4)
elif estimador=='ASTM':
df_1['P_s'] = 1 - (df_1['Rank, i']-0.5)/(N)
df_1[r'lnln(1/Ps)'] = np.log(np.log(1/df_1['P_s']))
df_1[r'ln(σ_f)'] = np.log(df_1['σ_f [MPa]'])
df_1 = df_1.round(decimals=2)
display(df_1)
def weibull_plot(N, m, sigma_0):
global weibull_m
global sigmas_f_m
global df_1
with np.errstate(divide='ignore'):
fig = plt.figure(figsize=(6, 6))
plt.plot(df_1[r'ln(σ_f)'], df_1[r'lnln(1/Ps)'], 'o', label='Dados Experimentais', c='navy', lw=3)
x = df_1[r'ln(σ_f)']
y = df_1[r'lnln(1/Ps)']
m_interp, sig_0_interp = np.polyfit(x, y, 1)
r_squared = 1 - (sum((y - (m_interp * x + sig_0_interp))**2) / ((len(y) - 1) * np.var(y, ddof=1)))
plt.plot(df_1[r'ln(σ_f)'], df_1[r'ln(σ_f)']*m_interp+sig_0_interp, '-', label='Interpolação', c='orangered', lw=3)
sigs = np.logspace(0, 6, 1000)
plt.plot(np.log(sigs), np.log(np.log(1/(1-weibull_m.cdf(sigs)))), '-k', label='Resultado Analítico', lw=3)
plt.text(0.1, 4.3, '$m_{exp}$ = ' + str(round(m_interp, 2)), fontsize=13)
plt.text(0.1, 3.6, '$\sigma_{0, exp}$ = ' + str(round(np.exp(-sig_0_interp/m_interp), 2)) + ' MPa', fontsize=13)
plt.text(0.1, 2.6, '$r^2$ = ' + str(round(r_squared, 2)), fontsize=13)
plt.ylabel(r'$\log{\log{\ \frac{1}{P_S}}}$', fontsize=15)
plt.xlabel(r'$\log{\sigma}$', fontsize=15)
plt.xticks(fontsize=16)
plt.yticks(fontsize=16)
plt.xlim(0, 6)
plt.ylim(-10, 5)
plt.legend(loc='lower right', fontsize=12)
plt.grid()
plt.show()
# select_est = widgets.Select(options=['Wayne', 'ASTM'])
slider_m = widgets.FloatSlider(value=5,
min=1,
max=100,
step=1,
description='$m$',
orientation='horizontal')
slider_sig_0 = widgets.FloatSlider(value=20,
min=2,
max=100,
step=5,
description='$\sigma_0$',
orientation='horizontal')
slider_N = widgets.FloatSlider(value=60,
min=5,
max=1000,
step=5,
description='# CP\'s, $N$',
orientation='horizontal')
pd_show_int = interactive(weibull_experiment, N=slider_N, m=slider_m, sigma_0=slider_sig_0)
weibull_plot_show_int = interactive(weibull_plot, N=slider_N, m=slider_m, sigma_0=slider_sig_0)
eqs_header = widgets.HTML(r'<p>Equação de Weibull:</p>')
equations_weib = widgets.HTMLMath(r'$\ln \ln \left(\frac{1}{1-P_{f}}\right) = m \ln\left( \frac{\sigma}{\sigma_{0}}\right)$')
eqs_prob = widgets.HTMLMath(r'$P_f(i) = \frac{i-0.3}{n+0.4}$')
eqs = widgets.VBox([eqs_header, equations_weib, eqs_prob])
# eqs.layout = widgets.Layout(display='flex',
# flex_flow='column',
# align_items='center',
# align_content='center')
controls = widgets.HBox([widgets.VBox(pd_show_int.children[:3]), eqs])
pd_out = pd_show_int.children[-1]
#HIDDEN
three_t = lambda m: (2*(m+1)**2)**(1/m)
four_t = lambda m: (4*(m+1)**2/(m+2))**(1/m)
space = widgets.HTMLMath('   ')
def vol_ensaios(m, X_pec, Y_pec, Z_pec, X_am, Y_am, Z_am, ensaio, N, sigma_0):
global sigmas_f_m
mean_sig = np.mean(sigmas_f_m)
V_amostra = X_am * 1e-3 * Y_am * 1e-3 * Z_am * 1e-3
V_peca = X_pec * Y_pec * Z_pec
if ensaio == 'Flexão 3 Pontos':
# disp_mean_sig = widgets.HTMLMath(f'σ_(f, 3pt, amostra): {round(mean_sig, 2)} MPa')
mean_sig_trac = mean_sig * (2*(m+1)**2)**(-1/m)
# disp_mean_sig_trac = widgets.HTMLMath(f'σ_(f, t, amostra): {round(mean_sig_trac, 2)} MPa')
elif ensaio == 'Flexão 4 Pontos':
# disp_mean_sig = widgets.HTMLMath(f'σ_(f, 4pt, amostra): {round(mean_sig, 2)} MPa')
mean_sig_trac = mean_sig * (4*(m+1)**2/(m+2))**(-1/m)
# disp_mean_sig_trac = widgets.HTMLMath(f'σ_(f, t, amostra): {round(mean_sig_trac, 2)} MPa')
mean_sig_peca = (V_amostra/V_peca)**(1/m) * mean_sig_trac
disp_mean_sig_peca = widgets.HTMLMath(f'σ_(f, t, peça): {round(mean_sig_peca, 6)} MPa')
# V_am = widgets.HTMLMath(f'Volume da CP: {round(V_amostra*1e6, 2)} cm³')
# V_pec = widgets.HTMLMath(f'Volume da Peça: {round(V_peca, 2)} m³')
# space_A = widgets.HTMLMath(' → ')
# display(widgets.VBox([widgets.HBox([V_am, space, V_pec]), widgets.HBox([disp_mean_sig, space_A, disp_mean_sig_trac]), disp_mean_sig_peca]))
print(f'Volume do CP: {round(V_amostra*1e6, 2)} cm³')
print(f'Volume da Peça: {round(V_peca, 2)} m³')
if ensaio == 'Flexão 3 Pontos':
print(f'σ_(f, 3pt, amostra): {round(mean_sig, 2)} MPa')
print(f'σ_(f, t, amostra): {round(mean_sig_trac, 2)} MPa')
if ensaio == 'Flexão 4 Pontos':
print(f'σ_(f, 4pt, amostra): {round(mean_sig, 2)} MPa')
print(f'σ_(f, t, amostra): {round(mean_sig_trac, 2)} MPa')
print(f'σ_(f, t, peça): {round(mean_sig_peca, 6)} MPa')
select_ensaio = widgets.Dropdown(description='Ensaio:', options=['Flexão 3 Pontos', 'Flexão 4 Pontos'])
slider_X_am = widgets.FloatSlider(value=125,
min=5,
max=1000,
step=5,
description='CP L [mm]:',
orientation='horizontal')
slider_Y_am = widgets.FloatSlider(value=25,
min=5,
max=1000,
step=5,
description='CP H [mm]:',
orientation='horizontal')
slider_Z_am = widgets.FloatSlider(value=25,
min=5,
max=1000,
step=5,
description='CP b [mm]:',
orientation='horizontal')
slider_X_pec = widgets.FloatSlider(value=10,
min=1,
max=10,
step=1,
description='Peça L [m]:',
orientation='horizontal')
slider_Y_pec = widgets.FloatSlider(value=1,
min=1,
max=10,
step=1,
description='Peça H [m]:',
orientation='horizontal')
slider_Z_pec = widgets.FloatSlider(value=1,
min=1,
max=10,
step=1,
description='Peça b [m]:',
orientation='horizontal')
weib_effects_int = interactive(vol_ensaios, m=slider_m, X_pec=slider_X_pec, Y_pec=slider_Y_pec, Z_pec=slider_Z_pec,
X_am=slider_X_am, Y_am=slider_Y_am, Z_am=slider_Z_am, ensaio=select_ensaio, N=slider_N, sigma_0=slider_sig_0)
weib_eff_controls = widgets.HBox([widgets.VBox(weib_effects_int.children[1:4]), widgets.VBox(weib_effects_int.children[4:7])])
aviso = widgets.HTMLMath('Observação: Altere o número de amostras, $N$, para plotar.')
weibull_output = widgets.VBox([aviso, controls, widgets.HBox([pd_out, weibull_plot_show_int.children[-1]]),
weib_eff_controls, weib_effects_int.children[7], weib_effects_int.children[10]])
weibull_output.layout = widgets.Layout(display='flex',
flex_flow='column',
align_items='center',
align_content='center')
out = widgets.Output()
with out:
df_1
weibull_output
Muito obrigado!
Muito obrigado por sua participação e boa sorte com seus estudos! Para quaisquer dúvidas estaremos disponíveis nas monitorias (toda Terça-Feira 14:00 - 15:00).
Esta aula, bem como a lista, e o material de apoio fizeram uso de recursos de ensino pouco comuns como Widgets Interativos e um site de suporte.
Caso tenha 5 minutos, por favor, responda o formulário abaixo.
Grato, Murilo.