diff --git a/publications/blog-posts/figures/rotated-boundary-conditions/generate-rotated-basis.py b/publications/blog-posts/figures/rotated-boundary-conditions/generate-rotated-basis.py new file mode 100644 index 000000000..b203819ba --- /dev/null +++ b/publications/blog-posts/figures/rotated-boundary-conditions/generate-rotated-basis.py @@ -0,0 +1,102 @@ +"""Geometry for the rotated-boundary-conditions figure. + +The cetz skill's rule: geometry computation happens in Python and Typst just +draws. Doing the triangulation here rather than in the figure is what stops +nodes being left out of the mesh -- an earlier version connected nodes by a +distance threshold and silently missed several. + +Emits the schema in the skill's `underworld-bridge.md`: + + {"vertices": [[x, y], ...], + "triangles": [[i, j, k], ...], + "surface": [i, ...], indices of the constrained nodes + "frames": [{"p": [x, y], "n": [nx, ny], "t": [tx, ty]}, ...], + "curve": [[x, y], ...]} the surface, finely sampled + +Run: python3 generate-rotated-basis.py +""" +import json +import pathlib + +import numpy as np +from scipy.spatial import Delaunay + +OUT = pathlib.Path(__file__).with_name("rotated-basis-data.json") + +X0, X1 = -3.5, 3.5 +BASE = -2.9 +NX = 9 # columns of nodes +NY = 4 # rows, surface included + + +def surface_y(x): + """A free surface: rises on the left, falls on the right, with an + inflection between, so the normal swings through a wide range and the + curvature changes sign. No global rotation straightens this out, which is + the reason the figure exists.""" + x = np.asarray(x, dtype=float) + return (0.95 * np.exp(-(((x + 1.75) / 1.30) ** 2)) + - 0.80 * np.exp(-(((x - 1.70) / 1.15) ** 2)) + + 0.75) + + +def surface_slope(x, eps=1.0e-4): + return (surface_y(x + eps) - surface_y(x - eps)) / (2 * eps) + + +# Nodes on a grid warped to sit under the surface. Columns are staggered on +# alternate rows so the Delaunay triangulation comes out as triangles rather +# than as near-degenerate right angles on a perfect lattice. +xs = np.linspace(X0, X1, NX) +pts = [] +surface_idx = [] +for row in range(NY): + frac = row / (NY - 1) # 0 at the base, 1 at the surface + offset = 0.0 if row % 2 == 0 else 0.5 * (xs[1] - xs[0]) + cols = xs + offset + if row == NY - 1: + cols = xs # surface row unstaggered + for x in cols: + if x < X0 - 1e-9 or x > X1 + 1e-9: + continue + top = float(surface_y(x)) + y = BASE + (top - BASE) * frac + if row == NY - 1: + surface_idx.append(len(pts)) + pts.append([float(x), float(y)]) + +pts = np.array(pts) +tri = Delaunay(pts) + +# Drop the slivers Delaunay leaves along a non-convex top edge: any triangle +# whose centroid sits above the surface is outside the domain. +keep = [] +for simplex in tri.simplices: + c = pts[simplex].mean(axis=0) + if c[1] <= float(surface_y(c[0])) + 1.0e-9: + keep.append([int(i) for i in simplex]) + +frames = [] +for i in surface_idx: + x = pts[i][0] + m = float(surface_slope(x)) + n = np.array([-m, 1.0]) + n /= np.linalg.norm(n) + t = np.array([1.0, m]) + t /= np.linalg.norm(t) + frames.append({"p": [float(pts[i][0]), float(pts[i][1])], + "n": [float(n[0]), float(n[1])], + "t": [float(t[0]), float(t[1])]}) + +curve_x = np.linspace(X0, X1, 121) +data = { + "vertices": [[float(a), float(b)] for a, b in pts], + "triangles": keep, + "surface": [int(i) for i in surface_idx], + "frames": frames, + "curve": [[float(a), float(b)] for a, b in zip(curve_x, surface_y(curve_x))], +} +OUT.write_text(json.dumps(data, indent=1)) +print("wrote %s: %d vertices, %d triangles, %d surface nodes" + % (OUT.name, len(data["vertices"]), len(data["triangles"]), + len(data["surface"]))) diff --git a/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis-data.json b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis-data.json new file mode 100644 index 000000000..76c066ed5 --- /dev/null +++ b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis-data.json @@ -0,0 +1,1011 @@ +{ + "vertices": [ + [ + -3.5, + -2.9 + ], + [ + -2.625, + -2.9 + ], + [ + -1.75, + -2.9 + ], + [ + -0.875, + -2.9 + ], + [ + 0.0, + -2.9 + ], + [ + 0.875, + -2.9 + ], + [ + 1.75, + -2.9 + ], + [ + 2.625, + -2.9 + ], + [ + 3.5, + -2.9 + ], + [ + -3.0625, + -1.5690676390839076 + ], + [ + -2.1875, + -1.4005781953894207 + ], + [ + -1.3125, + -1.4008544079276464 + ], + [ + -0.4375, + -1.5774930208125677 + ], + [ + 0.4375, + -1.744573187635388 + ], + [ + 1.3125, + -1.9201468750186008 + ], + [ + 2.1875, + -1.9061055410007908 + ], + [ + 3.0625, + -1.7488486663272358 + ], + [ + -3.5, + -0.3632396002330607 + ], + [ + -2.625, + -0.06405866884136113 + ], + [ + -1.75, + 0.1666008481044865 + ], + [ + -0.875, + -0.0676029423016371 + ], + [ + 0.0, + -0.4232124172195797 + ], + [ + 0.875, + -0.7747080243853817 + ], + [ + 1.75, + -0.9985423178094506 + ], + [ + 2.625, + -0.7459268346817169 + ], + [ + 3.5, + -0.5126942107107557 + ], + [ + -3.5, + 0.9051405996504092 + ], + [ + -2.625, + 1.3539119967379585 + ], + [ + -1.75, + 1.6999012721567301 + ], + [ + -0.875, + 1.3485955865475447 + ], + [ + 0.0, + 0.8151813741706304 + ], + [ + 0.875, + 0.2879379634219279 + ], + [ + 1.75, + -0.047813476714175795 + ], + [ + 2.625, + 0.331109747977425 + ], + [ + 3.5, + 0.6809586839338668 + ] + ], + "triangles": [ + [ + 31, + 30, + 21 + ], + [ + 30, + 20, + 21 + ], + [ + 26, + 17, + 18 + ], + [ + 27, + 26, + 18 + ], + [ + 24, + 25, + 33 + ], + [ + 25, + 34, + 33 + ], + [ + 20, + 12, + 21 + ], + [ + 12, + 20, + 11 + ], + [ + 16, + 7, + 8 + ], + [ + 25, + 16, + 8 + ], + [ + 16, + 25, + 24 + ], + [ + 32, + 24, + 33 + ], + [ + 20, + 19, + 11 + ], + [ + 27, + 19, + 28 + ], + [ + 19, + 27, + 18 + ], + [ + 29, + 20, + 30 + ], + [ + 19, + 29, + 28 + ], + [ + 29, + 19, + 20 + ], + [ + 17, + 9, + 18 + ], + [ + 9, + 17, + 0 + ], + [ + 13, + 12, + 4 + ], + [ + 5, + 13, + 4 + ], + [ + 13, + 5, + 14 + ], + [ + 12, + 13, + 21 + ], + [ + 12, + 3, + 4 + ], + [ + 3, + 12, + 11 + ], + [ + 2, + 3, + 11 + ], + [ + 16, + 15, + 7 + ], + [ + 15, + 16, + 24 + ], + [ + 5, + 6, + 14 + ], + [ + 15, + 6, + 7 + ], + [ + 6, + 15, + 14 + ], + [ + 32, + 23, + 24 + ], + [ + 23, + 15, + 24 + ], + [ + 15, + 23, + 14 + ], + [ + 9, + 10, + 18 + ], + [ + 10, + 19, + 18 + ], + [ + 19, + 10, + 11 + ], + [ + 10, + 2, + 11 + ], + [ + 1, + 9, + 0 + ], + [ + 10, + 1, + 2 + ], + [ + 1, + 10, + 9 + ], + [ + 22, + 32, + 31 + ], + [ + 22, + 23, + 32 + ], + [ + 22, + 31, + 21 + ], + [ + 13, + 22, + 21 + ], + [ + 22, + 13, + 14 + ], + [ + 23, + 22, + 14 + ] + ], + "surface": [ + 26, + 27, + 28, + 29, + 30, + 31, + 32, + 33, + 34 + ], + "frames": [ + { + "p": [ + -3.5, + 0.9051405996504092 + ], + "n": [ + -0.3058957249037569, + 0.9520650216690062 + ], + "t": [ + 0.9520650216690062, + 0.3058957249037569 + ] + }, + { + "p": [ + -2.625, + 1.3539119967379585 + ], + "n": [ + -0.5302120052673807, + 0.8478651009862024 + ], + "t": [ + 0.8478651009862024, + 0.5302120052673807 + ] + }, + { + "p": [ + -1.75, + 1.6999012721567301 + ], + "n": [ + 0.0005151017428689071, + 0.9999998673350886 + ], + "t": [ + 0.9999998673350886, + -0.0005151017428689071 + ] + }, + { + "p": [ + -0.875, + 1.3485955865475447 + ], + "n": [ + 0.5426587927926065, + 0.8399532335820077 + ], + "t": [ + 0.8399532335820077, + -0.5426587927926065 + ] + }, + { + "p": [ + 0.0, + 0.8151813741706304 + ], + "n": [ + 0.4836463257494078, + 0.8752635212260921 + ], + "t": [ + 0.8752635212260921, + -0.4836463257494078 + ] + }, + { + "p": [ + 0.875, + 0.2879379634219279 + ], + "n": [ + 0.5429867639410579, + 0.8397412543068359 + ], + "t": [ + 0.8397412543068358, + -0.5429867639410579 + ] + }, + { + "p": [ + 1.75, + -0.047813476714175795 + ], + "n": [ + -0.0574834315887112, + 0.9983464604498711 + ], + "t": [ + 0.9983464604498711, + 0.0574834315887112 + ] + }, + { + "p": [ + 2.625, + 0.331109747977425 + ], + "n": [ + -0.5055403695426544, + 0.8628029524536158 + ], + "t": [ + 0.8628029524536158, + 0.5055403695426544 + ] + }, + { + "p": [ + 3.5, + 0.6809586839338668 + ], + "n": [ + -0.1847046296744634, + 0.982794078012693 + ], + "t": [ + 0.982794078012693, + 0.1847046296744634 + ] + } + ], + "curve": [ + [ + -3.5, + 0.9051405996504089 + ], + [ + -3.441666666666667, + 0.924709919910546 + ], + [ + -3.3833333333333333, + 0.9459569983210258 + ], + [ + -3.325, + 0.9689047102584147 + ], + [ + -3.2666666666666666, + 0.9935569625066398 + ], + [ + -3.2083333333333335, + 1.019896408691709 + ], + [ + -3.15, + 1.0478823531124308 + ], + [ + -3.091666666666667, + 1.077448917340561 + ], + [ + -3.033333333333333, + 1.1085035443807112 + ], + [ + -2.975, + 1.1409259131181626 + ], + [ + -2.9166666666666665, + 1.1745673310574358 + ], + [ + -2.8583333333333334, + 1.2092506658746796 + ], + [ + -2.8, + 1.2447708660876644 + ], + [ + -2.7416666666666667, + 1.280896108319149 + ], + [ + -2.6833333333333336, + 1.31736959344396 + ], + [ + -2.625, + 1.353911996737958 + ], + [ + -2.5666666666666664, + 1.3902245584714161 + ], + [ + -2.5083333333333333, + 1.425992781792527 + ], + [ + -2.45, + 1.4608906848912475 + ], + [ + -2.3916666666666666, + 1.4945855350373591 + ], + [ + -2.333333333333333, + 1.526742973893888 + ], + [ + -2.275, + 1.5570324272567517 + ], + [ + -2.216666666666667, + 1.5851326787631834 + ], + [ + -2.158333333333333, + 1.6107374767717508 + ], + [ + -2.1, + 1.633561037067616 + ], + [ + -2.041666666666667, + 1.6533433016773964 + ], + [ + -1.9833333333333334, + 1.6698548161227493 + ], + [ + -1.925, + 1.6829010939636901 + ], + [ + -1.8666666666666667, + 1.6923263483659823 + ], + [ + -1.8083333333333333, + 1.6980164853777744 + ], + [ + -1.75, + 1.6999012721567306 + ], + [ + -1.6916666666666667, + 1.6979556149374022 + ], + [ + -1.6333333333333333, + 1.692199905332323 + ], + [ + -1.575, + 1.682699418789918 + ], + [ + -1.5166666666666666, + 1.6695627748036432 + ], + [ + -1.4583333333333335, + 1.6529394938808053 + ], + [ + -1.4, + 1.6330167104634032 + ], + [ + -1.3416666666666668, + 1.61001512313947 + ], + [ + -1.2833333333333332, + 1.5841842828851125 + ], + [ + -1.225, + 1.5557973361588586 + ], + [ + -1.1666666666666665, + 1.5251453520084297 + ], + [ + -1.1083333333333334, + 1.492531370690367 + ], + [ + -1.0499999999999998, + 1.4582643155616355 + ], + [ + -0.9916666666666667, + 1.4226529102631265 + ], + [ + -0.9333333333333331, + 1.3859997397146477 + ], + [ + -0.875, + 1.348595586547545 + ], + [ + -0.8166666666666664, + 1.3107141647846698 + ], + [ + -0.7583333333333333, + 1.272607360377266 + ], + [ + -0.7000000000000002, + 1.2345010741961506 + ], + [ + -0.6416666666666666, + 1.1965917478187014 + ], + [ + -0.5833333333333335, + 1.1590436364841392 + ], + [ + -0.5249999999999999, + 1.1219868773708934 + ], + [ + -0.4666666666666668, + 1.0855163852548482 + ], + [ + -0.4083333333333332, + 1.0496915919046683 + ], + [ + -0.3500000000000001, + 1.0145370304190584 + ], + [ + -0.2916666666666665, + 0.9800437511642502 + ], + [ + -0.2333333333333334, + 0.9461715419910942 + ], + [ + -0.17499999999999982, + 0.9128519118949112 + ], + [ + -0.1166666666666667, + 0.8799917840856701 + ], + [ + -0.058333333333333126, + 0.8474778314169195 + ], + [ + 0.0, + 0.8151813741706305 + ], + [ + 0.05833333333333357, + 0.7829637472749286 + ], + [ + 0.1166666666666667, + 0.7506820312100971 + ], + [ + 0.17500000000000027, + 0.7181950283314438 + ], + [ + 0.2333333333333334, + 0.6853693544451201 + ], + [ + 0.2916666666666665, + 0.6520855046989421 + ], + [ + 0.3500000000000001, + 0.6182437438110678 + ], + [ + 0.4083333333333332, + 0.5837696640761547 + ], + [ + 0.4666666666666668, + 0.548619251245789 + ], + [ + 0.5250000000000004, + 0.5127832990715608 + ], + [ + 0.583333333333333, + 0.4762910187656371 + ], + [ + 0.6416666666666666, + 0.4392127004940097 + ], + [ + 0.7000000000000002, + 0.4016613007004925 + ], + [ + 0.7583333333333337, + 0.363792851740313 + ], + [ + 0.8166666666666664, + 0.325805618850044 + ], + [ + 0.875, + 0.28793796342192796 + ], + [ + 0.9333333333333336, + 0.2504649100526243 + ], + [ + 0.9916666666666671, + 0.2136934567143398 + ], + [ + 1.0499999999999998, + 0.17795671114661105 + ], + [ + 1.1083333333333334, + 0.14360698042499975 + ], + [ + 1.166666666666667, + 0.11100798268452128 + ], + [ + 1.2249999999999996, + 0.08052638813756408 + ], + [ + 1.2833333333333332, + 0.05252292884049681 + ], + [ + 1.3416666666666668, + 0.02734334129772975 + ], + [ + 1.4000000000000004, + 0.005309421388709379 + ], + [ + 1.458333333333333, + -0.013289523913108559 + ], + [ + 1.5166666666666666, + -0.028204549692437686 + ], + [ + 1.5750000000000002, + -0.03923401019190231 + ], + [ + 1.6333333333333337, + -0.04622932722506001 + ], + [ + 1.6916666666666664, + -0.049099395459559214 + ], + [ + 1.75, + -0.047813476714175684 + ], + [ + 1.8083333333333336, + -0.042402488429535534 + ], + [ + 1.8666666666666671, + -0.032958648592165685 + ], + [ + 1.9249999999999998, + -0.01963349805658543 + ], + [ + 1.9833333333333334, + -0.002634378770072665 + ], + [ + 2.041666666666667, + 0.017780499734667887 + ], + [ + 2.0999999999999996, + 0.04130822544601709 + ], + [ + 2.158333333333333, + 0.06760836083270272 + ], + [ + 2.216666666666667, + 0.09631088206901939 + ], + [ + 2.2750000000000004, + 0.12702461196713855 + ], + [ + 2.333333333333333, + 0.15934587047687743 + ], + [ + 2.3916666666666666, + 0.19286707036237616 + ], + [ + 2.45, + 0.22718500012101506 + ], + [ + 2.5083333333333337, + 0.26190856024390546 + ], + [ + 2.5666666666666664, + 0.2966657509787571 + ], + [ + 2.625, + 0.3311097479774251 + ], + [ + 2.6833333333333336, + 0.36492394452614507 + ], + [ + 2.741666666666667, + 0.39782588332598795 + ], + [ + 2.8, + 0.429570044919817 + ], + [ + 2.8583333333333334, + 0.4599495019114082 + ], + [ + 2.916666666666667, + 0.4887964864025283 + ], + [ + 2.9750000000000005, + 0.5159819512151704 + ], + [ + 3.033333333333333, + 0.5414142324664324 + ], + [ + 3.091666666666667, + 0.5650369413091443 + ], + [ + 3.1500000000000004, + 0.5868262259143289 + ], + [ + 3.208333333333333, + 0.6067875511844959 + ], + [ + 3.2666666666666666, + 0.6249521436990915 + ], + [ + 3.325, + 0.6413732437141221 + ], + [ + 3.3833333333333337, + 0.6561222955694307 + ], + [ + 3.4416666666666664, + 0.6692851936251155 + ], + [ + 3.5, + 0.680958683933867 + ] + ] +} \ No newline at end of file diff --git a/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.png b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.png new file mode 100644 index 000000000..9d9ca0089 Binary files /dev/null and b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.png differ diff --git a/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.svg b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.svg new file mode 100644 index 000000000..85ba6b2ab --- /dev/null +++ b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.typ b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.typ new file mode 100644 index 000000000..0480ba4c8 --- /dev/null +++ b/publications/blog-posts/figures/rotated-boundary-conditions/rotated-basis.typ @@ -0,0 +1,180 @@ +// Rotated boundary conditions: where the basis changes, and where in the +// solver the rotation lives. +// +// Thesis: rotating the degrees of freedom leaves the discrete problem in a +// MIXED basis, but the obligation is CONTAINED -- it lives in the velocity +// block and its multigrid substructure, and the Schur/pressure machinery +// wrapping it never handles a rotated vector. +// +// House style follows publications/blog-posts/figures/arrays-sync-flow.typ +// (cetz 0.3.4, hex colours, helper-function pattern). + +#import "@preview/cetz:0.3.4" + +#set page(width: auto, height: auto, margin: 16pt) +#set text(size: 10pt) + +#cetz.canvas({ + import cetz.draw: * + + // Colours + let cart-fg = rgb("#4a7bf7") + let cart-bg = rgb("#dce8fc") + let rot-fg = rgb("#e57373") + let rot-bg = rgb("#fce4ec") + let plain-bg = rgb("#f2f2f0") + let mg-bg = rgb("#fdf1f4") + let schur-fg = rgb("#49a87c") // the un-rotated half + let schur-bg = rgb("#e8f4e8") + let ink = rgb("#1a1a1a") + let muted = rgb("#7a7a7a") + let hair = rgb("#d0d0d0") + + // ====================================================================== + // LEFT: a free surface, meshed. Geometry (nodes, Delaunay triangles, + // surface frames) comes from generate-rotated-basis.py via JSON -- the + // skill's rule, and the reason every node is now in the triangulation. + // An earlier version joined nodes by a distance threshold and left some out. + // ====================================================================== + let data = json("rotated-basis-data.json") + let vtx = data.vertices + let at(i) = (vtx.at(i).at(0), vtx.at(i).at(1)) + let is-surface(i) = data.surface.contains(i) + + // triangles first, then the surface, then the nodes: painter's algorithm + for tri in data.triangles { + line(at(tri.at(0)), at(tri.at(1)), at(tri.at(2)), close: true, + stroke: (paint: hair, thickness: 0.55pt)) + } + + line(..data.curve.map(q => (q.at(0), q.at(1))), + stroke: (paint: ink, thickness: 1.2pt)) + + // A rotated frame at every constrained node: n the outward surface normal, + // t the tangent. Both follow the node, which is the whole point. + for (k, f) in data.frames.enumerate() { + let p = (f.p.at(0), f.p.at(1)) + let n = f.n + let tg = f.t + let ln = 0.72 + let lt = 0.46 + line(p, (p.at(0) + ln * n.at(0), p.at(1) + ln * n.at(1)), + mark: (end: ">", scale: 0.4), stroke: (paint: rot-fg, thickness: 1.1pt)) + line(p, (p.at(0) + lt * tg.at(0), p.at(1) + lt * tg.at(1)), + mark: (end: ">", scale: 0.4), stroke: (paint: rot-fg, thickness: 1.1pt)) + if k == 4 { + content((p.at(0) + 1.20 * ln * n.at(0), p.at(1) + 1.20 * ln * n.at(1)), + text(fill: rot-fg, size: 9pt, $n$)) + content((p.at(0) + 1.55 * lt * tg.at(0) + 0.24 * n.at(0), + p.at(1) + 1.55 * lt * tg.at(1) + 0.24 * n.at(1)), + text(fill: rot-fg, size: 9pt, $t$)) + } + } + + // the Cartesian frame at one interior node -- identical at every other one, + // which is what makes the surface the odd one out + let ip = at(12) // an interior node with room around it + line(ip, (ip.at(0) + 0.66, ip.at(1)), mark: (end: ">", scale: 0.4), + stroke: (paint: cart-fg, thickness: 1.1pt)) + line(ip, (ip.at(0), ip.at(1) + 0.66), mark: (end: ">", scale: 0.4), + stroke: (paint: cart-fg, thickness: 1.1pt)) + content((ip.at(0) + 0.88, ip.at(1)), text(fill: cart-fg, size: 9pt, $x$)) + content((ip.at(0), ip.at(1) + 0.88), text(fill: cart-fg, size: 9pt, $y$)) + + for i in range(vtx.len()) { + if is-surface(i) { + circle(at(i), radius: 0.13, fill: rot-fg, + stroke: (paint: rot-fg, thickness: 1pt)) + } else { + circle(at(i), radius: 0.11, fill: cart-bg, + stroke: (paint: cart-fg, thickness: 1pt)) + } + } + + content((0, 3.05), text(weight: "bold", size: 10pt, fill: ink, + "A free surface has no preferred direction")) + + circle((-3.3, -3.55), radius: 0.13, fill: rot-fg, stroke: (paint: rot-fg)) + content((-3.05, -3.55), anchor: "west", + text(size: 9pt, fill: ink, [surface node --- solve for $(v_n, v_t)$, hold $v_n$])) + circle((-3.3, -4.1), radius: 0.11, fill: cart-bg, + stroke: (paint: cart-fg, thickness: 1pt)) + content((-3.05, -4.1), anchor: "west", + text(size: 9pt, fill: ink, [interior node --- solve for $(v_x, v_y)$])) + + line((4.35, -4.9), (4.35, 3.3), stroke: (paint: hair, thickness: 0.8pt)) + + // ====================================================================== + // RIGHT: where the rotation lives. Two blocks ABUTTING, not nested: the + // velocity solve is rotated, the Schur/pressure solve is not, and the + // single un-rotation on the boundary between them feeds both the Schur + // solve and everything outside. + // ====================================================================== + let panel(tl, br, title, subtitle, bg, edge) = { + import cetz.draw: * + rect(tl, br, fill: bg, stroke: (paint: edge, thickness: 1pt), radius: 5pt) + content((tl.at(0) + 0.30, tl.at(1) - 0.36), anchor: "west", + text(weight: "bold", size: 9.5pt, fill: edge, title)) + if subtitle != none { + // cetz content() lays out at natural width and spills over the rect. + // Box it to the panel's own width so the text wraps inside the border. + // The canvas default is 1cm per unit, so the arithmetic is direct. + let tw = (br.at(0) - tl.at(0) - 0.60) * 1cm + content((tl.at(0) + 0.30, tl.at(1) - 0.90), anchor: "north-west", + box(width: tw, text(size: 8.5pt, fill: muted, subtitle))) + } + } + + content((10.3, 3.05), text(weight: "bold", size: 10pt, fill: ink, + "Where the rotation lives")) + + // -- the rotated half --------------------------------------------------- + panel((4.9, 1.65), (11.75, -3.70), [Velocity solve --- rotated], + [$hat(A) = Q^T A Q$, #h(0.7em) $hat(b) = Q^T b$, #h(0.7em) + $v_n$ held strongly at surface nodes], rot-bg, rot-fg) + + panel((5.35, -0.15), (11.50, -3.25), "Multigrid", + [the transfers are the only further obligation], mg-bg, rot-fg) + + let mgrow(y, lhs, rhs) = { + import cetz.draw: * + content((5.65, y), anchor: "west", text(size: 9pt, fill: ink, lhs)) + content((8.05, y), anchor: "west", text(size: 9pt, fill: rot-fg, rhs)) + } + mgrow(-1.85, [prolongation], [$P -> Q^T P$]) + mgrow(-2.40, [coarse operators], [inherit $Q$ via $R A P$]) + mgrow(-2.93, [coarse solve], [SVD (rigid rotations)]) + + // -- the un-rotated half, abutting -------------------------------------- + panel((13.35, 1.65), (17.85, -1.35), [Fieldsplit / Schur solve], + [pressure and constraints. Isotropic, and carrying no boundary + condition of this kind, so it never sees a rotated vector.], + schur-bg, schur-fg) + + // -- one un-rotation on the boundary, feeding both ---------------------- + // The junction sits in the gap BETWEEN the two blocks, low enough to clear + // the velocity block's own subtitle -- one un-rotation, two consumers. + let jx = 12.55 + let jy = -0.40 + line((11.75, jy), (jx, jy), stroke: (paint: rot-fg, thickness: 1pt)) + circle((jx, jy), radius: 0.075, fill: rot-fg, stroke: (paint: rot-fg)) + content((jx, jy + 0.45), text(size: 9pt, fill: rot-fg, $v = Q hat(v)$)) + + // branch 1: into the Schur solve, which needs it un-rotated + line((jx, jy), (13.29, jy), mark: (end: ">", scale: 0.4), + stroke: (paint: rot-fg, thickness: 1pt)) + + // branch 2: out to everything else + line((jx, jy), (jx, -2.70), stroke: (paint: rot-fg, thickness: 1pt)) + line((jx, -2.70), (13.20, -2.70), mark: (end: ">", scale: 0.4), + stroke: (paint: rot-fg, thickness: 1pt)) + content((13.33, -2.70), anchor: "west", text(size: 8.5pt, fill: rot-fg, + [and to everything outside ---\ output, advection, the surface update])) + + // convention, stated -- a reader who assumes the transpose reads the + // whole figure backwards + line((4.75, -5.05), (18.1, -5.05), stroke: (paint: hair, thickness: 0.8pt)) + content((4.75, -5.55), anchor: "west", text(size: 9pt, fill: ink, + [Convention: the columns of $Q$ are the nodal frame, so $hat(v) = Q^T v$, + and $Q = I$ at every unconstrained node.])) +})