Merge pull request #81 from LemurPwned/streamlit/ui

Streamlit/UI

Author: LemurPwned

README.md

@@ -6,7 +6,7 @@
 [![pages-build-deployment](https://github.com/LemurPwned/cmtj/actions/workflows/pages/pages-build-deployment/badge.svg?branch=gh-pages)](https://github.com/LemurPwned/cmtj/actions/workflows/pages/pages-build-deployment)
 [![Version](https://img.shields.io/pypi/v/cmtj)](https://pypi.org/project/cmtj/)
 [![License](https://img.shields.io/pypi/l/cmtj.svg)](https://github.com/LemurPwned/cmtj/blob/master/LICENSE)
-[![Streamlit](https://static.streamlit.io/badges/streamlit_badge_black_white.svg)](http://cmtj-simulations.streamlit.app/)
+[![Streamlit](https://static.streamlit.io/badges/streamlit_badge_black_white.svg)](https://cmtj-app.streamlit.app/spectrum)
 ![Downloads](https://img.shields.io/pypi/dm/cmtj.svg)
 
 ## Table of contents
@@ -36,7 +36,14 @@ It is also possible to connect devices in parallel or in series to have electric
 
 ## Web GUI
 
-Check out the [streamlit hosted demo here](http://cmtj-simulations.streamlit.app/). You can simulate PIMM spectra and Spin-Diode spectra there. Let us know if you have any issues with the demo.
+Check out the [streamlit hosted demo here](https://cmtj-app.streamlit.app/spectrum). 
+You can simulate:
+
+* PIMM spectra and Spin-Diode spectra
+* Try some optimization fitting
+* Fit multi-domain or multi-level M(H) or R(H) loops in [Domain mode](https://cmtj-app.streamlit.app)
+
+Let us know if you have any issues with the demo.
 
 ## Quickstart
 

view/domain.py

@@ -5,22 +5,33 @@
 import streamlit as st
 from helpers import read_mh_data
 from simulation_fns import create_single_domain, get_axis_angles
-from utils import GENERIC_BOUNDS, GENERIC_UNITS
+from utils import GENERIC_BOUNDS, GENERIC_UNITS, get_init_kval
 
 from cmtj import AxialDriver, Junction
 from cmtj.utils import FieldScan
 
 _lock = RLock()
 
+
 with st.container():
     st.write("# Domain fitting")
-
-    st.write(
-        "Fit M(H) to multidomain model. "
-        "Upload your file with experimental data: with columns H, mx, my, mz.\n"
-        "First column is always H in A/m, the rest are the components of the magnetisation in range (-1, 1)."
-        "If you upload just two columns, the script will assume that the data is in the form H, mx."
-    )
+    with st.expander("Explanation"):
+        st.write(
+            "#### Fit M(H) to multidomain model. \n"
+            "##### Experimental data\n"
+            "Upload your file with experimental data: with columns H, mx, my, mz.\n"
+            "First column is always H in A/m, the rest are the components of the magnetisation in range (-1, 1)."
+            "If you upload just two columns, the script will assume that the data is in the form H, mx.\n"
+            "##### Simulation\n"
+            "Add new domains in the sidebar panel. Area of each domain is a weight. The resulting $m_\mathrm{mixed}(h)$ plot is produces with\n\n"
+            r"$m_\mathrm{mixed} = (1/{\sum_i a_i})\sum_i a_i m_i$"
+            "\n\n where $a_i$ is the area of $i$th domain and $m_i$ is $m(h)$ for that domain.\n"
+            "##### Coordinate system\n"
+            r"$\theta$ is the polar angle and $\phi$ is the azimuthal angle."
+            r" $\theta=90^\circ$ denotes fully in-plane magnetisation, "
+            r"$\theta=0^\circ$ denotes out-of-plane magnetisation (positive z-direction)"
+            r"$\phi=0^\circ$ denotes magnetisation pointing in the positive x direction."
+        )
 
     progress_bar = st.progress(0)
 
@@ -46,7 +57,8 @@ def simulate_domains():
     weights = [st.session_state[f"area{i}"] * 1e-18 for i in range(st.session_state.N)]
     wsum = sum(weights)
     weights = [w / wsum for w in weights]
-    Mmixed = np.zeros((len(Hscan), 3))
+    Mmixed = np.zeros((len(Hscan), 3), dtype=np.float16)
+    Mdomains = np.zeros((st.session_state.N, len(Hscan), 3), np.float16)
     for i, H in enumerate(Hvecs):
         j.clearLog()
         j.setLayerExternalFieldDriver("all", AxialDriver(*H))
@@ -60,11 +72,13 @@ def simulate_domains():
             mx = np.mean(log[f"domain_{l}_mx"][-st.session_state.last_N :])
             my = np.mean(log[f"domain_{l}_my"][-st.session_state.last_N :])
             mz = np.mean(log[f"domain_{l}_mz"][-st.session_state.last_N :])
-            Mmixed[i] += np.array([mx, my, mz]) * weights[l]
+            Mdomains[l, i] = np.array([mx, my, mz])
+            Mmixed[i] += Mdomains[l, i] * weights[l]
         # for each layer take last N values
         progress_bar.progress(int((i / maxlen) * 100))
     st.session_state["Mmixed"] = Mmixed
     st.session_state["Hscan"] = Hscan
+    st.session_state["Mdomains"] = Mdomains
 
 
 with st.sidebar:
@@ -77,6 +91,52 @@ def simulate_domains():
         accept_multiple_files=False,
         key="upload",
     )
+    st.checkbox("Show domains in polar angles", key="domain_use_angle", value=False)
+    domain_params = st.expander("Domain parameters", expanded=True)
+    with domain_params:
+        N = st.number_input(
+            "Number of domains", min_value=1, max_value=10, value=1, key="N"
+        )
+        for i in range(N):
+            st.markdown(f"#### Domain {i+1}")
+            unit_k = GENERIC_UNITS["K"]
+            st.number_input(
+                f"K ({i+1}) " r" ($\mathrm{" f"{unit_k}" "}$)",
+                min_value=GENERIC_BOUNDS["K"][0],
+                max_value=GENERIC_BOUNDS["K"][1],
+                value=1.2e3,
+                step=10.0,
+                key=f"K{i}",
+                help="Uniaxial anisotropy constant of the domain. Does not include shape anisotropy.",
+            )
+
+            st.number_input(
+                f"area ({i+1}) " r"($\mathrm{nm^2}$)",
+                min_value=1.0,
+                max_value=500.0,
+                value=10.0,
+                key=f"area{i}",
+                help="Area of the domain. In fact, this is a weight of the domain.",
+            )
+
+            st.number_input(
+                r"$\theta_\mathrm{K}$ (deg.)",
+                min_value=0.0,
+                max_value=180.0,
+                value=get_init_kval(i),
+                key=f"theta_K{i}",
+                help="Polar angle of the anisotropy axis",
+            )
+            st.number_input(
+                r"$\phi_\mathrm{K}$ (deg.)",
+                min_value=0.0,
+                max_value=180.0,
+                value=get_init_kval(i),
+                key=f"phi_K{i}",
+                help="Azimuthal angle of the anisotropy axis",
+            )
+            st.markdown("-----\n")
+
     with st.expander("Simulation parameters", expanded=True):
         st.selectbox(
             "H axis", options=["x", "y", "z", "xy", "xz", "yz"], key="H_axis", index=0
@@ -125,6 +185,7 @@ def simulate_domains():
             value=1.8,
             step=0.01,
             key="Ms_shared",
+            help="Saturation magnetisation of the system. Affects the shape anisotropy",
         )
         st.number_input(
             "thickness (nm)",
@@ -141,40 +202,9 @@ def simulate_domains():
             key="alpha_shared",
             format="%.3f",
         )
-    domain_params = st.expander("Domain parameters", expanded=True)
-    with domain_params:
-        N = st.number_input(
-            "Number of domains", min_value=1, max_value=10, value=1, key="N"
-        )
-        for i in range(N):
-            st.markdown(f"#### Domain {i+1}")
-            unit_k = GENERIC_UNITS["K"]
-            st.number_input(
-                f"K ({i+1}) " r" ($\mathrm{" f"{unit_k}" "}$)",
-                min_value=GENERIC_BOUNDS["K"][0],
-                max_value=GENERIC_BOUNDS["K"][1],
-                value=1.2e3,
-                step=10.0,
-                key=f"K{i}",
-            )
 
-            st.number_input(
-                f"area ({i+1}) " r"($\mathrm{nm^2}$)",
-                min_value=1.0,
-                max_value=500.0,
-                value=10.0,
-                key=f"area{i}",
-            )
-            st.radio(
-                f"anisotropy axis ({i+1})",
-                options=["x", "y", "z"],
-                key=f"anisotropy_axis{i}",
-                index=2,
-            )
-            st.markdown("-----\n")
 
-
-def render(Hscan, Mmixed):
+def render(Hscan, Mmixed, Mdomains=None):
     if len(Hscan) <= 0 or len(Mmixed) <= 0:
         return
     with _lock:
@@ -187,21 +217,35 @@ def render(Hscan, Mmixed):
                 fields, mh = read_mh_data()
                 render_from_exp(ax, fields=fields, mh=mh)
             render_from_sim(ax, Hscan, Mmixed)
-
+            fig.suptitle("Mixed Domain model")
+            st.pyplot(fig)
+            n = 3
+            if st.session_state["domain_use_angle"]:
+                n = 2
+            w, h = plt.figaspect(st.session_state.N / n)
+            fig, ax = plt.subplots(
+                st.session_state.N, n, dpi=400, figsize=(w, h), sharex=True, sharey=True
+            )
+            if st.session_state.N <= 1:
+                ax = np.array([ax])
+            try:
+                render_individual_domains(fig, ax, Hscan, Mdomains)
+            except IndexError:
+                pass
             st.pyplot(fig)
 
 
 def render_from_exp(ax, fields, mh):
     if len(fields) <= 0 or len(mh) <= 0:
         st.toast("No data to render")
-
+    fields = np.asarray(fields)
     m = np.asarray(mh)
     for i, l in zip(range(m.shape[1]), "xyz"):
         ax[i].plot(
-            fields,
+            fields / 1e3,  # A/m --> kA/m
             m[:, i],
             label=f"$m_{l}$",
-            color="black",
+            color="yellow",
             marker="x",
             alpha=0.5,
             linestyle="--",
@@ -221,11 +265,55 @@ def render_from_sim(ax, Hscan, Mmixed):
     ax[2].set_xlabel("H (kA/m)")
 
 
+def render_individual_domains(fig, ax, Hscan, M: list[list[float]]):
+    n = 3
+    if st.session_state["domain_use_angle"]:
+        n = 2
+        ax[0, 0].set_title(r"$\theta$")
+        ax[0, 1].set_title(r"$\phi$")
+        for i, domain_m in enumerate(M):
+            # convert to polar
+            theta_m, phi_m, _ = FieldScan.vector2angle(
+                domain_m[:, 0], domain_m[:, 1], domain_m[:, 2]
+            )
+            ax[i, 0].plot(
+                Hscan / 1e3, theta_m, label=rf"$\theta_{i+1}$", color="crimson"
+            )
+            ax[i, 1].plot(
+                Hscan / 1e3, phi_m, label=rf"$\phi_{i+1}$", color="forestgreen"
+            )
+            ax[i, 0].set_ylabel(rf"$\theta_{i+1}$ (deg.)")
+            ax[i, 1].set_ylabel(rf"$\phi_{i+1}$ (deg.)")
+    else:
+        n = 3
+        ax[0, 0].set_title(r"$m_\mathrm{x}$")
+        ax[0, 1].set_title(r"$m_\mathrm{y}$")
+        ax[0, 2].set_title(r"$m_\mathrm{z}$")
+        for i, domain_m in enumerate(M):
+            ax[i, 0].plot(
+                Hscan / 1e3, domain_m[:, 0], label=rf"$m_{i+1}$", color="crimson"
+            )
+            ax[i, 1].plot(
+                Hscan / 1e3, domain_m[:, 1], label=rf"$m_{i+1}$", color="forestgreen"
+            )
+            ax[i, 2].plot(
+                Hscan / 1e3, domain_m[:, 2], label=rf"$m_{i+1}$", color="royalblue"
+            )
+            ax[i, 0].set_ylabel(rf"$m_{i+1}$")
+    for j in range(n):
+        ax[-1, j].set_xlabel("H (kA/m)")
+    fig.suptitle("Individual Domains")  # fig.legend()
+
+
 simulate_btn = st.button(
     "Simulate",
     key="simulate",
     on_click=simulate_domains,
     type="primary",
 )
 
-render(st.session_state.get("Hscan", []), st.session_state.get("Mmixed", []))
+render(
+    st.session_state.get("Hscan", []),
+    st.session_state.get("Mmixed", []),
+    Mdomains=st.session_state.get("Mdomains", []),
+)

view/simulation_fns.py

@@ -14,11 +14,14 @@
 
 def create_single_domain(id_: str) -> Layer:
     demag = [CVector(0, 0, 0), CVector(0, 0, 0), CVector(0, 0, 1)]
-    Kdir1 = get_axis_cvector(st.session_state[f"anisotropy_axis{id_}"])
+    # Kdir1 = get_axis_cvector(st.session_state[f"anisotropy_axis{id_}"])
+    Kdir = FieldScan.angle2vector(
+        theta=st.session_state[f"theta_K{id_}"], phi=st.session_state[f"phi_K{id_}"]
+    )
     layer = Layer(
         id=f"domain_{id_}",
-        mag=Kdir1,
-        anis=Kdir1,
+        mag=Kdir,
+        anis=Kdir,
         Ms=st.session_state["Ms_shared"],
         thickness=st.session_state["thickness_shared"] * 1e-9,
         cellSurface=1e-16,
@@ -34,11 +37,13 @@ def create_single_domain(id_: str) -> Layer:
 def create_single_layer(id_: str) -> tuple:
     """Do not forget to rescale the units!"""
     demag = [CVector(0, 0, 0), CVector(0, 0, 0), CVector(0, 0, 1)]
-    Kdir1 = get_axis_cvector(st.session_state[f"anisotropy_axis{id_}"])
+    Kdir = FieldScan.angle2vector(
+        theta=st.session_state[f"theta_K{id_}"], phi=st.session_state[f"phi_K{id_}"]
+    )
     layer = Layer(
         id=f"layer_{id_}",
-        mag=Kdir1,
-        anis=Kdir1,
+        mag=Kdir,
+        anis=Kdir,
         Ms=st.session_state[f"Ms{id_}"],
         thickness=st.session_state[f"thickness{id_}"] * 1e-9,
         cellSurface=1e-16,