Skip to content

Commit 13f2bc7

Browse files
div-tygdiv-tyg
committed
modified parameter.json: units in keys and isogeometric elements
1 parent 246cc3b commit 13f2bc7

9 files changed

Lines changed: 98 additions & 57 deletions

File tree

‎examples/linear-elastic-plate-with-hole/fenics/Snakefile‎

Lines changed: 33 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,33 @@
1+
import json
2+
3+
rule all:
4+
input:
5+
"solution_metrics.json",
6+
"solution_field_data.zip"
7+
8+
rule create_mesh:
9+
input:
10+
script = "create_mesh.py",
11+
parameters = "parameters.json",
12+
output:
13+
mesh = "mesh.msh",
14+
conda: "environment_mesh.yml"
15+
shell:
16+
"""
17+
python3 {input.script} --input_parameter_file {input.parameters} --output_mesh_file {output.mesh}
18+
"""
19+
20+
rule run_simulation:
21+
input:
22+
script = "run_simulation.py",
23+
parameters = "parameters.json",
24+
mesh = "mesh.msh",
25+
output:
26+
zip = "solution_field_data.zip",
27+
metrics = "solution_metrics.json",
28+
conda:
29+
"environment_simulation.yml"
30+
shell:
31+
"""
32+
python3 {input.script} --input_parameter_file {input.parameters} --input_mesh_file {input.mesh} --output_solution_file_zip {output.zip} --output_metrics_file {output.metrics}
33+
"""

‎examples/linear-elastic-plate-with-hole/fenics/run_benchmark.py‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -43,8 +43,9 @@
4343
for file in root_unzipped_benchmark_dir.glob("parameters_*.json"):
4444
with open(file, "r") as f:
4545
data = json.load(f)
46-
if data.get("element-size").get("value") >= 0.025:
47-
46+
#if data.get("element_size[m]") >= 0.025:
47+
#if data.get("configuration") == "1":
48+
if 1==1:
4849
# Create output directory for the configuration
4950
output_dir = root_unzipped_benchmark_dir / "results" / data.get("configuration")
5051
output_dir.mkdir(parents=True, exist_ok=True)

‎examples/linear-elastic-plate-with-hole/fenics/run_simulation.py‎

Lines changed: 22 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@
1010
from dolfinx.fem.petsc import LinearProblem
1111
from mpi4py import MPI
1212
from pint import UnitRegistry
13+
from ufl.algorithms import estimate_total_polynomial_degree
1314

1415
# Add parent directory to sys.path
1516
sys.path.insert(0, str(Path(__file__).resolve().parent.parent))
@@ -29,7 +30,7 @@ def run_fenics_simulation(
2930
gdim=2,
3031
)
3132

32-
V = df.fem.functionspace(mesh, ("CG", parameters["element-degree"], (2,)))
33+
V = df.fem.functionspace(mesh, ("CG", parameters["isoparametric_element_degree"], (2,)))
3334

3435
tags_left = facet_tags.find(1)
3536
tags_bottom = facet_tags.find(2)
@@ -45,46 +46,46 @@ def run_fenics_simulation(
4546

4647
E = (
4748
ureg.Quantity(
48-
parameters["young-modulus"]["value"], parameters["young-modulus"]["unit"]
49+
parameters["young_modulus[Pa]"], "Pa"
4950
)
5051
.to_base_units()
5152
.magnitude
5253
)
5354
nu = (
5455
ureg.Quantity(
55-
parameters["poisson-ratio"]["value"], parameters["poisson-ratio"]["unit"]
56+
parameters["poisson_ratio"], ""
5657
)
5758
.to_base_units()
5859
.magnitude
5960
)
6061
radius = (
61-
ureg.Quantity(parameters["radius"]["value"], parameters["radius"]["unit"])
62+
ureg.Quantity(parameters["radius[m]"], "m")
6263
.to_base_units()
6364
.magnitude
6465
)
6566
L = (
66-
ureg.Quantity(parameters["length"]["value"], parameters["length"]["unit"])
67+
ureg.Quantity(parameters["length[m]"], "m")
6768
.to_base_units()
6869
.magnitude
6970
)
7071
load = (
71-
ureg.Quantity(parameters["load"]["value"], parameters["load"]["unit"])
72+
ureg.Quantity(parameters["load[MPa]"], "MPa")
7273
.to_base_units()
7374
.magnitude
7475
)
75-
displacement_evaluation_point = parameters["displacement-evaluation-point"]
76+
displacement_evaluation_point = parameters["displacement_evaluation_point[m]"]
7677
displacement_evaluation_x = (
7778
ureg.Quantity(
78-
displacement_evaluation_point["x"]["value"],
79-
displacement_evaluation_point["x"]["unit"],
79+
displacement_evaluation_point[0],
80+
"m",
8081
)
8182
.to_base_units()
8283
.magnitude
8384
)
8485
displacement_evaluation_y = (
8586
ureg.Quantity(
86-
displacement_evaluation_point["y"]["value"],
87-
displacement_evaluation_point["y"]["unit"],
87+
displacement_evaluation_point[1],
88+
"m",
8889
)
8990
.to_base_units()
9091
.magnitude
@@ -115,10 +116,6 @@ def as_tensor(v):
115116

116117
dx = ufl.Measure(
117118
"dx",
118-
metadata={
119-
"quadrature_degree": parameters["quadrature-degree"],
120-
"quadrature_scheme": parameters["quadrature-rule"],
121-
},
122119
)
123120
ds = ufl.Measure(
124121
"ds",
@@ -238,10 +235,10 @@ def project(
238235
return uh
239236

240237
plot_space_stress = df.fem.functionspace(
241-
mesh, ("DG", parameters["element-degree"] - 1, (2, 2))
238+
mesh, ("DG", parameters["isoparametric_element_degree"] - 1, (2, 2))
242239
)
243240
plot_space_mises = df.fem.functionspace(
244-
mesh, ("DG", parameters["element-degree"] - 1, (1,))
241+
mesh, ("DG", parameters["isoparametric_element_degree"] - 1, (1,))
245242
)
246243
stress_nodes_red = project(sigma(u), plot_space_stress, dx)
247244
stress_nodes_red.name = "stress"
@@ -291,7 +288,8 @@ def mises_stress(u):
291288
quad_element = basix.ufl.quadrature_element(
292289
mesh.topology.cell_name(),
293290
value_shape=(1,),
294-
degree=parameters["quadrature-degree"],
291+
scheme="default",
292+
degree=estimate_total_polynomial_degree(u_),
295293
)
296294

297295
Q_mises = df.fem.functionspace(mesh, quad_element)
@@ -321,12 +319,12 @@ def mises_stress(u):
321319

322320
# Save metrics
323321
metrics = {
324-
"max_von_mises_stress_nodes": max_mises_stress_nodes,
325-
"max_von_mises_stress_gauss_points": max_mises_stress_gauss_points,
326-
"l2_error_displacement": l2_error_displacement,
327-
"reaction_force_left_boundary_x": reaction_left_x,
328-
"reaction_force_left_boundary_y": reaction_left_y,
329-
f"displacement_x_at_evaluation_point (x={displacement_evaluation_x}, y={displacement_evaluation_y})": displacement_x_at_evaluation_point,
322+
"max_von_mises_stress_nodes[Pa]": max_mises_stress_nodes,
323+
"max_von_mises_stress_gauss_points[Pa]": max_mises_stress_gauss_points,
324+
"l2_error_displacement[m]": l2_error_displacement,
325+
"reaction_force_left_boundary_x[N]": reaction_left_x,
326+
"reaction_force_left_boundary_y[N]": reaction_left_y,
327+
f"displacement_x_at_evaluation_point (x={displacement_evaluation_x}, y={displacement_evaluation_y})[m]": displacement_x_at_evaluation_point,
330328
}
331329

332330
if MPI.COMM_WORLD.rank == 0:

examples/linear-elastic-plate-with-hole/kratos/Snakefile_kratos.smk renamed to examples/linear-elastic-plate-with-hole/kratos/Snakefile.smk

Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,18 @@ rule all:
99
"solution_metrics.json",
1010
"solution_field_data.zip"
1111

12+
rule create_mesh:
13+
input:
14+
script = "create_mesh.py",
15+
parameters = "parameters.json",
16+
output:
17+
mesh = "mesh.msh",
18+
conda: "environment_mesh.yml"
19+
shell:
20+
"""
21+
python3 {input.script} --input_parameter_file {input.parameters} --output_mesh_file {output.mesh}
22+
"""
23+
1224
rule mesh_to_mdpa:
1325
input:
1426
parameters = "parameters.json",

‎examples/linear-elastic-plate-with-hole/kratos/create_kratos_input.py‎

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -37,30 +37,30 @@ def create_kratos_input(
3737

3838
E = (
3939
ureg.Quantity(
40-
parameters["young-modulus"]["value"], parameters["young-modulus"]["unit"]
40+
parameters["young_modulus[Pa]"], "Pa"
4141
)
4242
.to_base_units()
4343
.magnitude
4444
)
4545
nu = (
4646
ureg.Quantity(
47-
parameters["poisson-ratio"]["value"], parameters["poisson-ratio"]["unit"]
47+
parameters["poisson_ratio"], ""
4848
)
4949
.to_base_units()
5050
.magnitude
5151
)
5252
radius = (
53-
ureg.Quantity(parameters["radius"]["value"], parameters["radius"]["unit"])
53+
ureg.Quantity(parameters["radius[m]"], "m")
5454
.to_base_units()
5555
.magnitude
5656
)
5757
L = (
58-
ureg.Quantity(parameters["length"]["value"], parameters["length"]["unit"])
58+
ureg.Quantity(parameters["length[m]"], "m")
5959
.to_base_units()
6060
.magnitude
6161
)
6262
load = (
63-
ureg.Quantity(parameters["load"]["value"], parameters["load"]["unit"])
63+
ureg.Quantity(parameters["load[MPa]"], "MPa")
6464
.to_base_units()
6565
.magnitude
6666
)

‎examples/linear-elastic-plate-with-hole/kratos/msh_to_mdpa.py‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -24,12 +24,12 @@ def msh_to_mdpa(parameter_file: str, mesh_file: str, mdpa_file: str):
2424
with open(parameter_file) as f:
2525
parameters = json.load(f)
2626
radius = (
27-
ureg.Quantity(parameters["radius"]["value"], parameters["radius"]["unit"])
27+
ureg.Quantity(parameters["radius[m]"], "m")
2828
.to_base_units()
2929
.magnitude
3030
)
3131
L = (
32-
ureg.Quantity(parameters["length"]["value"], parameters["length"]["unit"])
32+
ureg.Quantity(parameters["length[m]"], "m")
3333
.to_base_units()
3434
.magnitude
3535
)

‎examples/linear-elastic-plate-with-hole/kratos/postprocess_results.py‎

Lines changed: 16 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -24,47 +24,47 @@ def postprocess_results(input_parameter_file, input_result_vtk, output_metrics_f
2424

2525
E = (
2626
ureg.Quantity(
27-
parameters["young-modulus"]["value"], parameters["young-modulus"]["unit"]
27+
parameters["young_modulus[Pa]"], "Pa"
2828
)
2929
.to_base_units()
3030
.magnitude
3131
)
3232
nu = (
3333
ureg.Quantity(
34-
parameters["poisson-ratio"]["value"], parameters["poisson-ratio"]["unit"]
34+
parameters["poisson_ratio"], ""
3535
)
3636
.to_base_units()
3737
.magnitude
3838
)
3939
radius = (
40-
ureg.Quantity(parameters["radius"]["value"], parameters["radius"]["unit"])
40+
ureg.Quantity(parameters["radius[m]"], "m")
4141
.to_base_units()
4242
.magnitude
4343
)
4444
L = (
45-
ureg.Quantity(parameters["length"]["value"], parameters["length"]["unit"])
45+
ureg.Quantity(parameters["length[m]"], "m")
4646
.to_base_units()
4747
.magnitude
4848
)
4949
load = (
50-
ureg.Quantity(parameters["load"]["value"], parameters["load"]["unit"])
50+
ureg.Quantity(parameters["load[MPa]"], "MPa")
5151
.to_base_units()
5252
.magnitude
5353
)
5454

55-
displacement_evaluation_point = parameters["displacement-evaluation-point"]
55+
displacement_evaluation_point = parameters["displacement_evaluation_point[m]"]
5656
displacement_evaluation_x = (
5757
ureg.Quantity(
58-
displacement_evaluation_point["x"]["value"],
59-
displacement_evaluation_point["x"]["unit"],
58+
displacement_evaluation_point[0],
59+
"m",
6060
)
6161
.to_base_units()
6262
.magnitude
6363
)
6464
displacement_evaluation_y = (
6565
ureg.Quantity(
66-
displacement_evaluation_point["y"]["value"],
67-
displacement_evaluation_point["y"]["unit"],
66+
displacement_evaluation_point[1],
67+
"m",
6868
)
6969
.to_base_units()
7070
.magnitude
@@ -124,12 +124,12 @@ def postprocess_results(input_parameter_file, input_result_vtk, output_metrics_f
124124
displacement_x_at_evaluation_point = float(displacement_sampled[0, 0])
125125

126126
metrics = {
127-
"max_von_mises_stress_nodes": max_von_mises_stress_nodes,
128-
"max_von_mises_stress_gauss_points": max_von_mises_stress_gauss_points,
129-
"l2_error_displacement": l2_error_displacement,
130-
"reaction_force_left_boundary_x": reaction_force_left_boundary_x,
131-
"reaction_force_left_boundary_y": reaction_force_left_boundary_y,
132-
f"displacement_x_at_evaluation_point (x={displacement_evaluation_x}, y={displacement_evaluation_y})": displacement_x_at_evaluation_point,
127+
"max_von_mises_stress_nodes[Pa]": max_von_mises_stress_nodes,
128+
"max_von_mises_stress_gauss_points[Pa]": max_von_mises_stress_gauss_points,
129+
"l2_error_displacement[m]": l2_error_displacement,
130+
"reaction_force_left_boundary_x[N]": reaction_force_left_boundary_x,
131+
"reaction_force_left_boundary_y[N]": reaction_force_left_boundary_y,
132+
f"displacement_x_at_evaluation_point (x={displacement_evaluation_x}, y={displacement_evaluation_y})[m]": displacement_x_at_evaluation_point,
133133
}
134134
with open(output_metrics_file, "w") as f:
135135
json.dump(metrics, f, indent=4)

‎examples/linear-elastic-plate-with-hole/kratos/run_benchmark.py‎

Lines changed: 5 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -43,8 +43,8 @@
4343
for file in root_unzipped_benchmark_dir.glob("parameters_*.json"):
4444
with open(file, "r") as f:
4545
data = json.load(f)
46-
if data.get("element-size").get("value") >= 0.025:
47-
46+
#if data.get("element_size[m]") >= 0.025:
47+
if 1==1:
4848
# Create output directory for the configuration
4949
output_dir = root_unzipped_benchmark_dir / "results" / data.get("configuration")
5050
output_dir.mkdir(parents=True, exist_ok=True)
@@ -60,11 +60,8 @@
6060
continue
6161
else:
6262
shutil.copy(item, output_dir / item.name)
63-
64-
# Run the Snakemake workflow from the benchmark to create the mesh for the configuration
65-
subprocess.run(["snakemake", "--snakefile", "Snakefile", "--use-conda", "--force", "--cores", "all", "--conda-prefix", str(shared_env_dir), "--until", "create_mesh"], check=True, cwd=output_dir)
66-
67-
# Run the Snakemake workflow to run the simulation and postprocess results for the configuration
68-
subprocess.run(["snakemake", "--snakefile", "Snakefile_kratos.smk", "--use-conda", "--force", "--cores", "all", "--conda-prefix", str(shared_env_dir)], check=True, cwd=output_dir)
63+
64+
# Run the Snakemake workflow for the configuration
65+
subprocess.run(["snakemake", "--snakefile", "Snakefile.smk", "--use-conda", "--force", "--cores", "all", "--conda-prefix", str(shared_env_dir)], check=True, cwd=output_dir)
6966
print("Workflow executed successfully.")
7067

‎examples/linear-elastic-plate-with-hole/linear-elastic-plate-with-hole.zip‎

-1.26 KB
Binary file not shown.

0 commit comments

Comments
 (0)