Skip to content

Commit 74c0368

Browse files
th-skamth-skam
committed
[COR] Allow choosing name for ForceField
1 parent a5ece5b commit 74c0368

9 files changed

Lines changed: 52 additions & 25 deletions

File tree

‎examples/Freefem/MMS/1D/bar.py‎

Lines changed: 11 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -53,11 +53,14 @@ def compute(nodes, topology):
5353
return compute
5454

5555

56-
def build_bar_scene(root, mms, E_eff, nx):
56+
def build_bar_scene(root, mms, E_eff, nx,
57+
force_field="LinearSmallStrainFEMForceField"):
5758
"""Populate root with a static 1D bar scene on the non-dimensional domain [0,1].
5859
5960
BodyForce is assembled after init by BodyForceAssembler. BCs: Dirichlet at
6061
x=0, Neumann `mms.traction_bc(E_eff)` at x=1.
62+
63+
force_field : name of the FEM force field to test
6164
"""
6265
root.addObject('RequiredPlugin', pluginName=[
6366
"Elasticity",
@@ -103,7 +106,7 @@ def build_bar_scene(root, mms, E_eff, nx):
103106
name="dofs",
104107
template="Vec1d")
105108

106-
Bar.addObject('LinearSmallStrainFEMForceField',
109+
Bar.addObject(force_field,
107110
name="FEM",
108111
template="Vec1d",
109112
youngModulus=E_eff,
@@ -130,15 +133,16 @@ def case_scene(mms):
130133
"""Return a `createScene(rootNode)` bound to this MMS case."""
131134
def createScene(rootNode):
132135
cfg = load_params()
133-
build_bar_scene(rootNode, mms, cfg["E_eff"], cfg["reference"]["nx"])
136+
build_bar_scene(rootNode, mms, cfg["E_eff"], cfg["reference"]["nx"],
137+
force_field=cfg["forceField"])
134138
return rootNode
135139
return createScene
136140

137141

138-
def solve_bar(mms, E_eff, nx):
142+
def solve_bar(mms, E_eff, nx, force_field="LinearSmallStrainFEMForceField"):
139143
"""Build, run one static step, and return a BarSolution1D snapshot."""
140144
root = Sofa.Core.Node("root")
141-
build_bar_scene(root, mms, E_eff, nx)
145+
build_bar_scene(root, mms, E_eff, nx, force_field=force_field)
142146
Sofa.Simulation.init(root)
143147
Sofa.Simulation.animate(root, root.dt.value)
144148
Bar = root.Bar
@@ -180,7 +184,8 @@ def plot_solution(case, x, u_h, u_ex, label_ex):
180184
def run_reference_scene(mms):
181185
"""Solve one MMS case at the reference mesh, write the solution table and plot."""
182186
cfg = load_params()
183-
sol = solve_bar(mms, cfg["E_eff"], cfg["reference"]["nx"])
187+
sol = solve_bar(mms, cfg["E_eff"], cfg["reference"]["nx"],
188+
force_field=cfg["forceField"])
184189
l2 = l2_error_1d(sol.x0, sol.edges, sol.u_h, mms.u_ex, L2_QUADRATURE_1D)
185190
h1 = h1_semi_error_1d(sol.x0, sol.edges, sol.u_h, mms.du_ex, H1_QUADRATURE_1D)
186191
write_solution_table(f"solution_{mms.name}", sol.x0, sol.u_h, mms.u_ex,

‎examples/Freefem/MMS/1D/params.json‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,7 @@
11
{
22
"length": 2.0,
33
"youngModulus": 1e6,
4+
"forceField": "LinearSmallStrainFEMForceField",
45
"reference": {
56
"nx": 10
67
},

‎examples/Freefem/MMS/1D/run_convergence.py‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -24,7 +24,8 @@
2424
stem = f"{mms.name}_convergence"
2525
hs, errors = run_convergence_series(
2626
nx_values = cfg["convergence"]["nx_values"][mms.name],
27-
run_fn = lambda nx, _m=mms: solve_bar(_m, cfg["E_eff"], nx),
27+
run_fn = lambda nx, _m=mms: solve_bar(
28+
_m, cfg["E_eff"], nx, force_field=cfg["forceField"]),
2829
h_fn = lambda nx: 1.0 / (nx - 1),
2930
error_fns = {
3031
"L2": lambda sol, _m=mms: l2_error_1d(

‎examples/Freefem/MMS/2D/beam.py‎

Lines changed: 12 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -67,11 +67,13 @@ def compute(nodes, topology):
6767
return compute
6868

6969
def build_beam_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
70-
nx=10, ny=10, with_visual=True, dim="2d"):
70+
nx=10, ny=10, with_visual=True, dim="2d",
71+
force_field="LinearSmallStrainFEMForceField"):
7172
"""Build a SOFA scene for `mms` on the 2D `element` strategy.
7273
7374
dim : "2d" → Vec2d template / plane stress
7475
"3d" → Vec3d template / plane strain (z coordinate fixed at 0)
76+
force_field : name of the FEM force field to test
7577
Returns (dofs, topology). Nodes and connectivity become available
7678
after `Sofa.Simulation.init(root)` runs, via
7779
`dofs.rest_position.array()` and `element.read_connectivity(topology)`.
@@ -125,7 +127,7 @@ def build_beam_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
125127

126128
topology = element.add_topology(Beam)
127129

128-
Beam.addObject("LinearSmallStrainFEMForceField", name="FEM", template=tmpl,
130+
Beam.addObject(force_field, name="FEM", template=tmpl,
129131
youngModulus=E, poissonRatio=nu, topology="@topology")
130132

131133
mms.apply_bcs(Beam, nodes_2d, L, dim)
@@ -150,11 +152,13 @@ def build_beam_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
150152
# Simulation runner
151153
# ─────────────────────────────────────────────────────────────────────────────
152154

153-
def solve_beam(elem, mms, L, E, nu, nx, ny, dim="2d"):
155+
def solve_beam(elem, mms, L, E, nu, nx, ny, dim="2d",
156+
force_field="LinearSmallStrainFEMForceField"):
154157
"""Build, init, and run one static step. Returns a BeamSolution2D snapshot."""
155158
root = Sofa.Core.Node("root")
156159
dofs, topology = build_beam_scene(
157-
root, mms, elem, L=L, E=E, nu=nu, nx=nx, ny=ny, with_visual=False, dim=dim
160+
root, mms, elem, L=L, E=E, nu=nu, nx=nx, ny=ny, with_visual=False, dim=dim,
161+
force_field=force_field
158162
)
159163
Sofa.Simulation.init(root)
160164
# Read topology back from SOFA now that init has populated it.
@@ -246,9 +250,10 @@ def run_reference_scene(elem, mms):
246250
L, E = cfg["length"], cfg["youngModulus"]
247251
nu, dim = ref["nu"], ref["dim"]
248252
nx = ny = ref["nx"]
253+
ff = cfg["forceField"]
249254
hyp = "plane strain" if dim == "3d" else "plane stress"
250255

251-
sol = solve_beam(elem, mms, L, E, nu, nx, ny, dim=dim)
256+
sol = solve_beam(elem, mms, L, E, nu, nx, ny, dim=dim, force_field=ff)
252257
l2 = elem.compute_l2(sol, mms, L)
253258
h1 = elem.compute_h1(sol, mms, L)
254259

@@ -281,6 +286,7 @@ def createScene(rootNode):
281286
build_beam_scene(rootNode, mms, element,
282287
L=cfg["length"], E=cfg["youngModulus"],
283288
nu=ref["nu"], nx=ref["nx"], ny=ref["nx"],
284-
with_visual=True, dim=ref["dim"])
289+
with_visual=True, dim=ref["dim"],
290+
force_field=cfg["forceField"])
285291
return rootNode
286292
return createScene

‎examples/Freefem/MMS/2D/params.json‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,7 @@
11
{
22
"length": 1.0,
33
"youngModulus": 1e6,
4+
"forceField": "LinearSmallStrainFEMForceField",
45
"reference": {
56
"nx": 20,
67
"nu": 0.0,

‎examples/Freefem/MMS/2D/run_convergence.py‎

Lines changed: 7 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -19,14 +19,16 @@
1919
from output import plot_convergence
2020

2121

22-
def convergence_study(elem_specs, mms, L, E, nu, nx_values, dim="2d"):
22+
def convergence_study(elem_specs, mms, L, E, nu, nx_values, dim="2d",
23+
force_field="LinearSmallStrainFEMForceField"):
2324
"""
2425
Run a convergence series for each element type in elem_specs, write a
2526
per-(element) text table, and one shared plot with L²/H¹ for every
2627
element on the same axes.
2728
2829
elem_specs : list of dicts with keys 'elem', 'label', 'l2_style', 'h1_style'
2930
dim : "2d" (plane stress) or "3d" (plane strain)
31+
force_field : name of the FEM force field to test
3032
"""
3133
hyp = "plane strain" if dim == "3d" else "plane stress"
3234
print(f"\n PoissonRatio = {nu} ({hyp})", flush=True)
@@ -40,7 +42,7 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values, dim="2d"):
4042
hs, errors = run_convergence_series(
4143
nx_values = nx_values,
4244
run_fn = lambda nx, _e=elem: solve_beam(
43-
_e, mms, L, E, nu, nx, nx, dim=dim),
45+
_e, mms, L, E, nu, nx, nx, dim=dim, force_field=force_field),
4446
h_fn = lambda nx: L / (nx - 1),
4547
error_fns = {
4648
"L2": lambda sol, _e=elem: _e.compute_l2(sol, mms, L),
@@ -66,6 +68,7 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values, dim="2d"):
6668
cfg = load_params()
6769
L = cfg["length"]
6870
E = cfg["youngModulus"]
71+
ff = cfg["forceField"]
6972
conv = cfg["convergence"]
7073

7174
specs = [
@@ -81,4 +84,5 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values, dim="2d"):
8184
print(f"\n== {mms.name} ==")
8285
for DIM in conv["dim_values"]:
8386
for nu in conv["nu_values"]:
84-
convergence_study(specs, mms, L, E, nu, nx_vals, dim=DIM)
87+
convergence_study(specs, mms, L, E, nu, nx_vals, dim=DIM,
88+
force_field=ff)

‎examples/Freefem/MMS/3D/params.json‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,7 @@
11
{
22
"length": 1.0,
33
"youngModulus": 1e6,
4+
"forceField": "LinearSmallStrainFEMForceField",
45
"reference": {
56
"nx": 6,
67
"nu": 0.3

‎examples/Freefem/MMS/3D/run_convergence.py‎

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -17,13 +17,15 @@
1717
from output import plot_convergence
1818

1919

20-
def convergence_study(elem_specs, mms, L, E, nu, nx_values):
20+
def convergence_study(elem_specs, mms, L, E, nu, nx_values,
21+
force_field="LinearSmallStrainFEMForceField"):
2122
"""
2223
Run a convergence series for each element type in elem_specs, write a
2324
per-(element) text table, and one shared plot with L²/H¹ for every
2425
element on the same axes.
2526
2627
elem_specs : list of dicts with keys 'elem', 'label', 'l2_style', 'h1_style'
28+
force_field : name of the FEM force field to test
2729
"""
2830
print(f"\n PoissonRatio = {nu}", flush=True)
2931

@@ -36,7 +38,7 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values):
3638
hs, errors = run_convergence_series(
3739
nx_values = nx_values,
3840
run_fn = lambda nx, _e=elem: solve_solid(
39-
_e, mms, L, E, nu, nx, nx, nx),
41+
_e, mms, L, E, nu, nx, nx, nx, force_field=force_field),
4042
h_fn = lambda nx: L / (nx - 1),
4143
error_fns = {
4244
"L2": lambda sol, _e=elem: _e.compute_l2(sol, mms, L),
@@ -62,6 +64,7 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values):
6264
cfg = load_params()
6365
L = cfg["length"]
6466
E = cfg["youngModulus"]
67+
ff = cfg["forceField"]
6568
conv = cfg["convergence"]
6669

6770
specs = [
@@ -73,4 +76,4 @@ def convergence_study(elem_specs, mms, L, E, nu, nx_values):
7376
nx_vals = conv["nx_values"][mms.name]
7477
print(f"\n== {mms.name} ==")
7578
for nu in conv["nu_values"]:
76-
convergence_study(specs, mms, L, E, nu, nx_vals)
79+
convergence_study(specs, mms, L, E, nu, nx_vals, force_field=ff)

‎examples/Freefem/MMS/3D/solid.py‎

Lines changed: 11 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -56,12 +56,15 @@ def compute(nodes, topology):
5656
return compute
5757

5858
def build_solid_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
59-
nx=6, ny=6, nz=6, with_visual=True):
59+
nx=6, ny=6, nz=6, with_visual=True,
60+
force_field="LinearSmallStrainFEMForceField"):
6061
"""Build a SOFA scene for `mms` on the 3D `element` strategy.
6162
6263
Returns (dofs, topology). Nodes and connectivity become available
6364
after `Sofa.Simulation.init(root)` runs, via
6465
`dofs.rest_position.array()` and `element.read_connectivity(topology)`.
66+
67+
force_field : name of the FEM force field to test
6568
"""
6669
rootNode.addObject("RequiredPlugin", pluginName=[
6770
"Elasticity",
@@ -106,7 +109,7 @@ def build_solid_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
106109

107110
topology = element.add_topology(Solid)
108111

109-
Solid.addObject("LinearSmallStrainFEMForceField", name="FEM", template="Vec3d",
112+
Solid.addObject(force_field, name="FEM", template="Vec3d",
110113
youngModulus=E, poissonRatio=nu, topology="@topology")
111114

112115
mms.apply_bcs(Solid, nodes_3d, L)
@@ -130,12 +133,13 @@ def build_solid_scene(rootNode, mms, element, L=1.0, E=1e6, nu=0.3,
130133
# Simulation runner
131134
# ─────────────────────────────────────────────────────────────────────────────
132135

133-
def solve_solid(elem, mms, L, E, nu, nx, ny, nz):
136+
def solve_solid(elem, mms, L, E, nu, nx, ny, nz,
137+
force_field="LinearSmallStrainFEMForceField"):
134138
"""Build, init, and run one static step. Returns a SolidSolution3D snapshot."""
135139
root = Sofa.Core.Node("root")
136140
dofs, topology = build_solid_scene(
137141
root, mms, elem, L=L, E=E, nu=nu,
138-
nx=nx, ny=ny, nz=nz, with_visual=False
142+
nx=nx, ny=ny, nz=nz, with_visual=False, force_field=force_field
139143
)
140144
Sofa.Simulation.init(root)
141145
nodes_3d = dofs.rest_position.array().copy()
@@ -265,8 +269,9 @@ def run_reference_scene(elem, mms):
265269
L, E = cfg["length"], cfg["youngModulus"]
266270
nu = ref["nu"]
267271
nx = ny = nz = ref["nx"]
272+
ff = cfg["forceField"]
268273

269-
sol = solve_solid(elem, mms, L, E, nu, nx, ny, nz)
274+
sol = solve_solid(elem, mms, L, E, nu, nx, ny, nz, force_field=ff)
270275
l2 = elem.compute_l2(sol, mms, L)
271276
h1 = elem.compute_h1(sol, mms, L)
272277

@@ -300,6 +305,6 @@ def createScene(rootNode):
300305
L=cfg["length"], E=cfg["youngModulus"],
301306
nu=ref["nu"],
302307
nx=ref["nx"], ny=ref["nx"], nz=ref["nx"],
303-
with_visual=True)
308+
with_visual=True, force_field=cfg["forceField"])
304309
return rootNode
305310
return createScene

0 commit comments

Comments
 (0)