diff --git a/qsymm/json_serialization/convert_func.py b/qsymm/json_serialization/convert_func.py
new file mode 100644
index 0000000000000000000000000000000000000000..326eb23f35186db4ba68bd2a24609f76cfba7584
--- /dev/null
+++ b/qsymm/json_serialization/convert_func.py
@@ -0,0 +1,74 @@
+import numpy as np
+import sympy
+import re
+import json
+from astropy.utils.misc import JsonCustomEncoder
+import qsymm
+from icecream import ic
+
+def bloch_model_to_json(model, name):
+    name = name + '.json'
+    model_dict = {str(key) : model[key].tolist() for key in model.keys()}
+    with open(name, 'w') as f:
+        json.dump(model_dict, f, cls = JsonCustomEncoder)
+
+def bloch_model_from_json(name):
+    with open("json_models/" + name + ".json", "r") as f:
+        model = json.load(f, cls=JsonCustomDecoderExp)
+    return qsymm.Model(model, momenta=["k_x", "k_y"])
+
+
+class JsonCustomDecoderExp(json.JSONDecoder):
+    """
+    A decoder function tailored specifically to our purposes: take a dictionary of strings, output a dictionary
+    of symbolic keys (Mul types etc.) and arrays.
+    """
+
+    def __init__(self, *args, **kwargs):
+        json.JSONDecoder.__init__(self, object_hook=self.object_hook, *args, **kwargs)
+
+    def object_hook(self, obj):
+        new_dict = dict()
+        for key in obj.keys():
+            # back-convert the key from string to symbols
+            symbolic = []
+            # find exponents, works for arbitrary number of exponents
+            exp_full = re.findall('\*?e\*\*\(-?I\*k_\w\)\*?', key) # gives list with exponential terms in key
+            if exp_full:
+                coef = key
+                for exp in exp_full:
+                    coef = coef.replace(exp, '') # constant coefficient 
+                    
+                exp_power_pl = re.findall('\*?e\*\*\(I\*(k_\w)\)\*?', key) # gives list of positive exponents
+                exp_power_mn = re.findall('\*?e\*\*\(-I\*(k_\w)\)\*?', key) # gives list of negative exponents
+                
+                if len(exp_power_pl)>=1:
+                    sym_exp_pls = []
+                    for exp in exp_power_pl:
+                        sym_exp_pls.append(sympy.Pow(sympy.Symbol('e'), sympy.I*sympy.Symbol(exp)))
+                    if len(exp_power_mn)==0: # no negative powers
+                        symbolic_key = sympy.Mul(*sym_exp_pls, sympy.Symbol(coef))
+                    
+                if len(exp_power_mn)>=1:
+                    sym_exp_mns = []
+                    for exp in exp_power_mn:
+                        sym_exp_mns.append(sympy.Pow(sympy.Symbol('e'), -sympy.I*sympy.Symbol(exp)))
+                    if len(exp_power_pl)==0: # no positive powers
+                        symbolic_key = sympy.Mul(*sym_exp_mns, sympy.Symbol(coef))
+                    else: # both positive and negative powers
+                        symbolic_key = sympy.Mul(*sym_exp_mns, *sym_exp_pls, sympy.Symbol(coef))
+            else:
+                symbolic_key = sympy.Symbol(key)
+                
+            assert str(key)==str(symbolic_key) # make sure converter works
+            array = []
+            for line in obj[key]:
+                lineski = []
+                for entry in line:
+                    lineski.append(entry[0] + 1j * entry[1])
+                array.append(lineski)
+            new_dict[symbolic_key] = np.array(array)
+        return new_dict
+
+
+
diff --git a/qsymm/json_serialization/minimal_model.ipynb b/qsymm/json_serialization/minimal_model.ipynb
new file mode 100644
index 0000000000000000000000000000000000000000..36dfa66aee317a8a6900f7f145f868f8932add59
--- /dev/null
+++ b/qsymm/json_serialization/minimal_model.ipynb
@@ -0,0 +1,151 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "id": "ec96970a-359f-4e6f-b2d9-529749253dab",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "import tinyarray\n",
+    "import qsymm\n",
+    "import sympy\n",
+    "from sympy import Matrix\n",
+    "import string\n",
+    "from convert_func import bloch_model_from_json, bloch_model_to_json"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "id": "d7832d0c-1b99-47aa-9a40-00091eda980e",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/latex": [
+       "$\\displaystyle \\left[\\begin{matrix}0.5 e^{i k_{x}} - 0.5 i e^{i k_{y}} + 0.5 i e^{- i k_{y}} + 0.5 e^{- i k_{x}} & 0\\\\0 & - 0.5 e^{i k_{x}} - 0.5 i e^{i k_{y}} + 0.5 i e^{- i k_{y}} - 0.5 e^{- i k_{x}}\\end{matrix}\\right]$"
+      ],
+      "text/plain": [
+       "Matrix([\n",
+       "[0.5*e**(I*k_x) - 0.5*I*e**(I*k_y) + 0.5*I*e**(-I*k_y) + 0.5*e**(-I*k_x),                                                                        0],\n",
+       "[                                                                      0, -0.5*e**(I*k_x) - 0.5*I*e**(I*k_y) + 0.5*I*e**(-I*k_y) - 0.5*e**(-I*k_x)]])"
+      ]
+     },
+     "execution_count": 3,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "sigma_x = tinyarray.array([[0, 1], [1, 0]])\n",
+    "sigma_y = tinyarray.array([[0, -1j], [1j, 0]])\n",
+    "sigma_z = tinyarray.array([[1, 0], [0, -1]])\n",
+    "\n",
+    "#create qsymm model\n",
+    "ham = \" cos(k_x) * sigma_z + sin(k_y) * eye(2) \"\n",
+    "H = qsymm.Model(ham, momenta=[\"k_x\", \"k_y\"])\n",
+    "H.tosympy()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "id": "309fb778-7db3-4a02-a163-ab6479aa7711",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "8 1\n"
+     ]
+    },
+    {
+     "data": {
+      "text/latex": [
+       "$\\displaystyle \\left[\\begin{matrix}d + e^{i k_{x}} i + \\frac{e^{i k_{x}}}{2} + e^{i k_{y}} m + i e^{i k_{y}} q - \\frac{i e^{i k_{y}}}{2} + e^{- i k_{y}} m - i e^{- i k_{y}} q + \\frac{i e^{- i k_{y}}}{2} + e^{- i k_{x}} i + \\frac{e^{- i k_{x}}}{2} & e^{i k_{x}} j + i e^{i k_{x}} l + e^{i k_{y}} n + i e^{i k_{y}} r + f + i h + e^{- i k_{y}} o - i e^{- i k_{y}} s + e^{- i k_{x}} j + i e^{- i k_{x}} l\\\\e^{i k_{x}} j - i e^{i k_{x}} l + e^{i k_{y}} o + i e^{i k_{y}} s + f - i h + e^{- i k_{y}} n - i e^{- i k_{y}} r + e^{- i k_{x}} j - i e^{- i k_{x}} l & e^{i k_{x}} k - \\frac{e^{i k_{x}}}{2} + e^{i k_{y}} p + i e^{i k_{y}} t - \\frac{i e^{i k_{y}}}{2} + g + e^{- i k_{y}} p - i e^{- i k_{y}} t + \\frac{i e^{- i k_{y}}}{2} + e^{- i k_{x}} k - \\frac{e^{- i k_{x}}}{2}\\end{matrix}\\right]$"
+      ],
+      "text/plain": [
+       "Matrix([\n",
+       "[d + e**(I*k_x)*i + e**(I*k_x)/2 + e**(I*k_y)*m + I*e**(I*k_y)*q - I*e**(I*k_y)/2 + e**(-I*k_y)*m - I*e**(-I*k_y)*q + I*e**(-I*k_y)/2 + e**(-I*k_x)*i + e**(-I*k_x)/2,                          e**(I*k_x)*j + I*e**(I*k_x)*l + e**(I*k_y)*n + I*e**(I*k_y)*r + f + I*h + e**(-I*k_y)*o - I*e**(-I*k_y)*s + e**(-I*k_x)*j + I*e**(-I*k_x)*l],\n",
+       "[                         e**(I*k_x)*j - I*e**(I*k_x)*l + e**(I*k_y)*o + I*e**(I*k_y)*s + f - I*h + e**(-I*k_y)*n - I*e**(-I*k_y)*r + e**(-I*k_x)*j - I*e**(-I*k_x)*l, e**(I*k_x)*k - e**(I*k_x)/2 + e**(I*k_y)*p + I*e**(I*k_y)*t - I*e**(I*k_y)/2 + g + e**(-I*k_y)*p - I*e**(-I*k_y)*t + I*e**(-I*k_y)/2 + e**(-I*k_x)*k - e**(-I*k_x)/2]])"
+      ]
+     },
+     "execution_count": 4,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "#find symmetries\n",
+    "candidates = qsymm.groups.square()\n",
+    "discrete_symm, continuous_symm = qsymm.symmetries(H, candidates)\n",
+    "print(len(discrete_symm), len(continuous_symm))\n",
+    "\n",
+    "abcd_sympy = list(sympy.symbols('a:z'))\n",
+    "abcd_list = list(string.ascii_lowercase)\n",
+    "del abcd_sympy[4]\n",
+    "del abcd_list[4]\n",
+    "\n",
+    "norbs = [('A', 2)]  # A atom per unit cell, 4 orbitals each\n",
+    "hopping_vectors = [('A', 'A', [0, 1]), ('A', 'A', [1, 0])] \n",
+    "\n",
+    "#find all temrs compatible with 1 specific symmetry\n",
+    "family= qsymm.bloch_family(hopping_vectors, [discrete_symm[1]], norbs)\n",
+    "\n",
+    "#add all the terms compatible with the symmetry\n",
+    "H = (H + sum([symbol * term for symbol, term in zip(abcd_sympy[3:], family)]))\n",
+    "H.tosympy(nsimplify=True)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "id": "387f8212-1861-4875-8fc5-798cb5d83220",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "bloch_model_to_json(H, 'model')"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "de61fbf0-f6a6-43db-b32c-300e408fcbd3",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "f5cb9c98-cc14-448d-8a03-41df556ce3d2",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "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.9.10"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/qsymm/json_serialization/model.json b/qsymm/json_serialization/model.json
new file mode 100644
index 0000000000000000000000000000000000000000..8f80422ce4459d12508738506308cc1ded303a96
--- /dev/null
+++ b/qsymm/json_serialization/model.json
@@ -0,0 +1 @@
+{"e**(-I*k_x)": [[[0.5, 0.0], [0.0, 0.0]], [[0.0, 0.0], [-0.5, 0.0]]], "e**(-I*k_y)": [[[0.0, 0.5], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.5]]], "e**(I*k_x)": [[[0.5, 0.0], [0.0, 0.0]], [[0.0, 0.0], [-0.5, 0.0]]], "e**(I*k_y)": [[[0.0, -0.5], [0.0, 0.0]], [[0.0, 0.0], [0.0, -0.5]]], "e**(-I*k_x)*l": [[[-0.0, 0.0], [0.0, 1.0]], [[0.0, -1.0], [-0.0, 0.0]]], "e**(-I*k_x)*j": [[[0.0, 0.0], [1.0, 0.0]], [[1.0, 0.0], [0.0, 0.0]]], "e**(I*k_y)*s": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 1.0], [0.0, 0.0]]], "e**(-I*k_y)*o": [[[0.0, 0.0], [1.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(-I*k_y)*t": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, -1.0]]], "e**(-I*k_y)*m": [[[1.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(I*k_y)*o": [[[0.0, 0.0], [0.0, 0.0]], [[1.0, 0.0], [0.0, 0.0]]], "e**(I*k_y)*t": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 1.0]]], "e**(I*k_y)*m": [[[1.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(-I*k_y)*r": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, -1.0], [0.0, 0.0]]], "e**(-I*k_y)*n": [[[0.0, 0.0], [0.0, 0.0]], [[1.0, 0.0], [0.0, 0.0]]], "e**(-I*k_x)*i": [[[1.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(I*k_x)*j": [[[0.0, 0.0], [1.0, 0.0]], [[1.0, 0.0], [0.0, 0.0]]], "e**(I*k_x)*l": [[[-0.0, 0.0], [0.0, 1.0]], [[0.0, -1.0], [-0.0, 0.0]]], "e**(-I*k_x)*k": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [1.0, 0.0]]], "e**(-I*k_y)*p": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [1.0, 0.0]]], "d": [[[1.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "f": [[[0.0, 0.0], [1.0, 0.0]], [[1.0, 0.0], [0.0, 0.0]]], "e**(I*k_y)*n": [[[0.0, 0.0], [1.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "h": [[[-0.0, 0.0], [0.0, 1.0]], [[0.0, -1.0], [-0.0, 0.0]]], "e**(I*k_y)*p": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [1.0, 0.0]]], "e**(-I*k_y)*q": [[[0.0, -1.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(I*k_y)*r": [[[0.0, 0.0], [0.0, 1.0]], [[0.0, 0.0], [0.0, 0.0]]], "g": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [1.0, 0.0]]], "e**(I*k_y)*q": [[[0.0, 1.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(I*k_x)*i": [[[1.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [0.0, 0.0]]], "e**(I*k_x)*k": [[[0.0, 0.0], [0.0, 0.0]], [[0.0, 0.0], [1.0, 0.0]]], "e**(-I*k_y)*s": [[[0.0, 0.0], [0.0, -1.0]], [[0.0, 0.0], [0.0, 0.0]]]}
\ No newline at end of file