The strata near the fault did not undergo dislocation.

[kaikai430]

When using GemPy for 3D geological modeling, I set up the fault relationship matrix, and the model loaded successfully, but the strata near the fault did not show any displacement. What is the reason for this? How can I solve this problem?This is my first time using GitHub, and my question format might not be standard, so I hope you seniors can understand.
#-------------------------------------------------------------------------------------------------------
This is my code

%% [markdown]

# 1.5: Fault relations

GemPy断层关系建模示例

%% [markdown]

## 1. 导入依赖库

%%

import gempy as gp
import gempy_viewer as gpv
import numpy as np
import pandas as pd
from gempy.core.data import (
StructuralFrame,
SurfacePointsTable,
OrientationsTable,
ColorsGenerator,
StructuralGroup,
StructuralElement
)
from gempy_engine.core.data.stack_relation_type import StackRelationType
from gempy.core.data import FaultsRelationSpecialCase

固定随机种子保证可重复性

np.random.seed(1515)

%% [markdown]

## 2. 数据准备

%%

读取表面点和方向数据

df = pd.read_csv("6个别数据.CSV")
df_orientations = pd.read_csv("倾向数据_归一化.csv")

初始化颜色生成器

color_gen = ColorsGenerator()

%% [markdown]

## 3. 构建结构框架

%%

初始化结构框架（含地层、断层、基底三组）

structural_frame = gp.data.StructuralFrame(
structural_groups=[
gp.data.StructuralGroup(
name="Stratigraphy",
elements=[],
structural_relation=StackRelationType.ERODE
),
gp.data.StructuralGroup(
name="Faults",
elements=[],
structural_relation=StackRelationType.FAULT
),
gp.data.StructuralGroup(
name="Basement",
elements=[],
structural_relation=StackRelationType.BASEMENT
)
],
color_gen=color_gen
)

%% [markdown]

### 3.1 创建结构元素

%%

为每个surface创建结构元素

unique_surfaces = df['surface'].unique()

for name in unique_surfaces:
# 处理表面点数据
mask = df['surface'] == name
subset = df[mask]

sp_array = np.zeros(len(subset), dtype=SurfacePointsTable.dt)
sp_array["X"] = subset["X"].values
sp_array["Y"] = subset["Y"].values
sp_array["Z"] = subset["Z"].values
sp_array["id"] = 0
sp_array["nugget"] = 0.0
surface_points = SurfacePointsTable(data=sp_array)

# 处理方向数据
orientations = None
if name in df_orientations['surface'].unique():
    mask_ori = df_orientations['surface'] == name
    subset_ori = df_orientations[mask_ori]
    
    ori_array = np.zeros(len(subset_ori), dtype=OrientationsTable.dt)
    ori_array["X"] = subset_ori["X"].values
    ori_array["Y"] = subset_ori["Y"].values
    ori_array["Z"] = subset_ori["Z"].values
    ori_array["G_x"] = subset_ori["G_x"].values
    ori_array["G_y"] = subset_ori["G_y"].values
    ori_array["G_z"] = subset_ori["G_z"].values
    ori_array["id"] = 0
    ori_array["nugget"] = 0.0
    orientations = OrientationsTable(data=ori_array)

# 创建结构元素
element = StructuralElement(
    name=name,
    surface_points=surface_points,
    orientations=orientations,
    color=next(color_gen)
)

# 分配组别
group_name = "Faults" if "Fault" in name else "Basement" if "Basement" in name else "Stratigraphy"
group = structural_frame.get_group_by_name(group_name)
group.elements.append(element)

%% [markdown]

## 4. 创建地质模型

%%

设置模型范围

extent = [
df['X'].min(), df['X'].max(),
df['Y'].min(), df['Y'].max(),
df['Z'].min(), df['Z'].max()
]

初始化地质模型

geo_model = gp.create_geomodel(
project_name='Faults_relations',
extent=extent,
resolution=[30, 30, 30],
structural_frame=structural_frame
)

%% [markdown]

## 5. 配置断层关系

%%

获取断层元素

fault1 = geo_model.structural_frame.get_element_by_name("Fault_1")
fault2 = geo_model.structural_frame.get_element_by_name("Fault_2")
fault3 = geo_model.structural_frame.get_element_by_name("Fault_3")

移除默认的Faults组

Faults = geo_model.structural_frame.get_group_by_name("Faults")
Faults.elements.remove(fault1)
Faults.elements.remove(fault2)
Faults.elements.remove(fault3)

%% [markdown]

### 5.1 定义断层系列

%%

添加新的断层结构组

gp.add_structural_group(
model=geo_model,
group_index=0,
structural_group_name="fault_series_1",
elements=[fault1],
structural_relation=StackRelationType.FAULT,
fault_relations=FaultsRelationSpecialCase.OFFSET_ALL
)

gp.add_structural_group(
model=geo_model,
group_index=1,
structural_group_name="fault_series_2",
elements=[fault2],
structural_relation=StackRelationType.FAULT,
fault_relations=FaultsRelationSpecialCase.OFFSET_ALL
)

gp.add_structural_group(
model=geo_model,
group_index=2,
structural_group_name="fault_series_3",
elements=[fault3],
structural_relation=StackRelationType.FAULT,
fault_relations=FaultsRelationSpecialCase.OFFSET_ALL
)

执行清除

geo_model.structural_frame.structural_groups = [
g for g in geo_model.structural_frame.structural_groups
if g.name != "Faults"
]

print("清除后结构组:")
for g in geo_model.structural_frame.structural_groups:
print(" -", g.name)

%% [markdown]

### 5.2 设置断层关系矩阵

%%

配置断层关系矩阵

gp.set_fault_relation(
frame=geo_model.structural_frame,
rel_matrix=np.array([
[0, 0, 0, 1, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
])
)
print(geo_model.structural_frame)
print("\n 当前结构组顺序（与 fault_relation 矩阵行列一致）:")
for i, group in enumerate(geo_model.structural_frame.structural_groups):
print(f"[{i}] {group.name}")

更新模型

geo_model.update_transform()
gp.compute_model(geo_model)

#%%
#查看所有模型结构组
print(geo_model.structural_frame)
print("\n 当前结构组顺序（与 fault_relation 矩阵行列一致）:")
for i, group in enumerate(geo_model.structural_frame.structural_groups):
print(f"[{i}] {group.name}")

%% [markdown]

## 6. 可视化结果

%%

2D剖面可视化

gpv.plot_2d(geo_model, show_data=True)

%%

3D交互式可视化

gpv.plot_3d(
geo_model,
show_data=True,
image=False,
kwargs_plot_structured_grid={'opacity': 0.2}
)
#-------------------------------------------------------------------------------------------------
And this is my output
Setting Backend To: AvailableBackends.numpy
StructuralFrame(
structural_groups=[
StructuralGroup(
name=fault_series_1,
structural_relation=StackRelationType.FAULT,
elements=[
Element(
name=Fault_1,
color=#43b665,
is_active=True
)
]
),
StructuralGroup(
name=fault_series_2,
structural_relation=StackRelationType.FAULT,
elements=[
Element(
name=Fault_2,
color=#d14897,
is_active=True
)
]
),
StructuralGroup(
name=fault_series_3,
structural_relation=StackRelationType.FAULT,
elements=[
Element(
name=Fault_3,
color=#59e5b8,
is_active=True
)
]
),
StructuralGroup(
name=Stratigraphy,
structural_relation=StackRelationType.ERODE,
elements=[
Element(
name=Roof_Coal_8,
color=#d746be,
is_active=True
),
Element(
name=Floor_Coal_8,
color=#60eb91,
is_active=True
),
Element(
name=Overburden5,
color=#62dedb,
is_active=True
),
Element(
name=Roof_Coal_9,
color=#5e63d8,
is_active=True
),
Element(
name=Floor_Coal_9,
color=#52a22e,
is_active=True
),
Element(
name=Roof_Coal_10,
color=#371970,
is_active=True
),
Element(
name=Floor_Coal_10,
color=#e5c339,
is_active=True
),
Element(
name=Floor_Coal_11,
color=#d3dc76,
is_active=True
),
Element(
name=Roof_Coal_11,
color=#4d478e,
is_active=True
),
Element(
name=Roof_Coal_7,
color=#6d2b9f,
is_active=True
),
Element(
name=Floor_Coal_7,
color=#95ce4b,
is_active=True
),
Element(
name=Roof_Coal_4,
color=#b271d0,
is_active=True
),
Element(
name=Floor_Coal_4,
color=#5DA629,
is_active=True
),
Element(
name=Overburden4,
color=#e5421d,
is_active=True
),
Element(
name=Floor_Coal_5,
color=#72e54a,
is_active=True
),
Element(
name=Floor_Coal_6,
color=#d0e63d,
is_active=True
),
Element(
name=Roof_Coal_5,
color=#583bd1,
is_active=True
),
Element(
name=Roof_Coal_6,
color=#b949e2,
is_active=True
),
Element(
name=Roof_Coal_1,
color=#9f0052,
is_active=True
),
Element(
name=Floor_Coal_1,
color=#015482,
is_active=True
),
Element(
name=Roof_Coal_2,
color=#728f02,
is_active=True
),
Element(
name=Floor_Coal_2,
color=#ffbe00,
is_active=True
),
Element(
name=Floor_Coal_3,
color=#443988,
is_active=True
),
Element(
name=Roof_Coal_3,
color=#ff3f20,
is_active=True
)
]
),
StructuralGroup(
name=Basement,
structural_relation=StackRelationType.BASEMENT,
elements=[
Element(
name=Basement,
color=#df344e,
is_active=True
)
]
)
],
fault_relations=
[[False, False, False, False, False],
[False, False, False, False, False],
[False, False, False, False, False],
[False, False, False, False, False],
[False, False, False, False, False]],

当前结构组顺序（与 fault_relation 矩阵行列一致）:
[0] fault_series_1
[1] fault_series_2
[2] fault_series_3
[3] Stratigraphy
[4] Basement

[javoha]

Hello @kaikai430,

thanks for your question. It is really hard to figure out the exact problem in your code just by looking at it here. From what I can see, in the last thing you output your fault relations contains all False, which means all your faults do not cause any offset. So my assumption is that the mistake happens whenyou sert your fault relations:

gp.set_fault_relation(
frame=geo_model.structural_frame,
rel_matrix=np.array([
[0, 0, 0, 1, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 1, 0],
[0, 0, 0, 0, 0],
[0, 0, 0, 0, 0],
])
)

Ma suggestion would be to use the alternative way to set faults in gempy and use the default behavior without setting the fault relation matrix actively:

gp.set_is_fault(geo_model, ["fault_series_1", "fault_series_2", "fault_series_3"])

Also make sure to watch these videos to make sure you understand all the implications if you havent done so already.

If you are able to provide your full code with the input data as a zip file I could also take a closer look.

Cheers,
Jan

[kaikai430]

Thank you very much for your response! Following your method, my fault relationship matrix is correct now! I really appreciate your help during your busy schedule, thank you

very much!

[kaikai430]

    fault_relations=

[[False, False, False, True, True],
[False, False, False, True, True],
[False, False, False, True, True],
[False, False, False, False, False],
[False, False, False, False, False]],