diff --git a/research/physical-gaussian-expansion/EVIDENCE.md b/research/physical-gaussian-expansion/EVIDENCE.md new file mode 100644 index 0000000000..a9cccf6c4f --- /dev/null +++ b/research/physical-gaussian-expansion/EVIDENCE.md @@ -0,0 +1,107 @@ +# Verification evidence and limits + +Jose Manuel (`jmbarrancoml`), with substantial OpenAI Codex assistance. +I ran the local checks on Darwin arm64 with Python 3.12.13. I used the +repository's `run_python_replay` launcher with a fresh bytecode namespace, +`-I -O -Werror -B`, explicit process-group timeouts and no external Python +dependencies. The checks use explicit exceptions, so `-O` does not disable +their assertions. [VERIFY.json](VERIFY.json) registers eight commands with +1180 seconds of total timeout allowance. + +## Mathematical and implementation coverage + +| Suite | Evidence | +|---|---| +| sources | Complete 43-file core export, correction commit/archive, 14 live input pins, seven local runtime modules, one physical Count/Tape cohort and four binding controls | +| budget | 1248 exact production guard cases, complete block/tile word invoice, separately charged setup and inherited Fourier comparison | +| geometry | 78 exact exponent/crop identities, negative origins, period seams, `u*theta=4/7,4/13`, equal right endpoints and three encoded production radii | +| production | Four complete production tile descriptors at `p=4800,alpha=4,d=300`, including both period seams and a one-record ceiling transition; ten rejection controls | +| block | Two complete nonzero producer/block joins at p4/P32 and p8/P48, separate staged and periodic references, premature-extraction controls, 34 and 38 actual word calls | +| tile_overlap | Two complete blocks with equal source starts, values `[(3,0),(3,1)]`, 88 actual word calls, six fresh records and five retained overlap records | +| tile_wrap | Two complete blocks across a period seam, values `[(2,1),(3,0),(3,1)]`, 130 actual word calls, nine fresh records, seven retained overlap records, seven padding records and a clipped final block | +| ownership | Nineteen rejection/ownership controls, successful consumer reuse after failures and exact producer re-preparation | + +The displayed output pairs are integer p-grid numerators. The explicit +fixtures use reduced radii and work grids; I do not apply the production tail +theorem to them. For the separate periodic reference I sum directed Gaussian +intervals through distance `8alpha`, round interval endpoints outward to +96 bits and add `2^-256` for the omitted tails. The tail estimate is +`(1/alpha+1/(16pi))*exp(-64pi)<2^-256`. I bound the maximum component error +by two target ulps. [PROOF.md](PROOF.md) gives the separate complex-norm +production argument. + +At `p=4800`, the production checker derives `L=2^98`, radius280, w72, +P157286400, sigma19660800 and halo452984832000000. It checks tile indices +0,1,8,15. Tile8 has one fewer source record than the other three. The final +short block has length144. I inspect these encoded descriptors without +allocating `2^P`, a coefficient catalogue or a production-sized payload. + +The actual finite consumers use the fixed 4096-bit crossover of their +inherited FftWord owner, so their recorded calls take its base/Karatsuba path. +PR3787 supplies separate actual FFT-path evidence and the retained all-size +word argument. I do not describe these expansion fixtures as new large-word +FFT execution. + +## Stable payloads and preparation change + +I compared exact page and output bytes before and after reusing the even +kernel pairs `k_d=k_-d`. I reproduced these hashes: + +| Fixture | Page SHA256 | Output SHA256 | +|---|---|---| +| block p4 | `ea522ece071a20e7e4f3f07156aae149c425f0b642d04357eae7df51a563bcfe` | `17b837bf4383456d726bdee2c757903686c77451951319d6cdbcb0f52826e95e` | +| block p8 | `e7f9e13638be5ab74f68c2b33e2acb094c258bba998447277ff848038cd00b53` | `ba3a3fc42249d623b680b171fe61696ded2dab685d8e3b67ef891ccc1475a3e2` | +| overlap tile | `70d33b2d88b08d0f18af9eed039409c038b7f2e6ee97b2287ee4bcf418de24fc` | `cc050d5d43d1ba27893fe5658568da2157f278c545e550d08610ca4849b682fb` | +| wrapped tile | `645041b9abd04e9cba1c3ef5a3fa360d33c39392f6a6f26bb968db37a0b39b9b` | `d431d73e30d0270d3ce59de62f49169922ed409ff4c0dad8cfccd810f69b1d3c` | +| ownership tile | `1b91aa71588fe0296a04c1f0061cd0060c6b7f22135fa1247c9a3aee4b836247` | `bb2a4f35db7dd2f37c30c52d6f722af8eca987b1d549b3cb6d55cbf1332a0fa0` | + +Reflection reduced the wrapped tile from32 to25 exponential evaluations and +from647 to495 Taylor steps. The block pair takes6 and9 exponentials; the +overlap tile takes15. I retained the complete fixtures and negative controls. +Local reflected-kernel debug runs took69.996 seconds for the block pair, +57.286 for overlap and142.403 for the wrapped tile. The complete ownership +run took38.950 seconds. These timings describe this environment, rather than +a portable runtime promise. The official publication receipt records the +final bounded replay and hosted outcome. + +## Failures that led to corrections + +I retained the following distinctions when preparing the final checks: + +- The initial geometry checker mishandled marker initialization between + cases. I corrected the checker and reran it; I did not change the geometry + to accommodate those failures. +- The zero-scale ownership fixture did not distinguish early extraction + from correct extraction. I restricted that negative control to the two + nonzero-scale block cases and both complete tile cases, where it changes + the result. The ownership fixture still checks its exact output. +- A short imaginary source record exposed a real scratch-lifetime bug in + the tile validator. I nested the real component's ownership scope before + reading the imaginary component, then reran the nineteen failure/reuse + controls. The corrected consumer restores scratch and borrowed ports. + +I count completed runs only. A prior source's timeout or successful receipt +does not count as execution of this package. + +## Source and publication boundary + +The loader authenticates the compressed core export SHA256 +`5aaa84523bfa22eb817023d3130e8b8879153e9207e8f836cc6c61435fd051b8`, +its uncompressed SHA256 +`2d97e093634a4119d9990dcc5d175917df3eabda5e8d1c204f640381e92ab4a2`, +and receipt SHA256 +`d387e4d0a1a43b9b24104f76f8992603d4265bbe447601bea476d492e92e5340`. +It also checks the full inventory, member hashes and source head. The PR371 +correction archive SHA256 is +`31a553fc4aa4443e60f285633ceb309edc84ede4f64ef7c2646f8806e61da1e8`. + +The canonical [Discussion3664](https://github.com/CrocSwap/integer-mult-bounds/discussions/3664) +records the publication head, official export/receipt hashes, unchanged-byte +replay, exact CI tree and final handoff. Keeping those transport receipts +there avoids a self-referential archive hash in this file. Source integrity +does not establish mathematical acceptance. + +I leave the complete tensor/layout/CRT/integer join, production-sized +execution and independent review open. I retain the analytic, numerical +word-recovery and supplier premises from the pinned predecessors. I make no +new selected-exponent claim. diff --git a/research/physical-gaussian-expansion/PROOF.md b/research/physical-gaussian-expansion/PROOF.md new file mode 100644 index 0000000000..89642f4740 --- /dev/null +++ b/research/physical-gaussian-expansion/PROOF.md @@ -0,0 +1,262 @@ +# Complete physical Gaussian expansion tile + +Jose Manuel (`jmbarrancoml`), with substantial OpenAI Codex assistance. +I give a scoped implementation argument under the retained word, interval and +analytic source premises. The finite checks support the listed interfaces; +they do not establish a complete integer-multiplication theorem or independent +human review. + +## 1. Quantified boundary and sources + +I use the original global integers `s,t`, with `1100`, +integer `2<=alpha100`, `b²/u<5p/4`, `b=22`, `|h|<=2b+1`, `0<=z=2sigma+p+ceilLog2(32(n+1))+12`. The production values satisfy it with +room to spare. Thus internal complex error is below `2^(-p-12)`. + +Every omitted source index lies beyond distance `b` from the target center. +A monotone integral estimate bounds the two normalized tails by + +```text +(1/alpha + alpha/(2*pi*b))*exp(-pi*b²/u) + < exp(-pi*p) < 2^(-3p) < 2^(-p-12). +``` + +I charge final component truncations separately, with complex error below +`sqrt(2)*2^-p`. For the complete tile, use the maximum over its disjoint +target blocks, rather than summing over the number of blocks. The result is + +```text +||returned-Sx/2||_infinity + < (sqrt(2)+2^-11)*2^-p < p²*2^-p. +``` + +The shifted Gaussian sum is at most `1+alpha`, by its maximum plus its +integral. Hence the exact S/2 operator norm is at most +`1/2+1/(2alpha)<=3/4`. Before final truncation, the computed norm is below +`3/4+2^(-p-11)<1`; component truncation toward zero cannot increase it. +The consumer also checks its returned disk-grid values. + +Reduced-radius fixtures exercise the complete implementation and rounding +order. They do not satisfy this production tail proof. Their independent +reference reports bound the **maximum component error**, not complex modulus. + +## 4. Production halo and complete stream + +The production wrapper retains the paired compression domain: + +```text +A=4096p³, L=2^ceilLog2(p^8), t dyadic, +u>=ceilLog2(8p), u*d<=p, +1/(4d)<(t-s)/s<1/(2d-1), +w=ceil(sqrt(ceil(16p/u)))+2, d>4w+2, +L>=64dA, s>2(L+2A+2w). +``` + +For tile `k`, I use target core `[kL,(k+1)L)` and target halo +`[kL-A,(k+1)L+A)`. I invoke the retained `ProductionPlan` on that halo and +the source core `[floor(kappa*kL),floor(kappa*(k+1)L))`. For target endpoint +`v`, its source-window endpoint is + +```text +ceil(s*(2v-1)/(2t)), +``` + +because `q_j=floor(t*j/s+1/2)` belongs to an integer target interval exactly +between those two ceilings. I retain the principal complement, locality +margin and paired precision checks. Both positive and negative lifts use +this same endpoint convention. + +The entire expansion stencil is +`[floor(kappa*kL)-b-1, ceil(kappa*((k+1)L-1))+b+2)`. +Each side has clearance greater than `kappa*(A+1/2)-b-3>A/2-b-3>0` +inside the source halo. The ratio and `L>=64dA` bounds imply source count +at most L, including endpoint discrepancies. The implementation checks +that inequality and pads the supplied source halo to capacity L. + +For consecutive short blocks, their source left and right endpoints are +nondecreasing. Since `kappa<1` and each target step is at most b, the next +source window overlaps the previous one. The consumer validates the entire +capacity once, then scans it once in increasing order. It copies the first +window; for each next block it copies the retained suffix from one short +bank into the other and appends only newly reached source records. Equal +left endpoints and equal right endpoints are legal. It skips the remaining +source tail and padding once. It checks the complete ordered page catalogue, +common axis/precision, target positions and clipped final block. + +## 5. Complete cost and ownership + +For a block, count two input norm products, two A products per source record, +two B products, two c products and two output norm products per target record, +plus eight packed products for the two signed convolutions. The exact +FftWord call count is `4n+6m+8`. The tile adds `2*input_capacity` products +for its initial validation. These calls use the actual shared word owner. + +With the retained `O(q log²(q+2))` q-bit word contract, the eight packed +products have q=`O(bp)` and cost `O(bp log²(p+2))`; the scalar products, +packing, cropping, carries and head restoration fit the same block bound. +Across all blocks, `sum n_i=O(L+b)`. Copying overlap banks and streaming +coefficient pages also costs `O((L+b)p)`. The complete online bounds are + +```text +time O((input_capacity+L+b)*p*log²(p+2)), +space O((L+b)*p), +prepared catalogue O((L+b)*p). +``` + +For production capacity L, these give `O(Lp log²p)` time and `O(Lp)` +space. The output itself has L records. Metadata operations and fixed +scratch restoration fit these bounds. A caller must retain the complete +catalogue and reuse it across incoming payload lines. + +The producer uses exact rational metadata, directed exponentials and +sequential page copies. For fixed global parameters it performs a +polynomial number of operations on polynomial-length values in t,p. +Reflection scans add at most `O(b²p)` work per block. I charge total setup +time and scratch as **polynomial(t,p)**, separately from the per-line bill. +Only an enclosing tensor schedule with the required spectator stock may +absorb that term. This package does not implement that schedule. + +Under the retained bulk scaling `d=p^.997`, the complete expansion +arithmetic contributes `p^(.997+delta)` per data volume for each fixed +`delta>0`. At `delta=.0001`, `.9971<.99827619697`, the authenticated full +Fourier power. This comparison covers expansion arithmetic only. A full +caller must also pay its global permutations, short mixed-radix moves, +periodic halo copies, page access and compression. + +Each public owner reserves a fixed number of physical tapes before a call. +Inside a call it seals physical constructors, borrows canonical parameters, +suspends its scope around actual child owners and consumes protected child +returns. It restores borrowed values/heads and empties owned scratch on +success or failure. For paired reads, it opens the first value's ownership +scope before reading the second. The short-imaginary-input control exercises +that failure boundary. The checker disables decoded host dimensions during +the complete calls; host fractions and arrays belong only to reference tests. + +## 6. Remaining joins + +I leave the selected frontier unchanged. A complete transform still needs +the physical global coordinate gathers/scatters, short mixed-radix moves, +periodic halo routing, catalogue reuse and source-core assembly. A gather +that preserves canonical coarse-field order can move Theta(p) coordinates; +the caller must account for its full support, rather than pricing only the +selected low bits. + +PR3763 must receive the complete natural target tensor: D source lengths, +D-1 flags and one final contiguous axis. Microbox coarse indices cannot +become independent chirp row batches. Count the actual D scalar axes in +`gamma=sum_i(2u_i+ceilLog2(32p)+1)`; rad20 uses D=d-1. Keep fine normalization +within each endpoint and apply the final global restoration once. The +retained outer norm at most 1/2 leaves room for expansion error before the +chirp's .6 bound. The CRT entrance/return and both final integer-recovery S +factors remain separate obligations. Reviewers must also assess the +inherited numerical recovery, supplier and source-correspondence premises. diff --git a/research/physical-gaussian-expansion/README.md b/research/physical-gaussian-expansion/README.md new file mode 100644 index 0000000000..f11fe3dc45 --- /dev/null +++ b/research/physical-gaussian-expansion/README.md @@ -0,0 +1,121 @@ +# Complete physical Gaussian expansion tile + +Jose Manuel (`jmbarrancoml`), with substantial OpenAI Codex assistance. + +I implement the source-bound normalized Gaussian expansion `S/2` for a +complete target tile. The producer generates directed real coefficients and +ordered factor pages. The consumer executes the actual signed convolution, +keeps the fine grid through normalization, restores scale during final-grid +extraction, and streams overlapping source windows through two short banks. + +This closes the coefficient-generation and complete normalized-output +handoff in [Hanson Sun's PR1243](https://github.com/CrocSwap/integer-mult-bounds/pull/1243). +I retain [Artus Krohn-Grimberghe's PR371](https://github.com/CrocSwap/integer-mult-bounds/pull/371) +final-grid correction and the verified physical cohort from +[PR3851](https://github.com/CrocSwap/integer-mult-bounds/pull/3851). +The complete tensor/CRT caller remains open. I leave the selected frontier +unchanged and claim no complete multiplication theorem. + +For integer `p>100`, integer `2<=alpha=1`; the older combined +S/E statement has a narrower domain. [PROOF.md](PROOF.md) gives the exact +assumptions, crop identity, tail bound, arithmetic errors and complete bill. + +You can use these public interfaces inside `expansion_sources.opened()`: + +```python +producer = TilePreparation() +page = producer.run_production(s, t, d, alpha, p, tile_index) +result = Tile().run(page, source_halo) +``` + +Parameters are canonical borrowed `Count` values. Supply complex p-grid +records in increasing unwrapped global source order, then zero-pad through +the declared capacity. Consume each protected return before reusing its +owner. The explicit-grid methods support bounded fixtures; they do not +silently inherit the production tail theorem. The consumer requires the +matching producer's Gaussian pages, not merely arbitrary valid framing. + +I retain these implementation boundaries: + +| Files | Responsibility | +|---|---| +| `expansion_sources.py` | Authenticate the complete prior export, correction archive and common physical cohort | +| `expansion_geometry.py`, `expansion_records.py` | Signed global origins, short stencils, exact factor identity and fine-grid guards | +| `expansion_prepare.py`, `expansion_block.py` | Actual real-coefficient producer and complete normalized block | +| `expansion_tile_geometry.py` | Paired microbox halo and production admission | +| `expansion_tile_prepare.py`, `expansion_tile.py` | Complete ordered catalogue and monotone overlapping-window consumer | + +I execute two nonzero normalized blocks and two complete tiles. The tiles +exercise an unchanged source start, a periodic seam, padding and a clipped +last block. I compare them with separate exact staged arithmetic and directed +periodic Gaussian sums, and check the complete word invoices. I also test +nineteen failure/ownership controls and positive reuse. Production checks +exercise encoded metadata at four tiles, including a one-record ceiling +transition; they do not allocate the production-sized P-grid payload. +[EVIDENCE.md](EVIDENCE.md) records the scope, hashes and earlier failures. + +Run the eight registered checks with Python 3.12 and the standard library: + +```sh +python3.12 -B scripts/research_package.py --package research/physical-gaussian-expansion --check-only +python3.12 -B scripts/research_package.py --package research/physical-gaussian-expansion +``` + +The manifest allows at most 1180 seconds in total and 450 seconds for its +longest command. For a committed export, follow +[the repository exporter](../../docs/research/package-export.md); the +verifiers work outside Git with the declared inputs. I run the official +launcher to isolate Python bytecode, even with `-I -O -Werror -B`. + +## Open next questions + +- Can a contributor implement the full tensor caller with canonical global + coordinate gathers, short mixed-radix moves, periodic halos and page reuse? + Start with one complete tensor and restore its natural target order before + invoking PR3763 once. Full-support routing and the physical field-move + primitive remain obligations. +- Can a contributor join the physical CRT entrance and return to the complete + Gaussian maps? Start with one interval and check both endpoints and both S + factors in integer recovery. The complete Gaussian/CRT/integer join remains + outside this package. +- Can a reviewer validate the quantified source and precision correspondence? + Inspect the crop, final extraction, small-u-theta domain and inherited word + recovery against their exact pins. These are AI-assisted author checks; + independent human or formal-kernel review remains open. + +## Sources and credit + +I pin canonical main `20ecf2891d80b5acf18d980b2164fb44ea7faa90`, priority +revision 6, in [TASK.json](TASK.json) and [SOURCES.json](SOURCES.json). +I retain the pre-implementation decision and status in +[Discussion3664](https://github.com/CrocSwap/integer-mult-bounds/discussions/3664). +The whole PR3851 export has head `7ab67b4a84d28cd562b9711aa315ca5b47822b37`; +the PR371 correction archive has head `988ef67986c313445936fc7103908b6163b9d47e`. +The loader authenticates their bytes, inventory and receipts before import. + +I credit eumemic and icekylinx for the corrected real-chirp work; Hanson Sun +for the packed-correlation audit and explicit handoff; Artus +Krohn-Grimberghe for final-grid accounting; Chafik Boukhalfa for the Gaussian, +chirp and locality predecessors; RaD/hipotures for the retained blocked, +microbox and bulk-transfer reports; Mason Hale for the numerical reference; +Dugongue for the word-ABI and inverse precedents; Harvey and van der Hoeven +for the classical multiplication framework; and Douglas Colkitt and the +OpenAI manuscript contributors for the repository sources and research +framework. The complete nested exports retain their additional authors, +licenses and notices. I do not present their work as a new invention here. + +OpenAI Codex assisted with mathematics, implementation, source-analysis +agents, tests and prose. I made no independent-human-review claim. This +contribution follows the repository's Apache-2.0 license; retained sources +keep their source-specific notices. diff --git a/research/physical-gaussian-expansion/SOURCES.json b/research/physical-gaussian-expansion/SOURCES.json new file mode 100644 index 0000000000..811f919a4c --- /dev/null +++ b/research/physical-gaussian-expansion/SOURCES.json @@ -0,0 +1,22 @@ +{ + "schema_version": 1, + "base_commit": "20ecf2891d80b5acf18d980b2164fb44ea7faa90", + "files": [ + {"path":"docs/research/agent-task-catalog.json","usage":"provenance","git_blob":"67f085acd86acd7d4e03922d9cf4b26bca66fed5","sha256":"47b42433f5a85809a3aadd29c2c339da2c6ed728e4f7762df6b88bccb7874cf0"}, + {"path":"certificates/selected-result.json","usage":"provenance","git_blob":"b1fea0c94436b236cc2b7d41a7be032811b94c81","sha256":"bb2eaa58647b353f88d0aec2eafb9d30df1be53589231b23e12ceef3e4bad129"}, + {"path":"research/physical-gaussian-expansion/inputs/core-export.tar.gz","usage":"input","pin":"head:e463a6530298e2b12a2f1d87367be20d3a6c3e60:research/physical-gaussian-expansion/inputs/core-export.tar.gz","git_blob":"9719f1e08e1e9d90a5e5cab7b9cc7875db5cb340","sha256":"5aaa84523bfa22eb817023d3130e8b8879153e9207e8f836cc6c61435fd051b8"}, + {"path":"research/physical-gaussian-expansion/inputs/rounding-source.tar","usage":"input","pin":"head:e463a6530298e2b12a2f1d87367be20d3a6c3e60:research/physical-gaussian-expansion/inputs/rounding-source.tar","git_blob":"72db8f8da3bbfd8f184fa87e3d336cb946414ad9","sha256":"31a553fc4aa4443e60f285633ceb309edc84ede4f64ef7c2646f8806e61da1e8"}, + {"path":"references/semantic-bulk/rad20/code/downstream_blocked_gaussian.py","usage":"input","git_blob":"4e1066f40ab6fb8c40a397f197749b6c41118429","sha256":"2ca6a460a7a9718eaf50f96c1602a1f21c65e2254acb58d3a1693c7e38f1ef95"}, + {"path":"references/semantic-bulk/rad20/reports/downstream-microbox-resampling.md","usage":"input","git_blob":"27f962de34b6943aee33cc0d40b8c7883ccf56f2","sha256":"9598c236eed7e44ac5eb41cc20f7c23249829e195247d15e2857a34e14222c49"}, + {"path":"references/semantic-bulk/rad20/reports/bulk-resampling-tape-transfer.md","usage":"input","git_blob":"7141f84808633e22760e6ca760302ec153f3770a","sha256":"cda26d5aca04b9c824eb6c4f50bee564fb7fa908b09a72e76852df5e0b9d2622"}, + {"path":"references/semantic-bulk/rad20/reports/review-bulk-resampling.md","usage":"input","git_blob":"9731a03fdb7a1cecaddab737545befa238cb9ad3","sha256":"12aefaa9a0eee367fb71e8f910ed5f0094e1203d2caa29836b0c37887c446c29"}, + {"path":"references/semantic-bulk/rad20/LICENSE","usage":"input","git_blob":"0e259d42c996742e9e3cba14c677129b2c1b6311","sha256":"a2010f343487d3f7618affe54f789f5487602331c0a8d03f49e9a7c547cf0499"}, + {"path":"references/semantic-bulk/rad20/NOTICE","usage":"input","git_blob":"544634f5f3001c50330413f22dc0070e07c46738","sha256":"82fc72aff49769ae82a8ce886846d9bd0601827a09c955042d26a61177d14369"}, + {"path":"references/semantic-bulk/rad20/SOURCE.json","usage":"input","git_blob":"f19aa02a1d430ab30133727ddee5f6686b49ac37","sha256":"63477ee5103e4c1d0d1665bf2d1cfb91af5bcdef59cac6e46f4d6a1579b8eb2e"}, + {"path":"upstream/build/sections/07-resampling.tex","usage":"input","git_blob":"256243bf9416baefb82bf9203f80215b9ebe28f2","sha256":"8a6b91670f9163afb54892489fa2210c77fe1b058f0325ec7d18b44b94be6e2c"}, + {"path":"upstream/build/sections/08-assembly.tex","usage":"input","git_blob":"337810a1d9a58a2ecf14123420923c039ba322a6","sha256":"763d7945b1e1ec4ebae9b799bcefcc45fa9bae6e2ffad799868a9c3d513dbd4a"}, + {"path":"notes/chirp-scaling-correction.tex","usage":"input","git_blob":"a7917c55348fb4da5b095c227839c055349e6daa","sha256":"ed1eb5a225dbcf3189cf8236eed556cbbc24b0ea994eb36ee35fc2442b228e0c"}, + {"path":"research/chirp-packed-correlation-audit/README.md","usage":"input","git_blob":"2894f84c3076536858d015018bb060dcdd80de07","sha256":"4106b2a93625498b3eacec0b7e4a5a63215269f1c3a0aad152c0bdda071a1dc0"}, + {"path":"research/chirp-packed-correlation-audit/TASK.json","usage":"input","git_blob":"3e7bdfa3daa581912c5f6e63373035a505991564","sha256":"4f0436519041977cd76a644840cabcce06cfcaf91afece500a9ed4d184aa2256"} + ] +} diff --git a/research/physical-gaussian-expansion/TASK.json b/research/physical-gaussian-expansion/TASK.json new file mode 100644 index 0000000000..3e6297ff8a --- /dev/null +++ b/research/physical-gaussian-expansion/TASK.json @@ -0,0 +1,41 @@ +{ + "schema_version": 1, + "task_id": "uniform-transfer-audit", + "priority_revision": 6, + "base_commit": "20ecf2891d80b5acf18d980b2164fb44ea7faa90", + "original_card": { + "id": "uniform-transfer-audit", + "question": "Can one exact all-size/compiler, precision, hosting or restoration implication be proved or falsified?", + "starting_points": ["certificates/selected-result.json", "docs/audit.md"], + "methods": ["Prove one missing implication with explicit quantified hypotheses.", "Find an exact counterexample at the claimed interface boundary."], + "acceptance_gate": "Identify the proposition and what would change downstream; do not equate a finite replay with a universal proof.", + "first_step": { + "command": ["python3", "-B", "scripts/research_baseline.py"], + "expected": "PASS_RESEARCH_BASELINE_BINDINGS_ONLY", + "inspect": "Read the selected review scope and choose one exact all-size or precision implication.", + "scope": "Positive control / source-bound orientation only; this is not a new construction or full cold replay." + }, + "first_experiment": "Write its quantified inputs and conclusion; find a missing hypothesis or construct a minimal exact adversarial instance.", + "bottleneck_link": "Finite supplier improvements are useful only under the retained transfer obligations.", + "allowed_variations": ["Choose one tape, precision, hosting or restoration lemma."], + "scope_change_triggers": ["Replacing a quantified obligation with another finite smoke test."], + "overlap_terms": ["uniform transfer", "all-size", "precision"] + }, + "scope_status": "refined", + "scope_change_reason": "I take the source-bound S/2 expansion handoff after PR3851 and the explicit missing tasks in PR1243. Before implementation I recorded the matching producer, complete tile consumer, real-chirp identity, fine-grid correction, full invoice and stop rules in Discussion3664.", + "actual_scope": "One complete physical normalized Gaussian expansion tile, including directed coefficient generation, ordered pages, actual signed convolution, normalization, final-grid extraction, periodic source-window streaming and fixed-tape restoration.", + "bottleneck_link": "The retained bulk transfer required a complete expansion implementation with affordable real factors and final-grid accounting; a supplied-factor block or a dense per-line Gaussian sum would leave that interface open.", + "budget_minutes": 60, + "stop_rule": "Complete the producer and whole tile or report the exact obstruction. Stop at mismatched global geometry, alpha/u drift, premature extraction, unpaid factor access, growing real chirps, whole-tile rescans per block or failed restoration. The donor requested continued research beyond the default session; registered checks retain the 600/1200-second limits.", + "outcome": "completed", + "evidence": "I implemented complete encoded block and tile production/consumption, derived the quantified production error and bill, and checked nonzero blocks, repeated source starts, a periodic seam, zero padding and a clipped final block against separate staged and periodic models. Ownership controls cover malformed public boundaries, protected outputs, sealed allocation and positive reuse. Encoded production checks remain metadata-only. EVIDENCE.md and the public receipt distinguish local checks, official export replay and hosted verification.", + "author": "Jose Manuel (jmbarrancoml)", + "ai_assistance": "Substantial OpenAI Codex assistance in mathematics, implementation, source-analysis agents, verification and prose. No independent human or formal-kernel review claimed.", + "canonical_discussion": "https://github.com/CrocSwap/integer-mult-bounds/discussions/3664", + "overlap": "I started from canonical main20ecf2891d80b5acf18d980b2164fb44ea7faa90, revision6. The broad search through2026-10-18 23:59:51UTC covered577open PRs,2964closed/merged PRs,2532Discussion records and4126local records with complete pagination and zero warnings. I inspected17leads and found no competing complete expansion in that coverage. At2026-10-19 00:55:29UTC I refreshed main with no guidance delta and checked the dependency heads, this thread,35hub comments with replies and100recent all-state PRs covering all3updates since the broad search. Dugongue's PR3872 supplies a complementary rational-operation sliding inverse and leaves physical/bit/tape integration open. The separate expired precision audit had no published outcome. The public thread retains subsequent checkpoints and the publication receipt.", + "open_questions": [ + {"question":"Can the complete physical tensor caller reuse these pages with its full paid movement schedule?","why_it_matters":"One expansion tile does not execute global coordinate gathers, short mixed-radix moves, halos or the complete chirp/compression join.","first_step":"Implement one complete tensor, pay full-support canonical gathers and restore natural target order before PR3763; count D actual axes in gamma.","blocked_by":"Physical mixed-radix moves, halo/catalogue routing and source-core assembly remain separate."}, + {"question":"Can the physical CRT entrance and return preserve the valid interval through the complete maps?","why_it_matters":"Coordinate movement alone does not implement arithmetic field conversion or integer recovery.","first_step":"Implement one source/target interval conversion and check both endpoints and both final S factors.","blocked_by":"The complete Gaussian/CRT/integer join remains open."}, + {"question":"Can a reviewer validate the complete numerical and source correspondence?","why_it_matters":"Finite execution does not prove the universal Gaussian tail, word-recovery or inherited transfer premises.","first_step":"Inspect the exact crop, small-u-theta bounds, fine-grid normalization and final extraction against the pinned sources.","blocked_by":"Independent human or formal-kernel review remains open."} + ] +} diff --git a/research/physical-gaussian-expansion/VERIFY.json b/research/physical-gaussian-expansion/VERIFY.json new file mode 100644 index 0000000000..eebfea8d6d --- /dev/null +++ b/research/physical-gaussian-expansion/VERIFY.json @@ -0,0 +1,32 @@ +{ + "schema_version": 1, + "status": "ready", + "python": "3.12", + "requirements": [], + "inputs": [ + "research/physical-gaussian-expansion/inputs/core-export.tar.gz", + "research/physical-gaussian-expansion/inputs/rounding-source.tar", + "references/semantic-bulk/rad20/code/downstream_blocked_gaussian.py", + "references/semantic-bulk/rad20/reports/downstream-microbox-resampling.md", + "references/semantic-bulk/rad20/reports/bulk-resampling-tape-transfer.md", + "references/semantic-bulk/rad20/reports/review-bulk-resampling.md", + "references/semantic-bulk/rad20/LICENSE", + "references/semantic-bulk/rad20/NOTICE", + "references/semantic-bulk/rad20/SOURCE.json", + "upstream/build/sections/07-resampling.tex", + "upstream/build/sections/08-assembly.tex", + "notes/chirp-scaling-correction.tex", + "research/chirp-packed-correlation-audit/README.md", + "research/chirp-packed-correlation-audit/TASK.json" + ], + "commands": [ + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","sources"],"timeout_seconds":15,"checks":"Complete source exports, live analytic and license pins, shared physical types, actual acyclic owners and four source-binding controls."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","budget"],"timeout_seconds":15,"checks":"Exact stencil, factor-growth, precision, tail and complete arithmetic-budget controls bound to the retained full Fourier power."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","geometry"],"timeout_seconds":20,"checks":"Signed origins, exact factor/crop identities, low u*theta, equal right endpoints and encoded production short-block metadata."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","production"],"timeout_seconds":30,"checks":"Encoded complete production tile metadata at both seams and an interior ceiling transition, paired compression halo and ten domain controls; no production-sized payload."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","block"],"timeout_seconds":230,"checks":"Actual coefficient producer and two complete nonzero normalized blocks against separate staged and periodic Gaussian models, premature-extraction controls and complete word invoices."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","tile_overlap"],"timeout_seconds":240,"checks":"Whole producer/tile join with an unchanged source start, retained overlap, exact outputs, separate periodic reference and complete stream/word invoices."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","tile_wrap"],"timeout_seconds":450,"checks":"Whole periodic-crossing tile, clipped final block, source padding, actual convolution and normalization, exact outputs and full stream/word invoices."}, + {"argv":["{python}","-I","-O","-Werror","-B","{package}/verify.py","--suite","ownership"],"timeout_seconds":180,"checks":"Complete public producer/consumer ownership with nineteen malformed-boundary, protected-output and allocation controls, failed-call restoration and deterministic positive reuse."} + ] +} diff --git a/research/physical-gaussian-expansion/expansion_block.py b/research/physical-gaussian-expansion/expansion_block.py new file mode 100644 index 0000000000..be35473b70 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_block.py @@ -0,0 +1,162 @@ +"""Complete normalized real-Gaussian block with the actual word convolution.""" +from expansion_sources import active +active() +from cardinal_counts import Count +from encoded import compare, subtract +from expansion_records import BlockHeader +from fft_arithmetic import Arithmetic +from fft_word_engine import FftWord +from gaussian_convolution import Convolution +from gaussian_page import skip_bits +from records import Tape +from ring_words import Signed, read_signed, scaled, truncate, write_signed +from router_metadata import indices, owned +from scratch_pool import ScratchPool +from session import (OuterScope, borrowed_tape, erase_marker, initialize_marker, + sealed_constructors) +from words import BLANK, need + + +class Block: + def __init__(self): + self.word = FftWord(base_bits=4096) + self.convolution = Convolution(self.word) + self.pool = ScratchPool(256) + self.scope = OuterScope(self.pool) + self.arithmetic = Arithmetic(self.word, self.scope) + self.output = Tape() + self.closed = True + self.calls = 0 + for tape in (*self.pool.tapes, self.output): + erase_marker(tape) + + def idle(self): + self.convolution.idle() + need(self.convolution.closed and not self.convolution.pool.live + and self.convolution.output.read() == BLANK, "consumed expansion convolution return") + + def norm(self, real, imaginary, header): + with owned(self.arithmetic.multiply(real, real, header.input_width)) as (first,): + with owned(self.arithmetic.multiply(imaginary, imaginary, header.input_width)) as (second,): + first.magnitude.add(second.magnitude) + need(compare(first.magnitude, header.norm_bound) <= 0, "expansion disk-grid packet") + + def factors(self, page, target, length, header): + for _ in indices(length): + with owned(read_signed(page, header.width)) as (value,): + need(not value.negative and compare(value.magnitude, header.unit) <= 0, + "nonnegative normalized expansion coefficient") + write_signed(value, target, header.width) + target.rewind() + + def scalar(self, source, a, kernel, b, normalization, destination, header, component): + weighted = returned = None + source.rewind() + a.rewind() + b.rewind() + try: + weighted = self.pool.acquire() + for _ in indices(header.n): + with owned(read_signed(source, header.input_width)) as (real,): + with owned(read_signed(source, header.input_width)) as (imaginary,): + if component == 0: + self.norm(real, imaginary, header) + with owned(scaled(real if component == 0 else imaginary, header.promote)) as (value,): + with owned(read_signed(a, header.width)) as (coefficient,): + with owned(self.arithmetic.multiply(value, coefficient, header.width)) as (product,): + with owned(truncate(product, header.P)) as (rounded,): + write_signed(rounded, weighted, header.width) + need(source.read() == a.read() == BLANK, "complete expansion source window and input factors") + with owned(subtract(header.nc, header.n), Signed(Count.literal(0))) as (padding, zero): + for _ in indices(padding): + write_signed(zero, weighted, header.width) + weighted.rewind() + with self.scope.suspended(): + returned = self.convolution.run(weighted, kernel, header.nc, header.width, header.P) + with owned(Count.literal(1)) as (one,): + with owned(subtract(header.n, one)) as (prefix,): + with owned(Count.product(prefix, header.width)) as (bits,): + skip_bits(returned, bits) + for _ in indices(header.m): + with owned(read_signed(returned, header.width)) as (value,): + with owned(read_signed(b, header.width)) as (coefficient,): + with owned(self.arithmetic.multiply(value, coefficient, header.width)) as (product,): + with owned(truncate(product, header.P)) as (post_chirp,): + with owned(self.arithmetic.multiply(post_chirp, normalization, header.width)) as (product,): + with owned(truncate(product, header.P)) as (fine_normalized,): + with owned(truncate(fine_normalized, header.extract)) as (final,): + write_signed(final, destination, header.input_width) + need(returned.read() == b.read() == BLANK, "complete cropped expansion output") + destination.rewind() + finally: + if returned is not None: + returned.clear() + erase_marker(returned) + if weighted is not None: + weighted.clear() + for value in (source, a, kernel, b): + value.rewind() + + def run(self, page, source): + page, source = borrowed_tape(page), borrowed_tape(source) + forbidden = (*self.pool.tapes, self.output, *self.convolution.pool.tapes, + self.convolution.output, *self.word.pool.tapes, self.word.output, + *self.word.canonical) + need(page is not source and page not in forbidden and source not in forbidden, + "external expansion block ports") + need(self.closed and not self.pool.live and self.output.read() == BLANK, + "consume prior expansion block output") + self.idle() + for tape in (*self.pool.tapes, self.output): + initialize_marker(tape) + self.closed = False + header = a = kernel = b = normalization = real = imaginary = None + succeeded = False + try: + with sealed_constructors(), self.scope.installed(): + try: + page.rewind() + header = BlockHeader(page) + a, kernel, b, real, imaginary = (self.pool.acquire() for _ in range(5)) + self.factors(page, a, header.n, header) + self.factors(page, kernel, header.nc, header) + self.factors(page, b, header.m, header) + normalization = read_signed(page, header.width) + with owned(Count(header.unit)) as (quarter,): + quarter.divide_small(4) + need(not normalization.negative and not normalization.magnitude.zero() + and compare(normalization.magnitude, quarter) <= 0, + "positive normalized expansion reciprocal alpha") + need(page.read() == BLANK, "complete expansion page extent") + self.scalar(source, a, kernel, b, normalization, real, header, 0) + self.scalar(source, a, kernel, b, normalization, imaginary, header, 1) + for _ in indices(header.m): + with owned(read_signed(real, header.input_width)) as (first,): + with owned(read_signed(imaginary, header.input_width)) as (second,): + self.norm(first, second, header) + write_signed(first, self.output, header.input_width) + write_signed(second, self.output, header.input_width) + need(real.read() == imaginary.read() == BLANK, "complete normalized complex expansion block") + self.idle() + self.calls += 1 + succeeded = True + finally: + for value in (a, kernel, b, real, imaginary): + if value is not None: + value.clear() + if normalization is not None: + normalization.close() + if header is not None: + header.close() + page.rewind() + source.rewind() + self.pool.empty() + finally: + for tape in self.pool.tapes: + erase_marker(tape) + self.output.rewind() + if not succeeded: + self.output.clear() + erase_marker(self.output) + self.closed = True + return self.output diff --git a/research/physical-gaussian-expansion/expansion_check_block.py b/research/physical-gaussian-expansion/expansion_check_block.py new file mode 100644 index 0000000000..670125bf0d --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_block.py @@ -0,0 +1,208 @@ +"""Whole producer/block checks against staged arithmetic and Gaussian intervals.""" +from fractions import Fraction as Q +from hashlib import sha256 +from expansion_sources import require + + +def packet_bytes(tape): + from words import BLANK + tape.rewind() + result = bytearray() + try: + while tape.read() != BLANK: + result.append(tape.read()) + tape.step(1) + return bytes(result) + finally: + tape.rewind() + + +def inputs(start, length, period): + pairs = ((Q(3, 4), Q(-1, 4)), (Q(1, 2), Q(1, 4)), (Q(-1, 8), Q(1, 2)), + (Q(5, 8), Q(-3, 8)), (Q(-1, 4), Q(1, 4)), (Q(3, 8), Q(-1, 8)), + (Q(1, 2), Q(3, 8))) + return [pairs[(j % period) % len(pairs)] for j in range(start, start+length)] + + +def expected_page(page, s, t, alpha, radius, start, length, p, bits, sigma): + from expansion_records import BlockHeader + from gaussian_check_exponentials import oracle + from ring_check_words import decode + from ring_words import read_signed + from router_metadata import owned + from words import BLANK + page.rewind() + with owned(BlockHeader(page)) as (header,): + signed = lambda value: (-1 if value.negative else 1)*decode(value.magnitude) + kappa = Q(s, t) + j0 = (kappa*start).__floor__() + gamma = j0-kappa*start + low = -radius-1 + high = (kappa*(length-1)-gamma).__ceil__()+radius+1 + n, nc = high-low+1, high-low+length + require((decode(header.n), decode(header.m), decode(header.nc), signed(header.source_start), + signed(header.target_start)) == (n, length, nc, j0+low, start), + "producer binds the full global source and target windows") + factors = {} + for kind, count in (("A", n), ("K", nc), ("B", length)): + factors[kind] = [] + for index in range(count): + with owned(read_signed(page, header.width)) as (value,): + actual = Q(signed(value), 1 << bits) + if kind == "A": + h = low+index + exponent = ((1-kappa)*h*h+2*gamma*h+gamma*gamma)/(alpha*alpha) + elif kind == "K": + d = -high+index + exponent = kappa*d*d/(alpha*alpha) + else: + z = index + exponent = ((kappa*kappa-kappa)*z*z-2*gamma*kappa*z)/(alpha*alpha) + lower, upper = oracle(abs(exponent), exponent < 0, 0 if kind == "K" else sigma) + require(0 <= actual <= lower <= upper <= 1 and upper-actual < Q(4, 1 << bits), + "complete normalized factor encloses its original real Gaussian coefficient") + factors[kind].append(actual) + with owned(read_signed(page, header.width)) as (value,): + c = Q(signed(value), 1 << bits) + require(c == Q((1 << bits)//(2*alpha), 1 << bits) and page.read() == BLANK, + "exact encoded reciprocal alpha and complete factor page") + page.rewind() + return j0+low, factors, c + + +def models(data, factors, c, s, t, alpha, source_start, target_start, m, p, bits, sigma, check_premature=True): + from gaussian_check_exponentials import oracle + n = len(data) + rounding = lambda value: Q(int(value*(1 << bits)), 1 << bits) + actual, premature = [], [] + errors = Q(0) + for z in range(m): + pair, wrong = [], [] + for component in (0, 1): + weighted = [rounding(factors["A"][i]*data[i][component]) for i in range(n)] + convolution = rounding(sum((weighted[i]*factors["K"][n-1+z-i] for i in range(n)), Q(0))) + post_chirp = rounding(factors["B"][z]*convolution) + fine = rounding(c*post_chirp) + value = int(fine*(1 << (2*sigma+p))) + pair.append(value) + early = Q(int(post_chirp*(1 << p)), 1 << p) + wrong.append(int(c*early*(1 << (2*sigma+p)))) + lower = upper = Q(0) + for i in range(n): + exponent = (source_start+i-Q(s, t)*(target_start+z))**2/(alpha*alpha) + lo, hi = oracle(exponent) + ends = sorted((data[i][component]*lo, data[i][component]*hi)) + lower += ends[0]/(2*alpha) + upper += ends[1]/(2*alpha) + errors = max(errors, abs(Q(value, 1 << p)-lower), abs(Q(value, 1 << p)-upper)) + actual.append(pair) + premature.append(wrong) + require(errors < Q(2, 1 << p), "whole normalized finite block agrees with direct Gaussian intervals") + if check_premature: + require(actual != premature, "premature p-grid extraction changes the complete normalized result") + return actual, premature, Q((errors*(1 << 96)).__ceil__(), 1 << 96) + + +def periodic_error(actual, s, t, alpha, target_start, p): + """Independent periodic S/2 enclosure, with a rigorous omitted tail.""" + from gaussian_check_exponentials import oracle + require(alpha >= 2, "periodic Gaussian oracle alpha domain") + result = Q(0) + scale = 1 << 96 + tail = Q(1, 1 << 256) + for z, pair in enumerate(actual): + center = Q(s, t)*(target_start+z) + origin = center.__floor__() + lower = [Q(0), Q(0)] + upper = [Q(0), Q(0)] + for j in range(origin-8*alpha, origin+8*alpha+1): + lo, hi = oracle((j-center)**2/(alpha*alpha)) + values = inputs(j, 1, s)[0] + for component, value in enumerate(values): + ends = sorted((value*lo/(2*alpha), value*hi/(2*alpha))) + lower[component] += Q((ends[0]*scale).__floor__(), scale) + upper[component] += Q((ends[1]*scale).__ceil__(), scale) + # Beyond distance8alpha the two normalized tails are bounded by + # (1/alpha+1/(16pi))*exp(-64pi)<2^-256 for alpha>=2. + for component, value in enumerate(pair): + result = max(result, abs(Q(value, 1 << p)-(lower[component]-tail)), + abs(Q(value, 1 << p)-(upper[component]+tail))) + require(result < Q(2, 1 << p), "finite complete output agrees with the periodic Gaussian map") + return Q((result*scale).__ceil__(), scale) + + +def run(): + from cardinal_counts import Count + from check_support import no_decode + from expansion_block import Block + from expansion_prepare import Preparation + from records import Tape + from ring_check_words import body, decode + from ring_words import Signed, read_signed, write_signed + from router_metadata import owned + from session import erase_marker + from words import BLANK + preparation, block = Preparation(), Block() + records = [] + for s, t, alpha, radius, start, length, p, bits, sigma in ( + (7, 8, 2, 1, 0, 1, 4, 32, 4), + (7, 8, 2, 1, -1, 1, 8, 48, 4)): + page = returned = None + source = Tape() + with owned(*(Count.literal(v) for v in (s, t, alpha, radius, abs(start), length, p, bits, sigma))) as args: + source_length, target_length, a, b, origin, m, target_p, P, scale = args + saved = [(bytes(value.tape.cells), value.tape.head) for value in args] + try: + with no_decode(): + page = preparation.run_explicit(source_length, target_length, a, b, + Signed(origin, start < 0), m, target_p, P, scale) + page_bytes = packet_bytes(page) + source_start, factors, c = expected_page(page, s, t, alpha, radius, start, length, p, bits, sigma) + data = inputs(source_start, len(factors["A"]), s) + expected, premature, error = models(data, factors, c, s, t, alpha, source_start, + start, length, p, bits, sigma) + global_error = periodic_error(expected, s, t, alpha, start, p) + with owned(Count.literal(p+2)) as (width,): + for pair in data: + for value in pair: + numerator = int(value*(1 << p)) + with owned(Signed(Count.literal(abs(numerator)), numerator < 0)) as (encoded,): + write_signed(encoded, source, width) + source.rewind() + source_bytes = body(source) + before = dict(block.word.calls) + with no_decode(): + returned = block.run(page, source) + raw = body(returned) + actual = [] + for _ in range(length): + with owned(read_signed(returned, width), read_signed(returned, width)) as (real, imaginary): + actual.append([(-1 if v.negative else 1)*decode(v.magnitude) for v in (real, imaginary)]) + require(actual == expected and returned.read() == BLANK, + "complete produced expansion block equals the exact staged model") + require(any(v for pair in actual for v in pair), "nonzero complete expansion return") + word_calls = {key: block.word.calls[key]-before[key] for key in before} + require(sum(word_calls.values()) == 4*len(data)+6*length+8, + "complete expansion word invoice includes norms and normalization") + require(packet_bytes(page) == page_bytes and body(source) == source_bytes, + "complete expansion preserves borrowed page and source") + require(saved == [(bytes(value.tape.cells), value.tape.head) for value in args], + "complete factor producer preserves borrowed parameters") + require(preparation.closed and not preparation.pool.live and block.closed and not block.pool.live, + "complete expansion owners close their scratch leases") + records.append({"s": s, "t": t, "alpha": alpha, "radius": radius, + "target": [start, length], "source": [source_start, len(data)], + "p": p, "P": bits, "sigma": sigma, "values": actual, + "premature_values": premature, "direct_error": str(error), + "periodic_error": str(global_error), + "page_sha256": sha256(page_bytes).hexdigest(), "core_sha256": sha256(raw).hexdigest(), + "word_calls": word_calls, "preparation_counts": preparation.last_counts}) + finally: + for tape in (returned, page, source): + if tape is not None: + tape.clear() + erase_marker(tape) + block.idle() + return {"blocks": records, "preparation_pool_peak": preparation.pool.peak, + "block_pool_peak": block.pool.peak, + "scope": "complete encoded coefficient producer and normalized finite block; reduced radii do not assert production Gaussian tail bounds"} diff --git a/research/physical-gaussian-expansion/expansion_check_budget.py b/research/physical-gaussian-expansion/expansion_check_budget.py new file mode 100644 index 0000000000..6c664bd727 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_budget.py @@ -0,0 +1,44 @@ +"""Exact controls for the quantified short-block proof and complete tile bill.""" +from fractions import Fraction as Q +import importlib.util +from math import isqrt +from expansion_sources import active, require + + +def run(): + cases = 0 + for p in (*range(101, 513), 1000, 4096, 4800, 16384): + root = isqrt(p)+(isqrt(p)**2 < p) + largest = isqrt(p-1) + for alpha in sorted({2, largest, max(2, largest//2)}): + b, u = alpha*root, alpha*alpha + require(Q(b*b, u) < Q(5*p, 4), "production radius squared bound") + require(3*b+3 < 4*p and 4*b+2 < 5*p, "whole short stencil and prefix convolution bound") + require(Q((2*b+2)**2, u) < 20*p and Q((b+1)**2, u) < 3*p + and Q((2*b+1)**2, u) < 12*p, "real Gaussian exponent magnitude controls") + P, sigma = 32768*p, 4096*p + guard = (32*(3*b+4)-1).bit_length() + require(P >= 2*sigma+p+guard+12, "production fine grid covers normalization before extraction") + require(Q(1, alpha)+Q(alpha, 6*b) < 1 and 3*p > p+12, + "Gaussian tail coefficient and final-grid error separation") + require(2*sigma > 2*80*p and 2*sigma < P-p, + "normalized chirp factors and retained fine bits") + cases += 1 + path = active().ring_package / "ring_check_budget.py" + spec = importlib.util.spec_from_file_location("expansion_inherited_fourier_budget", path) + module = importlib.util.module_from_spec(spec) + spec.loader.exec_module(module) + inherited = module.run() + phi = Q(inherited["fourier_power"]) + epsilon, delta = Q(997, 1000), Q(1, 10000) + expansion = epsilon+delta + require(expansion == Q(9971, 10000) and phi-expansion == Q(117619697, 100000000000), + "complete expansion arithmetic fits the inherited Fourier allowance") + return {"production_guard_cases": cases, + "block_word_calls": "4*n+6*m+8", + "tile_word_calls": "2*input_capacity+sum_i(4*n_i+6*m_i+8)", + "tile_online_time": "O((input_capacity+L+b)*p*log2(p+2)^2)", + "tile_online_space": "O((L+b)*p)", "coefficient_catalogue": "O((L+b)*p)", + "setup": "polynomial(t,p), retained separately from online payload work", + "expansion_arithmetic_power": str(expansion), "fourier_power": str(phi), + "gap": str(phi-expansion), "scope": "finite algebraic controls; quantified claims and inherited premises in PROOF.md"} diff --git a/research/physical-gaussian-expansion/expansion_check_geometry.py b/research/physical-gaussian-expansion/expansion_check_geometry.py new file mode 100644 index 0000000000..d199887cc8 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_geometry.py @@ -0,0 +1,128 @@ +"""Exact global origins, signed stencils, crop indices and production metadata.""" +from fractions import Fraction as Q +from math import isqrt +from expansion_sources import require + + +def run(): + from cardinal_counts import Count + from check_support import no_decode + from expansion_geometry import BlockGeometry, ProductionBlock + from expansion_records import BlockHeader + from ring_check_words import decode + from ring_words import Signed + from router_metadata import owned + from scratch_pool import ScratchPool + from session import OuterScope, erase_marker, initialize_marker, sealed_constructors + from words import BLANK + pool = ScratchPool(256) + scope = OuterScope(pool) + for tape in pool.tapes: + erase_marker(tape) + records = [] + identities = 0 + cases = ((7, 8, 2, 1, 0, 1, 4, 32, 4), + (7, 8, 2, 1, 1, 1, 4, 32, 4), + (7, 8, 2, 2, 7, 2, 4, 40, 5), + (13, 14, 2, 2, -2, 2, 8, 48, 8), + (13, 14, 2, 2, 11, 1, 8, 48, 8), + (3, 4, 2, 2, 2, 2, 4, 32, 3), + (3, 4, 2, 2, 4, 1, 4, 32, 3)) + for s, t, alpha, radius, start, length, p, bits, sigma in cases: + with owned(*(Count.literal(v) for v in (s, t, alpha, radius, abs(start), length, p, bits, sigma))) as args: + source_length, target_length, a, b, origin, m, target_p, P, scale = args + saved = [(bytes(value.tape.cells), value.tape.head) for value in args] + plan = header = page = None + for tape in pool.tapes: + initialize_marker(tape) + with sealed_constructors(), scope.installed(): + try: + with no_decode(): + plan = BlockGeometry(source_length, target_length, a, b, + Signed(origin, start < 0), m, target_p, P, scale) + def signed(value): + return (-1 if value.negative else 1)*decode(value.magnitude) + def rational(value): + return Q(signed(value.numerator), decode(value.denominator)) + kappa = Q(s, t) + j0 = (kappa*start).__floor__() + gamma = j0-kappa*start + low = -radius-1 + high = (kappa*(length-1)-gamma).__ceil__()+radius+1 + n, nc = high-low+1, high-low+length + require(signed(plan.j0) == j0 and rational(plan.gamma) == gamma + and signed(plan.source_start) == j0+low + and (decode(plan.n), decode(plan.nc)) == (n, nc), + "exact lifted expansion source window and convolution length") + arrays = {} + for kind, count in (("A", n), ("K", nc), ("B", length)): + values = [] + for index in range(count): + with owned(Count.literal(index)) as (position,): + with no_decode(): + value = plan.exponent(kind, position) + with owned(value): + values.append(rational(value)) + arrays[kind] = values + for index, h in enumerate(range(low, high+1)): + for z in range(length): + kernel_index = n-1+z-index + require(0 <= kernel_index < nc, "complete prefix crop contains every term") + direct = Q((j0+h-kappa*(start+z))**2, alpha*alpha) + require(arrays["A"][index]+arrays["K"][kernel_index]+arrays["B"][z] == direct, + "exact original real Gaussian exponent after crop") + identities += 1 + page = pool.acquire() + with no_decode(): + BlockHeader.dump(plan, page) + page.rewind() + header = BlockHeader(page) + require(page.read() == BLANK, "complete expansion metadata frame") + header.repeat(page) + require(page.read() == BLANK and decode(header.extract) == bits-p-2*sigma, + "normalization remains before final-grid extraction") + records.append({"s": s, "t": t, "alpha": alpha, "radius": radius, + "target": [start, length], "source": [j0+low, n], + "nc": nc, "p": p, "P": bits, "sigma": sigma, + "u_theta": str(Q(alpha*alpha*(t-s), s))}) + require(saved == [(bytes(value.tape.cells), value.tape.head) for value in args], + "geometry preserves borrowed parameters") + finally: + for value in (header, plan): + if value is not None: + value.close() + if page is not None: + page.clear() + pool.empty() + for tape in pool.tapes: + erase_marker(tape) + require(records[-2]["source"] == [-2, 9] and records[-1]["source"] == [0, 7], + "positive left advance with an equal source right endpoint") + production = [] + for p, alpha in ((101, 2), (127, 3), (16384, 6)): + with owned(Count.literal(p), Count.literal(alpha)) as args: + for tape in pool.tapes: + initialize_marker(tape) + with sealed_constructors(), scope.installed(): + plan = None + try: + with no_decode(): + plan = ProductionBlock(*args) + root = isqrt(p) + root += root*root < p + actual = {name: decode(getattr(plan, name)) for name in ("u", "radius", "P", "sigma")} + require(actual == {"u": alpha*alpha, "radius": alpha*root, + "P": 32768*p, "sigma": 4096*p}, + "encoded production expansion metadata") + production.append({"p": p, "alpha": alpha, **actual}) + finally: + if plan is not None: + plan.close() + pool.empty() + for tape in pool.tapes: + erase_marker(tape) + for tape in pool.tapes: + require(tape.head == 1 and all(value == BLANK for value in tape.cells), "geometry scratch restoration") + return {"blocks": records, "factor_identities": identities, "production": production, + "pool_peak": pool.peak, + "scope": "exact encoded geometry and metadata only; coefficient pages and data consumers remain separate"} diff --git a/research/physical-gaussian-expansion/expansion_check_ownership.py b/research/physical-gaussian-expansion/expansion_check_ownership.py new file mode 100644 index 0000000000..8690c15122 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_ownership.py @@ -0,0 +1,189 @@ +"""Complete tile ownership, malformed boundaries, sealed allocation and reuse.""" +from hashlib import sha256 +from unittest.mock import patch +from expansion_check_block import expected_page, inputs, models, packet_bytes +from expansion_sources import require + + +def run(): + from cardinal_counts import Count + from check_support import no_decode + from expansion_tile import Tile + from expansion_tile_geometry import TileHeader + from expansion_tile_prepare import TilePreparation, copy_cells + import expansion_prepare + from gaussian_rationals import read_count + from records import Tape + from ring_check_words import decode + from ring_words import Signed, read_signed, write_signed + from router_metadata import owned + from session import erase_marker + from words import BLANK, RIGHT + preparation, tile = TilePreparation(), Tile() + + def packet(data): + result = Tape() + for value in data: + result.append(value) + result.rewind() + return result + + def count_bytes(value): + return bytes([(value >> bit) & 1 for bit in range(max(1, value.bit_length()))]+[RIGHT]) + + def position(data, fields): + start = 0 + for _ in range(fields): + start = data.index(RIGHT, start)+1 + return start + + def changed(data, field, value): + start = position(data, field) + end = data.index(RIGHT, start)+1 + return data[:start]+count_bytes(value)+data[end:] + + def clear(tape): + tape.clear() + erase_marker(tape) + + def restored(owner): + require(owner.closed and not owner.pool.live + and all(t.head == 1 and all(v == BLANK for v in t.cells) for t in owner.pool.tapes), + "complete tile fixed scratch restoration") + if hasattr(owner, "parameter_pool"): + require(not owner.parameter_pool.live + and all(t.head == 1 and all(v == BLANK for v in t.cells) for t in owner.parameter_pool.tapes), + "complete tile parameter scratch restoration") + + def rejected(call, message): + try: + with no_decode(): + call() + except ValueError as error: + require(message in str(error), "specific tile rejection: " + message + ": " + str(error)) + else: + raise ValueError("missing tile rejection: " + message) + + source = Tape() + controls = [] + with owned(*(Count.literal(v) for v in (3, 4, 2, 1, 2, 5, 8, 0, 1, 4, 24, 0))) as args: + s, t, alpha, radius, start, length, capacity, target, target_length, p, bits, sigma = args + arguments = (s, t, alpha, radius, Signed(start, True), length, capacity, + Signed(target), target_length, p, bits, sigma) + saved = [(bytes(v.tape.cells), v.tape.head) for v in args] + with no_decode(): + page = preparation.run_explicit(*arguments) + page_bytes = packet_bytes(page) + restored(preparation) + restored(preparation.block) + rejected(lambda: preparation.run_explicit(*arguments), "consume prior expansion tile page") + controls.append("protected_producer_page") + require(packet_bytes(page) == page_bytes, "protected producer page survives rejected reuse") + block_page = Tape() + with owned(TileHeader(page)) as (header,): + with owned(read_count(page)) as (cells,): + copy_cells(page, block_page, cells) + block_page.rewind() + origin, factors, c = expected_page(block_page, 3, 4, 2, 1, 0, 1, 4, 24, 0) + data = inputs(origin, 5, 3) + expected, premature, error = models(data, factors, c, 3, 4, 2, origin, 0, 1, 4, 24, 0, + check_premature=False) + clear(block_page) + page.rewind() + with owned(Count.literal(6)) as (width,): + for pair in data+[(0, 0)]*3: + for value in pair: + integer = int(value*16) + with owned(Signed(Count.literal(abs(integer)), integer < 0)) as (encoded,): + write_signed(encoded, source, width) + source.rewind() + source_bytes = packet_bytes(source) + with no_decode(): + result = tile.run(page, source) + with owned(read_signed(result, width), read_signed(result, width)) as pair: + actual = [[(-1 if v.negative else 1)*decode(v.magnitude) for v in pair]] + require(actual == expected and result.read() == BLANK, "whole ownership tile matches separate staged oracle") + result_bytes = packet_bytes(result) + rejected(lambda: tile.run(page, source), "consume prior expansion tile output") + require(packet_bytes(result) == result_bytes, "protected full tile output survives rejected reuse") + controls.append("protected_tile_output") + clear(result) + child_start = position(page_bytes, 15) + child = page_bytes[child_start:] + factor_start = position(child, 14) + factor_width = 30 + negative = bytes([1]*factor_width) + bad_factor = child[:factor_start]+negative+child[factor_start+factor_width:] + # Header edits keep the nested cell count truthful, so a semantic + # rejection cannot hide behind an unrelated short-page failure. + prefix = page_bytes[:position(page_bytes, 14)] + def child_packet(value): + return prefix+count_bytes(len(value))+value + target_position = position(child, 13) + bad_target = child[:target_position]+bytes([0])+count_bytes(1)+child[position(child, 14):] + source_position = position(child, 12) + bad_origin = child[:source_position]+bytes([1])+count_bytes(3)+child[position(child, 13):] + padding = source_bytes[:60]+bytes([1])+source_bytes[61:] + non_disk = bytes((24 >> i) & 1 for i in range(6))+source_bytes[6:] + variants = ( + ("tile_version", changed(page_bytes, 0, 1), source_bytes, "tile version"), + ("precision_guard", changed(page_bytes, 2, 23), source_bytes, "precision guard"), + ("tile_block_count", changed(page_bytes, 11, 2), source_bytes, "block count"), + ("source_capacity", changed(page_bytes, 8, 9), source_bytes, "padded capacity"), + ("short_page", page_bytes[:-1], source_bytes, "catalogue cell"), + ("extra_page", page_bytes+b"\x00", source_bytes, "streamed expansion tile"), + ("short_input", page_bytes, source_bytes[:-1], "complete"), + ("extra_input", page_bytes, source_bytes+b"\x00", "source extent"), + ("nonzero_padding", page_bytes, padding, "zero expansion source padding"), + ("input_outside_disk", page_bytes, non_disk, "disk-grid input"), + ("child_precision", child_packet(changed(child, 1, 5)), source_bytes, "precision guard"), + ("target_order", child_packet(bad_target), source_bytes, "ordered expansion target"), + ("source_origin", child_packet(bad_origin), source_bytes, "begins in supplied halo"), + ("negative_factor", child_packet(bad_factor), source_bytes, "nonnegative normalized"), + ("zero_normalization", child_packet(child[:-factor_width]+bytes(factor_width)), + source_bytes, "positive normalized expansion reciprocal"), + ) + for name, bad_page, bad_source, message in variants: + first, second = packet(bad_page), packet(bad_source) + try: + rejected(lambda: tile.run(first, second), message) + require(packet_bytes(first) == bad_page and packet_bytes(second) == bad_source, + "failed tile preserves both malformed borrowed ports") + tile.idle() + restored(tile) + restored(tile.block) + require(tile.output.read() == BLANK, "failed tile clears partial output") + controls.append(name) + finally: + clear(first) + clear(second) + def allocate(*args, **kwargs): + Tape() + with patch.object(tile, "copy_records", allocate): + rejected(lambda: tile.run(page, source), "physical Tape allocation inside") + controls.append("dynamic_physical_tape") + with no_decode(): + result = tile.run(page, source) + require(packet_bytes(result) == result_bytes and packet_bytes(page) == page_bytes + and packet_bytes(source) == source_bytes, "positive tile reuse after failed calls") + clear(result) + clear(page) + with patch.object(expansion_prepare, "prepare_page", allocate): + rejected(lambda: preparation.run_explicit(*arguments), "physical Tape allocation inside") + controls.append("producer_dynamic_physical_tape") + restored(preparation) + restored(preparation.block) + require(preparation.output.read() == BLANK, "failed producer clears partial catalogue") + with no_decode(): + page = preparation.run_explicit(*arguments) + require(packet_bytes(page) == page_bytes, "positive producer reuse preserves exact generated coefficients") + require(saved == [(bytes(v.tape.cells), v.tape.head) for v in args], "complete producer preserves parameter ports") + clear(page) + tile.idle() + restored(tile) + restored(preparation) + clear(source) + return {"values": actual, "premature_values": premature, "finite_error": str(error), + "controls": controls, "producer_calls": preparation.calls, "tile_calls": tile.calls, + "page_sha256": sha256(page_bytes).hexdigest(), "core_sha256": sha256(result_bytes).hexdigest(), + "producer_pool_peak": preparation.pool.peak, "tile_pool_peak": tile.pool.peak} diff --git a/research/physical-gaussian-expansion/expansion_check_production.py b/research/physical-gaussian-expansion/expansion_check_production.py new file mode 100644 index 0000000000..29fbc44dd9 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_production.py @@ -0,0 +1,101 @@ +"""Encoded whole-tile production metadata without production-sized coefficients.""" +from fractions import Fraction as Q +from math import isqrt +from expansion_sources import require + + +def run(): + from cardinal_counts import Count + from check_support import no_decode + from expansion_tile_geometry import ProductionTile, stencil_bounds + from ring_check_words import decode + from router_metadata import owned + from scratch_pool import ScratchPool + from session import OuterScope, erase_marker, initialize_marker, sealed_constructors + from words import BLANK + pool = ScratchPool(256) + scope = OuterScope(pool) + for tape in pool.tapes: + erase_marker(tape) + signed = lambda v: (-1 if v.negative else 1)*decode(v.magnitude) + p, alpha, d = 4800, 4, 300 + L = 1 << ((p**8-1).bit_length()) + t = 16*L + s = t-(1 << 92)-1 + halo = 4096*p**3 + records, controls = [], [] + cases = [("first", s, t, d, alpha, p, 0, None), + ("middle", s, t, d, alpha, p, 1, None), + ("rounding_transition", s, t, d, alpha, p, 8, None), + ("last", s, t, d, alpha, p, 15, None), + ("target_precision", s, t, d, alpha, 100, "production expansion precision"), + ("alpha_one", s, t, d, 1, p, "integer alpha"), + ("alpha_square", s, t, d, 70, p, "alpha squared"), + ("concentration", s, t, d, 2, p, "logarithmic Gaussian concentration"), + ("similarity", s, t, 301, alpha, p, "similarity guard"), + ("cell_separation", s, t, 290, alpha, p, "cell/band separation"), + ("nondyadic_target", s, t-1, d, alpha, p, "dyadic target"), + ("incomplete_tile", (1 << 64)-(1 << 54)-1, 1 << 64, d, alpha, p, "complete target tile"), + ("tile_out_of_range", s, t, d, alpha, p, "complete target tile"), + ("global_period", 2*L-(1 << 89)-1, 2*L, d, alpha, p, "global period")] + for case in cases: + if len(case) == 8: + name, S, T, D, a, target_p, index, rejection = case + else: + name, S, T, D, a, target_p, rejection = case + index = 16 if name == "tile_out_of_range" else 0 + with owned(*(Count.literal(v) for v in (S, T, D, a, target_p, index))) as parameters: + saved = [(bytes(v.tape.cells), v.tape.head) for v in parameters] + for tape in pool.tapes: + initialize_marker(tape) + plan = None + try: + with sealed_constructors(), scope.installed(): + try: + with no_decode(): + plan = ProductionTile(*parameters) + require(rejection is None, "missing production domain rejection: " + name) + root = isqrt(p)+(isqrt(p)**2 < p) + expected = {"p": p, "alpha": alpha, "u": alpha*alpha, "d": d, + "radius": alpha*root, "P": 32768*p, "sigma": 4096*p, + "L": L, "halo": halo, "target_length": L, "input_capacity": L, + "blocks": (L+alpha*root-1)//(alpha*root)} + actual = {key: decode(getattr(plan, key)) for key in expected} + require(actual == expected and signed(plan.target_start) == index*L, + "encoded production expansion tile parameters") + # q_j=floor(t*j/s+1/2) lies in [a,b) precisely on + # [ceil(s*(2a-1)/(2t)), ceil(s*(2b-1)/(2t))). + left = -((-s*(2*(index*L-halo)-1))//(2*t)) + right = -((-s*(2*((index+1)*L+halo)-1))//(2*t)) + require(signed(plan.source_start) == left and decode(plan.source_count) == right-left, + "same original global compression window at both period seams") + with no_decode(): + low, high = stencil_bounds(plan.s, plan.t, plan.radius, + plan.target_start, plan.target_length) + with owned(low, high): + require(left <= signed(low) <= signed(high) <= right, + "production halo contains complete expansion stencils") + require(decode(plan.inverse_plan.w) == 72 and right-left <= L + and L % (alpha*root) == 144, + "paired band, source padding and clipped last production block") + records.append({"tile_index": index, **actual, "source_start": left, + "source_count": right-left, "source_padding": L-(right-left), + "w": decode(plan.inverse_plan.w), "u_theta": str(Q(alpha*alpha*(t-s), s))}) + except ValueError as error: + if rejection is None: + raise + require(rejection in str(error), "specific production rejection: " + name + ": " + str(error)) + controls.append(name) + finally: + if plan is not None: + plan.close() + pool.empty() + finally: + for tape in pool.tapes: + erase_marker(tape) + require(saved == [(bytes(v.tape.cells), v.tape.head) for v in parameters], + "production metadata preserves borrowed parameters") + require(all(tape.head == 1 and all(v == BLANK for v in tape.cells) for tape in pool.tapes), + "production metadata restores owned tapes after success and rejection") + return {"tiles": records, "controls": controls, "pool_peak": pool.peak, + "scope": "encoded production metadata only; no production-sized coefficient page or payload execution"} diff --git a/research/physical-gaussian-expansion/expansion_check_sources.py b/research/physical-gaussian-expansion/expansion_check_sources.py new file mode 100644 index 0000000000..9ece37daad --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_sources.py @@ -0,0 +1,66 @@ +"""Complete source exports, live analytic pins and one actual physical cohort.""" +import importlib +import json +from pathlib import Path +from unittest.mock import patch +import expansion_sources +from expansion_sources import HERE, active, digest, require + + +def references(rows): + root = HERE.parents[1] + for row in rows: + if row["usage"] == "input": + path = root / row["path"] + require(path.is_file() and not path.is_symlink() and digest(path) == row["sha256"], + "pinned expansion live source: " + row["path"]) + + +def run(): + from cardinal_counts import Count + from records import Tape + from word_product import WordProduct + from expansion_tile import Tile + from expansion_tile_prepare import TilePreparation + rows = json.loads((HERE / "SOURCES.json").read_text())["files"] + references(rows) + modules = ("expansion_records", "expansion_geometry", "expansion_prepare", "expansion_block", + "expansion_tile_geometry", "expansion_tile_prepare", "expansion_tile") + for name in modules: + module = importlib.import_module(name) + require(Path(module.__file__).resolve() == HERE / (name + ".py"), "package-local expansion runtime") + for key, value in (("Count", Count), ("Tape", Tape)): + if hasattr(module, key): + require(getattr(module, key) is value, "one expansion Count/Tape cohort") + tile, producer = Tile(), TilePreparation() + require(tile.word is tile.block.word is tile.block.convolution.word + and type(tile.word.base) is WordProduct and type(tile.word.product.word) is WordProduct + and tile.word.base_bits == 4096, "actual acyclic word and convolution owners") + tile.idle() + producer.idle() + require(tile.closed and producer.closed, "idle complete expansion owners") + controls = [] + def reject(call, fragment, name): + try: + call() + except ValueError as error: + require(fragment in str(error), "specific expansion source rejection") + controls.append(name) + else: + raise ValueError("missing expansion source rejection") + altered = [dict(row) for row in rows] + next(row for row in altered if row["usage"] == "input")["sha256"] = "0"*64 + reject(lambda: references(altered), "pinned expansion live source", "live_pin") + with patch.object(expansion_sources, "CORE_DIGEST", "0"*64): + core_tar = active().core_package.parents[2] / "core.tar" + reject(lambda: expansion_sources.authenticate_core(core_tar), "complete uncompressed", "core_archive") + with patch.object(expansion_sources, "ROUNDING_HEAD", "0"*40): + reject(lambda: expansion_sources.authenticate_rounding(HERE / "inputs/rounding-source.tar"), + "exact final-grid correction commit", "rounding_head") + def second(): + with expansion_sources.opened(): + raise ValueError("unexpected second expansion cohort") + reject(second, "cold expansion namespace", "second_cohort") + return {"core_head": expansion_sources.CORE_HEAD, "rounding_head": expansion_sources.ROUNDING_HEAD, + "core_export_files": 43, "live_inputs": sum(row["usage"] == "input" for row in rows), + "runtime_modules": list(modules), "controls": controls} diff --git a/research/physical-gaussian-expansion/expansion_check_tile.py b/research/physical-gaussian-expansion/expansion_check_tile.py new file mode 100644 index 0000000000..c76a8a7ade --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_check_tile.py @@ -0,0 +1,135 @@ +"""Whole prepared tiles: periodic seams, repeated starts and a clipped last block.""" +from hashlib import sha256 +from expansion_sources import require +from expansion_check_block import expected_page, inputs, models, packet_bytes, periodic_error + + +def run(case): + from cardinal_counts import Count + from check_support import no_decode + from expansion_tile import Tile + from expansion_tile_geometry import TileHeader + from expansion_tile_prepare import TilePreparation, copy_cells + from gaussian_rationals import read_count + from records import Tape + from ring_check_words import body, decode + from ring_words import Signed, read_signed, write_signed + from router_metadata import owned + from session import erase_marker + from words import BLANK + require(case in ("overlap", "wrap"), "declared complete tile fixture") + if case == "overlap": + s, t, alpha, radius, source_start, source_count, capacity, start, length = 3, 4, 2, 1, -3, 8, 8, 0, 2 + else: + s, t, alpha, radius, source_start, source_count, capacity, start, length = 3, 4, 2, 2, -1, 9, 16, 3, 3 + p, bits, sigma = 4, 32, 3 + preparation, tile = TilePreparation(), Tile() + page = returned = None + source = Tape() + with owned(*(Count.literal(v) for v in + (s, t, alpha, radius, abs(source_start), source_count, capacity, abs(start), length, p, bits, sigma))) as args: + S, T, a, b, source_origin, count, padded, target_origin, target_count, target_p, P, scale = args + saved = [(bytes(value.tape.cells), value.tape.head) for value in args] + try: + with no_decode(): + page = preparation.run_explicit(S, T, a, b, Signed(source_origin, source_start < 0), + count, padded, Signed(target_origin, start < 0), target_count, target_p, P, scale) + page_bytes = packet_bytes(page) + expected, descriptions = [], [] + word_count = 2*capacity + with owned(TileHeader(page)) as (header,): + require(decode(header.blocks) == (length+radius-1)//radius + and decode(header.input_capacity) == capacity and decode(header.source_count) == source_count, + "complete expansion catalogue metadata") + offset = 0 + for _ in range(decode(header.blocks)): + block_page = Tape() + try: + with owned(read_count(page)) as (cells,): + copy_cells(page, block_page, cells) + block_page.rewind() + m = min(radius, length-offset) + block_start, factors, c = expected_page(block_page, s, t, alpha, radius, + start+offset, m, p, bits, sigma) + data = inputs(block_start, len(factors["A"]), s) + values, premature, error = models(data, factors, c, s, t, alpha, block_start, + start+offset, m, p, bits, sigma) + expected.extend(values) + descriptions.append({"source_start": block_start, "source_length": len(data), + "target_start": start+offset, "target_length": m, + "values": values, "premature_values": premature, + "finite_error": str(error)}) + word_count += 4*len(data)+6*m+8 + offset += m + finally: + block_page.clear() + erase_marker(block_page) + require(page.read() == BLANK and offset == length, "complete expansion block catalogue") + page.rewind() + global_error = periodic_error(expected, s, t, alpha, start, p) + with owned(Count.literal(p+2)) as (width,): + data = inputs(source_start, source_count, s)+[(0, 0)]*(capacity-source_count) + for pair in data: + for value in pair: + numerator = int(value*(1 << p)) + with owned(Signed(Count.literal(abs(numerator)), numerator < 0)) as (encoded,): + write_signed(encoded, source, width) + source.rewind() + source_bytes = body(source) + before = dict(tile.word.calls) + with no_decode(): + returned = tile.run(page, source) + raw = body(returned) + actual = [] + for _ in range(length): + with owned(read_signed(returned, width), read_signed(returned, width)) as (real, imaginary): + actual.append([(-1 if v.negative else 1)*decode(v.magnitude) for v in (real, imaginary)]) + require(actual == expected and returned.read() == BLANK, + "complete streamed expansion equals independently staged blocks") + require(any(v for pair in actual for v in pair), "nonzero complete expansion tile") + calls = {key: tile.word.calls[key]-before[key] for key in before} + require(sum(calls.values()) == word_count, "whole tile word invoice includes padding validation") + first = descriptions[0]["source_start"]-source_start + end = first+descriptions[0]["source_length"] + fresh = descriptions[0]["source_length"] + overlap = 0 + for record in descriptions[1:]: + begin = record["source_start"]-source_start + overlap += end-begin + following = begin+record["source_length"] + fresh += following-end + end = following + expected_counts = {"blocks": len(descriptions), "validation_records": capacity, + "fresh_records": fresh, "overlap_records": overlap, + "initial_skipped": first, "tail_skipped": capacity-end} + require(tile.last_counts == expected_counts + and first+fresh+capacity-end == capacity, + "one validation scan and one monotone source scan with bounded overlaps") + if case == "overlap": + require(descriptions[0]["source_start"] == descriptions[1]["source_start"] + and overlap > 0, "actual zero-advance source-window reuse") + else: + require(start < t < start+length and descriptions[-1]["target_length"] < radius + and capacity > source_count, "actual period crossing, clipped last block and zero padding") + require(packet_bytes(page) == page_bytes and body(source) == source_bytes, + "tile preserves borrowed complete catalogue and input") + require(saved == [(bytes(value.tape.cells), value.tape.head) for value in args], + "tile producer preserves borrowed global parameters") + require(preparation.closed and not preparation.pool.live and tile.closed and not tile.pool.live, + "complete tile owner restoration") + result = {"case": case, "s": s, "t": t, "alpha": alpha, "radius": radius, + "source": [source_start, source_count, capacity], "target": [start, length], + "p": p, "P": bits, "sigma": sigma, "values": actual, + "page_sha256": sha256(page_bytes).hexdigest(), "core_sha256": sha256(raw).hexdigest(), + "periodic_error": str(global_error), "blocks": descriptions, + "word_calls": calls, "stream_counts": tile.last_counts, + "preparation_counts": preparation.last_counts, + "preparation_pool_peak": preparation.pool.peak, "tile_pool_peak": tile.pool.peak, + "scope": "complete encoded tile production and streaming at explicit finite grids; no production-size tile execution"} + finally: + for tape in (returned, page, source): + if tape is not None: + tape.clear() + erase_marker(tape) + tile.idle() + return result diff --git a/research/physical-gaussian-expansion/expansion_geometry.py b/research/physical-gaussian-expansion/expansion_geometry.py new file mode 100644 index 0000000000..816a371a21 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_geometry.py @@ -0,0 +1,155 @@ +"""Global lifted geometry and rational real-chirp exponents for one block.""" +from expansion_sources import active +active() +from cardinal_counts import Count, positive +from encoded import compare, subtract +from expansion_records import Fields, compile_frame, work_width +from gaussian_rationals import Rational +from ring_words import Signed, add, fits +from router_metadata import owned +from words import need + + +def integer(value): + return Rational(Signed(Count(value)), Count.literal(1), False) + + +def absolute(value): + return Rational(Signed(Count(value.numerator.magnitude)), Count(value.denominator), False) + + +class BlockGeometry(Fields): + def __init__(self, s, t, alpha, radius, target_start, target_length, p, precision, sigma): + super().__init__() + try: + for name, value in (("s", s), ("t", t), ("alpha", alpha), ("radius", radius), + ("m", target_length), ("p", p), ("P", precision), ("sigma", sigma)): + self.keep(name, Count(value)) + if name != "sigma": + positive(value) + self.keep("version", Count.literal(1)) + self.keep("target_start", target_start.copy()) + with owned(Count(s)) as (twice_s,): + twice_s.multiply_small(2) + need(compare(s, t) < 0 and compare(t, twice_s) < 0, "expansion source ratio") + need(compare(target_length, radius) <= 0, "expansion target block within radius") + fits(s, p) + fits(t, p) + fits(alpha, p) + with owned(Count(p)) as (origin_width,): + origin_width.increment() + origin_width.increment() + fits(target_start.magnitude, origin_width) + self.keep("u", Count.product(alpha, alpha)) + self.keep("kappa", Rational(Signed(Count(s)), Count(t))) + self.keep("u_value", integer(self.u)) + with owned(Rational.integer(target_start)) as (target,): + with owned(self.kappa.times(target)) as (center,): + self.keep("j0", center.integer_part()) + with owned(Rational.integer(self.j0)) as (origin,): + self.keep("gamma", origin.minus(center)) + with owned(Rational.literal(-1), Rational.literal(0)) as (minus_one, zero): + need(self.gamma.compare(minus_one) > 0 and self.gamma.compare(zero) <= 0, + "global expansion fractional origin") + low_magnitude = Count(radius) + low_magnitude.increment() + self.keep("h_low", Signed(low_magnitude, True)) + with owned(Count.literal(1)) as (one,): + with owned(subtract(target_length, one)) as (last,): + with owned(integer(last)) as (position,): + with owned(self.kappa.times(position)) as (scaled,): + with owned(scaled.minus(self.gamma)) as (relative,): + high = relative.integer_part(True) + try: + need(not high.negative, "nonnegative expansion upper origin") + high.magnitude.add(radius) + high.magnitude.increment() + self.keep("h_high", high) + high = None + finally: + if high is not None: + high.close() + self.keep("n", Count(self.h_high.magnitude)) + self.n.add(radius) + self.n.increment() + self.n.increment() + with owned(Count(self.n), Count.literal(1)) as (total, one): + total.add(target_length) + self.keep("nc", subtract(total, one)) + self.keep("source_start", add(self.j0, self.h_low)) + self.keep("kernel_start", self.h_high.negated()) + self.keep("width", work_width(precision, self.n)) + compile_frame(self) + self.prepare_quadratics() + except BaseException: + self.close() + raise + + def prepare_quadratics(self): + with owned(Rational.literal(1), Rational.literal(2)) as (one, two): + with owned(one.minus(self.kappa)) as (difference,): + self.keep("a_square", difference.over(self.u_value)) + with owned(two.times(self.gamma)) as (twice_gamma,): + self.keep("a_linear", twice_gamma.over(self.u_value)) + with owned(twice_gamma.times(self.kappa)) as (product,): + with owned(product.negated()) as (negative,): + self.keep("b_linear", negative.over(self.u_value)) + with owned(self.gamma.times(self.gamma)) as (square,): + self.keep("a_constant", square.over(self.u_value)) + with owned(self.kappa.times(self.kappa)) as (square,): + with owned(square.minus(self.kappa)) as (difference,): + self.keep("b_square", difference.over(self.u_value)) + self.keep("kernel_square", self.kappa.over(self.u_value)) + + def exponent(self, kind, index): + need(kind in ("A", "K", "B"), "fixed expansion factor family") + length = self.n if kind == "A" else self.nc if kind == "K" else self.m + need(compare(index, length) < 0, "expansion factor index") + with owned(Signed(Count(index))) as (ordinal,): + if kind == "A": + position = add(self.h_low, ordinal) + elif kind == "K": + position = add(self.kernel_start, ordinal) + else: + position = ordinal.copy() + with owned(position): + with owned(Rational.integer(position)) as (x,): + with owned(x.times(x)) as (square,): + if kind == "K": + return self.kernel_square.times(square) + square_factor = self.a_square if kind == "A" else self.b_square + linear_factor = self.a_linear if kind == "A" else self.b_linear + with owned(square_factor.times(square), linear_factor.times(x)) as (quadratic, linear): + total = quadratic.plus(linear) + if kind == "B": + return total + with owned(total): + return total.plus(self.a_constant) + + +class ProductionBlock(Fields): + """Derive the retained all-size radius, work grid and normalization scale.""" + def __init__(self, p, alpha): + super().__init__() + try: + positive(p) + positive(alpha) + with owned(Count.literal(101), Count.literal(2)) as (minimum, two): + need(compare(p, minimum) >= 0 and compare(alpha, two) >= 0, + "production expansion precision and integer alpha") + self.keep("u", Count.product(alpha, alpha)) + need(compare(self.u, p) < 0, "production alpha squared below precision") + self.keep("root", Count.literal(1)) + while True: + with owned(Count.product(self.root, self.root)) as (square,): + if compare(square, p) >= 0: + break + self.root.increment() + self.keep("radius", Count.product(alpha, self.root)) + self.keep("P", Count(p)) + self.P.multiply_small(32768) + self.keep("sigma", Count(p)) + self.sigma.multiply_small(4096) + except BaseException: + self.close() + raise diff --git a/research/physical-gaussian-expansion/expansion_prepare.py b/research/physical-gaussian-expansion/expansion_prepare.py new file mode 100644 index 0000000000..69f862e7a1 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_prepare.py @@ -0,0 +1,193 @@ +"""Encoded real-Gaussian factor production and protected complete block pages.""" +from expansion_sources import active +active() +from cardinal_counts import Count +from encoded import compare, subtract +from expansion_geometry import BlockGeometry, ProductionBlock, absolute +from expansion_records import BlockHeader +from gaussian_exponentials import Exponentials +from gaussian_page import skip_bits +from gaussian_rationals import read_rational, write_rational +from records import Tape +from ring_words import Signed, read_signed, right_shift, write_signed +from router_metadata import indices, owned +from scratch_pool import ScratchPool +from session import (OuterScope, borrowed_count, erase_marker, initialize_marker, + sealed_constructors) +from twiddle_setup import quotient +from words import BLANK, need + + +def reflected_kernel(plan, index, cache): + if compare(index, plan.h_high.magnitude) <= 0: + return None + with owned(Count(plan.h_high.magnitude)) as (twice,): + twice.multiply_small(2) + if compare(index, twice) > 0: + return None + with owned(subtract(twice, index)) as (earlier,): + cache.rewind() + with owned(Count.product(earlier, plan.width)) as (prefix,): + skip_bits(cache, prefix) + value = read_signed(cache, plan.width) + try: + with owned(subtract(index, earlier), Count.literal(1)) as (distance, one): + with owned(subtract(distance, one)) as (remaining,): + with owned(Count.product(remaining, plan.width)) as (suffix,): + skip_bits(cache, suffix) + need(cache.read() == BLANK, "complete reflected expansion kernel prefix") + return value + except BaseException: + value.close() + raise + + +def prepare_page(plan, target, pool): + streams = {} + bound = exponentials = kernel_cache = None + exact_units = reflected = 0 + try: + bound = Count.literal(1) + for kind, length in (("A", plan.n), ("K", plan.nc), ("B", plan.m)): + stream = pool.acquire() + streams[kind] = stream + for index in indices(length): + with owned(plan.exponent(kind, index)) as (exponent,): + write_rational(exponent, stream) + with owned(absolute(exponent)) as (size,): + with owned(size.integer_part(True)) as (ceiling,): + if compare(ceiling.magnitude, bound) > 0: + following = Count(ceiling.magnitude) + bound.close() + bound = following + stream.rewind() + exponentials = Exponentials(plan.P, bound) + kernel_cache = pool.acquire() + BlockHeader.dump(plan, target) + for kind, length in (("A", plan.n), ("K", plan.nc), ("B", plan.m)): + stream = streams[kind] + scale = None if kind == "K" else plan.sigma + for index in indices(length): + with owned(read_rational(stream)) as (exponent,): + if kind == "K": + need(not exponent.numerator.negative, "nonnegative expansion Gaussian kernel") + reused = reflected_kernel(plan, index, kernel_cache) if kind == "K" else None + if reused is not None: + with owned(reused): + write_signed(reused, target, plan.width) + write_signed(reused, kernel_cache, plan.width) + reflected += 1 + elif exponent.numerator.magnitude.zero(): + value = Count(plan.unit) if scale is None else right_shift(plan.unit, scale) + with owned(Signed(value)) as (coefficient,): + write_signed(coefficient, target, plan.width) + if kind == "K": + write_signed(coefficient, kernel_cache, plan.width) + exact_units += 1 + else: + with owned(absolute(exponent)) as (argument,): + lower, upper = exponentials.grid_interval(argument, exponent.numerator.negative, scale) + with owned(Signed(lower), upper) as (coefficient, bound_upper): + need(compare(bound_upper, plan.unit) <= 0, + "certified expansion factor within normalized unit interval") + write_signed(coefficient, target, plan.width) + if kind == "K": + write_signed(coefficient, kernel_cache, plan.width) + need(stream.read() == BLANK, "complete prepared expansion factor family") + with owned(Count(plan.alpha)) as (denominator,): + denominator.multiply_small(2) + with owned(Signed(quotient(plan.unit, denominator))) as (normalization,): + need(not normalization.magnitude.zero(), "positive expansion normalization") + write_signed(normalization, target, plan.width) + return {"exponentials": exponentials.calls, "exact_unit_factors": exact_units, + "reflected_kernel_factors": reflected, + "taylor_steps": exponentials.taylor_steps, "early_stops": exponentials.early_stops} + finally: + if exponentials is not None: + exponentials.close() + if bound is not None: + bound.close() + for stream in streams.values(): + stream.clear() + if kernel_cache is not None: + kernel_cache.clear() + + +class Preparation: + def __init__(self): + self.pool = ScratchPool(256) + self.scope = OuterScope(self.pool) + self.parameter_pool = ScratchPool(64) + self.parameter_scope = OuterScope(self.parameter_pool) + self.output = Tape() + self.closed = True + self.calls = 0 + self.last_counts = None + for tape in (*self.pool.tapes, *self.parameter_pool.tapes, self.output): + erase_marker(tape) + + def run_explicit(self, s, t, alpha, radius, target_start, target_length, p, precision, sigma): + s, t, alpha, radius, target_length, p, precision, sigma = (borrowed_count(value) for value in + (s, t, alpha, radius, target_length, p, precision, sigma)) + target_start = Signed(borrowed_count(target_start.magnitude), target_start.negative) + parameters = (s, t, alpha, radius, target_length, p, precision, sigma, target_start.magnitude) + need(all(value.tape not in (*self.pool.tapes, self.output) for value in parameters), + "external expansion preparation parameters") + need(self.closed and not self.pool.live and self.output.read() == BLANK, + "consume prior expansion prepared page") + for tape in (*self.pool.tapes, self.output): + initialize_marker(tape) + self.closed = False + self.last_counts = None + plan = None + succeeded = False + try: + with sealed_constructors(), self.scope.installed(): + try: + plan = BlockGeometry(s, t, alpha, radius, target_start, target_length, p, precision, sigma) + self.last_counts = prepare_page(plan, self.output, self.pool) + self.calls += 1 + succeeded = True + finally: + if plan is not None: + plan.close() + for value in parameters: + value.rewind() + self.pool.empty() + finally: + for tape in self.pool.tapes: + erase_marker(tape) + self.output.rewind() + if not succeeded: + self.output.clear() + erase_marker(self.output) + self.closed = True + return self.output + + def run_production(self, s, t, alpha, target_start, target_length, p): + need(self.closed and not self.pool.live and not self.parameter_pool.live and self.output.read() == BLANK, + "consume prior production expansion page") + s, t, alpha, target_length, p = (borrowed_count(value) for value in (s, t, alpha, target_length, p)) + target_start = Signed(borrowed_count(target_start.magnitude), target_start.negative) + parameters = (s, t, alpha, target_length, p, target_start.magnitude) + need(all(value.tape not in (*self.pool.tapes, *self.parameter_pool.tapes, self.output) + for value in parameters), "external production expansion parameters") + for tape in self.parameter_pool.tapes: + initialize_marker(tape) + plan = None + try: + with sealed_constructors(), self.parameter_scope.installed(): + try: + plan = ProductionBlock(p, alpha) + with self.parameter_scope.suspended(): + return self.run_explicit(s, t, alpha, plan.radius, target_start, + target_length, p, plan.P, plan.sigma) + finally: + if plan is not None: + plan.close() + for value in parameters: + value.rewind() + self.parameter_pool.empty() + finally: + for tape in self.parameter_pool.tapes: + erase_marker(tape) diff --git a/research/physical-gaussian-expansion/expansion_records.py b/research/physical-gaussian-expansion/expansion_records.py new file mode 100644 index 0000000000..0a90dafb5c --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_records.py @@ -0,0 +1,152 @@ +"""Canonical expansion page metadata and complete numerical frame checks.""" +from expansion_sources import active +active() +from cardinal_counts import Count, positive +from encoded import compare, subtract +from gaussian_rationals import read_count, write_count +from ring_words import Signed, power_two +from router_metadata import ceiling_log, owned +from words import need + + +class Fields: + def __init__(self): + self.fields = {} + + def keep(self, name, value): + need(name not in self.fields, "one expansion metadata role") + self.fields[name] = value + return value + + def __getattr__(self, name): + return self.fields[name] + + def close(self): + for value in reversed(tuple(self.fields.values())): + value.close() + self.fields.clear() + + +def write_integer(value, target): + target.append(int(value.negative)) + write_count(value.magnitude, target) + + +def read_integer(source): + sign = source.read() + need(sign in (0, 1), "expansion integer sign") + source.step(1) + magnitude = read_count(source) + try: + need(not sign or not magnitude.zero(), "canonical expansion integer zero") + return Signed(magnitude, bool(sign)) + except BaseException: + magnitude.close() + raise + + +def same_integer(first, second): + return first.negative == second.negative and compare(first.magnitude, second.magnitude) == 0 + + +def work_width(precision, length): + result = Count(precision) + try: + with owned(Count(length), Count.literal(3)) as (more, three): + more.increment() + with owned(ceiling_log(more)) as (guard,): + result.add(guard) + result.add(three) + return result + except BaseException: + result.close() + raise + + +def compile_frame(owner): + for name in BlockHeader.names: + if name != "sigma": + positive(getattr(owner, name)) + with owned(Count.literal(1), Count.literal(2)) as (one, two): + need(compare(owner.version, one) == 0, "expansion page version") + need(compare(owner.alpha, two) >= 0, "expansion integer alpha at least two") + with owned(Count(owner.s)) as (twice_s,): + twice_s.multiply_small(2) + need(compare(owner.s, owner.t) < 0 and compare(owner.t, twice_s) < 0, + "near-one expansion ratio") + need(compare(owner.m, owner.radius) <= 0 and compare(owner.m, owner.n) <= 0, + "complete short expansion block shape") + with owned(Count(owner.radius), Count.literal(3)) as (minimum, three): + minimum.multiply_small(2) + minimum.add(three) + with owned(Count(minimum)) as (maximum,): + maximum.add(owner.m) + need(compare(minimum, owner.n) <= 0 and compare(owner.n, maximum) <= 0, + "complete bounded expansion stencil span") + with owned(Count(owner.n), Count.literal(1)) as (total, one): + total.add(owner.m) + with owned(subtract(total, one)) as (expected,): + need(compare(expected, owner.nc) == 0, "complete prefix convolution extent") + with owned(Count(owner.p), Count(owner.p)) as (max_n, max_nc): + max_n.multiply_small(4) + max_nc.multiply_small(5) + need(compare(owner.n, max_n) <= 0 and compare(owner.nc, max_nc) <= 0, + "expansion block lengths within precision allowance") + with owned(work_width(owner.P, owner.n)) as (expected,): + need(compare(owner.width, expected) == 0, "expansion work width") + with owned(Count(owner.n), Count(owner.sigma), Count(owner.p), + Count.literal(12)) as (scaled_length, twice_sigma, minimum, guard): + scaled_length.increment() + scaled_length.multiply_small(32) + twice_sigma.multiply_small(2) + minimum.add(twice_sigma) + minimum.add(guard) + with owned(ceiling_log(scaled_length)) as (logarithm,): + minimum.add(logarithm) + need(compare(owner.P, minimum) >= 0, "expansion pre-extraction precision guard") + with owned(subtract(owner.P, owner.p)) as (shift,): + owner.keep("promote", Count(shift)) + owner.keep("extract", subtract(shift, twice_sigma)) + owner.keep("input_width", Count(owner.p)) + owner.input_width.increment() + owner.input_width.increment() + owner.keep("unit", power_two(owner.P)) + with owned(Count(owner.p)) as (twice_p,): + twice_p.multiply_small(2) + owner.keep("norm_bound", power_two(twice_p)) + for name, length in (("input_bits", owner.n), ("kernel_bits", owner.nc), + ("output_bits", owner.m)): + owner.keep(name, Count.product(length, owner.width)) + + +class BlockHeader(Fields): + names = ("version", "p", "P", "width", "sigma", "s", "t", "alpha", "radius", "n", "m", "nc") + origins = ("source_start", "target_start") + + def __init__(self, page): + super().__init__() + try: + for name in self.names: + self.keep(name, read_count(page)) + for name in self.origins: + self.keep(name, read_integer(page)) + compile_frame(self) + except BaseException: + self.close() + raise + + @classmethod + def dump(cls, plan, page): + for name in cls.names: + write_count(getattr(plan, name), page) + for name in cls.origins: + write_integer(getattr(plan, name), page) + + def repeat(self, page): + page.rewind() + for name in self.names: + with owned(read_count(page)) as (value,): + need(compare(value, getattr(self, name)) == 0, "unchanged expansion block header") + for name in self.origins: + with owned(read_integer(page)) as (value,): + need(same_integer(value, getattr(self, name)), "unchanged expansion block origin") diff --git a/research/physical-gaussian-expansion/expansion_sources.py b/research/physical-gaussian-expansion/expansion_sources.py new file mode 100644 index 0000000000..f8073c8214 --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_sources.py @@ -0,0 +1,129 @@ +"""Authenticate the complete Gaussian cohort and the final-grid correction.""" +from contextlib import contextmanager +import gzip +from hashlib import sha256 +import importlib.util +import json +from pathlib import Path, PurePosixPath +import shutil +import sys +import tarfile +import tempfile + +HERE = Path(__file__).resolve().parent +CORE_PACKAGE = "research/physical-gaussian-core" +CORE_HEAD = "7ab67b4a84d28cd562b9711aa315ca5b47822b37" +COMPRESSED_DIGEST = "5aaa84523bfa22eb817023d3130e8b8879153e9207e8f836cc6c61435fd051b8" +CORE_DIGEST = "2d97e093634a4119d9990dcc5d175917df3eabda5e8d1c204f640381e92ab4a2" +CORE_RECEIPT = "d387e4d0a1a43b9b24104f76f8992603d4265bbe447601bea476d492e92e5340" +ROUNDING_HEAD = "988ef67986c313445936fc7103908b6163b9d47e" +ROUNDING_DIGEST = "31a553fc4aa4443e60f285633ceb309edc84ede4f64ef7c2646f8806e61da1e8" +_active = None + + +def require(condition, message): + if not condition: + raise ValueError(message) + + +def digest(path): + value = sha256() + with path.open("rb") as source: + for block in iter(lambda: source.read(65536), b""): + value.update(block) + return value.hexdigest() + + +def regular(item): + path = PurePosixPath(item.name) + require((item.isfile() or item.isdir()) and not path.is_absolute() + and ".." not in path.parts and "\\" not in item.name + and path.as_posix() == item.name and item.size <= 2**27, + "bounded regular expansion source member") + + +def authenticate_core(path): + require(digest(path) == CORE_DIGEST, "complete uncompressed Gaussian core export") + with tarfile.open(path, "r:") as archive: + members = archive.getmembers() + names = [item.name for item in members] + require(len(names) == len(set(names)) == 43, "complete unique core export inventory") + for item in members: + regular(item) + require(item.isfile(), "regular core export member") + raw = archive.extractfile(".research-export.json").read() + require(sha256(raw).hexdigest() == CORE_RECEIPT, "pinned complete core receipt") + receipt = json.loads(raw) + require((receipt["commit"], receipt["package"], len(receipt["files"])) + == (CORE_HEAD, CORE_PACKAGE, 42), "Gaussian core source identity") + require(set(names) == {row["path"] for row in receipt["files"]} | {".research-export.json"}, + "exact core source member set") + for row in receipt["files"]: + data = archive.extractfile(row["path"]).read() + require(len(data) == row["size"] and sha256(data).hexdigest() == row["sha256"], + "core source member: " + row["path"]) + + +def authenticate_rounding(path): + require(path.is_file() and not path.is_symlink() and digest(path) == ROUNDING_DIGEST, + "pinned final-grid correction source") + with tarfile.open(path, "r:") as archive: + require(archive.pax_headers.get("comment") == ROUNDING_HEAD, + "exact final-grid correction commit") + members = archive.getmembers() + require(len({item.name for item in members}) == len(members), "unique correction source members") + expected = {"LICENSE", "NOTICE", "notes/fast-gaussian-resampling.tex", + "scripts/check_fast_gaussian_rounding.py"} + require({item.name for item in members if item.isfile()} == expected, + "complete correction files and notices") + for item in members: + regular(item) + require(item.name in expected | {"notes", "scripts"}, "correction source roots") + + +def active(): + require(_active is not None, "authenticated expansion physical cohort") + return _active + + +@contextmanager +def opened(): + global _active + require(_active is None and "gaussian_sources" not in sys.modules + and "cardinal_counts" not in sys.modules, "cold expansion namespace") + compressed = HERE / "inputs/core-export.tar.gz" + require(compressed.is_file() and not compressed.is_symlink() + and digest(compressed) == COMPRESSED_DIGEST, "pinned compressed core source") + rounding = HERE / "inputs/rounding-source.tar" + authenticate_rounding(rounding) + old_path = sys.path[:] + with tempfile.TemporaryDirectory(prefix="gaussian-expansion-source-") as directory: + root = Path(directory) + source_tar = root / "core.tar" + with gzip.open(compressed, "rb") as source, source_tar.open("xb") as destination: + shutil.copyfileobj(source, destination, 65536) + require(source_tar.stat().st_size <= 2**27, "bounded decompressed core source") + authenticate_core(source_tar) + core_root, rounding_root = root / "core", root / "rounding" + core_root.mkdir() + rounding_root.mkdir() + for path, target in ((source_tar, core_root), (rounding, rounding_root)): + with tarfile.open(path, "r:") as archive: + archive.extractall(target, filter="data") + package = core_root / CORE_PACKAGE + spec = importlib.util.spec_from_file_location("gaussian_sources", package / "gaussian_sources.py") + source = importlib.util.module_from_spec(spec) + sys.modules["gaussian_sources"] = source + try: + spec.loader.exec_module(source) + sys.path.insert(0, str(package)) + with source.opened() as state: + _active = state + state.core_package = package + state.rounding_root = rounding_root + yield state + finally: + _active = None + sys.path[:] = old_path + if sys.modules.get("gaussian_sources") is source: + del sys.modules["gaussian_sources"] diff --git a/research/physical-gaussian-expansion/expansion_tile.py b/research/physical-gaussian-expansion/expansion_tile.py new file mode 100644 index 0000000000..d1b4bebb2e --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_tile.py @@ -0,0 +1,184 @@ +"""Whole target tiles with a monotone source stream and two short overlap banks.""" +from expansion_sources import active +active() +from cardinal_counts import Count, positive +from compact_counts import copy_bits +from encoded import compare, subtract +from expansion_block import Block +from expansion_records import BlockHeader, same_integer +from expansion_tile_geometry import TileHeader, advance +from expansion_tile_prepare import copy_cells +from fft_arithmetic import Arithmetic +from gaussian_page import skip_bits +from gaussian_rationals import read_count +from records import Tape +from ring_words import difference, read_signed +from router_metadata import indices, owned +from scratch_pool import ScratchPool +from session import (OuterScope, borrowed_tape, erase_marker, initialize_marker, + sealed_constructors) +from words import BLANK, need + + +class Tile: + def __init__(self): + self.block = Block() + self.word = self.block.word + self.pool = ScratchPool(256) + self.scope = OuterScope(self.pool) + self.arithmetic = Arithmetic(self.word, self.scope) + self.output = Tape() + self.closed = True + self.calls = 0 + self.last_counts = None + for tape in (*self.pool.tapes, self.output): + erase_marker(tape) + + def idle(self): + self.block.idle() + need(self.block.closed and not self.block.pool.live and self.block.output.read() == BLANK, + "consumed expansion tile child return") + + def norm(self, real, imaginary, header): + with owned(self.arithmetic.multiply(real, real, header.input_width)) as (first,): + with owned(self.arithmetic.multiply(imaginary, imaginary, header.input_width)) as (second,): + first.magnitude.add(second.magnitude) + need(compare(first.magnitude, header.norm_bound) <= 0, "expansion tile disk-grid input") + + def validate(self, source, header): + source.rewind() + try: + for index in indices(header.input_capacity): + with owned(read_signed(source, header.input_width)) as (real,): + with owned(read_signed(source, header.input_width)) as (imaginary,): + self.norm(real, imaginary, header) + if compare(index, header.source_count) >= 0: + need(real.magnitude.zero() and imaginary.magnitude.zero(), "zero expansion source padding") + self.last_counts["validation_records"] += 1 + need(source.read() == BLANK, "complete padded expansion source extent") + finally: + source.rewind() + + def copy_records(self, source, target, length, width, counter): + for _ in indices(length): + copy_bits(source, target, width) + self.last_counts[counter] += 1 + + def skip_records(self, source, length, width, counter): + for _ in indices(length): + skip_bits(source, width) + self.last_counts[counter] += 1 + + def run(self, page, source): + page, source = borrowed_tape(page), borrowed_tape(source) + forbidden = (*self.pool.tapes, self.output, *self.block.pool.tapes, self.block.output, + *self.block.convolution.pool.tapes, self.block.convolution.output, + *self.word.pool.tapes, self.word.output, *self.word.canonical) + need(page is not source and page not in forbidden and source not in forbidden, + "external expansion tile ports") + need(self.closed and not self.pool.live and self.output.read() == BLANK, + "consume prior expansion tile output") + self.idle() + for tape in (*self.pool.tapes, self.output): + initialize_marker(tape) + self.closed = False + self.last_counts = dict(blocks=0, validation_records=0, fresh_records=0, + overlap_records=0, initial_skipped=0, tail_skipped=0) + plan = block_page = previous = current = previous_offset = consumed = emitted = next_target = None + succeeded = False + try: + with sealed_constructors(), self.scope.installed(): + try: + page.rewind() + plan = TileHeader(page) + self.validate(source, plan) + block_page, previous, current = (self.pool.acquire() for _ in range(3)) + consumed, emitted, next_target = Count.literal(0), Count.literal(0), plan.target_start.copy() + for _ in indices(plan.blocks): + borrowed_tape(block_page).clear() + with owned(read_count(page)) as (cells,): + positive(cells) + copy_cells(page, block_page, cells) + block_page.rewind() + with owned(BlockHeader(block_page)) as (header,): + for name in ("s", "t", "alpha", "radius", "p", "P", "sigma"): + need(compare(getattr(header, name), getattr(plan, name)) == 0, + "one global expansion tile axis and precision") + need(same_integer(header.target_start, next_target), "ordered expansion target blocks") + with owned(subtract(plan.target_length, emitted)) as (remaining,): + wanted = plan.radius if compare(plan.radius, remaining) <= 0 else remaining + need(compare(header.m, wanted) == 0, "complete or clipped final expansion block") + with owned(difference(header.source_start, plan.source_start)) as (offset,): + need(not offset.negative, "expansion source block begins in supplied halo") + with owned(Count(offset.magnitude)) as (end,): + end.add(header.n) + need(compare(end, plan.source_count) <= 0, + "expansion source block ends in supplied halo") + borrowed_tape(current).clear() + if previous_offset is None: + self.skip_records(source, offset.magnitude, plan.record_width, "initial_skipped") + consumed.add(offset.magnitude) + self.copy_records(source, current, header.n, plan.record_width, "fresh_records") + consumed.add(header.n) + else: + need(compare(previous_offset, offset.magnitude) <= 0 + and compare(offset.magnitude, consumed) <= 0, + "monotone overlapping expansion source windows") + previous.rewind() + with owned(subtract(offset.magnitude, previous_offset)) as (skip,): + with owned(Count.product(skip, plan.record_width)) as (bits,): + skip_bits(previous, bits) + with owned(subtract(consumed, offset.magnitude)) as (overlap,): + need(compare(overlap, header.n) <= 0, + "monotone expansion source right endpoints") + self.copy_records(previous, current, overlap, plan.record_width, "overlap_records") + need(previous.read() == BLANK, "complete retained expansion overlap") + with owned(subtract(header.n, overlap)) as (fresh,): + self.copy_records(source, current, fresh, plan.record_width, "fresh_records") + consumed.add(fresh) + current.rewind() + block_page.rewind() + with self.scope.suspended(): + returned = self.block.run(block_page, current) + try: + with owned(Count.product(header.m, plan.record_width)) as (bits,): + copy_bits(returned, self.output, bits) + need(returned.read() == BLANK, "complete expansion block return in tile") + finally: + returned.clear() + erase_marker(returned) + if previous_offset is not None: + previous_offset.close() + previous_offset = Count(offset.magnitude) + previous, current = current, previous + emitted.add(header.m) + following = advance(next_target, header.m) + next_target.close() + next_target = following + self.last_counts["blocks"] += 1 + with owned(subtract(plan.input_capacity, consumed)) as (tail,): + self.skip_records(source, tail, plan.record_width, "tail_skipped") + need(source.read() == page.read() == BLANK and compare(emitted, plan.target_length) == 0, + "complete streamed expansion tile and catalogue") + self.idle() + self.calls += 1 + succeeded = True + finally: + for value in (block_page, previous, current): + if value is not None: + value.clear() + for value in (previous_offset, consumed, emitted, next_target, plan): + if value is not None: + value.close() + page.rewind() + source.rewind() + self.pool.empty() + finally: + for tape in self.pool.tapes: + erase_marker(tape) + self.output.rewind() + if not succeeded: + self.output.clear() + erase_marker(self.output) + self.closed = True + return self.output diff --git a/research/physical-gaussian-expansion/expansion_tile_geometry.py b/research/physical-gaussian-expansion/expansion_tile_geometry.py new file mode 100644 index 0000000000..30ff87e1ff --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_tile_geometry.py @@ -0,0 +1,197 @@ +"""Complete tile frames, contained short windows and retained bulk metadata.""" +from expansion_sources import active +active() +from cardinal_counts import Count, positive +from encoded import compare, subtract +from expansion_geometry import ProductionBlock +from expansion_records import Fields, read_integer, write_integer +from gaussian_prepare import ProductionPlan +from gaussian_rationals import Rational, read_count, write_count +from ring_words import Signed, add, difference, fits, power_two +from router_metadata import ceiling_log, owned +from twiddle_setup import quotient +from words import need + + +def advance(origin, length): + with owned(Signed(Count(length))) as (distance,): + return add(origin, distance) + + +def source_index(s, t, target, ceiling=False): + with owned(Rational(Signed(Count.product(s, target.magnitude), target.negative), Count(t))) as (ratio,): + return ratio.integer_part(ceiling) + + +def stencil_bounds(s, t, radius, target_start, target_length): + with owned(source_index(s, t, target_start), Count(radius)) as (first, distance): + distance.increment() + with owned(Signed(Count(distance))) as (margin,): + start = difference(first, margin) + try: + with owned(advance(target_start, target_length), Signed(Count.literal(1))) as (end, one): + with owned(difference(end, one)) as (last_target,): + with owned(source_index(s, t, last_target, True), Count(radius)) as (last, distance): + distance.increment() + distance.increment() + finish = advance(last, distance) + return start, finish + except BaseException: + start.close() + raise + + +def compile_tile(owner): + for name in TileHeader.names: + if name != "sigma": + positive(getattr(owner, name)) + with owned(Count.literal(2)) as (two,): + need(compare(owner.version, two) == 0 and compare(owner.alpha, two) >= 0, + "expansion tile version and integer alpha") + with owned(Count(owner.s)) as (twice,): + twice.multiply_small(2) + need(compare(owner.s, owner.t) < 0 and compare(owner.t, twice) < 0, + "expansion tile near-one ratio") + for value in (owner.s, owner.t, owner.alpha): + fits(value, owner.p) + need(compare(owner.source_count, owner.input_capacity) <= 0, "expansion halo fits its padded capacity") + with owned(quotient(owner.target_length, owner.radius, True)) as (blocks,): + need(compare(blocks, owner.blocks) == 0, "complete expansion tile block count") + with owned(Count(owner.radius), Count(owner.radius), Count.literal(3), Count.literal(2), + Count(owner.p), Count(owner.p)) as (max_n, max_nc, three, two, p4, p5): + max_n.multiply_small(3) + max_n.add(three) + max_nc.multiply_small(4) + max_nc.add(two) + p4.multiply_small(4) + p5.multiply_small(5) + need(compare(max_n, p4) <= 0 and compare(max_nc, p5) <= 0, + "uniform short expansion tile blocks") + max_n.increment() + max_n.multiply_small(32) + with owned(ceiling_log(max_n), Count(owner.p), Count(owner.sigma), Count.literal(12)) as (logarithm, minimum, scale, guard): + scale.multiply_small(2) + minimum.add(scale) + minimum.add(logarithm) + minimum.add(guard) + need(compare(owner.P, minimum) >= 0, "complete expansion tile precision guard") + owner.keep("input_width", Count(owner.p)) + owner.input_width.increment() + owner.input_width.increment() + owner.keep("record_width", Count(owner.input_width)) + owner.record_width.multiply_small(2) + with owned(Count(owner.p)) as (twice_p,): + twice_p.multiply_small(2) + owner.keep("norm_bound", power_two(twice_p)) + + +class TileHeader(Fields): + names = ("version", "p", "P", "sigma", "s", "t", "alpha", "radius", + "source_count", "input_capacity", "target_length", "blocks") + origins = ("source_start", "target_start") + + def __init__(self, page): + super().__init__() + try: + for name in self.names: + self.keep(name, read_count(page)) + for name in self.origins: + self.keep(name, read_integer(page)) + compile_tile(self) + except BaseException: + self.close() + raise + + @classmethod + def dump(cls, plan, target): + for name in cls.names: + write_count(getattr(plan, name), target) + for name in cls.origins: + write_integer(getattr(plan, name), target) + + +class TileGeometry(Fields): + def __init__(self, s, t, alpha, radius, source_start, source_count, input_capacity, + target_start, target_length, p, precision, sigma): + super().__init__() + try: + for name, value in (("s", s), ("t", t), ("alpha", alpha), ("radius", radius), + ("source_count", source_count), ("input_capacity", input_capacity), + ("target_length", target_length), ("p", p), ("P", precision), ("sigma", sigma)): + self.keep(name, Count(value)) + if name != "sigma": + positive(value) + self.keep("version", Count.literal(2)) + self.keep("source_start", source_start.copy()) + self.keep("target_start", target_start.copy()) + self.keep("blocks", quotient(target_length, radius, True)) + compile_tile(self) + first, last = stencil_bounds(s, t, radius, target_start, target_length) + with owned(first, last, advance(source_start, source_count)) as (first, last, source_end): + with owned(difference(first, source_start), difference(source_end, last)) as (left, right): + need(not left.negative and not right.negative, + "complete expansion source halo contains every short stencil") + except BaseException: + self.close() + raise + + +class ProductionTile(Fields): + """Retain the microbox halo and the matching compression metadata.""" + def __init__(self, s, t, d, alpha, p, tile_index): + super().__init__() + try: + self.keep("block_plan", ProductionBlock(p, alpha)) + for name, value in (("s", s), ("t", t), ("d", d), ("alpha", alpha), ("p", p), + ("radius", self.block_plan.radius), ("P", self.block_plan.P), + ("sigma", self.block_plan.sigma), ("u", self.block_plan.u)): + self.keep(name, Count(value)) + with owned(Count.product(p, p)) as (square,): + self.keep("halo", Count.product(square, p)) + self.halo.multiply_small(4096) + with owned(Count.product(square, square)) as (fourth,): + with owned(Count.product(fourth, fourth)) as (eighth,): + with owned(ceiling_log(eighth)) as (logarithm,): + self.keep("L", power_two(logarithm)) + with owned(ceiling_log(t)) as (logarithm,): + with owned(power_two(logarithm)) as (dyadic,): + need(compare(dyadic, t) == 0, "production expansion dyadic target") + self.keep("target_start", Signed(Count.product(tile_index, self.L))) + self.keep("target_length", Count(self.L)) + self.keep("input_capacity", Count(self.L)) + with owned(advance(self.target_start, self.L)) as (end,): + need(not end.negative and compare(end.magnitude, t) <= 0, + "production complete target tile") + core_start = source_index(s, t, self.target_start) + try: + core_end = source_index(s, t, end) + except BaseException: + core_start.close() + raise + with owned(core_start, core_end) as (core_start, core_end): + with owned(difference(core_end, core_start), Signed(Count(self.halo))) as (core_length, halo): + need(not core_length.negative and not core_length.magnitude.zero(), "positive source core") + with owned(difference(self.target_start, halo), Count(self.halo)) as (halo_start, halo_length): + halo_length.multiply_small(2) + halo_length.add(self.L) + self.keep("inverse_plan", ProductionPlan(s, t, d, self.u, p, halo_start, + halo_length, core_start, core_length.magnitude)) + self.keep("source_start", self.inverse_plan.start.copy()) + self.keep("source_count", Count(self.inverse_plan.m)) + with owned(Count(self.L), Count(self.halo), Count(self.inverse_plan.w)) as (covered, halos, bands): + halos.multiply_small(2) + bands.multiply_small(2) + covered.add(halos) + covered.add(bands) + covered.multiply_small(2) + need(compare(s, covered) > 0, "production global period exceeds twice the padded halo") + with owned(Count.product(d, self.halo)) as (required,): + required.multiply_small(64) + need(compare(self.L, required) >= 0, "production source-padding headroom") + need(compare(self.P, self.inverse_plan.P) == 0, "paired Gaussian work grid") + self.keep("version", Count.literal(2)) + self.keep("blocks", quotient(self.L, self.radius, True)) + compile_tile(self) + except BaseException: + self.close() + raise diff --git a/research/physical-gaussian-expansion/expansion_tile_prepare.py b/research/physical-gaussian-expansion/expansion_tile_prepare.py new file mode 100644 index 0000000000..501fca979c --- /dev/null +++ b/research/physical-gaussian-expansion/expansion_tile_prepare.py @@ -0,0 +1,164 @@ +"""Sequential catalogue of complete real-Gaussian expansion block pages.""" +from expansion_sources import active +active() +from cardinal_counts import Count +from encoded import compare, subtract +from expansion_prepare import Preparation +from expansion_records import BlockHeader, same_integer +from expansion_tile_geometry import ProductionTile, TileGeometry, TileHeader, advance +from gaussian_rationals import write_count +from records import Tape +from ring_words import Signed, difference +from router_metadata import indices, owned +from scratch_pool import ScratchPool +from session import (OuterScope, borrowed_count, erase_marker, initialize_marker, + sealed_constructors) +from words import BLANK, RIGHT, need + + +def copy_cells(source, target, length): + for _ in indices(length): + value = source.read() + need(value in (0, 1, RIGHT), "complete expansion catalogue cell") + target.append(value) + source.step(1) + + +def measure_cells(source): + result = Count.literal(0) + source.rewind() + try: + while source.read() != BLANK: + need(source.read() in (0, 1, RIGHT), "canonical expansion catalogue alphabet") + result.increment() + source.step(1) + return result + except BaseException: + result.close() + raise + finally: + source.rewind() + + +class TilePreparation: + def __init__(self): + self.block = Preparation() + self.pool = ScratchPool(256) + self.scope = OuterScope(self.pool) + self.parameter_pool = ScratchPool(256) + self.parameter_scope = OuterScope(self.parameter_pool) + self.output = Tape() + self.closed = True + self.calls = 0 + self.last_counts = None + for tape in (*self.pool.tapes, *self.parameter_pool.tapes, self.output): + erase_marker(tape) + + def idle(self): + need(self.block.closed and not self.block.pool.live and not self.block.parameter_pool.live + and self.block.output.read() == BLANK, "consumed expansion catalogue child page") + + def run_explicit(self, s, t, alpha, radius, source_start, source_count, input_capacity, + target_start, target_length, p, precision, sigma): + values = tuple(borrowed_count(value) for value in + (s, t, alpha, radius, source_count, input_capacity, target_length, p, precision, sigma)) + s, t, alpha, radius, source_count, input_capacity, target_length, p, precision, sigma = values + source_start = Signed(borrowed_count(source_start.magnitude), source_start.negative) + target_start = Signed(borrowed_count(target_start.magnitude), target_start.negative) + parameters = (*values, source_start.magnitude, target_start.magnitude) + need(all(value.tape not in (*self.pool.tapes, self.output) for value in parameters), + "external expansion tile preparation parameters") + need(self.closed and not self.pool.live and self.output.read() == BLANK, + "consume prior expansion tile page") + self.idle() + for tape in (*self.pool.tapes, self.output): + initialize_marker(tape) + self.closed = False + self.last_counts = dict(blocks=0, exponentials=0, exact_unit_factors=0, + reflected_kernel_factors=0, taylor_steps=0, early_stops=0) + plan = next_target = emitted = None + succeeded = False + try: + with sealed_constructors(), self.scope.installed(): + try: + plan = TileGeometry(s, t, alpha, radius, source_start, source_count, input_capacity, + target_start, target_length, p, precision, sigma) + TileHeader.dump(plan, self.output) + next_target, emitted = target_start.copy(), Count.literal(0) + for _ in indices(plan.blocks): + with owned(subtract(target_length, emitted)) as (remaining,): + with owned(Count(radius if compare(radius, remaining) <= 0 else remaining)) as (length,): + with self.scope.suspended(): + returned = self.block.run_explicit(s, t, alpha, radius, next_target, + length, p, precision, sigma) + try: + with owned(BlockHeader(returned)) as (header,): + need(same_integer(header.target_start, next_target) + and compare(header.m, length) == 0, "complete ordered expansion catalogue block") + with owned(difference(header.source_start, source_start)) as (offset,): + need(not offset.negative, "expansion block begins within source halo") + with owned(Count(offset.magnitude)) as (end,): + end.add(header.n) + need(compare(end, source_count) <= 0, + "expansion block ends within source halo") + with owned(measure_cells(returned)) as (cells,): + write_count(cells, self.output) + copy_cells(returned, self.output, cells) + need(returned.read() == BLANK, "complete catalogue page copy") + self.last_counts["blocks"] += 1 + for key, value in self.block.last_counts.items(): + self.last_counts[key] += value + finally: + returned.clear() + erase_marker(returned) + emitted.add(length) + following = advance(next_target, length) + next_target.close() + next_target = following + need(compare(emitted, target_length) == 0, "complete expansion catalogue target extent") + self.idle() + self.calls += 1 + succeeded = True + finally: + for value in (next_target, emitted, plan): + if value is not None: + value.close() + for value in parameters: + value.rewind() + self.pool.empty() + finally: + for tape in self.pool.tapes: + erase_marker(tape) + self.output.rewind() + if not succeeded: + self.output.clear() + erase_marker(self.output) + self.closed = True + return self.output + + def run_production(self, s, t, d, alpha, p, tile_index): + need(self.closed and not self.pool.live and not self.parameter_pool.live and self.output.read() == BLANK, + "consume prior production expansion tile page") + parameters = tuple(borrowed_count(value) for value in (s, t, d, alpha, p, tile_index)) + need(all(value.tape not in (*self.pool.tapes, *self.parameter_pool.tapes, self.output) + for value in parameters), "external production expansion tile parameters") + for tape in self.parameter_pool.tapes: + initialize_marker(tape) + plan = None + try: + with sealed_constructors(), self.parameter_scope.installed(): + try: + plan = ProductionTile(*parameters) + with self.parameter_scope.suspended(): + return self.run_explicit(plan.s, plan.t, plan.alpha, plan.radius, + plan.source_start, plan.source_count, plan.input_capacity, plan.target_start, + plan.target_length, plan.p, plan.P, plan.sigma) + finally: + if plan is not None: + plan.close() + for value in parameters: + value.rewind() + self.parameter_pool.empty() + finally: + for tape in self.parameter_pool.tapes: + erase_marker(tape) diff --git a/research/physical-gaussian-expansion/inputs/core-export.tar.gz b/research/physical-gaussian-expansion/inputs/core-export.tar.gz new file mode 100644 index 0000000000..9719f1e08e Binary files /dev/null and b/research/physical-gaussian-expansion/inputs/core-export.tar.gz differ diff --git a/research/physical-gaussian-expansion/inputs/rounding-source.tar b/research/physical-gaussian-expansion/inputs/rounding-source.tar new file mode 100644 index 0000000000..72db8f8da3 Binary files /dev/null and b/research/physical-gaussian-expansion/inputs/rounding-source.tar differ diff --git a/research/physical-gaussian-expansion/verify.py b/research/physical-gaussian-expansion/verify.py new file mode 100644 index 0000000000..58f11e3565 --- /dev/null +++ b/research/physical-gaussian-expansion/verify.py @@ -0,0 +1,28 @@ +"""Cold entry point for physical Gaussian expansion checks.""" +import argparse +import importlib +import json +from pathlib import Path +import sys + +sys.path.insert(0, str(Path(__file__).resolve().parent)) +from expansion_sources import opened, require + + +def main(): + parser = argparse.ArgumentParser() + parser.add_argument("--suite", required=True) + args = parser.parse_args() + require(args.suite in ("sources", "budget", "geometry", "block", "tile_overlap", "tile_wrap", "ownership", "production"), "declared expansion check") + with opened(): + name = "expansion_check_tile" if args.suite.startswith("tile_") else "expansion_check_" + args.suite + module = importlib.import_module(name) + require(Path(module.__file__).resolve() == Path(__file__).resolve().parent / (name + ".py"), + "package-local expansion checker") + evidence = module.run(args.suite.removeprefix("tile_")) if args.suite.startswith("tile_") else module.run() + print(json.dumps({"status": "PASS_PHYSICAL_GAUSSIAN_EXPANSION", "suite": args.suite, + "evidence": evidence}, sort_keys=True)) + + +if __name__ == "__main__": + main()