Scribd
Cargar
10 vistas8 páginas

Análisis de deslizamientos en terremotos

Este documento presenta el análisis de datos obtenidos de la simulación dinámica de un terremoto. Se cargan datos de deslizamiento, velocidad de deslizamiento, esfuerzo de corte y otros parámetros. Se grafican estos datos en 3D y se calculan propiedades como la energía de fractura en diferentes puntos de la falla. Finalmente, se analizan las relaciones entre la caída promedio y ponderada del esfuerzo de corte.

Cargado por

William Segundo Araujo Navarro
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd
10 vistas8 páginas

Análisis de deslizamientos en terremotos

Este documento presenta el análisis de datos obtenidos de la simulación dinámica de un terremoto. Se cargan datos de deslizamiento, velocidad de deslizamiento, esfuerzo de corte y otros parámetros. Se grafican estos datos en 3D y se calculan propiedades como la energía de fractura en diferentes puntos de la falla. Finalmente, se analizan las relaciones entre la caída promedio y ponderada del esfuerzo de corte.

Cargado por

William Segundo Araujo Navarro
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd

# -*- coding: utf-8 -*-

"""
Created on Mon Nov 5 15:21:14 2018

@author: Alejandro
"""
#%%
import numpy as np
import [Link] as plt
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D

#datos = [Link]("3/09141_step3171.dat",skiprows=9)
#datos = [Link]("3/09141_step3187.dat",skiprows=9)
#datos = [Link]("3/09141_step3229.dat",skiprows=9)
#datos = [Link]("3/09141_step3247.dat",skiprows=9)
#datos = [Link]("3/09141_step3264.dat",skiprows=9)

datos = [Link]("537/00537_step2104.dat",skiprows=9)

datos2 = [Link]("537/[Link]",skiprows=9)

grid = 2 ; ds = 4*4**(grid)
dt = ds/2/6000

x = datos[:,0]; y = datos[:,1];
slip_vel = datos[:,2]; # m/s
slip = datos[:,3] # mm
shear_stres = datos[:,4] # MPa
t = datos2[:,0] # s
moment_rate = datos2[:,1] # Nm/s
moment = datos2[:,2] # Nm
Mw = datos2[:,3]

#%%

#% Scatter 3D
function = slip*4**(grid)

fig = [Link](figsize=(10,6))
ax = [Link](projection='3d')

#tt = list(range(-1,8))
sct = [Link](x*ds, y*ds, function,vmin=0, vmax=[Link](function), s = 50, alpha
= 1 , c=function, cmap='jet', linewidth=0.1, marker = "s"); #velocidad
deslizamientos
#sct = [Link](x*ds, y*ds, function,vmin=0, vmax=4, s = 40, alpha = 1 ,
c=function, cmap='jet', linewidth=0.1, marker = "s"); #deslizamientos
#sct = [Link](x*4*4**(grid), y*4*4**(grid), function, s = 50, c=function,
cmap='jet', linewidth=0.1);

ax.set_xlabel('x [m]')
ax.set_ylabel('y [m]')
ax.set_zticklabels([]); #ax.set_zlabel('MPa')
#ax.set_zlim(-1, 7);
ax.set_xlim([Link](x*4*4**(grid)), [Link](x*4*4**(grid)));
ax.set_ylim([Link](y*4*4**(grid)), [Link](y*4*4**(grid)))
[Link](sct, shrink=0.5, aspect=5, label='Deslizamiento [mm]');#
fig.set_label('label')

ax.view_init(elev=90., azim=0) # Vista en planta


#ax.view_init(elev=20., azim=30) # Vista 3D

[Link]('Deslizamiento terremoto Mw = 3.77943, t = 0.554667 [s]')

[Link]()

#%%
"""
LEER TODOS LOS ARCHIVOS PRODUCIDOS POR HIDEO AOCHI
"""

import glob
import numpy as np
import [Link] as plt
import obspy as ob

filenames = sorted([Link]('537/00537_step2*.dat'))

# 0: x, 1: y, 2: slip vel [m/s], 3: slip [mm], 4: shear stress [MPa]


col = 2

grid = 2 ; ds = 4*4**(grid)

ds = 4*4**(grid)
dt = ds/2/6000

fun_0 = [Link]('537/00537_step0001.dat')
X = fun_0[:,0]
Y = fun_0[:,1]

fun = fun_0[:,col]
fun2 = fun_0[:,4]
fun3 = fun_0[:,3]
for f in filenames:
data = [Link](fname=f)
fun = np.c_[fun, data[:,col]]
fun2 = np.c_[fun2, data[:,4]]
fun3 = np.c_[fun3, data[:,3]]

t2 = dt*[Link](0,len(fun[0,:]),len(fun[0,:]))

#%%
Datos_finales = [Link]("537/[Link]",skiprows=9)

time_step = Datos_finales[:,0]
time = Datos_finales[:,1] # [s]
moment2 = Datos_finales[:,2] # [Nm/s]
moment_rate2 = Datos_finales[:,3] # [Nm]

#%%
from [Link] import pearsonr
ss = [Link](fun == [Link](fun[:,:]))

int = [Link](ss[0])

[Link]()
[Link](221)
[Link](t2,fun[int,:],'k',linewidth=1.5)
[Link]('Tiempo [s]'); [Link]('Velocidad de deslizamiento [m/s]')
[Link](t2[0], t2[len(t2)-1])

#[Link]()
[Link](222)
[Link](t2,fun2[int,:],'k',linewidth=1.5)
[Link]('Tiempo [s]'); [Link]('Esfuerzo de Corte [MPa]')
[Link](t2[0], t2[len(t2)-1])

#[Link]()
[Link](223)
[Link](t2,fun3[int,:]*4**(grid),'k',linewidth=1.5)
[Link]('Tiempo [s]'); [Link]('Deslizamientos [mm]')
[Link](t2[0], t2[len(t2)-1])

[Link](224)
[Link](time, moment_rate2,'k',linewidth=1.5)
[Link]('Tiempo [s]'); [Link]('Tasa de momento sísmico [Nm/s] \n de toda la
ruptura')
[Link](time[0], time[len(time)-1])

[Link]('Parámetros de ruptura \n simulación ruptura dinámica terremoto Mw =


3.77943 ',size=16)

p = pearsonr(fun[int,:],fun2[int,:]) # Correlation and p-value


print('La correlación es', p[0])

print('El punto está en', (x[int])*4*4**(grid),(y[int])*4*4**(grid))


#%% realizar filtro pasabajo 5 Hz

#%% calcular Gc (energía de fractura)

## fun2[int,:] # Esfuerzo de corte


## fun3[int,:] # Deslizamientos
##import [Link] as plt
##
##fun2[int,0:65] = fun2[int,0:65]
##fun2[int,65:] = 0

[Link]()
[Link]((fun3[int,:])*(10**(-3)),fun2[int,:],'o-k')
[Link]('Deslizamientos [m]'); [Link]('Esfuerzo de corete [MPa]')

#[Link](fun3[500,:],fun2[500,:],'k')
index = [Link]( fun2[int,:] == 0 ) # indices para encontrar el valor cuando la
tensión cae a cero.
index = [Link](index)

dc1 = (fun3[int,index[0,0]])*(10**(-3)) # Distancia a la cual se alcanza por


primera vez el esfuerzo cero.
tau1 = [Link](fun2[int,:])
G_c1 = (tau1*dc1)/2

print('La energía de fractura para ese punto es:', G_c1,'M J/m/m')

#%% Energía de Fractura para toda la falla

n = len(fun2[:,0])

dc = ([Link]((n)))
tau = ([Link]((n)))
G_c = ([Link]((n)))

for k in range(n):
if (( max(fun3[k,:]) > 0 )):
if (min(fun2[k,:]) == 0 ):
index = [Link]( fun2[k,:] == 0 ) # indices para encontrar el valor
cuando la tensión cae a cero.
index = [Link](index)
dc[k] = (fun3[k,index[0,0]])*(10**(-3)) # en metros
tau[k] = ([Link](fun2[k,:])) # MPa
G_c[k] = (((tau[k])*(dc[k]))/2 )*(10**(6)) # J/m/m
# else:
# dc[k] = 0
# tau[k] = 0
# G_c[k] = (0.16)
else:
dc[k] = 0
tau[k] = 0
G_c[k] = 3* (1.25*(10**(3))) # J/m/m

# pass
#%% Graficación de energía de fractura
fig = [Link](figsize=(10,6))
ax = [Link](projection='3d')

sct = [Link](X*ds, Y*ds, G_c, s = 50, alpha = 1 , c=G_c, cmap='jet',


linewidth=0.1, marker = "s"); #velocidad deslizamientos

ax.set_xlabel('x [m]')
ax.set_ylabel('y [m]')
ax.set_zticklabels([]); #ax.set_zlabel('MPa')
#ax.set_zlim(-1, 7);
ax.set_xlim([Link](X*4*4**(grid)), [Link](X*4*4**(grid)));
ax.set_ylim([Link](Y*4*4**(grid)), [Link](Y*4*4**(grid)))
[Link](sct, shrink=0.5, aspect=5, label='Energía de Fractura [J/m**2]');#
fig.set_label('label')

ax.view_init(elev=90., azim=0) # Vista en planta


#ax.view_init(elev=20., azim=30) # Vista 3D

[Link]('Deslizamiento terremoto Mw = 3.77943, t = 0.554667 [s]')

[Link]()

#%% Relación de caída de esfuerzo propemdio y ponderado

# slip : fun3
# shear stres : fun2

I_1 = [Link]((len(fun3)))
I_2 = [Link]((len(fun3)))
I_3 = [Link]((len(fun3)))
#dsigma = [Link]((len(fun3)))
du = [Link]((len(fun3)))

for i in range(len(I_1)):
dsigma = max(fun2[i,:]) - fun2[i,len(fun2[0,:])-1]
du[i] = fun3[i,len(fun3[0,:])-1]
I_1[i] = dsigma*(du[i])*ds
I_2[i] = (du[i])*ds
I_3[i] = dsigma*ds

ddddd = [Link](du > 0)


A = (len(ddddd[0]))*(4*4**(grid))

dsigmaE = (([Link](I_1)) / ([Link](I_2)))


dsigmaA = [Link](I_3) /A
% Programa de correlaciones TwoPoint Statistic
clear all; clc;% close all;

datos = load('1534/01534_step3066.dat');
Z = 3 ;
grid = 4*(4^(Z)) ;
x = (datos(:,1))*grid ; % metros
y = (datos(:,2))*grid ; % metros
slip_vel = datos(:,3) ; % m/s
slip = datos(:,4)*grid/4 ; % mm
shear_stress = datos(:,5) ; % MPa

H = vec2mat(slip,64);
%

X = 1:64; X = X*grid;
Y = 1:64; Y = Y*grid;
[X,Y] = meshgrid(X,Y);

figure()

s = surf(X,Y,H);colormap jet; axis xy;

[Link] = 'none';
c = colorbar;
[Link] = 'Esfuerzo de corte [MPa]';
xlim( [x(1) , x(end)] ); ylim( [y(1) , y(end)] );
title('Esfuerzo de corte final Terremoto Mw4.62181')
xlabel({'Distancia [m]'});ylabel({'Distancia [m]'})
view(2)

%% Two point statistic


%[ corrfun r rw] = twopointcorr( x,y,dr,blksize,verbose)

[kx , ky ] = find( 0 < H & H <= 2000 ) ;


%k = (0 < H) & (H < 200);
%ky = (0< H) & (H < 220);

%dr = grid*2
dr =550

[ corrfun, r, rw] = twopointcorr( kx*grid,ky*grid,dr,1000);


%[ corrfun, r, rw] = tw opointcorr( x(k),y(k),dr,1000);
%[ corrfun, r, rw] = twopointcorr( kx,ky,dr);
%figure(); plot(r, corrfun/max(corrfun), '-k'); xlabel('Distancia [m]');
ylabel('Funci�n de correlaci�n')
figure()
hold on
plot(r, corrfun); xlabel('Distancia [m]'); ylabel('Funci�n de correlaci�n')
# -*- coding: utf-8 -*-
"""
Created on Tue Dec 4 15:53:58 2018

@author: Alejandro
"""

#%% Datos de caida de esfuerzos

import numpy as np
import [Link] as plt
from matplotlib import cm
#from mpl_toolkits.mplot3d import Axes3D

datos = [Link]("[Link]")
numero = datos[:,0]
Mw = datos[:,1]
S_A = datos[:,2]
S_E = datos[:,3]

function = Mw
#[Link](S_A, S_E, function,vmin=[Link](function), vmax=[Link](function), s =
50, alpha = 1 , c=function, cmap='jet', linewidth=0.1, marker = "s"); #velocidad
deslizamientos
sc = [Link](S_E, S_A,s = 100 ,vmin=[Link](function), vmax=[Link](function),
alpha = 1 , c=function, cmap='jet', linewidth=0.1);
[Link]('Caída de esfuerzos')
[Link]('$\Delta \sigma_E$ [MPa]'); [Link]('$\Delta \sigma_A$ [MPa]')
#[Link](0, n);[Link](0, n)
[Link](sc,label='Magnitud de momento Mw')

[Link]()

También podría gustarte