Задачи

Задача 1

У данной задачи есть 2 варианта (1.1, 1.2).

Ниже переставлен код выполняющий расчет скорости поперечных и продольных волн в породе состоящий из двух компонентов - минеральной (Кальцит) и пор наполненных водой. Расчет выполняется с использованием методов теории эффективных сред в зависимости от объемных долей компонентов, т.е. от пористости породы.

В данном коде представлено два метода: метод Фойгта-Ройсса-Хилла - Voigt_Reuss_Hill_method и метод Хашина-Штрикмана для двух компонентов - Hashin_Shtrikman_method_2comp. Оба метода основаны на усреднении некоторых граничных значений. В первом случае границы Фойгта - Voigt_bound и Ройсса - Reuss_bound, во втором Hashin_Shtrikman_bound_2comp.

Задача

  • задача 1.1 для метода Фойгта-Ройсса-Хилла
  • задача 1.2 для метода Хашина-Штрикмана

Построить на одном рисунке результаты работы одного из методов (два или четыре графика на одном figure):

  • Модуль упругости сжатия K, ГПа
  • Модуль упругости сдвига Mu, ГПа
  • Скорость продольных волн Vp, км/с
  • Скорость поперечных волн Vs, км/с

в виде графиков и записать полученные результаты в CSV файл. На графиках должен быть заголовок, подписи по осям, координатная сетка.

Код выполняющий задание свойств компонентов среды и выполняющий построение графиков и запись результатов должен быть размещен в отдельном файле от методов расчета. Код не относящийся к используемому методу должен быть удален.

Базовый код 

import numpy as np
import matplotlib.pyplot as plt

#%% Расчет

# Число компонентов

N_components = 2

# Число шагов

N_steps = 100

# Инициализируем массивы исходных данных и результатов

var_c = ['Vp', 'Vs', 'Rho', 'K', 'Mu', 'Vol']

c = np.zeros([N_components,N_steps], dtype=[(x, np.float64) for x in var_c])

var_r = ['K_V', 'Mu_V', 'Vp_V', 'Vs_V',
         'K_R', 'Mu_R', 'Vp_R', 'Vs_R',
         'K_VRH', 'Mu_VRH', 'Vp_VRH', 'Vs_VRH',
         'K_HS_UP', 'Mu_HS_UP', 'K_HS_LO', 'Mu_HS_LO',
         'K_HS', 'Mu_HS','Vp_HS','Vs_HS']

r = np.zeros([N_steps], dtype=[(x, np.float64) for x in var_r])

# Задаем исходные данные

# Calcite
c[0]['Vp'] = 6.54
c[0]['Vs'] = 3.35
c[0]['Rho'] = 2.7

# Water
c[1]['Vp'] = 1.6
c[1]['Vs'] = 1e-6
c[1]['Rho'] = 1.1

# Задаем пористость

poros = np.linspace(0.0, 1.0, num=N_steps)

# Определяем содержание компонентов

c[0]['Vol'] = 1 - poros
c[1]['Vol'] = poros

# Расчет

# Расчет K и Mu

c['K'] = c['Rho']*(c['Vp']**2-4./3.*c['Vs']**2)
c['Mu'] = c['Rho']*c['Vs']**2

# axis = 0 : Суммирование по компонентам N_components
# результат - одномерный массив из N_steps элементов

Rho_eff = np.sum(c['Vol']*c['Rho'], axis=0)

# Метод Фойгта

def Voigt_bound(c):
    K = np.sum(c['Vol']*c['K'], axis=0)
    Mu = np.sum(c['Vol']*c['Mu'], axis=0)
    return K, Mu


r['K_V'], r['Mu_V'] = Voigt_bound(c)

# Метод Ройсса

def Reuss_bound(c):
    K = 1.0 / np.sum(c['Vol']/c['K'], axis=0)
    Mu = 1.0 / np.sum(c['Vol']/c['Mu'], axis=0)
    return K, Mu

r['K_R'], r['Mu_R'] = Reuss_bound(c)

# Метод Фойгта-Ройсса-Хилла

def Voigt_Reuss_Hill_method(c):
    K_V, Mu_V = Voigt_bound(c)
    K_R, Mu_R = Reuss_bound(c)
    K = (K_V+K_R) / 2.0
    Mu = (Mu_V+Mu_R) / 2.0
    return K, Mu

r['K_VRH'] , r['Mu_VRH'] = Voigt_Reuss_Hill_method(c)

# Метод Хашина-Штрикмана для двух компонентов

def Hashin_Shtrikman_bound_2comp(c0,c1):
    K1, Mu1, f1 = c0['K'], c0['Mu'], c0['Vol']
    K2, Mu2, f2 = c1['K'], c1['Mu'], c1['Vol']
    
    Z0 = K1+4./3.*Mu1
    YK = 1./(K2-K1)
    YMu = 1./(Mu2-Mu1)
    
    K_HS = K1+f2/(YK+f1/Z0)
    
    Z1 = 2.*f1*(K1+2.*Mu1)/(5.*Mu1*Z0)
    Mu_HS = Mu1+f2/(YMu+Z1)
    
    return K_HS, Mu_HS

r['K_HS_UP'], r['Mu_HS_UP'] = Hashin_Shtrikman_bound_2comp(c[0],c[1])
r['K_HS_LO'], r['Mu_HS_LO'] = Hashin_Shtrikman_bound_2comp(c[1],c[0])

def Hashin_Shtrikman_method_2comp(c0,c1):
    K_HS_UP, Mu_HS_UP = Hashin_Shtrikman_bound_2comp(c0,c1)
    K_HS_LO, Mu_HS_LO = Hashin_Shtrikman_bound_2comp(c1,c0)
    K = (K_HS_UP+K_HS_LO) / 2.0
    Mu = (Mu_HS_UP+Mu_HS_LO) / 2.0
    return K, Mu

r['K_HS'], r['Mu_HS'] = Hashin_Shtrikman_method_2comp(c[0],c[1])

# Для каждого метода рассчитываем Vp и Vs

for m in ['V', 'R', 'VRH', 'HS']:
    r[f'Vp_{m}'] = np.sqrt((r[f'K_{m}']+4.0/3.0*r[f'Mu_{m}'])/Rho_eff)
    r[f'Vs_{m}'] = np.sqrt(r[f'Mu_{m}']/Rho_eff)

Задача 2

У данной задачи есть 4 варианта (2.1, 2.2, 2.3, 2.4).

2.1 Сейсмические события на удалении от станции

Задача

Зададимся координатами любого города. Предположим, что в этом городе находится сейсмостанция. Отметить на карте сам город и выбрать из каталога сейсмических событий события, отстоящие от него по координатам на не менее чем 15 градусов и не более чем на 60 градусов. Карту следует построить в стереографической проекции с центром в выбранном городе.

Данные

Каталог сейсмических событий с магнитудой больше 7 USGS

Базовый код 

import numpy as np
import matplotlib.pyplot as plt
from cartopy import crs as ccrs

fig = plt.figure(figsize=(10,10))

central_point = [60,160] # сейсмостанция

ax = fig.add_subplot(projection=ccrs.Stereographic(
     central_point[0],central_point[1],scale_factor=3))

ax.set_global()

ax.coastlines(resolution='110m')

2.2 Количество сейсмических событий на сетке

Построить на карте мира теплокарту числа сейсмических событий из каталога с шагом 3 на 3 градуса. Использовать проекцию Робинсона с центром по линии перемены дат.

Данные

Каталог сейсмических событий с магнитудой больше 7 USGS

Базовый код 

import numpy as np
import matplotlib.pyplot as plt
from cartopy import crs as ccrs
import matplotlib.ticker as mticker

fig = plt.figure(figsize=(10,10))

ax = fig.add_subplot(projection=ccrs.Robinson(central_longitude=180))

ax.set_global()

ax.coastlines(resolution='110m')

2.3 Карта напряженности магнитного поля

Построить глобальную карту напряженности магнитного поля за 1980, 1990 и 2000 годы. Использовать проекцию Робинсона. Пример карты https://ngdc.noaa.gov/geomag/magfield-wist/ (Total Intensity: Main field) Использовать контуры с подписями.

Данные

Данные NGDC NOAA 1980 год 1980 год 1980 год

Базовый код 

import numpy as np
import matplotlib.pyplot as plt
from cartopy import crs as ccrs
import matplotlib.ticker as mticker

fig = plt.figure(figsize=(10,10))

ax = fig.add_subplot(projection=ccrs.Robinson(central_longitude=0))

ax.set_global()

ax.coastlines(resolution='110m')

2.4 Построение смещений на карте

В данном Shape-файле описаны точки трех разных видов (kind): base: опорные точки, profile1 и profile2 - результаты измерений смещений вдоль профилей. Построить карту используя подложку tile сервера OpenStreetMap на которую нанести: опорные точки, точки профилей 1 и 2 соединенные линиями различных цветов. Для каждой из точек профиля указать направление (az в градусах) и величину смещения (value) в виде стрелки.

Данные

Синтетические данные

Базовый код 

import numpy as np
import matplotlib.pyplot as plt
from cartopy import crs as ccrs
import  shapefile
from cartopy.io.img_tiles import OSM

shpfile = shapefile.Reader("/home/f0ma/tdata/tdata")

fields = [x[0] for x in shpfile.fields[1:] ]

for item in shpfile.iterShapeRecords():
    print("Position:", item.shape.points[0])
    print("Data:", dict(zip(fields, item.record)))


tiles_provider = OSM()


plt.figure(figsize=(16, 16))

ax = plt.axes(projection=tiles_provider.crs)

ax.set_extent([shpfile.bbox[0]-0.1, shpfile.bbox[2]+0.1,
               shpfile.bbox[1]-0.1, shpfile.bbox[3]+0.1], ccrs.PlateCarree())

ax.add_image(tiles_provider, 11)
tasks.txt · Последнее изменение: 2023/04/08 20:18 — root
Если не указано иное, содержимое этой вики предоставляется на условиях следующей лицензии: CC Attribution-Share Alike 4.0 International


  • Показать исходный текст
  • История страницы
  • Ссылки сюда
  • Недавние изменения
  • Управление медиафайлами
  • Все страницы