replace old example and update performance docs

py: format using black
poetry: add flake8-black
flake: add buildInputs for flake8-black

README.md

@@ -29,12 +29,7 @@ nix develop
 
 - RAM usage and execution time estimate for simulation **ONLY**
 
-    1. 100,000,000 time steps: 6GiB RAM and 12 seconds 
+    - 500,000,000 time steps: 8GiB RAM and 55 seconds
-    
-    2. 200,000,000 time steps: 11GiB RAM and 20 seconds
-
-    3. 300,000,000 time steps: 16GiB RAM and 35 seconds
-
 
 <br>
 <details><summary><b>WRPLL formulas</b></summary>

flake.lock

@@ -66,11 +66,11 @@
         "treefmt-nix": "treefmt-nix"
       },
       "locked": {
-        "lastModified": 1702365004,
+        "lastModified": 1705060653,
-        "narHash": "sha256-IRFvmyP1uk1hchRVxaXTqu6YoZCvMM/NVtUf2hD2Tag=",
+        "narHash": "sha256-puYyylgrBS4AFAHeyVRTjTUVD8DZdecJfymWJe7H438=",
         "owner": "nix-community",
         "repo": "poetry2nix",
-        "rev": "c12ac880114d52a3cad5fa02b00f2e2090e89982",
+        "rev": "e0b44e9e2d3aa855d1dd77b06f067cd0e0c3860d",
         "type": "github"
       },
       "original": {

flake.nix

@@ -32,8 +32,13 @@
             old: {
               buildInputs = (old.buildInputs or [ ]) ++ [ super.poetry-core ];
             });
-          scipy = super.scipy.override {
+          flake8-black = super.flake8-black.overridePythonAttrs (
-            # fail to build
+            # ModuleNotFoundError: No module named 'setuptools'
+            old: {
+              buildInputs = (old.buildInputs or [ ]) ++ [ super.setuptools ];
+            });
+          jupytext = super.jupytext.override {
+            # bypass error: metadata-generation-failed
             preferWheel = true;
           };
         });

poetry.lock

@@ -204,6 +204,50 @@ soupsieve = ">1.2"
 html5lib = ["html5lib"]
 lxml = ["lxml"]
 
+[[package]]
+name = "black"
+version = "23.12.1"
+description = "The uncompromising code formatter."
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "black-23.12.1-cp310-cp310-macosx_10_9_x86_64.whl", hash = "sha256:e0aaf6041986767a5e0ce663c7a2f0e9eaf21e6ff87a5f95cbf3675bfd4c41d2"},
+    {file = "black-23.12.1-cp310-cp310-macosx_11_0_arm64.whl", hash = "sha256:c88b3711d12905b74206227109272673edce0cb29f27e1385f33b0163c414bba"},
+    {file = "black-23.12.1-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:a920b569dc6b3472513ba6ddea21f440d4b4c699494d2e972a1753cdc25df7b0"},
+    {file = "black-23.12.1-cp310-cp310-win_amd64.whl", hash = "sha256:3fa4be75ef2a6b96ea8d92b1587dd8cb3a35c7e3d51f0738ced0781c3aa3a5a3"},
+    {file = "black-23.12.1-cp311-cp311-macosx_10_9_x86_64.whl", hash = "sha256:8d4df77958a622f9b5a4c96edb4b8c0034f8434032ab11077ec6c56ae9f384ba"},
+    {file = "black-23.12.1-cp311-cp311-macosx_11_0_arm64.whl", hash = "sha256:602cfb1196dc692424c70b6507593a2b29aac0547c1be9a1d1365f0d964c353b"},
+    {file = "black-23.12.1-cp311-cp311-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:9c4352800f14be5b4864016882cdba10755bd50805c95f728011bcb47a4afd59"},
+    {file = "black-23.12.1-cp311-cp311-win_amd64.whl", hash = "sha256:0808494f2b2df923ffc5723ed3c7b096bd76341f6213989759287611e9837d50"},
+    {file = "black-23.12.1-cp312-cp312-macosx_10_9_x86_64.whl", hash = "sha256:25e57fd232a6d6ff3f4478a6fd0580838e47c93c83eaf1ccc92d4faf27112c4e"},
+    {file = "black-23.12.1-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:2d9e13db441c509a3763a7a3d9a49ccc1b4e974a47be4e08ade2a228876500ec"},
+    {file = "black-23.12.1-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:6d1bd9c210f8b109b1762ec9fd36592fdd528485aadb3f5849b2740ef17e674e"},
+    {file = "black-23.12.1-cp312-cp312-win_amd64.whl", hash = "sha256:ae76c22bde5cbb6bfd211ec343ded2163bba7883c7bc77f6b756a1049436fbb9"},
+    {file = "black-23.12.1-cp38-cp38-macosx_10_9_x86_64.whl", hash = "sha256:1fa88a0f74e50e4487477bc0bb900c6781dbddfdfa32691e780bf854c3b4a47f"},
+    {file = "black-23.12.1-cp38-cp38-macosx_11_0_arm64.whl", hash = "sha256:a4d6a9668e45ad99d2f8ec70d5c8c04ef4f32f648ef39048d010b0689832ec6d"},
+    {file = "black-23.12.1-cp38-cp38-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:b18fb2ae6c4bb63eebe5be6bd869ba2f14fd0259bda7d18a46b764d8fb86298a"},
+    {file = "black-23.12.1-cp38-cp38-win_amd64.whl", hash = "sha256:c04b6d9d20e9c13f43eee8ea87d44156b8505ca8a3c878773f68b4e4812a421e"},
+    {file = "black-23.12.1-cp39-cp39-macosx_10_9_x86_64.whl", hash = "sha256:3e1b38b3135fd4c025c28c55ddfc236b05af657828a8a6abe5deec419a0b7055"},
+    {file = "black-23.12.1-cp39-cp39-macosx_11_0_arm64.whl", hash = "sha256:4f0031eaa7b921db76decd73636ef3a12c942ed367d8c3841a0739412b260a54"},
+    {file = "black-23.12.1-cp39-cp39-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:97e56155c6b737854e60a9ab1c598ff2533d57e7506d97af5481141671abf3ea"},
+    {file = "black-23.12.1-cp39-cp39-win_amd64.whl", hash = "sha256:dd15245c8b68fe2b6bd0f32c1556509d11bb33aec9b5d0866dd8e2ed3dba09c2"},
+    {file = "black-23.12.1-py3-none-any.whl", hash = "sha256:78baad24af0f033958cad29731e27363183e140962595def56423e626f4bee3e"},
+    {file = "black-23.12.1.tar.gz", hash = "sha256:4ce3ef14ebe8d9509188014d96af1c456a910d5b5cbf434a09fef7e024b3d0d5"},
+]
+
+[package.dependencies]
+click = ">=8.0.0"
+mypy-extensions = ">=0.4.3"
+packaging = ">=22.0"
+pathspec = ">=0.9.0"
+platformdirs = ">=2"
+
+[package.extras]
+colorama = ["colorama (>=0.4.3)"]
+d = ["aiohttp (>=3.7.4)", "aiohttp (>=3.7.4,!=3.9.0)"]
+jupyter = ["ipython (>=7.8.0)", "tokenize-rt (>=3.2.0)"]
+uvloop = ["uvloop (>=0.15.2)"]
+
 [[package]]
 name = "bleach"
 version = "6.1.0"
@@ -593,6 +637,40 @@ files = [
 [package.extras]
 devel = ["colorama", "json-spec", "jsonschema", "pylint", "pytest", "pytest-benchmark", "pytest-cache", "validictory"]
 
+[[package]]
+name = "flake8"
+version = "7.0.0"
+description = "the modular source code checker: pep8 pyflakes and co"
+optional = false
+python-versions = ">=3.8.1"
+files = [
+    {file = "flake8-7.0.0-py2.py3-none-any.whl", hash = "sha256:a6dfbb75e03252917f2473ea9653f7cd799c3064e54d4c8140044c5c065f53c3"},
+    {file = "flake8-7.0.0.tar.gz", hash = "sha256:33f96621059e65eec474169085dc92bf26e7b2d47366b70be2f67ab80dc25132"},
+]
+
+[package.dependencies]
+mccabe = ">=0.7.0,<0.8.0"
+pycodestyle = ">=2.11.0,<2.12.0"
+pyflakes = ">=3.2.0,<3.3.0"
+
+[[package]]
+name = "flake8-black"
+version = "0.3.6"
+description = "flake8 plugin to call black as a code style validator"
+optional = false
+python-versions = ">=3.7"
+files = [
+    {file = "flake8-black-0.3.6.tar.gz", hash = "sha256:0dfbca3274777792a5bcb2af887a4cad72c72d0e86c94e08e3a3de151bb41c34"},
+    {file = "flake8_black-0.3.6-py3-none-any.whl", hash = "sha256:fe8ea2eca98d8a504f22040d9117347f6b367458366952862ac3586e7d4eeaca"},
+]
+
+[package.dependencies]
+black = ">=22.1.0"
+flake8 = ">=3"
+
+[package.extras]
+develop = ["build", "twine"]
+
 [[package]]
 name = "flask"
 version = "3.0.0"
@@ -1131,6 +1209,35 @@ files = [
     {file = "jupyterlab_widgets-3.0.9.tar.gz", hash = "sha256:6005a4e974c7beee84060fdfba341a3218495046de8ae3ec64888e5fe19fdb4c"},
 ]
 
+[[package]]
+name = "jupytext"
+version = "1.16.1"
+description = "Jupyter notebooks as Markdown documents, Julia, Python or R scripts"
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "jupytext-1.16.1-py3-none-any.whl", hash = "sha256:796ec4f68ada663569e5d38d4ef03738a01284bfe21c943c485bc36433898bd0"},
+    {file = "jupytext-1.16.1.tar.gz", hash = "sha256:68c7b68685e870e80e60fda8286fbd6269e9c74dc1df4316df6fe46eabc94c99"},
+]
+
+[package.dependencies]
+markdown-it-py = ">=1.0"
+mdit-py-plugins = "*"
+nbformat = "*"
+packaging = "*"
+pyyaml = "*"
+toml = "*"
+
+[package.extras]
+dev = ["jupytext[test-cov,test-external]"]
+docs = ["myst-parser", "sphinx", "sphinx-copybutton", "sphinx-rtd-theme"]
+test = ["pytest", "pytest-randomly", "pytest-xdist"]
+test-cov = ["jupytext[test-integration]", "pytest-cov (>=2.6.1)"]
+test-external = ["autopep8", "black", "flake8", "gitpython", "isort", "jupyter-fs (<0.4.0)", "jupytext[test-integration]", "pre-commit", "sphinx-gallery (<0.8)"]
+test-functional = ["jupytext[test]"]
+test-integration = ["ipykernel", "jupyter-server (!=2.11)", "jupytext[test-functional]", "nbconvert"]
+test-ui = ["calysto-bash"]
+
 [[package]]
 name = "llvmlite"
 version = "0.41.1"
@@ -1164,6 +1271,30 @@ files = [
     {file = "llvmlite-0.41.1.tar.gz", hash = "sha256:f19f767a018e6ec89608e1f6b13348fa2fcde657151137cb64e56d48598a92db"},
 ]
 
+[[package]]
+name = "markdown-it-py"
+version = "3.0.0"
+description = "Python port of markdown-it. Markdown parsing, done right!"
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "markdown-it-py-3.0.0.tar.gz", hash = "sha256:e3f60a94fa066dc52ec76661e37c851cb232d92f9886b15cb560aaada2df8feb"},
+    {file = "markdown_it_py-3.0.0-py3-none-any.whl", hash = "sha256:355216845c60bd96232cd8d8c40e8f9765cc86f46880e43a8fd22dc1a1a8cab1"},
+]
+
+[package.dependencies]
+mdurl = ">=0.1,<1.0"
+
+[package.extras]
+benchmarking = ["psutil", "pytest", "pytest-benchmark"]
+code-style = ["pre-commit (>=3.0,<4.0)"]
+compare = ["commonmark (>=0.9,<1.0)", "markdown (>=3.4,<4.0)", "mistletoe (>=1.0,<2.0)", "mistune (>=2.0,<3.0)", "panflute (>=2.3,<3.0)"]
+linkify = ["linkify-it-py (>=1,<3)"]
+plugins = ["mdit-py-plugins"]
+profiling = ["gprof2dot"]
+rtd = ["jupyter_sphinx", "mdit-py-plugins", "myst-parser", "pyyaml", "sphinx", "sphinx-copybutton", "sphinx-design", "sphinx_book_theme"]
+testing = ["coverage", "pytest", "pytest-cov", "pytest-regressions"]
+
 [[package]]
 name = "markupsafe"
 version = "2.1.3"
@@ -1247,6 +1378,47 @@ files = [
 [package.dependencies]
 traitlets = "*"
 
+[[package]]
+name = "mccabe"
+version = "0.7.0"
+description = "McCabe checker, plugin for flake8"
+optional = false
+python-versions = ">=3.6"
+files = [
+    {file = "mccabe-0.7.0-py2.py3-none-any.whl", hash = "sha256:6c2d30ab6be0e4a46919781807b4f0d834ebdd6c6e3dca0bda5a15f863427b6e"},
+    {file = "mccabe-0.7.0.tar.gz", hash = "sha256:348e0240c33b60bbdf4e523192ef919f28cb2c3d7d5c7794f74009290f236325"},
+]
+
+[[package]]
+name = "mdit-py-plugins"
+version = "0.4.0"
+description = "Collection of plugins for markdown-it-py"
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "mdit_py_plugins-0.4.0-py3-none-any.whl", hash = "sha256:b51b3bb70691f57f974e257e367107857a93b36f322a9e6d44ca5bf28ec2def9"},
+    {file = "mdit_py_plugins-0.4.0.tar.gz", hash = "sha256:d8ab27e9aed6c38aa716819fedfde15ca275715955f8a185a8e1cf90fb1d2c1b"},
+]
+
+[package.dependencies]
+markdown-it-py = ">=1.0.0,<4.0.0"
+
+[package.extras]
+code-style = ["pre-commit"]
+rtd = ["myst-parser", "sphinx-book-theme"]
+testing = ["coverage", "pytest", "pytest-cov", "pytest-regressions"]
+
+[[package]]
+name = "mdurl"
+version = "0.1.2"
+description = "Markdown URL utilities"
+optional = false
+python-versions = ">=3.7"
+files = [
+    {file = "mdurl-0.1.2-py3-none-any.whl", hash = "sha256:84008a41e51615a49fc9966191ff91509e3c40b939176e643fd50a5c2196b8f8"},
+    {file = "mdurl-0.1.2.tar.gz", hash = "sha256:bb413d29f5eea38f31dd4754dd7377d4465116fb207585f97bf925588687c1ba"},
+]
+
 [[package]]
 name = "mistune"
 version = "3.0.2"
@@ -1258,6 +1430,17 @@ files = [
     {file = "mistune-3.0.2.tar.gz", hash = "sha256:fc7f93ded930c92394ef2cb6f04a8aabab4117a91449e72dcc8dfa646a508be8"},
 ]
 
+[[package]]
+name = "mypy-extensions"
+version = "1.0.0"
+description = "Type system extensions for programs checked with the mypy type checker."
+optional = false
+python-versions = ">=3.5"
+files = [
+    {file = "mypy_extensions-1.0.0-py3-none-any.whl", hash = "sha256:4392f6c0eb8a5668a69e23d168ffa70f0be9ccfd32b5cc2d26a34ae5b844552d"},
+    {file = "mypy_extensions-1.0.0.tar.gz", hash = "sha256:75dbf8955dc00442a438fc4d0666508a9a97b6bd41aa2f0ffe9d2f2725af0782"},
+]
+
 [[package]]
 name = "nbclient"
 version = "0.9.0"
@@ -1642,6 +1825,17 @@ files = [
 qa = ["flake8 (==3.8.3)", "mypy (==0.782)"]
 testing = ["docopt", "pytest (<6.0.0)"]
 
+[[package]]
+name = "pathspec"
+version = "0.12.1"
+description = "Utility library for gitignore style pattern matching of file paths."
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "pathspec-0.12.1-py3-none-any.whl", hash = "sha256:a0d503e138a4c123b27490a4f7beda6a01c6f288df0e4a8b79c7eb0dc7b4cc08"},
+    {file = "pathspec-0.12.1.tar.gz", hash = "sha256:a482d51503a1ab33b1c67a6c3813a26953dbdc71c31dacaef9a838c4e29f5712"},
+]
+
 [[package]]
 name = "pexpect"
 version = "4.9.0"
@@ -1790,6 +1984,17 @@ files = [
 [package.extras]
 tests = ["pytest"]
 
+[[package]]
+name = "pycodestyle"
+version = "2.11.1"
+description = "Python style guide checker"
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "pycodestyle-2.11.1-py2.py3-none-any.whl", hash = "sha256:44fe31000b2d866f2e41841b18528a505fbd7fef9017b04eff4e2648a0fadc67"},
+    {file = "pycodestyle-2.11.1.tar.gz", hash = "sha256:41ba0e7afc9752dfb53ced5489e89f8186be00e599e712660695b7a75ff2663f"},
+]
+
 [[package]]
 name = "pycparser"
 version = "2.21"
@@ -1801,6 +2006,17 @@ files = [
     {file = "pycparser-2.21.tar.gz", hash = "sha256:e644fdec12f7872f86c58ff790da456218b10f863970249516d60a5eaca77206"},
 ]
 
+[[package]]
+name = "pyflakes"
+version = "3.2.0"
+description = "passive checker of Python programs"
+optional = false
+python-versions = ">=3.8"
+files = [
+    {file = "pyflakes-3.2.0-py2.py3-none-any.whl", hash = "sha256:84b5be138a2dfbb40689ca07e2152deb896a65c3a3e24c251c5c62489568074a"},
+    {file = "pyflakes-3.2.0.tar.gz", hash = "sha256:1c61603ff154621fb2a9172037d84dca3500def8c8b630657d1701f026f8af3f"},
+]
+
 [[package]]
 name = "pygments"
 version = "2.17.2"
@@ -2458,6 +2674,17 @@ webencodings = ">=0.4"
 doc = ["sphinx", "sphinx_rtd_theme"]
 test = ["flake8", "isort", "pytest"]
 
+[[package]]
+name = "toml"
+version = "0.10.2"
+description = "Python Library for Tom's Obvious, Minimal Language"
+optional = false
+python-versions = ">=2.6, !=3.0.*, !=3.1.*, !=3.2.*"
+files = [
+    {file = "toml-0.10.2-py2.py3-none-any.whl", hash = "sha256:806143ae5bfb6a3c6e736a764057db0e6a0e05e338b5630894a5f779cabb4f9b"},
+    {file = "toml-0.10.2.tar.gz", hash = "sha256:b3bda1d108d5dd99f4a20d24d9c348e91c4db7ab1b749200bded2f839ccbe68f"},
+]
+
 [[package]]
 name = "tornado"
 version = "6.4"
@@ -2739,4 +2966,4 @@ testing = ["big-O", "jaraco.functools", "jaraco.itertools", "more-itertools", "p
 [metadata]
 lock-version = "2.0"
 python-versions = "^3.11"
-content-hash = "121a69b8179bfa06910de0f2fd41aa135a5a4822061538b5c021ba0641e5fa1a"
+content-hash = "5c5b288075581046408471099dc252467fad16a37215d5ad9f1ccb1b5f079821"

pyproject.toml

@@ -16,6 +16,8 @@ numba = "^0.58.1"
 pandas = "^2.1.4"
 plotly = "^5.18.0"
 plotly-resampler = "^0.9.1"
+flake8-black = "^0.3.6"
+jupytext = "^1.16.1"
 
 
 [build-system]

src/WRPLL_timesim_example.py

@@ -0,0 +1,68 @@
+# %%
+from plotly_resampler import FigureWidgetResampler
+from plotly.subplots import make_subplots
+import plotly.graph_objects as go
+import numpy as np
+
+from wrpll_simulation.config import Timesim_Config, PI_Config
+from wrpll_simulation.timesim import WRPLL_Timesim
+
+# %%
+KP, KI = 6, 2
+sim_config = Timesim_Config(
+    timestep_size=1e-11,
+    sim_length=500_000_000,
+    helper_PI=PI_Config(KP=KP, KI=KI),
+    main_PI=PI_Config(KP=KP, KI=KI),
+    has_jitter=False,
+    step_input_time=10e-6,
+    step_frequency=6.25,
+    step_phase=0,
+)
+sim = WRPLL_Timesim(sim_config, np.random.default_rng(1))
+
+fig = FigureWidgetResampler(make_subplots(rows=2))
+# ROW 1
+fig.add_trace(
+    go.Scattergl(name="main - gtx freq diff"),
+    hf_x=sim.time,
+    hf_y=sim.freq_diff,
+    row=1,
+    col=1,
+)
+
+# ROW 2
+fig.add_trace(
+    go.Scattergl(name="main - gtx phase diff"),
+    hf_x=sim.time,
+    hf_y=sim.phase_diff,
+    row=2,
+    col=1,
+)
+
+
+title = (
+    f'step freq = {sim_config.step_frequency}Hz, '
+    f'step phase = {sim_config.step_phase}° | '
+    f'KP = {KP}, KI = {KI}'
+)
+
+fig.update_layout(
+    xaxis2=dict(title="Time (sec)", exponentformat="SI"),
+    yaxis1=dict(title="Frequency (Hz)"),
+    yaxis2=dict(title="Phase (degree)"),
+    height=950,
+    showlegend=True,
+    title_text=title,
+    legend=dict(
+        orientation="h",
+        yanchor="bottom",
+        y=1.02,
+        xanchor="right",
+        x=1,
+    ),
+)
+
+fig
+
+# %%

src/both_PLL_example.ipynb

@@ -1,158 +0,0 @@
-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Both PLL mode example\n",
-    "\n",
-    "Time domain simulation with helper and main PLL mode. \n",
-    "\n",
-    "- Period error (**Helper PLL**)\n",
-    "    - $\\Delta{period} = N - (tag_{gtx}[n] - tag_{gtx}[n-1]) $\n",
-    "    - where ideally $f_{helper} = \\dfrac{f_{in} * (N-1)}{N}$ \n",
-    "\n",
-    "- Phase error (**Main PLL**)\n",
-    "    - $\\text{Let } \\Delta tag[n] = tag_{main}[n] - tag_{gtx}[n] \\text{ mod N}$\n",
-    "\n",
-    "    - $\\Delta\\phi[n] = \\begin{cases}\n",
-    "                \\Delta tag[n] - N, & \\text{if } \\Delta tag > N/2  \\\\\n",
-    "                \\Delta tag[n], & otherwise\n",
-    "        \\end{cases} \\quad$\n",
-    " \n",
-    "\n",
-    "- ADPLL PID (common for main and helper PLL)\n",
-    "    - $P[n] = err[n] * K_P$\n",
-    "    - $I[n] = I[n-1] + err[n] * K_I$\n",
-    "    - $D[n] = (err[n] - err[n-1]) * K_D$\n",
-    "    - $adpll[n] = \\text{base adpll} + P[i] + I[n] + D[n] $\n",
-    "    - where $\\text{base adpll}$ is constant and obtain from frequency counter in HW"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "from plotly_resampler import FigureResampler, FigureWidgetResampler\n",
-    "from plotly.subplots import make_subplots\n",
-    "import plotly.graph_objects as go\n",
-    "import numpy as np\n",
-    "from wrpll_simulation.wrpll import WRPLL_simulator"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# settings\n",
-    "timestep = 1e-10\n",
-    "total_steps = 200_000_000\n",
-    "sim_mode = \"both\"\n",
-    "adpll_period = 200e-6     # in seconds, the period that pll will trigger, (minimum > total DCXO frequency change delay)\n",
-    "start_up_delay = 100e-6   # in seconds, the frequency adjustment is DISABLE until time > start_up_delay\n",
-    "\n",
-    "gtx_freq = 125_001_519\n",
-    "\n",
-    "\n",
-    "helper_filter = {\n",
-    "    \"KP\": 2,\n",
-    "    \"KI\": 0.5,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "main_filter = {\n",
-    "    \"KP\": 12,\n",
-    "    \"KI\": 0,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "\n",
-    "# simulation have RNG for\n",
-    "# - gtx, main and helper jitter\n",
-    "# - starting phase for main and helper\n",
-    "# - base_adpll error\n",
-    "\n",
-    "wrpll_sim = WRPLL_simulator(\n",
-    "    timestep=timestep,\n",
-    "    total_steps=total_steps,\n",
-    "    sim_mode=sim_mode,\n",
-    "    helper_filter=helper_filter,\n",
-    "    main_filter=main_filter,\n",
-    "    gtx_freq=gtx_freq,\n",
-    "    adpll_write_period=adpll_period,\n",
-    "    start_up_delay=start_up_delay,\n",
-    ")\n",
-    "wrpll_sim.run()"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# faster than pyplot with resampling feature\n",
-    "# see https://github.com/predict-idlab/plotly-resampler\n",
-    "\n",
-    "fig = FigureWidgetResampler(make_subplots(rows=4, shared_xaxes=True))\n",
-    "fig.add_trace(go.Scattergl(name='phase error'),  hf_x=wrpll_sim.time, hf_y=wrpll_sim.phase_err, row=1, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='freq error (ppm)'),  hf_x=wrpll_sim.time, hf_y=(\n",
-    "    wrpll_sim.mainfreq-gtx_freq) * (1e6/gtx_freq), row=2, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='period error'),  hf_x=wrpll_sim.time, hf_y=wrpll_sim.period_err, row=3, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='gtx'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.gtx+1, row=4, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='main'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.main, row=4, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='helper'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.helper-1, row=4, col=1)\n",
-    "\n",
-    "\n",
-    "fig.update_layout(\n",
-    "    xaxis4=dict(title=\"time (sec)\"),\n",
-    "    yaxis1=dict(title=\"phase error\"),\n",
-    "    yaxis2=dict(title=\"freq error (ppm)\"),\n",
-    "    yaxis3=dict(title=\"beating period error\"),\n",
-    "    yaxis4=dict(title=\"Signal\"),\n",
-    "\n",
-    "    height=1000,\n",
-    "    showlegend=True,\n",
-    "    title_text=\"PLL example\",\n",
-    "    legend=dict(\n",
-    "        orientation=\"h\",\n",
-    "        yanchor=\"bottom\",\n",
-    "        y=1.02,\n",
-    "        xanchor=\"right\",\n",
-    "        x=1,\n",
-    "    ),\n",
-    ")\n",
-    "\n",
-    "fig"
-   ]
-  }
- ],
- "metadata": {
-  "kernelspec": {
-   "display_name": "Python 3 (ipykernel)",
-   "language": "python",
-   "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.11.6"
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 4
-}

src/helper_PLL_example.ipynb

@@ -1,149 +0,0 @@
-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Helper PLL mode example\n",
-    "\n",
-    "Time domain simulation with helper PLL mode\n",
-    "\n",
-    "- Period error (**Helper PLL**)\n",
-    "    - $\\Delta{period} = N - (tag_{gtx}[n] - tag_{gtx}[n-1]) $\n",
-    "    - where ideally $f_{helper} = \\dfrac{f_{in} * (N-1)}{N}$ \n",
-    "\n",
-    "- Phase error (**Main PLL**)\n",
-    "    - $\\text{Let } \\Delta tag[n] = tag_{main}[n] - tag_{gtx}[n] \\text{ mod N}$\n",
-    "\n",
-    "    - $\\Delta\\phi[n] = \\begin{cases}\n",
-    "                \\Delta tag[n] - N, & \\text{if } \\Delta tag > N/2  \\\\\n",
-    "                \\Delta tag[n], & otherwise\n",
-    "        \\end{cases} \\quad$\n",
-    " \n",
-    "\n",
-    "- ADPLL PID (common for main and helper PLL)\n",
-    "    - $P[n] = err[n] * K_P$\n",
-    "    - $I[n] = I[n-1] + err[n] * K_I$\n",
-    "    - $D[n] = (err[n] - err[n-1]) * K_D$\n",
-    "    - $adpll[n] = \\text{base adpll} + P[i] + I[n] + D[n] $\n",
-    "    - where $\\text{base adpll}$ is constant and obtain from frequency counter in HW"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "from plotly_resampler import FigureResampler, FigureWidgetResampler\n",
-    "from plotly.subplots import make_subplots\n",
-    "import plotly.graph_objects as go\n",
-    "import numpy as np\n",
-    "from wrpll_simulation.wrpll import WRPLL_simulator"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# settings\n",
-    "timestep = 1e-10\n",
-    "total_steps = 100_000_000\n",
-    "sim_mode = \"helper_pll\"\n",
-    "adpll_period = 200e-6     # in seconds, the period that pll will trigger, (minimum > total DCXO frequency change delay)\n",
-    "start_up_delay = 100e-6   # in seconds, the frequency adjustment is DISABLE until time > start_up_delay\n",
-    "\n",
-    "gtx_freq = 125_001_519\n",
-    "\n",
-    "helper_filter = {\n",
-    "    \"KP\": 2,\n",
-    "    \"KI\": 0.5,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "main_filter = {  # unused\n",
-    "    \"KP\": 12,\n",
-    "    \"KI\": 0,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "# simulation have RNG for\n",
-    "# - gtx, main and helper jitter\n",
-    "# - starting phase for main and helper\n",
-    "# - base_adpll error\n",
-    "\n",
-    "wrpll_sim = WRPLL_simulator(\n",
-    "    timestep=timestep,\n",
-    "    total_steps=total_steps,\n",
-    "    sim_mode=sim_mode,\n",
-    "    helper_filter=helper_filter,\n",
-    "    main_filter=main_filter,\n",
-    "    gtx_freq=gtx_freq,\n",
-    "    adpll_write_period=adpll_period,\n",
-    "    start_up_delay=start_up_delay,\n",
-    ")\n",
-    "wrpll_sim.run()"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# faster than pyplot with resampling feature\n",
-    "# see https://github.com/predict-idlab/plotly-resampler\n",
-    "\n",
-    "fig = FigureWidgetResampler(make_subplots(rows=2, shared_xaxes=True))\n",
-    "fig.add_trace(go.Scattergl(name='period error'),  hf_x=wrpll_sim.time, hf_y=wrpll_sim.period_err, row=1, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='gtx'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.gtx+1, row=2, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='main'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.main, row=2, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='helper'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.helper-1, row=2, col=1)\n",
-    "\n",
-    "\n",
-    "fig.update_layout(\n",
-    "    xaxis2=dict(title=\"time (sec)\"),\n",
-    "\n",
-    "    yaxis1=dict(title=\"beating period error\"),\n",
-    "    yaxis2=dict(title=\"Signal\"),\n",
-    "    height=1000,\n",
-    "    showlegend=True,\n",
-    "    title_text=\"PLL example\",\n",
-    "    legend=dict(\n",
-    "        orientation=\"h\",\n",
-    "        yanchor=\"bottom\",\n",
-    "        y=1.02,\n",
-    "        xanchor=\"right\",\n",
-    "        x=1,\n",
-    "    ),\n",
-    ")\n",
-    "\n",
-    "fig"
-   ]
-  }
- ],
- "metadata": {
-  "kernelspec": {
-   "display_name": "Python 3 (ipykernel)",
-   "language": "python",
-   "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.11.6"
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 4
-}

src/main_PLL_example.ipynb

@@ -1,164 +0,0 @@
-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "metadata": {},
-   "source": [
-    "## Main PLL mode example\n",
-    "\n",
-    "Time domain simulation with main PLL mode, and helper PLL is assumed to be locked \n",
-    "\n",
-    "(`helper_init_freq` variable need to be set)\n",
-    "\n",
-    "\n",
-    "- Period error (**Helper PLL**)\n",
-    "    - $\\Delta{period} = N - (tag_{gtx}[n] - tag_{gtx}[n-1]) $\n",
-    "    - where ideally $f_{helper} = \\dfrac{f_{in} * (N-1)}{N}$ \n",
-    "\n",
-    "- Phase error (**Main PLL**)\n",
-    "    - $\\text{Let } \\Delta tag[n] = tag_{main}[n] - tag_{gtx}[n] \\text{ mod N}$\n",
-    "\n",
-    "    - $\\Delta\\phi[n] = \\begin{cases}\n",
-    "                \\Delta tag[n] - N, & \\text{if } \\Delta tag > N/2  \\\\\n",
-    "                \\Delta tag[n], & otherwise\n",
-    "        \\end{cases} \\quad$\n",
-    " \n",
-    "\n",
-    "- ADPLL PID (common for main and helper PLL)\n",
-    "    - $P[n] = err[n] * K_P$\n",
-    "    - $I[n] = I[n-1] + err[n] * K_I$\n",
-    "    - $D[n] = (err[n] - err[n-1]) * K_D$\n",
-    "    - $adpll[n] = \\text{base adpll} + P[i] + I[n] + D[n] $\n",
-    "    - where $\\text{base adpll}$ is constant and obtain from frequency counter in HW"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "from plotly_resampler import FigureResampler, FigureWidgetResampler\n",
-    "from plotly.subplots import make_subplots\n",
-    "import plotly.graph_objects as go\n",
-    "import numpy as np\n",
-    "from wrpll_simulation.wrpll import WRPLL_simulator"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# settings\n",
-    "timestep = 1e-10\n",
-    "total_steps = 100_000_000\n",
-    "sim_mode = \"main_pll\"\n",
-    "adpll_period = 200e-6     # in seconds, the period that pll will trigger, (minimum > total DCXO frequency change delay)\n",
-    "start_up_delay = 100e-6   # in seconds, the frequency adjustment is DISABLE until time > start_up_delay\n",
-    "\n",
-    "gtx_freq = 125_001_519\n",
-    "\n",
-    "helper_init_freq = gtx_freq * (4096-1)/4096\n",
-    "\n",
-    "helper_filter = {  # unused\n",
-    "    \"KP\": 2,\n",
-    "    \"KI\": 0.5,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "main_filter = {\n",
-    "    \"KP\": 12,\n",
-    "    \"KI\": 0,\n",
-    "    \"KD\": 0,\n",
-    "}\n",
-    "\n",
-    "\n",
-    "# simulation have RNG for\n",
-    "# - gtx, main and helper jitter\n",
-    "# - starting phase for main and helper\n",
-    "# - base_adpll error\n",
-    "\n",
-    "\n",
-    "wrpll_sim = WRPLL_simulator(\n",
-    "    timestep=timestep,\n",
-    "    total_steps=total_steps,\n",
-    "    sim_mode=sim_mode,\n",
-    "    helper_filter=helper_filter,\n",
-    "    main_filter=main_filter,\n",
-    "    gtx_freq=gtx_freq,\n",
-    "    adpll_write_period=adpll_period,\n",
-    "    start_up_delay=start_up_delay,\n",
-    "    helper_init_freq=helper_init_freq\n",
-    ")\n",
-    "wrpll_sim.run()"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# faster than pyplot with resampling feature\n",
-    "# see https://github.com/predict-idlab/plotly-resampler\n",
-    "\n",
-    "fig = FigureWidgetResampler(make_subplots(rows=4, shared_xaxes=True))\n",
-    "fig.add_trace(go.Scattergl(name='phase error'),  hf_x=wrpll_sim.time, hf_y=wrpll_sim.phase_err, row=1, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='freq error (ppm)'),  hf_x=wrpll_sim.time, hf_y=(\n",
-    "    wrpll_sim.mainfreq-gtx_freq) * (1e6/gtx_freq), row=2, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='period error'),  hf_x=wrpll_sim.time, hf_y=wrpll_sim.period_err, row=3, col=1)\n",
-    "\n",
-    "fig.add_trace(go.Scattergl(name='gtx'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.gtx+1, row=4, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='main'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.main, row=4, col=1)\n",
-    "fig.add_trace(go.Scattergl(name='helper'), hf_x=wrpll_sim.time, hf_y=wrpll_sim.helper-1, row=4, col=1)\n",
-    "\n",
-    "\n",
-    "fig.update_layout(\n",
-    "    xaxis4=dict(title=\"time (sec)\"),\n",
-    "    yaxis1=dict(title=\"phase error\"),\n",
-    "    yaxis2=dict(title=\"freq error (ppm)\"),\n",
-    "    yaxis3=dict(title=\"beating period error\"),\n",
-    "    yaxis4=dict(title=\"Signal\"),\n",
-    "\n",
-    "    height=1000,\n",
-    "    showlegend=True,\n",
-    "    title_text=\"PLL example\",\n",
-    "    legend=dict(\n",
-    "        orientation=\"h\",\n",
-    "        yanchor=\"bottom\",\n",
-    "        y=1.02,\n",
-    "        xanchor=\"right\",\n",
-    "        x=1,\n",
-    "    ),\n",
-    ")\n",
-    "\n",
-    "fig"
-   ]
-  }
- ],
- "metadata": {
-  "kernelspec": {
-   "display_name": "Python 3 (ipykernel)",
-   "language": "python",
-   "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.11.6"
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 4
-}

src/wrpll_simulation/config.py

@@ -0,0 +1,89 @@
+from numba.core.types import *
+from numba.experimental import jitclass
+
+
+@jitclass
+class PI_Config(object):
+    KP: int64
+    KI: int64
+
+    def __init__(self, KP, KI):
+        self.KP = KP
+        self.KI = KI
+
+
+@jitclass
+class Timesim_Config(object):
+    timestep_size: float32
+    sim_length: int64
+    helper_PI: PI_Config
+    main_PI: PI_Config
+    has_jitter: bool
+
+    step_input_time: float64
+    step_frequency: float64
+    step_phase: float64
+
+    # preset
+    beating_period: int64
+    blind_period: int16
+    adpll_limit: int32
+
+    gtx_init_phase: float64
+    gtx_init_freq: float64
+    gtx_jitter: float64
+
+    helper_init_phase: float64
+    main_init_phase: float64
+    helper_init_freq: float64
+    main_init_freq: float64
+    dcxo_jitter: float64
+
+    irq_delay: float64
+    i2c_comm_delay: float64
+    dcxo_settling_delay: float64
+
+    # jitclass does not support **kwargs
+    def __init__(
+        self,
+        timestep_size: float32,
+        sim_length: int64,
+        helper_PI: PI_Config,
+        main_PI: PI_Config,
+        has_jitter: bool,
+        step_input_time: float64,
+        step_frequency: float64,
+        step_phase: float64,
+    ):
+        # jitter << timestep_size
+        # otherwise simulation will add negative phase to Phase_Accumlator
+        self.timestep_size = timestep_size
+        self.sim_length = sim_length
+        self.helper_PI = helper_PI
+        self.main_PI = main_PI
+        self.has_jitter = has_jitter
+
+        # step input
+        self.step_input_time = step_input_time
+        self.step_frequency = step_frequency
+        self.step_phase = step_phase
+
+        # preset
+        self.beating_period = 32768
+        self.blind_period = 200
+        self.adpll_limit = 8161512
+
+        self.gtx_init_phase = 0.0
+        self.gtx_init_freq = 125_000_000
+        self.gtx_jitter = 200e-17
+
+        self.helper_init_phase = 0.0
+        self.main_init_phase = 0.0
+        self.helper_init_freq = 125_000_000 * (1 - (1 / self.beating_period))
+        self.main_init_freq = 125_000_000
+        self.dcxo_jitter = 95e-17
+
+        # hardware delay
+        self.irq_delay = 2e-6  # ~2us for the interrupt handling
+        self.i2c_comm_delay = 85.6e-6
+        self.dcxo_settling_delay = 100e-6

src/wrpll_simulation/sim.py

@@ -1,305 +0,0 @@
-import numpy as np
-from numba import njit
-from wrpll_simulation.wave_gen import square
-
-
-@njit
-def simulation_jit(
-    time,
-    gtx_freq,
-    gtx_jitter,
-    helper_pll,
-    main_pll,
-    h_KP,
-    h_KI,
-    h_KD,
-    m_KP,
-    m_KI,
-    m_KD,
-    dcxo_freq,
-    h_jitter,
-    m_jitter,
-    base_adpll,
-    N,
-    adpll_write_period,
-    i2c_comm_delay,
-    dcxo_settling_delay,
-    blind_period,
-    start_up_delay,
-    helper_init_freq=0
-):
-
-    arr_len = len(time)
-
-    main = np.zeros(arr_len, dtype=np.int8)
-    helper = np.zeros(arr_len, dtype=np.int8)
-    gtx = np.zeros(arr_len, dtype=np.int8)
-
-    gtx_tag = 0
-    gtx_ready = 0
-    gtx_FF = 0
-    gtx_beating = np.zeros(arr_len, dtype=np.int8)
-
-    main_tag = 0
-    main_ready = 0
-    main_FF = 0
-    main_beating = np.zeros(arr_len, dtype=np.int8)
-
-    phase_err_arr = np.zeros(arr_len, dtype=np.int16)
-    period_err_arr = np.zeros(arr_len, dtype=np.int16)
-    helper_adpll_arr = np.zeros(arr_len, dtype=np.int32)
-    main_adpll_arr = np.zeros(arr_len, dtype=np.int32)
-
-    helperfreq = np.zeros(arr_len, dtype=np.int32)
-    mainfreq = np.zeros(arr_len, dtype=np.int32)
-
-    phase_collector_r = 0
-    colr_gtx_tag = colr_main_tag = 0
-    phase_collector_state = 0
-
-    period_collector_r = 0
-    colr_last_gtx_tag = 0
-    beating_period = 0
-    period_collector_state = 0
-
-    # initial condition
-    timestep = time[1] - time[0]
-    gtx_phase = 0
-    h_phase = np.random.uniform(0, 360)
-    m_phase = np.random.uniform(0, 360)
-
-    if main_pll and not helper_pll:
-        helper_freq = helper_init_freq
-    else:
-        helper_freq = dcxo_freq * (1 + base_adpll * 0.0001164 / 1_000_000) * ((N-1) / N)
-    main_freq = dcxo_freq * (1 + base_adpll * 0.0001164 / 1_000_000)
-
-    last_gtx_tag = last_gtx_FF = last_gtx_beat = 0
-    last_main_tag = last_main_FF = last_main_beat = 0
-
-    last_helper = 0
-
-    counter = 0
-    gtx_blind_counter = 0
-    gtx_blinded = False
-
-    main_blind_counter = 0
-    main_blinded = False
-
-    wait_gtx = True
-    wait_main = True
-
-    # firmware values
-    FW_gtx_tag = 0
-    FW_main_tag = 0
-
-    period_err = last_period_err = 0
-    h_prop = h_integrator = h_derivative = 0
-    h_adpll = base_adpll
-    h_i2c_active_index = 0
-    period_colr_arm = h_i2c_active = False
-
-    phase_err = last_phase_err = 0
-    m_prop = m_integrator = m_derivative = 0
-    m_adpll = base_adpll
-    m_i2c_active_index = 0
-    phase_colr_arm = m_i2c_active = False
-
-    adpll_max = 8161512
-    def clip(n, minn, maxn): return max(min(maxn, n), minn)
-
-    for i, t in enumerate(time):
-
-        h_phase += 360 * helper_freq * (timestep + h_jitter[i])
-        helper[i] = square(h_phase)
-
-        m_phase += 360 * main_freq * (timestep + m_jitter[i])
-        main[i] = square(m_phase)
-
-        gtx_phase += 360 * gtx_freq * (timestep + gtx_jitter[i])
-        gtx[i] = square(gtx_phase)
-
-        if not last_helper and helper[i]:
-
-            gtx_FF, gtx_beating[i] = DDMTD(gtx[i], last_gtx_FF)
-            main_FF, main_beating[i] = DDMTD(main[i], last_main_FF)
-
-            gtx_tag, gtx_ready, gtx_blind_counter, gtx_blinded = Deglitcher(gtx_beating[i], gtx_tag, gtx_ready,
-                                                                            gtx_blind_counter, gtx_blinded, blind_period,
-                                                                            last_gtx_beat, last_gtx_tag, counter)
-
-            main_tag, main_ready, main_blind_counter, main_blinded = Deglitcher(main_beating[i], main_tag, main_ready,
-                                                                                main_blind_counter, main_blinded, blind_period,
-                                                                                last_main_beat, last_main_tag, counter)
-
-            phase_collector_r, wait_gtx, wait_main, colr_gtx_tag, colr_main_tag, phase_collector_state = phase_collector_FSM(gtx_ready, main_ready, gtx_tag, main_tag,
-                                                                                                                             wait_gtx, wait_main, colr_gtx_tag, colr_main_tag, phase_collector_state)
-
-            if phase_collector_r:
-                FW_gtx_tag = colr_gtx_tag
-                FW_main_tag = colr_main_tag
-
-            period_collector_r, colr_last_gtx_tag, beating_period, period_collector_state = period_collector_FSM(gtx_ready, gtx_tag,
-                                                                                                                 colr_last_gtx_tag, period_collector_state)
-
-            counter += 1
-
-            last_gtx_beat = gtx_beating[i]
-            last_gtx_FF = gtx_FF
-            last_gtx_tag = gtx_tag
-            last_main_beat = main_beating[i]
-            last_main_FF = main_FF
-            last_main_tag = main_tag
-
-        else:
-            gtx_beating[i] = last_gtx_beat
-            gtx_FF = last_gtx_FF
-            gtx_tag = last_gtx_tag
-
-            main_beating[i] = last_main_beat
-            main_FF = last_main_FF
-            main_tag = last_main_tag
-
-        if i > start_up_delay:
-
-            if i % adpll_write_period == 0:
-                period_colr_arm = phase_colr_arm = True
-
-            # Firmware filters
-
-            if period_colr_arm and period_collector_r:
-
-                period_colr_arm = False
-                period_err = N - beating_period
-
-                if helper_pll:
-                    h_prop = period_err * h_KP
-                    h_integrator += period_err * h_KI
-                    h_derivative = (period_err - last_period_err) * h_KD
-
-                    h_adpll = clip(int(base_adpll + h_prop + h_integrator + h_derivative), -adpll_max, adpll_max)
-                    last_period_err = period_err
-                    h_i2c_active_index = i
-                    h_i2c_active = True
-
-            if phase_colr_arm and phase_collector_r:
-
-                phase_colr_arm = False
-                tag_diff = ((FW_main_tag - FW_gtx_tag) % N)
-                if tag_diff > N/2:
-                    phase_err = tag_diff - N
-                else:
-                    phase_err = tag_diff
-
-                if main_pll:
-                    m_prop = phase_err * m_KP
-                    m_integrator += phase_err * m_KI
-                    m_derivative = (phase_err - last_phase_err) * m_KD
-
-                    m_adpll = clip(int(base_adpll + m_prop + m_integrator + m_derivative), -adpll_max, adpll_max)
-                    last_phase_err = phase_err
-
-                    m_i2c_active_index = i
-                    m_i2c_active = True
-
-            # i2c communication delay
-
-            if h_i2c_active and i >= i2c_comm_delay + dcxo_settling_delay + h_i2c_active_index:
-                helper_freq = dcxo_freq * (1 + h_adpll * 0.0001164 / 1_000_000) * ((N-1) / N)
-                h_i2c_active = False
-
-            if m_i2c_active and i >= i2c_comm_delay + dcxo_settling_delay + m_i2c_active_index:
-                main_freq = dcxo_freq * (1 + m_adpll * 0.0001164 / 1_000_000)
-                m_i2c_active = False
-
-        last_helper = helper[i]
-
-        # Data
-        period_err_arr[i] = period_err
-        phase_err_arr[i] = phase_err
-        helper_adpll_arr[i] = h_adpll
-        main_adpll_arr[i] = m_adpll
-        helperfreq[i] = helper_freq
-        mainfreq[i] = main_freq
-
-    return period_err_arr, phase_err_arr, helper_adpll_arr, main_adpll_arr, gtx_beating, main_beating, gtx, helper, main, helperfreq, mainfreq
-
-
-@njit
-def DDMTD(sig_in, last_FF):
-    return sig_in, last_FF
-
-
-@njit
-def Deglitcher(beating, t_out, t_ready, blind_counter, blinded, blind_period, last_beat, last_tag, counter):
-
-    if blind_counter == 0 and beating and not last_beat:  # rising
-        t_out = counter
-        t_ready = 1
-        blinded = True
-    else:
-        t_out = last_tag
-        t_ready = 0
-
-    if beating:
-        blind_counter = blind_period - 1
-
-    if blind_counter != 0:
-        blind_counter -= 1
-
-    return t_out, t_ready, blind_counter, blinded
-
-
-@njit
-def phase_collector_FSM(g_tag_r, m_tag_r, gtx_tag, main_tag, wait_gtx, wait_main, colr_gtx_tag, colr_main_tag, FSM_state):
-
-    collector_r = 0
-
-    match FSM_state:
-        case 0:  # IDEL
-            if g_tag_r and m_tag_r:
-                colr_gtx_tag = gtx_tag
-                colr_main_tag = main_tag
-                FSM_state = 3  # OUTPUT
-
-            elif g_tag_r:
-                colr_gtx_tag = gtx_tag
-                wait_main = True
-                FSM_state = 2  # WAITMAIN
-
-            elif m_tag_r:
-                colr_main_tag = main_tag
-                wait_gtx = True
-                FSM_state = 1  # WAITGTX
-        case 1:  # WAITGTX
-            if g_tag_r:
-                colr_gtx_tag = gtx_tag
-                FSM_state = 3  # OUTPUT
-        case 2:  # WAITMAIN
-            if m_tag_r:
-                colr_main_tag = main_tag
-                FSM_state = 3  # OUTPUT
-        case 3:  # OUTPUT
-            wait_gtx = wait_main = False
-            collector_r = 1
-            FSM_state = 0
-
-    return collector_r, wait_gtx, wait_main, colr_gtx_tag, colr_main_tag, FSM_state
-
-
-@njit
-def period_collector_FSM(g_tag_r, gtx_tag, colr_last_gtx_tag, FSM_state):
-    collector_r = 0
-    beating_period = 0
-    match FSM_state:
-        case 0:  # IDEL
-            if g_tag_r:
-                colr_last_gtx_tag = gtx_tag
-                FSM_state = 1
-        case 1:  # OUTPUT
-            if g_tag_r:
-                beating_period = gtx_tag - colr_last_gtx_tag
-                collector_r = 1
-                FSM_state = 0  # IDEL
-
-    return collector_r, colr_last_gtx_tag, beating_period, FSM_state

src/wrpll_simulation/timesim.py

@@ -0,0 +1,132 @@
+import typing as tp
+import numpy as np
+from numba.core.types import *
+from numba.experimental import jitclass
+
+from wrpll_simulation.config import Timesim_Config
+from wrpll_simulation.timesim_node import *
+
+
+@jitclass
+class WRPLL_Timesim(object):
+    cfg: Timesim_Config
+
+    time: float32[:]
+
+    freq_diff: float32[:]
+    phase_diff: float32[:]
+    helper_error: int16[:]
+    main_error: int16[:]
+
+    def __init__(self, cfg: Timesim_Config, rng: tp.Generator):
+        # subclass/inheritance is not supported by numba jitclass
+        # https://github.com/numba/numba/issues/1694
+        self.cfg = cfg
+
+        sim_length = cfg.sim_length
+        stop_time = cfg.timestep_size * sim_length
+        self.time = np.linspace(0, stop_time, sim_length).astype(np.float32)
+
+        self.freq_diff = np.zeros(sim_length, dtype=np.float32)
+        self.phase_diff = np.zeros(sim_length, dtype=np.float32)
+        self.helper_error = np.zeros(sim_length, dtype=np.int16)
+        self.main_error = np.zeros(sim_length, dtype=np.int16)
+
+        # __post_init__ is not supported by numba jitclass
+        # https://github.com/numba/numba/issues/4037
+        self.simulate(rng)
+
+    def simulate(self, rng: tp.Generator):
+        cfg = self.cfg
+
+        timestep_size = cfg.timestep_size
+        step_input_time = self.seconds_to_step(cfg.step_input_time)
+        irq_delay = self.seconds_to_step(cfg.irq_delay)
+        i2c_comm_delay = self.seconds_to_step(cfg.i2c_comm_delay)
+
+        # simulation node
+        gtx = Phase_Accumlator(cfg.gtx_init_freq, cfg.gtx_init_phase)
+        helper = Phase_Accumlator(cfg.helper_init_freq, cfg.helper_init_phase)
+        main = Phase_Accumlator(cfg.main_init_freq, cfg.main_init_phase)
+
+        ddmtd_gtx = DDMTD(cfg.blind_period)
+        ddmtd_main = DDMTD(cfg.blind_period)
+        gtx_tag_irq = EventManager_IRQ()
+        main_tag_irq = EventManager_IRQ()
+
+        tag_collector = Tag_Collector(cfg.beating_period)
+        helper_PLL = PI_loop(cfg.helper_PI, cfg.helper_init_freq, 0, cfg.adpll_limit)
+        main_PLL = PI_loop(cfg.main_PI, cfg.main_init_freq, 0, cfg.adpll_limit)
+
+        counter = 0
+
+        print("Running...")
+        for i in range(cfg.sim_length):
+            if i == step_input_time:
+                gtx.set_freq(gtx.freq + cfg.step_frequency)
+                gtx.set_phase(gtx.phase + cfg.step_phase)
+
+            if cfg.has_jitter:
+                gtx.update(timestep_size + rng.normal(0, cfg.gtx_jitter))
+                helper.update(timestep_size + rng.normal(0, cfg.dcxo_jitter))
+                main.update(timestep_size + rng.normal(0, cfg.dcxo_jitter))
+            else:
+                gtx.update(timestep_size)
+                helper.update(timestep_size)
+                main.update(timestep_size)
+
+            # GATEWARE
+            if helper.is_rising():
+                ddmtd_gtx.sync_update(gtx.o, counter)
+                ddmtd_main.sync_update(main.o, counter)
+
+                # for clock domain crossing
+                gtx_tag_irq.multireg(ddmtd_gtx.tag_ready, ddmtd_gtx.tag)
+                main_tag_irq.multireg(ddmtd_main.tag_ready, ddmtd_main.tag)
+
+                counter += 1
+
+            if main.is_rising():
+                # Generate interrupt request
+                gtx_tag_irq.sync_update(i)
+                main_tag_irq.sync_update(i)
+
+            # FIRMWARE
+            if gtx_tag_irq.is_due(i, irq_delay):
+                tag_collector.collect_gtx_tag(gtx_tag_irq.tag_csr)
+                helper_PLL.update(i, tag_collector.get_period_error())
+
+                if tag_collector.is_phase_error_ready():
+                    tag_collector.set_phase_error_ready(False)
+                    main_PLL.update(i, tag_collector.get_phase_error())
+
+            if main_tag_irq.is_due(i, irq_delay):
+                tag_collector.collect_main_tag(main_tag_irq.tag_csr)
+
+                if tag_collector.is_phase_error_ready():
+                    tag_collector.set_phase_error_ready(False)
+                    main_PLL.update(i, tag_collector.get_phase_error())
+
+            if helper_PLL.i2c_is_due(i, i2c_comm_delay):
+                helper.set_freq(helper_PLL.get_new_freq())
+
+            if main_PLL.i2c_is_due(i, i2c_comm_delay):
+                main.set_freq(main_PLL.get_new_freq())
+
+            # Data Logging
+            self.freq_diff[i] = np.float32(main.freq - gtx.freq)
+            self.phase_diff[i] = np.float32((main.phase - gtx.phase) % 360)
+            if self.phase_diff[i] > 180:
+                self.phase_diff[i] -= 360
+
+            if helper_PLL.i2c_is_due(i, i2c_comm_delay):
+                self.helper_error[i] = tag_collector.get_period_error()
+            elif i > 0:
+                self.helper_error[i] = self.helper_error[i - 1]
+            if main_PLL.i2c_is_due(i, i2c_comm_delay):
+                self.main_error[i] = tag_collector.get_phase_error()
+            elif i > 0:
+                self.main_error[i] = self.main_error[i - 1]
+
+    def seconds_to_step(self, seconds: float64):
+        return int(seconds / self.cfg.timestep_size)

src/wrpll_simulation/timesim_node.py

@@ -0,0 +1,223 @@
+from numba.core.types import *
+from numba.experimental import jitclass
+
+from wrpll_simulation.config import PI_Config
+
+
+@jitclass
+class Phase_Accumlator(object):
+    last_o: int8
+    o: int8
+    freq: float64
+    phase: float64
+
+    def __init__(self, freq: float64, phase: float64):
+        self.last_o = 0
+        self.o = 0
+        self.freq = freq
+        self.phase = phase
+
+    def update(self, time_increment: float64):
+        self.last_o = self.o
+
+        self.phase = (self.phase + 360 * self.freq * time_increment) % 360
+
+        # square wave function
+        if self.phase < 180:
+            self.o = 0
+        else:
+            self.o = 1
+
+    def set_freq(self, freq: float64):
+        self.freq = freq
+
+    def set_phase(self, phase: float64):
+        self.phase = phase
+
+    def is_rising(self) -> bool:
+        return not self.last_o and self.o
+
+
+# GATEWARE
+@jitclass
+class DDMTD(object):
+    FF: int8
+    beating: int8
+    last_beating: int8
+
+    tag: int32
+    tag_ready: int8
+    blind_counter: int16
+    blind_period: int16
+
+    def __init__(self, blind_period: int16):
+        # back to back Flip Flop
+        self.FF = 0
+        self.beating = 0
+        self.last_beating = 0
+
+        # deglticher
+        self.tag = 0
+        self.tag_ready = 0
+        self.blind_counter = 0
+        self.blind_period = blind_period
+
+    def sync_update(self, D_in: int8, counter: int32):
+        self.last_beating = self.beating
+
+        # FF shifting
+        self.beating = self.FF
+        self.FF = D_in
+
+        self.deglitcher_first_edge(counter)
+
+    def deglitcher_first_edge(self, counter: int32):
+        if self.blind_counter == 0 and self.beating_is_rising():
+            self.tag = counter
+            self.tag_ready = 1
+        else:
+            self.tag_ready = 0
+
+        if self.beating:
+            self.blind_counter = self.blind_period - 1
+
+        if self.blind_counter != 0:
+            self.blind_counter -= 1
+
+    def beating_is_rising(self) -> bool:
+        return not self.last_beating and self.beating
+
+
+@jitclass
+class EventManager_IRQ(object):
+    trigger: int8
+    tag_csr: int32
+
+    trigger_index: int64
+
+    def __init__(self):
+        self.trigger = 0
+        self.tag_csr = 0
+
+        # for simulating delay
+        self.trigger_index = 0
+
+    def sync_update(self, index):
+        if self.trigger:
+            self.trigger = 0
+            self.trigger_index = index
+
+    def multireg(self, trigger: int8, tag: int32):
+        if trigger:
+            self.trigger = 1
+            self.tag_csr = tag
+
+    def is_due(self, index, delay):
+        return index - self.trigger_index == delay
+
+
+# FIRMWARE
+@jitclass
+class Tag_Collector(object):
+    setpt_beating_period: int64
+
+    gtx_tag_ready: bool
+    gtx_tag: int32
+    main_tag_ready: bool
+    main_tag: int32
+
+    def __init__(self, setpt_beating_period: int64):
+        self.setpt_beating_period = setpt_beating_period
+
+        # for main PLL
+        self.gtx_tag_ready = False
+        self.gtx_tag = 0
+        self.main_tag_ready = False
+        self.main_tag = 0
+
+    def collect_gtx_tag(self, tag: int32):
+        self.gtx_tag = tag
+        self.gtx_tag_ready = True
+
+    def collect_main_tag(self, tag: int32):
+        self.main_tag = tag
+        self.main_tag_ready = True
+
+    def get_period_error(self) -> int32:
+        # period_error = n * ideal_beating - GTX_tag(n)
+        period_error = (-self.gtx_tag) % self.setpt_beating_period
+
+        if period_error > self.setpt_beating_period / 2:
+            return period_error - self.setpt_beating_period
+
+        return period_error
+
+    def get_phase_error(self) -> int32:
+        # tag_diff = main_tag(n) - gtx_tag(n)
+        tag_diff = (self.main_tag - self.gtx_tag) % self.setpt_beating_period
+
+        # mapping tags from [0, 2π] -> [-π, π]
+        if tag_diff > self.setpt_beating_period / 2:
+            return tag_diff - self.setpt_beating_period
+
+        return tag_diff
+
+    def set_phase_error_ready(self, ready: bool):
+        self.main_tag_ready = ready
+        self.gtx_tag_ready = ready
+
+    def is_phase_error_ready(self) -> bool:
+        return self.main_tag_ready and self.gtx_tag_ready
+
+
+@jitclass
+class PI_loop(object):
+    KP: int64
+    KI: int64
+    integrator: int64
+
+    center_freq: float64
+    adpll: int32
+    base_adpll: int32
+    adpll_limit: int32
+
+    i2c_transfer_index: int64
+
+    def __init__(
+        self,
+        PI_Conifg: PI_Config,
+        center_freq: float64,
+        base_adpll: int32,
+        adpll_limit: int32,
+    ):
+        # PI controller
+        self.KP = PI_Conifg.KP
+        self.KI = PI_Conifg.KI
+        self.integrator = 0
+
+        # ADPLL calcuation
+        self.center_freq = center_freq
+        self.adpll = base_adpll
+        self.base_adpll = base_adpll
+        self.adpll_limit = adpll_limit
+
+        # for simulating delay
+        self.i2c_transfer_index = 0
+
+    def update(self, index, tag_error: int32):
+        self.i2c_transfer_index = index
+        self.adpll = self.get_adpll(tag_error)
+
+    def get_adpll(self, tag_error: int32) -> int32:
+        prop = tag_error * self.KP
+        self.integrator += tag_error * self.KI
+        return self.cramp_adpll(self.base_adpll + int(prop + self.integrator))
+
+    def get_new_freq(self):
+        return self.center_freq * (1 + self.adpll * 0.0001164 / 1_000_000)
+
+    def cramp_adpll(self, adpll: int32) -> int32:
+        return max(min(self.adpll_limit, adpll), -self.adpll_limit)
+
+    def i2c_is_due(self, index, delay):
+        return index - self.i2c_transfer_index == delay

src/wrpll_simulation/wave_gen.py

@@ -1,29 +0,0 @@
-import numpy as np
-import numba as nb
-from numba import njit
-
-
-
-def gussian_jitter(RMS_jitter, size, seed=None):
-    return np.random.default_rng(seed).normal(0, RMS_jitter/2, size)
-
-
-@njit(fastmath=True)
-def square_with_jitter(time, freq, jitter):
-    n = len(time)
-    wave = np.empty(n)
-    timestep = time[1] - time[0]
-
-    phase = 0.
-    for i in range(n):
-        phase += 360 * freq * (timestep + jitter[i])
-        wave[i] = square(phase)
-    return wave
-
-
-@njit
-def square(x):
-    if np.mod(x, 360) < 180:
-        return 0
-    else:
-        return 1

src/wrpll_simulation/wrpll.py

@@ -1,102 +0,0 @@
-import numpy as np
-from wrpll_simulation.sim import simulation_jit
-from wrpll_simulation.wave_gen import gussian_jitter
-
-
-class WRPLL_simulator():
-
-    def __init__(
-        self,
-        timestep,
-        total_steps,
-        sim_mode,
-        helper_filter,
-        main_filter,
-        gtx_freq,
-        adpll_write_period,
-        start_up_delay,
-        i2c_comm_delay=85.6e-6,
-        dcxo_settling_delay=100e-6,
-        gtx_jitter=200e-15,
-        dcxo_freq=125_000_000,
-        dcxo_jitter=95e-15,
-        freq_acquisition_error=100,
-        N=4069,
-        blind_period=128,
-        helper_init_freq=None,
-        seed=None
-    ):
-
-        self.time = np.linspace(0, timestep*total_steps, total_steps)
-        self.sim_mode = sim_mode
-        self.h_KP = helper_filter["KP"]
-        self.h_KI = helper_filter["KI"]
-        self.h_KD = helper_filter["KD"]
-        self.m_KP = main_filter["KP"]
-        self.m_KI = main_filter["KI"]
-        self.m_KD = main_filter["KD"]
-
-        # init condition
-        self.gtx_freq = gtx_freq
-        self.gtx_jitter = gussian_jitter(gtx_jitter, len(self.time), seed)
-
-        self.dcxo_freq = dcxo_freq
-        self.h_jitter = gussian_jitter(dcxo_jitter, len(self.time), seed)
-        self.m_jitter = gussian_jitter(dcxo_jitter, len(self.time), seed)
-        self.N = N
-        self.helper_init_freq = helper_init_freq
-
-        # freq_acquisition() error
-        freq_diff = gtx_freq - dcxo_freq + \
-            np.random.default_rng(seed).uniform(-freq_acquisition_error, freq_acquisition_error)
-        self.base_adpll = int(freq_diff * (1 / dcxo_freq) * (1e6 / 0.0001164))
-
-        # sim config
-        self.i2c_comm_delay = int(i2c_comm_delay/timestep)
-        self.dcxo_settling_delay = int(dcxo_settling_delay/timestep)
-        self.blind_period = blind_period
-        self.adpll_write_period = int(adpll_write_period/timestep)
-        self.start_up_delay = int(start_up_delay/timestep)
-
-        if type(self.sim_mode) is not str:
-            raise ValueError(f"pll_type {type(self.sim_mode)} is not a string")
-
-        self.helper_pll = self.main_pll = False
-        if self.sim_mode.lower() == "both":
-            self.helper_pll = self.main_pll = True
-        elif self.sim_mode.lower() == "helper_pll":
-            self.helper_pll = True
-        elif self.sim_mode.lower() == "main_pll":
-            if self.helper_init_freq == None:
-                raise ValueError("main pll mode need to set a helper frequency")
-            self.main_pll = True
-        else:
-            raise ValueError("sim_mode is not helper_pll nor main_pll")
-
-    def run(self):
-        print("running simulation...")
-        self.period_err, self.phase_err, self.helper_adpll, self.main_adpll, self.gtx_beating, self.main_beating, self.gtx, self.helper, self.main, self.helperfreq, self.mainfreq = simulation_jit(
-            self.time,
-            self.gtx_freq,
-            self.gtx_jitter,
-            self.helper_pll,
-            self.main_pll,
-            self.h_KP,
-            self.h_KI,
-            self.h_KD,
-            self.m_KP,
-            self.m_KI,
-            self.m_KD,
-            self.dcxo_freq,
-            self.h_jitter,
-            self.m_jitter,
-            self.base_adpll,
-            self.N,
-            self.adpll_write_period,
-            self.i2c_comm_delay,
-            self.dcxo_settling_delay,
-            self.blind_period,
-            self.start_up_delay,
-            self.helper_init_freq
-        )
-        print("Done!")