Coverage for pySDC/projects/FEM_with_FEniCS/run_examples.py: 98%
87 statements
« prev ^ index » next coverage.py v7.17.0, created at 2026-11-02 12:30 +0000
1"""
2Run every declared combination and write the tables the README quotes to
3``data/fem_with_fenics_out.txt``, and the element-order study as a plot to
4``data/fem_with_fenics_orders.png``.
5"""
7from pathlib import Path
9import numpy as np
11from pySDC.helpers.stats_helper import get_sorted
12from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
13from pySDC.projects.FEM_with_FEniCS.setups import (
14 EXAMPLES,
15 get_coarsenings,
16 get_description,
17 get_families,
18 get_order_study,
19 get_pfasst_procs,
20)
22#: Cost of one solution restriction as a fraction of one fine-level sweep. Measured across the
23#: examples as 0.006-0.014 for the L2 projection and 0.023-0.037 for point sampling, so it is a
24#: 1-10%% correction either way -- small, but not zero, and the metric used to assume it was zero.
25RESTRICTION_COST = 0.01
28def run(example, nlevels=1, num_procs=1, **kwargs):
29 """
30 Run one configuration.
32 Returns
33 -------
34 dict
35 ``niter`` (mean over the steps), ``uend``, ``dofs`` per level and ``work``, the number of
36 fine-level sweep equivalents: iterations times the summed dof ratio of the hierarchy, plus
37 the solution restrictions, which are not free -- see :data:`RESTRICTION_COST`.
38 """
39 description, controller_params, t0, Tend = get_description(example, nlevels=nlevels, **kwargs)
40 controller = controller_nonMPI(num_procs=num_procs, controller_params=controller_params, description=description)
42 step = controller.MS[0]
43 dofs = [level.prob.init.dim() for level in step.levels]
44 prob = step.levels[0].prob
46 uend, stats = controller.run(u0=prob.u_exact(t0), t0=t0, Tend=Tend)
47 niter = np.mean([item[1] for item in get_sorted(stats, type='niter', sortby='time')])
49 # one restriction per coarse level per iteration, projecting the M node values and u0
50 num_nodes = step.levels[0].sweep.coll.num_nodes
51 restrictions = (num_nodes + 1) * (len(dofs) - 1) * RESTRICTION_COST
53 return {
54 'niter': niter,
55 'uend': uend,
56 'dofs': dofs,
57 'work': niter * (sum(n / dofs[0] for n in dofs) + restrictions),
58 }
61def compare_mlsdc(example, family='CG', coarsening='h', out=print, **kwargs):
62 """SDC against MLSDC with 2 and 3 levels. Returns the per-level results keyed by nlevels."""
63 results = {
64 nlevels: run(example, nlevels=nlevels, family=family, coarsening=coarsening, **kwargs) for nlevels in (1, 2, 3)
65 }
66 ref = results[1]
68 ladder = 'mesh coarsening' if coarsening == 'h' else 'order coarsening'
69 out(f'\n{example} [{family}, {coarsening}]: SDC vs MLSDC ({ladder}, collocation nodes kept)')
70 out(f' {"levels":>7s} {"dofs":>22s} {"niter":>7s} {"work":>7s} {"speed-up":>9s} {"|u - u_SDC|":>12s}')
71 for nlevels, res in results.items():
72 speedup = ref['work'] / res['work']
73 diff = abs(res['uend'] - ref['uend'])
74 out(
75 f' {nlevels:7d} {str(res["dofs"]):>22s} {res["niter"]:7.2f} {res["work"]:7.2f} '
76 f'{speedup:8.2f}x {diff:12.2e}'
77 )
78 return results
81def check_pfasst(example, nlevels=2, family='CG', coarsening='h', procs=None, out=print, **kwargs):
82 """PFASST over a growing number of parallel steps. Returns the results keyed by num_procs."""
83 procs = get_pfasst_procs(example) if procs is None else procs
84 results = {
85 p: run(example, nlevels=nlevels, num_procs=p, family=family, coarsening=coarsening, **kwargs) for p in procs
86 }
87 ref = results[procs[0]]
89 out(f'\n{example} [{family}, {coarsening}]: PFASST with {nlevels} levels, iterations over parallel steps')
90 out(f' {"procs":>6s} {"niter":>7s} {"|u - u_serial|":>15s}')
91 for p, res in results.items():
92 out(f' {p:6d} {res["niter"]:7.2f} {abs(res["uend"] - ref["uend"]):15.2e}')
93 return results
96def compare_orders(example, orders, family='CG', out=print, **kwargs):
97 """
98 MLSDC speed-up against the element order, at identical dof counts on every level.
100 This is the table behind "use high-order elements": the refinement ladder is shifted so that
101 ``CG1``, ``CG2`` and ``CG4`` all give the same dofs per level, which leaves the *quality* of the
102 coarse space as the only thing that varies. Returns the speed-ups keyed by (order, nlevels).
103 """
104 speedups = {}
105 for order in orders:
106 results = {
107 nlevels: run(example, nlevels=nlevels, family=family, coarsening='h', order=order, **kwargs)
108 for nlevels in (1, 2, 3)
109 }
110 for nlevels in (2, 3):
111 speedups[(order, nlevels)] = results[1]['work'] / results[nlevels]['work']
113 out(f'\n{example} [{family}]: MLSDC speed-up against element order, at equal dofs per level')
114 out(f' {"levels":>7s}' + ''.join(f'{"CG%d" % o:>12s}' for o in orders))
115 for nlevels in (2, 3):
116 out(f' {nlevels:7d}' + ''.join(f'{speedups[(o, nlevels)]:11.2f}x' for o in orders))
117 return speedups
120def plot_orders(speedups, fname='data/fem_with_fenics_orders.png', nlevels=3):
121 """
122 Plot the MLSDC speed-up against the element order, one line per example.
124 Parameters
125 ----------
126 speedups : dict
127 The results of :func:`compare_orders`, keyed by example.
128 fname : str
129 Where to save the figure.
130 nlevels : int
131 Which hierarchy to show, 2 or 3 levels.
132 """
133 import matplotlib.pyplot as plt
135 colors = ['#2a78d6', '#eb6834', '#1baf7a', '#eda100']
136 markers = ['o', 's', 'D', '^']
138 fig, ax = plt.subplots(figsize=(4.8, 2.6), constrained_layout=True)
139 ax.axhline(1.0, color='#8a8a85', lw=1, ls='--', zorder=1)
140 ax.text(0.01, 0.985, 'as fast as SDC', color='#5f5e5a', fontsize=9, va='top', transform=ax.get_yaxis_transform())
141 for i, (example, study) in enumerate(speedups.items()):
142 orders = sorted({order for order, _ in study})
143 ax.plot(
144 range(len(orders)),
145 [study[(order, nlevels)] for order in orders],
146 color=colors[i],
147 marker=markers[i],
148 ms=8,
149 lw=2.5,
150 mec='white',
151 mew=1.5,
152 label=example,
153 zorder=2,
154 )
155 ax.set_xticks(range(len(orders)), [f'CG{order}' for order in orders])
156 ax.set_xlabel('element order, at equal dofs per level')
157 ax.set_ylabel(f'speed-up over SDC,\n{nlevels}-level MLSDC')
158 ax.grid(axis='y', color='#e6e5df', lw=0.8)
159 ax.set_axisbelow(True)
160 for side in ('top', 'right'):
161 ax.spines[side].set_visible(False)
162 ax.legend(frameon=False, loc='upper left')
163 fig.savefig(fname, dpi=150)
164 plt.close(fig)
167def main():
168 Path('data').mkdir(parents=True, exist_ok=True)
169 with open('data/fem_with_fenics_out.txt', 'w') as f:
171 def out(line=''):
172 print(line)
173 f.write(str(line) + '\n')
175 out('FEniCS + pySDC, mass-matrix formulation throughout (no mass inverse anywhere).')
176 out(
177 f'work = iterations x [sum(dofs_l / dofs_0) + (M+1)(nlevels-1) x {RESTRICTION_COST}], '
178 'i.e. fine-level sweep equivalents including the solution restrictions.'
179 )
180 out('[family, coarsening]: CG/DG elements, h = coarser mesh, p = lower element order.')
181 speedups = {}
182 for example in EXAMPLES:
183 for family in get_families(example):
184 for coarsening in get_coarsenings(example):
185 compare_mlsdc(example, family=family, coarsening=coarsening, out=out)
186 check_pfasst(example, family=family, coarsening=coarsening, out=out)
187 if get_order_study(example):
188 speedups[example] = compare_orders(example, orders=get_order_study(example), out=out)
189 if speedups:
190 plot_orders(speedups)
193if __name__ == '__main__':
194 main()