Skip to content

Commit

Permalink
More tests, updated stubs and updating vexp_bounce.ipynb.
Browse files Browse the repository at this point in the history
  • Loading branch information
vitenti committed Apr 19, 2024
1 parent 026f3f1 commit 4ef487b
Show file tree
Hide file tree
Showing 7 changed files with 283 additions and 116 deletions.
134 changes: 50 additions & 84 deletions notebooks/primordial_perturbations/vexp_bounce.ipynb
Expand Up @@ -70,13 +70,10 @@
"outputs": [],
"source": [
"import sys\n",
"import math\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"from matplotlib.pyplot import cm\n",
"\n",
"from scipy.integrate import odeint\n",
"\n",
"from numcosmo_py import Ncm, Nc\n",
"from numcosmo_py.plotting.tools import set_rc_params_article, latex_float"
]
Expand Down Expand Up @@ -146,18 +143,10 @@
"Vexp_gauss = Nc.HICosmoVexp.new()\n",
"Vexp_gauss.set_em_coupling(Nc.HICosmoVexpEMCoupling.GAUSS)\n",
"\n",
"### SET 1\n",
"# current_set = \"set1\"\n",
"# Vexp.props.alphab = 8.3163e-2 # Alpha (# of e-folds) at the bounce\n",
"# Vexp.props.sigmaphi = 9.0 # Width of the Gaussian solution for the WdW equation\n",
"# Vexp.props.xb = 2.0e36 # Inverse of the scale factor at the bounce (Initial condition)\n",
"\n",
"\n",
"### SET 2\n",
"current_set = {\n",
" \"alphab\": 7.4847e-3, # Alpha (# of e-fold\\s) at the bounce\n",
" \"sigmaphi\": 100.0, # Width of the Gaussian solution for the WdW equation\n",
" \"xb\": 4.0e36, # Inverse of the scale factor at the bounce (Initial condition)\n",
" \"xb\": 2.0e36, # Inverse of the scale factor at the bounce (Initial condition)\n",
" \"dphi\": -9.0e-4, # Deviation of the Gaussian solution for the WdW equation\n",
" \"OmegaL\": 1.0, # H²(a when w=-1)/H²(a0). Basically gives the DE-dominated phase\n",
" \"Omegac\": 1.0, # Omega_d???\n",
Expand Down Expand Up @@ -203,8 +192,8 @@
"tau_a = np.linspace(tau_min, tau_max, npoints)\n",
"tau_q_a = np.linspace(tau_qt_c, tau_qt_e, npoints)\n",
"\n",
"print(f\"Vexp model prepared in the interval ({tau_min},{tau_max})\")\n",
"print(f\"Vexp quantum dominated interval ({tau_qt_c},{tau_qt_e})\")"
"print(f\"Vexp model prepared in the interval ({tau_min}, {tau_max})\")\n",
"print(f\"Vexp quantum dominated interval ({tau_qt_c}, {tau_qt_e})\")"
]
},
{
Expand Down Expand Up @@ -500,19 +489,36 @@
"metadata": {},
"outputs": [],
"source": [
"tau_f = 10.0\n",
"\n",
"pem_none = Nc.HIPertEM.new()\n",
"pem_none.set_initial_condition_type(Ncm.CSQ1DInitialStateType.ADIABATIC2)\n",
"pem_none.set_abstol(\n",
" 1.0e-200\n",
") # The vacuum is a fixed point at zero, we need an absolute tolerance to evolve the system\n",
"pem_none.set_k(1.0e3)\n",
"pem_none.set_ti(Vexp_none.tau_min())\n",
"pem_none.set_tf(tau_f) # We do not want to evolve through the singularity\n",
"pem_none.set_vacuum_max_time(-12.0)\n",
"pem_none.set_vacuum_reltol(1.0e-8)\n",
"\n",
"\n",
"pem_cauchy = Nc.HIPertEM.new()\n",
"pem_cauchy.set_k(1.0e3)\n",
"pem_cauchy.set_initial_condition_type(Ncm.CSQ1DInitialStateType.ADIABATIC4)\n",
"pem_cauchy.set_ti(Vexp_cauchy.tau_min())\n",
"pem_cauchy.set_tf(tau_f) # We do not want to evolve through the singularity\n",
"pem_cauchy.set_vacuum_max_time(-1.0e-1)\n",
"pem_cauchy.set_vacuum_reltol(1.0e-8)\n",
"\n",
"\n",
"pem_gauss = Nc.HIPertEM.new()\n",
"pem_gauss.set_initial_condition_type(Ncm.CSQ1DInitialStateType.ADIABATIC2)\n",
"\n",
"pem_none.set_k(1.0e3)\n",
"pem_cauchy.set_k(1.0e3)\n",
"pem_gauss.set_k(1.0e3)"
"pem_gauss.set_k(1.0e3)\n",
"pem_gauss.set_ti(Vexp_gauss.tau_min())\n",
"pem_gauss.set_tf(tau_f) # We do not want to evolve through the singularity\n",
"pem_gauss.set_vacuum_max_time(-1.0e-1)\n",
"pem_gauss.set_vacuum_reltol(1.0e-8)"
]
},
{
Expand Down Expand Up @@ -544,15 +550,15 @@
"source": [
"print(f\"Searching for WKB scale in ({tau_min}, {tau_qt_c})\")\n",
"\n",
"(Found_tau_i_none, tau_i_none) = pem_none.find_adiab_time_limit(\n",
" Vexp_none, tau_min, -12.0, 1.0e-9\n",
")\n",
"(Found_tau_i_cauchy, tau_i_cauchy) = pem_cauchy.find_adiab_time_limit(\n",
" Vexp_cauchy, tau_min, -1.0e-2, 1.0e-9\n",
" Vexp_cauchy, tau_min, -1.0e-5, 1.0e-9\n",
")\n",
"(Found_tau_i_gauss, tau_i_gauss) = pem_gauss.find_adiab_time_limit(\n",
" Vexp_gauss, tau_min, -1.0e-2, 1.0e-9\n",
")\n",
"\n",
"# The vacuum for no coupling is a fixed point, all points are adiabatic, using the same interval as Cauchy\n",
"Found_tau_i_none, tau_i_none = Found_tau_i_cauchy, tau_i_cauchy"
" Vexp_gauss, tau_min, -1.0e-5, 1.0e-9\n",
")"
]
},
{
Expand All @@ -570,8 +576,6 @@
"metadata": {},
"outputs": [],
"source": [
"tau_f = 10.0\n",
"\n",
"print(\n",
" f\"Found ({Found_tau_i_none }) WKB scale for no-potential at {tau_i_none:7.3f}\"\n",
")\n",
Expand Down Expand Up @@ -610,76 +614,38 @@
"tau_ad_a = np.linspace(-2.35, -1.0, 100)\n",
"tau_ad_a = tau_a\n",
"\n",
"F1_a_none = np.array([pem_none.eval_F1(Vexp_none, tau) for tau in tau_ad_a])\n",
"F2_a_none = np.array([pem_none.eval_F2(Vexp_none, tau) for tau in tau_ad_a])\n",
"nu_a_none = np.array([pem_none.eval_nu(Vexp_none, tau) for tau in tau_ad_a])\n",
"\n",
"F1_a_cauchy = np.array([pem_cauchy.eval_F1(Vexp_cauchy, tau) for tau in tau_ad_a])\n",
"F2_a_cauchy = np.array([pem_cauchy.eval_F2(Vexp_cauchy, tau) for tau in tau_ad_a])\n",
"nu_a_cauchy = np.array([pem_cauchy.eval_nu(Vexp_cauchy, tau) for tau in tau_ad_a])\n",
"F1_a_cauchy = np.abs([pem_cauchy.eval_F1(Vexp_cauchy, tau) for tau in tau_ad_a])\n",
"F2_a_cauchy = np.abs([pem_cauchy.eval_F2(Vexp_cauchy, tau) for tau in tau_ad_a])\n",
"\n",
"F1_a_gauss = np.array([pem_gauss.eval_F1(Vexp_gauss, tau) for tau in tau_ad_a])\n",
"F2_a_gauss = np.array([pem_gauss.eval_F2(Vexp_gauss, tau) for tau in tau_ad_a])\n",
"nu_a_gauss = np.array([pem_gauss.eval_nu(Vexp_gauss, tau) for tau in tau_ad_a])\n",
"F1_a_gauss = np.abs([pem_gauss.eval_F1(Vexp_gauss, tau) for tau in tau_ad_a])\n",
"F2_a_gauss = np.abs([pem_gauss.eval_F2(Vexp_gauss, tau) for tau in tau_ad_a])\n",
"\n",
"plt.plot(tau_ad_a, F1_a_none, label=r\"no-pot $F_1$\")\n",
"plt.plot(tau_ad_a, F1_a_cauchy, label=r\"cauchy $F_1$\")\n",
"plt.plot(tau_ad_a, F1_a_gauss, label=r\"gauss $F_1$\")\n",
"plt.plot(tau_ad_a, F1_a_cauchy, label=r\"cauchy $|F_1|$\")\n",
"plt.plot(tau_ad_a, F1_a_gauss, label=r\"gauss $|F_1|$\")\n",
"\n",
"plt.plot(tau_ad_a, F2_a_none, label=r\"no-pot $F_2$\")\n",
"plt.plot(tau_ad_a, F2_a_cauchy, label=r\"cauchy $F_2$\")\n",
"plt.plot(tau_ad_a, F2_a_gauss, label=r\"gauss $F_2$\")\n",
"plt.plot(tau_ad_a, F2_a_cauchy, label=r\"cauchy $|F_2|$\")\n",
"plt.plot(tau_ad_a, F2_a_gauss, label=r\"gauss $|F_2|$\")\n",
"\n",
"plt.grid(visible=True, which=\"both\", linestyle=\":\", color=\"0.75\", linewidth=0.5)\n",
"leg = plt.legend(loc=\"best\", ncol=2, fontsize=13)\n",
"\n",
"plt.ylim(1.0e-20, None)\n",
"plt.xlabel(r\"$\\tau$\")\n",
"plt.xscale(\"symlog\")\n",
"plt.yscale(\"symlog\", linthresh=1.0e-20)\n",
"plt.yscale(\"log\")\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "6da0d3d9",
"metadata": {},
"source": [
"tau_i = -0.94\n",
"tau_f = 5.0\n",
"print(csq1d0.eval_adiab_at (Vexp, tau_i))\n",
"print(csq1d1.eval_adiab_at (Vexp, tau_i))"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "825f9b1d",
"metadata": {},
"outputs": [],
"source": [
"pem_none.set_init_cond_adiab(Vexp_none, tau_i_none)\n",
"pem_cauchy.set_init_cond_adiab(Vexp_cauchy, tau_i_cauchy)\n",
"pem_gauss.set_init_cond_adiab(Vexp_gauss, tau_i_gauss)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "1603d89c",
"metadata": {},
"outputs": [],
"source": [
"pem_none.set_ti(tau_i_none)\n",
"pem_none.set_tf(tau_f)\n",
"pem_none.prepare(Vexp_none)\n",
"\n",
"pem_cauchy.set_ti(tau_i_cauchy)\n",
"pem_cauchy.set_tf(tau_f)\n",
"pem_cauchy.prepare(Vexp_cauchy)\n",
"\n",
"pem_gauss.set_ti(tau_i_gauss)\n",
"pem_gauss.set_tf(tau_f)\n",
"pem_gauss.prepare(Vexp_gauss)\n",
"\n",
"tau_evol_a_none, tau_s_none = pem_none.get_time_array()\n",
Expand Down Expand Up @@ -716,23 +682,23 @@
"\n",
"for tau in tau_evol_a_none:\n",
" PE, PB = pem_none.eval_PE_PB(Vexp_none, tau)\n",
" J11, _, J22 = pem_none.eval_at(tau, state_none).get_J()\n",
" J11, _, J22 = pem_none.eval_at(Vexp_none, tau, state_none).get_J()\n",
" PE_a_none.append(PE)\n",
" PB_a_none.append(PB)\n",
" J11_a_none.append(J11)\n",
" J22_a_none.append(J22)\n",
"\n",
"for tau in tau_evol_a_cauchy:\n",
" PE, PB = pem_cauchy.eval_PE_PB(Vexp_cauchy, tau)\n",
" J11, _, J22 = pem_cauchy.eval_at(tau, state_cauchy).get_J()\n",
" J11, _, J22 = pem_cauchy.eval_at(Vexp_cauchy, tau, state_cauchy).get_J()\n",
" PE_a_cauchy.append(PE)\n",
" PB_a_cauchy.append(PB)\n",
" J11_a_cauchy.append(J11)\n",
" J22_a_cauchy.append(J22)\n",
"\n",
"for tau in tau_evol_a_gauss:\n",
" PE, PB = pem_gauss.eval_PE_PB(Vexp_gauss, tau)\n",
" J11, _, J22 = pem_gauss.eval_at(tau, state_gauss).get_J()\n",
" J11, _, J22 = pem_gauss.eval_at(Vexp_gauss, tau, state_gauss).get_J()\n",
" PE_a_gauss.append(PE)\n",
" PB_a_gauss.append(PB)\n",
" J11_a_gauss.append(J11)\n",
Expand Down Expand Up @@ -763,7 +729,7 @@
"source": [
"mylw = 1\n",
"\n",
"plt.figure(figsize=(16, 8))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_none,\n",
Expand All @@ -782,7 +748,7 @@
"plt.yscale(\"log\")\n",
"plt.show()\n",
"\n",
"plt.figure(figsize=(16, 8))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_cauchy,\n",
Expand Down Expand Up @@ -817,7 +783,7 @@
"metadata": {},
"outputs": [],
"source": [
"plt.figure(figsize=(14, 7))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_none,\n",
Expand All @@ -836,7 +802,7 @@
"plt.yscale(\"log\")\n",
"plt.show()\n",
"\n",
"plt.figure(figsize=(14, 7))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_cauchy,\n",
Expand All @@ -857,7 +823,7 @@
"metadata": {},
"outputs": [],
"source": [
"plt.figure(figsize=(14, 7))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_none,\n",
Expand All @@ -884,7 +850,7 @@
"# plt.ylim(1.0e-60,None)\n",
"plt.show()\n",
"\n",
"plt.figure(figsize=(14, 7))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_none,\n",
Expand Down Expand Up @@ -920,7 +886,7 @@
"metadata": {},
"outputs": [],
"source": [
"plt.figure(figsize=(14, 7))\n",
"plt.figure()\n",
"\n",
"plt.plot(\n",
" tau_evol_a_none,\n",
Expand Down

0 comments on commit 4ef487b

Please sign in to comment.