Last active
May 19, 2019 19:55
-
-
Save PedroBern/77e800cbbe4dd69d8111b483834d5334 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import numpy as np | |
from mpl_toolkits.mplot3d import Axes3D | |
import matplotlib.pyplot as plt | |
from matplotlib import cm | |
from matplotlib.ticker import LinearLocator, FormatStrFormatter | |
# Dados do problema | |
gama = 19 | |
m = 0.5 | |
k = 0.5 | |
P = 1500 | |
ti = 0 | |
r_arr = z_arr = np.arange(0.5, 3.5, 0.5) | |
# configuracoes para o araquivo csv | |
file_name = "tensoes_geostaticas.csv" | |
delimiter = ";" | |
fmt = "%.2f" | |
# Parte 0: definicao de utilidades | |
def get_R(r, z): | |
return np.sqrt(np.power(r,2) + np.power(z,2)) | |
def get_phi(k): | |
return np.degrees(np.arcsin(1 - k)) | |
def get_vi(gama, z): | |
return gama * z | |
def get_hi(k, vi): | |
return k * vi | |
def get_delta_v(P, R, z): | |
return (3 * P * np.power(z,3)) / (2 * np.pi * np.power(R, 5)) | |
def get_delta_h(P, R, z, r, m): | |
return P / (2 * np.pi) * ((3 * np.power(r, 2) * z) / np.power(R, 5) - (1 - 2 * m) / (R * (R + z))) | |
def get_delta_t(P, R, z, r): | |
return (3 * P * r * np.power(z, 2)) / (2 * np.pi * np.power(R, 5)) | |
def get_vf(vi, delta_v): | |
return vi + delta_v | |
def get_hf(hi, delta_h): | |
return hi + delta_h | |
def get_tf(ti, delta_t): | |
return ti + delta_t | |
# Parte 1: Gráficos | |
fig = plt.figure(figsize=plt.figaspect(0.5)) | |
fig.subplots_adjust(hspace=0.5) | |
X = r_arr | |
Y = z_arr | |
X, Y = np.meshgrid(X, Y) | |
R = get_R(X,Y) | |
strains = [ | |
get_vf(get_vi(gama, Y), get_delta_v(P,R,Y)), | |
get_hf(get_hi(k, get_vi(gama, Y)), get_delta_h(P, R, Y, X, m)), | |
get_tf(ti, get_delta_t(P, R, Y, X)), | |
] | |
titles = [ | |
"Tensões Verticais", | |
"Tensões Horizontais", | |
"Tensões Cisalhantes", | |
] | |
for index, Z in enumerate(strains): | |
ax = fig.add_subplot(2, 2, index + 1, projection='3d') | |
ax.view_init(30, 45) | |
surf = ax.plot_surface(X, Y, Z, cmap=cm.coolwarm, linewidth=0, antialiased=False) | |
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f')) | |
ax.set_xlabel('Distancia horizontal (r)') | |
ax.set_ylabel('Profundidade (z)') | |
ax.set_title(titles[index]) | |
plt.show() | |
# Parte 2: Arquvio pro excel | |
result = [] | |
for r in r_arr: | |
for z in z_arr: | |
R = get_R(r, z) | |
vi = get_vi(gama, z) | |
hi = get_hi(k, vi) | |
delta_v = get_delta_v(P, R, z) | |
delta_h = get_delta_h(P, R, z, r, m) | |
delta_t = get_delta_t(P, R, z, r) | |
vf = get_vf(vi, delta_v) | |
hf = get_hf(hi, delta_h) | |
tf = get_tf(ti, delta_t) | |
result.append((r,z,vi,hi,delta_v, delta_h, delta_t,vf,hf,tf)) | |
file = np.asarray(result) | |
np.savetxt( | |
file_name, | |
file, | |
delimiter=delimiter, | |
fmt=fmt, | |
header="r;z;vi;hi;delta_v; delta_h; delta_t;vf;hf;tf" | |
) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment