first debug on metapeg. At the moment there is only one frequency being optimized. Although we are not using the epigraph yet, it would still be convenient to still consider a small range around the central frequency and optimize the average objective between them. This would still be better than considering only one frequency, but for the moment I need to understand if this model can work or not. Before running a full batch, more debugging is required to better understand how everything works. In this release, the code has not been tested in the results part, but the main optimization function run a couple of rounds without throwing errors. We still need to assure the normalization we make on the far fields is reasonable or not. The next step is to run a dummy simulation until the end and observe the results generation works, then the true test would be to let it run a true optimization with the right resolution.

metapeg.ipynb

@@ -2,17 +2,16 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "id": "8658f8cd",
    "metadata": {},
    "source": [
-    "This code implements an inverse designed metasurface that is used to recreate a given intenstity\n",
-    "pattern. It is based on meep and tidy3d tutorials."
+    " This code implements an inverse designed metasurface that is used to recreate a given intenstity\n",
+    "\n",
+    " pattern. It is based on meep and tidy3d tutorials."
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "4c8450d1",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -27,17 +26,16 @@
     "import scipy as sp\n",
     "import os\n",
     "\n",
-    "mp.verbosity(0)\n",
+    "mp.verbosity(2)\n",
     "\n",
     "#setup a results folder\n",
     "if not os.path.exists(\"results\"):\n",
-    "    os.makedirs(\"results\")"
+    "    os.makedirs(\"results\")\n"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "c728065b",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -57,13 +55,12 @@
     "\n",
     "######PML######\n",
     "pml_size = 1.0\n",
-    "pml_layers = [mp.PML(pml_size)]"
+    "pml_layers = [mp.PML(pml_size)]\n"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "ca9ab89e",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -74,19 +71,18 @@
     ")\n",
     "eta_e = 0.55  # erosion design field thresholding point (between 0 and 1)\n",
     "eta_d = 1 - eta_e  # dilation design field thresholding point (between 0 and 1)\n",
-    "filter_radius = mpa.get_conic_radius_from_eta_e(minimum_length, eta_e)"
+    "filter_radius = mpa.get_conic_radius_from_eta_e(minimum_length, eta_e)\n"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "5fc04b92",
    "metadata": {},
    "outputs": [],
    "source": [
     "######DESIGN REGION######\n",
-    "design_region_width = 3.0\n",
-    "design_region_length = 3.0\n",
+    "design_region_width = 4.0\n",
+    "design_region_length = 4.0\n",
     "design_region_height = (wavelength/np.sqrt(subs_index))/2 #half-wavelength in the substrate\n",
     "\n",
     "resolution = 20 # pixels/um\n",
@@ -143,7 +139,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "3163836f",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -160,7 +155,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "d2f566e1",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -179,7 +173,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "41b6bdaa",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -201,7 +194,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "3c05d8e1",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -246,22 +238,20 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "456d598e",
    "metadata": {},
    "outputs": [],
    "source": [
     "xs = ys = np.linspace(-design_region_width / 2, design_region_width / 2, Nx)\n",
-    "intensity_desired = intensity_desired_fn_logo(xs,ys)\n",
-    "print(intensity_desired.shape)"
+    "intensity_desired = intensity_desired_fn_logo(xs,ys)"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "5c6ed0bc",
    "metadata": {},
    "outputs": [],
    "source": [
+    "plt.figure()\n",
     "plt.pcolormesh(xs, ys, intensity_desired.T, cmap=\"magma\")\n",
     "plt.gca().set_aspect(\"equal\")\n",
     "plt.xlabel(\"x\")\n",
@@ -274,7 +264,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "8d1d7dc0",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -306,7 +295,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "c725fa9c",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -333,7 +321,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "362842c5",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -354,18 +341,20 @@
     "norm_sim.plot2D(output_plane=output_plane,show_sources=True, show_boundaries=True, show_geometry=True)\n",
     "plt.savefig(\"results/normalization_simulation_geometry.pdf\")\n",
     "\n",
-    "norm_sim.run(until=200)\n",
+    "norm_sim.run(until=1500)\n",
     "ref_fields = np.array([norm_sim.get_farfield(norm_near2far, point) for point in far_z_plane])\n",
     "#TODO this todo is linked to the one above. Need to check if the far-field need to be flipped or something\n",
     "#at this time, the two definition are consistent with each other.\n",
-    "intensity_norm_mean = npa.mean(npa.abs(ref_fields[:,:,0])**2)\n",
+    "#NOTE: read carefully the docs for get_farfield.\n",
+    "#it returns a list of the type [Ex1,Ey1,Ez1,Hx1,Hy1,Hz1,Ex2,Ey2,Ez2,Hx2,Hy2,Hz2....]\n",
+    "#for frequencies 1,2...\n",
+    "intensity_norm_mean = npa.mean(npa.abs(ref_fields[:,0])**2) #get the Ex component for all elements in the the list\n",
     "print(\"intensity_norm_mean shape (should be 1)\",intensity_norm_mean.shape)"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "710abef1",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -402,7 +391,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "aa694ac2",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -430,7 +418,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "c182593d",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -483,7 +470,7 @@
     "    if gradient.size > 0:\n",
     "        # NOTE: gradient is assigned in-place! No need to return.\n",
     "        gradient[...] = tensor_jacobian_product(mapping, 0)(\n",
-    "            v, eta_i, cur_beta, np.sum(dJ_du, axis=1)\n",
+    "            v, eta_i, cur_beta, dJ_du #NOTE: here we are using a single frequency!!\n",
     "        )\n",
     "    \n",
     "    evaluation_history.append(np.real(f0))\n",
@@ -559,7 +546,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "5225b3d4",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -609,7 +595,6 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "5c0ec83f",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -641,16 +626,11 @@
     "\n",
     "opt.sim.plot2D(output_plane=full_output_plane,fields=mp.Ex)\n",
     "plt.savefig(\"results/full_simulation_geometry.pdf\")\n",
-    "opt.sim.run(until=100)\n"
+    "opt.sim.run(until=1500)"
    ]
   }
  ],
  "metadata": {
-  "kernelspec": {
-   "display_name": "mp",
-   "language": "python",
-   "name": "python3"
-  },
   "language_info": {
    "codemirror_mode": {
     "name": "ipython",
@@ -661,9 +641,9 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.13.12"
+   "version": 3
   }
  },
  "nbformat": 4,
- "nbformat_minor": 5
+ "nbformat_minor": 2
 }