diff --git a/fleck/jax.py b/fleck/jax.py index be13d26..5fba6ce 100644 --- a/fleck/jax.py +++ b/fleck/jax.py @@ -123,6 +123,10 @@ def rotation_model(self, f0=0, t0_rot=0, u1=0, u2=0): spot_model : array Flux as a function of time and wavelength """ + u1 = jnp.atleast_1d(u1) + u2 = jnp.atleast_1d(u2) + u_ld = jnp.column_stack([u1, u2]) + ( spot_position_x, spot_position_y, spot_position_z, major_axis, minor_axis, angle, rad, contrast @@ -135,22 +139,32 @@ def rotation_model(self, f0=0, t0_rot=0, u1=0, u2=0): ) radial_coord = 1 - jnp.geomspace(1e-5, 1, 100)[::-1] + unspotted_total_flux = trapezoid( y=( - 2 * np.pi * radial_coord * - self.limb_darkening(radial_coord, u1, u2) - ), + 2 * np.pi * radial_coord[:, None] * + self.limb_darkening( + radial_coord[:, None], *u_ld.T + ) + ).T, x=radial_coord ) + limb_dark = self.limb_darkening( + mu, + u1=u1[None, None, :, None], + u2=u2[None, None, :, None] + ) + # Morris 2020 Eqn 6-7 spot_model = f0 - jnp.sum( np.pi * rad ** 2 * (1 - contrast) * - self.limb_darkening(mu, u1, u2) * - mask_behind_star, + limb_dark * + mask_behind_star / + unspotted_total_flux[None, None, :, None], axis=1 - ) / unspotted_total_flux + ) f_S = rad ** 2 * mu * (spot_position_z < 0).astype(int) return spot_model, f_S @@ -272,24 +286,22 @@ def add_spot(self, lon, lat, rad, contrast=None, temperature=None, spectrum=None grid is ``ActiveStar.phot`` """ if contrast is None and spectrum is None and temperature is not None: - self.phot = self._blackbody(self.wavelength, self.T_eff) + if self.phot is None: + self.phot = self._blackbody(self.wavelength, self.T_eff) spectrum = self._blackbody(self.wavelength, temperature) for attr, new_value in zip("lon, lat, rad, spectrum, temperature".split(', '), [lon, lat, rad, spectrum, temperature]): - prop = getattr(self, attr) - if not hasattr(new_value, 'ndim'): - new_value = jnp.array([new_value]) + new_value = jnp.atleast_1d(new_value) - if prop is not None: - if prop.ndim > 1 or (len(prop) > 1 and len(prop) == len(new_value)): + if attr == 'spectrum': + if len(prop): new_value = jnp.vstack([prop, new_value]) - else: - new_value = jnp.concatenate([prop, new_value]) - - setattr(self, attr, new_value) + elif len(prop.shape) and len(new_value.shape): + new_value = jnp.concatenate([prop, new_value]) + setattr(self, attr, new_value) @jit def _blackbody(self, wavelength_meters, temperature): @@ -439,11 +451,17 @@ def transit_model(self, t0, period, rp, a, inclination, jnp.sum(spot_coverages * spot_spectra, axis=1) ) + # if rp is scalar, turn it into a vector with length == N_wavelengths: + rp = rp * jnp.ones(u1.shape[0]) + transit = vmap( - lambda u_ld: jaxoplanet.core.light_curve( - u1=u_ld[0], u2=u_ld[1], b=jnp.hypot(X, Y), r=rp + lambda u_ld, rp: jaxoplanet.core.light_curve( + u1=u_ld[0], + u2=u_ld[1], + b=jnp.hypot(X, Y), + r=rp ), in_axes=0, out_axes=1 - )(u_ld) + )(u_ld, rp) contaminated_transit = ( time_series_spectrum - jnp.abs(transit) * self.phot[None, :] @@ -463,7 +481,7 @@ def transit_model(self, t0, period, rp, a, inclination, spot_position_x - Y[:, None, None, None] ) occultation_possible = jnp.squeeze( - (planet_spot_distance < (major_axis + rp)) & + (planet_spot_distance < (major_axis + rp.mean())) & (spot_position_z < 0) ) @@ -471,7 +489,7 @@ def transit_model(self, t0, period, rp, a, inclination, def time_step( carry, j, X=X, Y=Y, spot_position_y=spot_position_y, spot_position_x=spot_position_x, major_axis=major_axis, - minor_axis=minor_axis, rp=rp, angle=angle, + minor_axis=minor_axis, rp=rp.mean(), angle=angle, occultation_possible=occultation_possible ): return carry, lax.cond( @@ -509,7 +527,7 @@ def time_step( return ( out_of_transit[..., 0] * (contaminated_transit + scaled_occultation), - apparent_rprs2, X, Y, + apparent_rprs2, out_of_transit[..., 0], (contaminated_transit + scaled_occultation), spectrum_at_transit ) diff --git a/notebooks/AU-Mic.ipynb b/notebooks/AU-Mic.ipynb new file mode 100644 index 0000000..754f376 --- /dev/null +++ b/notebooks/AU-Mic.ipynb @@ -0,0 +1,613 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "4b3bdeee-cd82-485f-8305-98b0c3491cc6", + "metadata": {}, + "source": [ + "# AU Mic with `fleck`\n", + "\n", + "Fit the TESS Sector 1 light curve using `fleck` with three cool active regions. Using the best-fit spot map, extrapolate to find the wavelength dependence of the rotational modulation at other wavelengths." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cd5f73ae-a46e-4d43-8899-215d88503f8a", + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib inline\n", + "import numpyro\n", + "from numpyro.infer import MCMC, NUTS, Predictive\n", + "from numpyro import distributions as dist\n", + "\n", + "# Set the number of cores on your machine for parallel computing:\n", + "cpu_cores = 4\n", + "numpyro.set_host_device_count(cpu_cores)\n", + "\n", + "import jax\n", + "from jax import jit, numpy as jnp, config\n", + "from jax.random import PRNGKey, split\n", + "\n", + "# we need float64 support:\n", + "config.update(\"jax_enable_x64\", True)\n", + "\n", + "import arviz\n", + "from corner import corner\n", + "\n", + "import matplotlib.pyplot as plt\n", + "from matplotlib.gridspec import GridSpec\n", + "\n", + "import numpy as np\n", + "import astropy.units as u\n", + "from astropy.time import Time\n", + "from expecto import get_spectrum\n", + "from lightkurve import search_lightcurve\n", + "\n", + "from fleck.jax import ActiveStar, bin_spectrum" + ] + }, + { + "cell_type": "markdown", + "id": "fdf04d91-cf00-4361-9e2e-94b672821107", + "metadata": {}, + "source": [ + "The commented code below finds the transmittance weighted mean wavelength of the TESS bandpass." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6302e7be-d9a6-4f39-8dc5-0f7a4716caaa", + "metadata": {}, + "outputs": [], + "source": [ + "# import numpy as np\n", + "# import astropy.units as u\n", + "# from tynt import FilterGenerator\n", + "\n", + "# f = FilterGenerator()\n", + "# tess_filt = f.reconstruct('TESS/TESS.Red')\n", + "# tess_mean_wavelength = np.average(tess_filt.wavelength, weights=tess_filt.transmittance).to_value(u.m)\n", + "# tess_mean_wavelength\n", + "\n", + "# this is the answer you'd get if you ran the above code:\n", + "tess_mean_wavelength = 8.004867649770393e-07 # [m]" + ] + }, + { + "cell_type": "markdown", + "id": "8d7dcdbf-e6f5-476c-ad88-d3b72ea46c7b", + "metadata": {}, + "source": [ + "The sigma clipping below helps remove flares:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "44082f4b-9009-4683-b9af-4266e79dc843", + "metadata": {}, + "outputs": [], + "source": [ + "tess_lc = search_lightcurve(\n", + " \"AU Mic\", mission=\"TESS\", author=\"SPOC\", sector=1\n", + ").download_all()[0].normalize().remove_nans().remove_outliers(sigma_upper=2.8)" + ] + }, + { + "cell_type": "markdown", + "id": "2edb6973-88cd-4711-8351-279975d79a31", + "metadata": {}, + "source": [ + "Assume spots are *very* cold. They're not likely *this* cold." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "016c4123-123d-4725-889a-ad4744f0adfe", + "metadata": {}, + "outputs": [], + "source": [ + "T_phot = 3700 # Plavchan (2009)\n", + "T_spot1 = T_spot2 = T_spot3 = 2300\n", + "\n", + "blackbody = lambda *args: ActiveStar()._blackbody(*args)\n", + "\n", + "phot = jnp.array([float(blackbody(tess_mean_wavelength, T_phot))])\n", + "spectrum = jnp.array(\n", + " [[blackbody(tess_mean_wavelength, T_spot1)],\n", + " [blackbody(tess_mean_wavelength, T_spot2)],\n", + " [blackbody(tess_mean_wavelength, T_spot3)]]\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "d5f65c19-1293-407d-b64a-e979fa767e5b", + "metadata": {}, + "source": [ + "Below we construct a model that you can tweak by hand, to see how the parameters change:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ae0dce21-044a-44d8-b278-b701303fee75", + "metadata": {}, + "outputs": [], + "source": [ + "u1 = 0.453\n", + "u2 = 0.207\n", + "\n", + "r = ActiveStar(\n", + " times=tess_lc.time.value, \n", + " inclination=np.pi/2,\n", + " T_eff=T_phot,\n", + " wavelength=jnp.array([tess_mean_wavelength]),\n", + " P_rot=4.863,\n", + " phot=phot,\n", + " spectrum=spectrum\n", + ")\n", + "\n", + "r.lon = jnp.array([2.2, -0.5, 2.2])\n", + "r.lat = jnp.array([np.pi/2, np.pi/3, 0.1])\n", + "r.rad = jnp.array([0.2, 0.2, 0.2])\n", + "r.temperature = jnp.array([T_spot1, T_spot2, T_spot3])\n", + "\n", + "lc, contam = r.rotation_model(f0=1.05, u1=u1, u2=u2, t0_rot=tess_lc.time.value[0])\n", + "lc = np.squeeze(lc)\n", + "\n", + "ax = plt.gca()\n", + "tess_lc.plot(ax=ax)\n", + "ax.plot(tess_lc.time.value, lc)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5092769d-633a-4a40-9625-2d1db08ab5a6", + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(10, 4))\n", + "\n", + "r.plot_star(0, 2, 10, np.pi/2, ax=ax[0])\n", + "ax[0].set(\n", + " title='$\\phi = 0$'\n", + ")\n", + "\n", + "r.plot_star(0, 2, 10, np.pi/2, t0_rot=r.P_rot/2, ax=ax[1])\n", + "ax[1].set(\n", + " title='$\\phi = \\pi$'\n", + ");" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0a832743-b143-424c-af37-c75fc8789b19", + "metadata": {}, + "outputs": [], + "source": [ + "def numpyro_model(\n", + " y=jnp.array(np.array(tess_lc.flux.value)), \n", + " y_err=jnp.array(np.array(tess_lc.flux_err.value)), \n", + " n_spots=3, save_model=False\n", + "): \n", + "\n", + " theta = numpyro.sample(\n", + " 'theta', dist.Uniform(-np.pi, np.pi), \n", + " sample_shape=(n_spots,)\n", + " )\n", + "\n", + " x1, x2 = jnp.sin(theta), jnp.cos(theta)\n", + " lon = numpyro.deterministic(\n", + " 'lon', jnp.arctan2(x1, x2)\n", + " ) + np.pi\n", + "\n", + " rad = numpyro.sample(\n", + " 'rad', dist.Uniform(low=0, high=0.3), sample_shape=(n_spots,)\n", + " )\n", + "\n", + " f0 = numpyro.sample(\n", + " 'f0', dist.Uniform(low=0.9, high=1.1)\n", + " )\n", + " \n", + " r = ActiveStar(\n", + " times=tess_lc.time.value, \n", + " inclination=np.pi/2,\n", + " T_eff=T_phot,\n", + " phot=phot,\n", + " spectrum=spectrum,\n", + " wavelength=jnp.array([tess_mean_wavelength]),\n", + " P_rot=4.863\n", + " )\n", + " \n", + " r.lon = lon\n", + " r.lat = jnp.ones(n_spots) * np.pi / 2\n", + " r.rad = rad\n", + "\n", + " lc, contam = r.rotation_model(f0=f0, u1=u1, u2=u2, t0_rot=tess_lc.time.value[0])\n", + " lc = jnp.squeeze(lc)\n", + "\n", + " log_beta = numpyro.sample('log_beta', dist.Uniform(-1, 2))\n", + " \n", + " if save_model:\n", + " # this gets used to produce posterior predictive samples later on:\n", + " numpyro.deterministic(\"_lc_model\", lc)\n", + " contaminated_depth = 1e6 * contam\n", + " numpyro.deterministic(\"_depth_model\", contaminated_depth)\n", + "\n", + " # Normally distributed likelihood\n", + " numpyro.sample(\n", + " \"obs\", dist.Normal(\n", + " loc=lc, \n", + " scale=y_err * jnp.exp(log_beta)\n", + " ), obs=y\n", + " )" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "faab30fd-4770-4833-89b5-3fb4e60650b9", + "metadata": {}, + "outputs": [], + "source": [ + "from numpyro.infer import (\n", + " SVI, autoguide, Trace_ELBO, \n", + ")\n", + "from numpyro import optim\n", + "\n", + "guide = autoguide.AutoMultivariateNormal(numpyro_model)\n", + "\n", + "# # some alternatives to this guide are:\n", + "# guide = autoguide.AutoDAIS(numpyro_model)\n", + "# guide = autoguide.AutoBNAFNormal(numpyro_model)\n", + "\n", + "svi = SVI(\n", + " model=numpyro_model, \n", + " guide=guide, \n", + " optim=optim.Adagrad(step_size=0.15),\n", + " loss=Trace_ELBO()\n", + ")\n", + "svi_result = svi.run(\n", + " rng_key=PRNGKey(1), \n", + " num_steps=500\n", + ")\n", + "plt.loglog(svi_result.losses + 1 - np.nanmin(svi_result.losses))\n", + "plt.show()\n", + "\n", + "params = svi_result.params\n", + "posteriors = guide.sample_posterior(PRNGKey(1), params, sample_shape=(2_000,))\n", + "labels = [k for k, v in posteriors.items() if not k.startswith(\"_\")]\n", + "\n", + "\n", + "samples = []\n", + "iter_labels = []\n", + "for key, label in zip(posteriors, labels):\n", + " \n", + " if label == 'x':\n", + " continue\n", + " \n", + " if posteriors[key].ndim > 1:\n", + " for i, col in enumerate(posteriors[key].T):\n", + " samples.append(posteriors[key][:, i])\n", + " iter_labels.append(f\"{label}_{i}\")\n", + " else: \n", + " samples.append(posteriors[key][None, :])\n", + " iter_labels.append(label)\n", + " \n", + "samples = np.vstack(samples).T\n", + "\n", + "\n", + "corner(samples, labels=iter_labels)\n", + "fig = plt.gcf()\n", + "fig.suptitle('SVI')\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "a092e7a6-773e-4c91-bf5a-6f3f4c5fed7c", + "metadata": {}, + "source": [ + "Below is what you'd run for HMC, but this particular model is highly degenerate, and the fits won't converge:" + ] + }, + { + "cell_type": "raw", + "id": "4b3dc31b-d883-4ed9-96a1-f8f3c051db42", + "metadata": {}, + "source": [ + "rng_seed = 0\n", + "rng_keys = split(\n", + " PRNGKey(rng_seed), \n", + " cpu_cores\n", + ")\n", + "\n", + "sampler = NUTS(\n", + " numpyro_model, \n", + " dense_mass=True\n", + ")\n", + "\n", + "mcmc = MCMC(\n", + " sampler, \n", + " num_warmup=300, \n", + " num_samples=500, \n", + " num_chains=cpu_cores\n", + ")\n", + "\n", + "# Run the MCMC\n", + "mcmc.run(rng_keys)\n", + "\n", + "result = arviz.from_numpyro(mcmc)\n", + "\n", + "corner(\n", + " result,\n", + " var_names=['f0', 'lon', 'rad'],\n", + " quiet=True,\n", + ");" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "60865082-07d1-4409-9aad-8afc3912f367", + "metadata": {}, + "outputs": [], + "source": [ + "posterior_predictive = Predictive(\n", + " model=numpyro_model, \n", + " posterior_samples=posteriors,\n", + " return_sites=['_lc_model', '_depth_model'],\n", + ")\n", + "\n", + "pred = posterior_predictive(\n", + " rng_key=PRNGKey(1), \n", + " save_model=True\n", + ")\n", + "\n", + "y_pred = pred['_lc_model'] # contaminated transit model over wavelength\n", + "depth_pred = pred['_depth_model'] # contaminated transit depth at mid-transit time" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "70265d51-68e9-4a8f-9d5e-70fedf20b34d", + "metadata": {}, + "outputs": [], + "source": [ + "low, mid, high = np.percentile(y_pred, [16, 50, 84], axis=0)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "df1a7a47-efed-4c6d-a992-822e85283308", + "metadata": {}, + "outputs": [], + "source": [ + "plt.scatter(tess_lc.time.value, tess_lc.flux.value, color='silver', s=2)\n", + "\n", + "plt.plot(tess_lc.time.value, mid, color='DodgerBlue', lw=2)\n", + "\n", + "plt.fill_between(tess_lc.time.value, low, high, color='DodgerBlue', alpha=0.4, lw=0)\n", + "\n", + "plt.gca().set(\n", + " xlabel='Time [d]',\n", + " ylabel='Flux'\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "5d892e24-3180-4e11-9b71-b6ffdd68a055", + "metadata": {}, + "source": [ + "### extrapolate wavelength dependence from best fit in TESS\n", + "\n", + "This time, use PHOENIX model spectra to predict the wavelength dependence of variability. The lowest temperature model available with this package is 2300 K." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "2931762d-04b6-4c08-9a20-2cf2eca78559", + "metadata": {}, + "outputs": [], + "source": [ + "from fleck.jax import bin_spectrum\n", + "\n", + "wavelengths = jnp.linspace(0.5, 5, 100) * u.um\n", + "times = jnp.linspace(tess_lc.time.value.min(), tess_lc.time.value.max(), 300)\n", + "\n", + "kwargs = dict(\n", + " bins=wavelengths, \n", + " min=wavelengths.min(), \n", + " max=wavelengths.max(), \n", + " log=False\n", + ")\n", + "\n", + "\n", + "phot, spot = [\n", + " bin_spectrum(\n", + " get_spectrum(temp, 5.0, cache=True), **kwargs\n", + " ) for temp in [T_phot, T_spot1]\n", + "]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "78e3577f-3342-4d68-9e42-e699b67a6f0e", + "metadata": {}, + "outputs": [], + "source": [ + "spectrum = jnp.vstack([\n", + " spot.flux.value, \n", + " spot.flux.value, \n", + " spot.flux.value\n", + "])\n", + "\n", + "r = ActiveStar(\n", + " times=times, \n", + " inclination=np.pi/2,\n", + " T_eff=T_phot,\n", + " wavelength=phot.wavelength.to_value(u.m),\n", + " phot=phot.flux.value,\n", + " spectrum=spectrum,\n", + " P_rot=4.863,\n", + " temperature=jnp.array([T_spot1, T_spot2, T_spot3])\n", + ")\n", + "\n", + "r.lon = posteriors['lon'].mean(axis=0)\n", + "r.lat = jnp.array([np.pi / 2] * len(r.lon))\n", + "r.rad = posteriors['rad'].mean(axis=0)\n", + "\n", + "lc, contam = r.rotation_model(\n", + " f0=posteriors['f0'].mean(), u1=u1, u2=u2, \n", + " t0_rot=tess_lc.time.value[0]\n", + ")\n", + "lc = np.squeeze(lc)" + ] + }, + { + "cell_type": "markdown", + "id": "c7470263-8da4-4150-8e00-73c5eec9de5f", + "metadata": {}, + "source": [ + "Compute rotational modulation without transits:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c567aaa9-3b23-405c-81d6-ce1d524e1a62", + "metadata": {}, + "outputs": [], + "source": [ + "# blue to red colormap over `wavelengths`\n", + "cmap = lambda x: plt.cm.Spectral_r((x - wavelengths.min()) / np.ptp(wavelengths))\n", + "\n", + "for i in range(0, len(wavelengths), 20):\n", + " plt.plot(\n", + " times, lc[:, i], \n", + " color=cmap(wavelengths[i]), \n", + " label=f\"{wavelengths[i].to_value(u.um):.1f} µm\"\n", + " )\n", + " \n", + "plt.legend(loc='upper right', framealpha=1)\n", + "plt.gca().set(\n", + " xlabel='Time [d]',\n", + " ylabel='Flux'\n", + ");" + ] + }, + { + "cell_type": "markdown", + "id": "fcb7b7d3-4f65-403c-bf56-c8074ca933e4", + "metadata": {}, + "source": [ + "Predict transmission contamination from this best-fit model for AU Mic b:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c2b28584-3ca6-4915-9eba-93fef738be29", + "metadata": {}, + "outputs": [], + "source": [ + "# Wittrock (2023)\n", + "period = 8.46308\n", + "rp = 0.0488\n", + "inclination = np.radians(89.57917)\n", + "t0 = 2458322.77 - 2457000\n", + "a = 18.79\n", + "\n", + "# compute model near a transit:\n", + "r.times = jnp.linspace(t0 - 0.25, t0 + 0.25, 250)\n", + "\n", + "transit_lc, contam = r.transit_model(\n", + " t0=t0, \n", + " period=period, \n", + " rp=rp,\n", + " a=a,\n", + " inclination=inclination, \n", + " u1=u1, \n", + " u2=u2,\n", + " t0_rot=tess_lc.time.value[0]\n", + ")[:2]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c65f0778-4ba3-4187-8154-a038a62e64fc", + "metadata": {}, + "outputs": [], + "source": [ + "# blue to red colormap over `wavelengths`\n", + "fig, ax = plt.subplots(1, 3, figsize=(15, 4))\n", + "\n", + "cmap = lambda x: plt.cm.Spectral_r((x - wavelengths.min()) / np.ptp(wavelengths))\n", + "\n", + "for i in range(0, len(wavelengths), 20):\n", + " ax[0].plot(\n", + " r.times, transit_lc[:, i], \n", + " color=cmap(wavelengths[i]), \n", + " label=f\"{wavelengths[i].to_value(u.um):.1f} µm\"\n", + " )\n", + "\n", + "ax[0].legend(loc='lower right', framealpha=1)\n", + "ax[0].set(\n", + " xlabel='Time [d]',\n", + " ylabel='Flux'\n", + ")\n", + "\n", + "r.plot_star(t0, rp, a, inclination, ax=ax[1], t0_rot=tess_lc.time.value[0])\n", + "\n", + "ax[2].plot(wavelengths.to_value(u.um)[:-1], 1e6 * contam)\n", + "ax[2].set(\n", + " xlabel='Wavelength [µm]',\n", + " ylabel='Contaminated transit\\ndepth [ppm]'\n", + ")\n", + "fig.tight_layout()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ef326c45-d46c-4a04-908b-5e2947cfa6fa", + "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.11.5" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}