{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "假设在红移为1处去观测一个现在红移为2的星系，其红移为多少？"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$z_{1} = 1$, $z_{2} = 2$, $z_{1,2}$为题目所求红移 <p>\n",
    "$1+z_{1} = \\frac{a(t_0)}{a(t_1)}$ <p>\n",
    "$1+z_{2} = \\frac{a(t_0)}{a(t_2)}$ <p>\n",
    "$1+z_{1,2} = \\frac{a(t_1)}{a(t_2)} = \\frac{1+z_{2}}{1+z_{1}} = \\frac{3}{2}$ <p>\n",
    "$z_{1,2} = 0.5$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "假设宇宙的膨胀是一个线性过程，请根据今天的哈勃常数估算宇宙的年龄。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "假设宇宙膨胀的速度为常数v，那么尺度因子$a(t)=vt$, $\\dot a(t) = v$ <p>\n",
    "$H(t) = \\frac{\\dot a(t)}{a(t)} = \\frac{1}{t}$ <p>\n",
    "那么$H_0 = \\frac{1}{t0}$，即$t_0 = \\frac{1}{H_0} = \\frac{3.086 \\times 10^{19}}{70} = 1.39 \\times 10^{10} {\\rm yr} = 13.9 {\\rm Gyr}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "一团1000个太阳质量的中性气体云（全部为H原子）, 温度为30K，当云的密度大于多少时，气体云将发生塌缩？塌缩（自由下落）时标为多少？<p>\n",
    "云塌缩：势能 > 热能"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "粒子数(都是H原子)$N=\\frac{M}{m_H}$ <p>\n",
    "势能$E_p = \\frac{3}{5} \\times \\frac{GM^2}{R}$ <p>\n",
    "热能$E_th = \\frac{3}{2}NkT = \\frac{3}{2} \\frac{M}{m_H} kT$ <p>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "当密度大于4.84e-23 g/cm**3时，气体云会塌缩\n",
      "自由下落时标为9.58 Myr\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "#势能和热能相当时\n",
    "G = 6.67408e-8\n",
    "mSun = 1.988e33\n",
    "mH = 1.67e-24\n",
    "k = 1.38e-16\n",
    "yr = 31536000\n",
    "\n",
    "M_cl = 1000 * mSun\n",
    "T_cl = 30\n",
    "\n",
    "E_th = 1.5 * (M_cl/mH) * k * T_cl\n",
    "\n",
    "R_cl = (3/5) * G * M_cl**2 / E_th\n",
    "\n",
    "den_cl = M_cl / (4/3 * np.pi * R_cl**3)\n",
    "print(\"当密度大于%.2e g/cm**3时，气体云会塌缩\" %den_cl)\n",
    "\n",
    "t_ff = np.sqrt(3 * np.pi / (32 * G * den_cl))\n",
    "t_ff_Myr = t_ff/(1e6*yr)\n",
    "print(\"自由下落时标为%.2f Myr\" %t_ff_Myr)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这团气体中能形成的最大质量的恒星的质量是多少？<p>\n",
    "• 假设分子云中的恒星形成遵循Salpter初始质量函数（最小质量的恒星为0.08太阳质量）<p>\n",
    "• 这道题不是很好解析求解（需要估算一个合理的ΔM），或者采用MC方法 <p>\n",
    "• http://cluster.shao.ac.cn/~shen/Lecture/IMF.pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "最大质量恒星的平均质量：41.90 ± 21.77 Msun\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "M_cl = 1000  \n",
    "M_min = 0.08  \n",
    "M_max_guess = 100  # 初始猜测的最大恒星质量\n",
    "alpha = 2.35\n",
    "\n",
    "def salpeter_sample(M_min, M_max, size=1):\n",
    "    U = np.random.uniform(0, 1, size=size)\n",
    "    exponent = -1 / (alpha - 1)\n",
    "    return (M_min**(1 - alpha) + U * (M_max**(1 - alpha) - M_min**(1 - alpha)))**exponent\n",
    "\n",
    "def MC_simulation(num_trials):\n",
    "    max_masses = []\n",
    "\n",
    "    for _ in range(num_trials):\n",
    "        total_mass = 0\n",
    "        masses = []\n",
    "\n",
    "        while total_mass < M_cl:\n",
    "            mass = salpeter_sample(M_min, M_max_guess)\n",
    "            total_mass += mass\n",
    "            masses.append(mass)\n",
    "\n",
    "        max_masses.append(max(masses))\n",
    "\n",
    "    return np.mean(max_masses), np.std(max_masses)\n",
    "\n",
    "mean_max_mass, std_max_mass = MC_simulation(10000)\n",
    "print(f\"最大质量恒星的平均质量：{mean_max_mass:.2f} ± {std_max_mass:.2f} Msun\")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# good"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "太阳中的氢大概有10%在主序阶段被燃烧，请估算太阳处于主序阶段的时间？太阳表面的温度大概是~5500K，请由此估算地球表面的温度。<p>\n",
    "• 地球的反射率~0.3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "太阳在主序阶段的时间约为7.47 Gyr\n",
      "假设地球为黑体时表面温度为358 K\n"
     ]
    }
   ],
   "source": [
    "mSun = 1.988e33\n",
    "m_H = 0.75 * mSun\n",
    "\n",
    "E_H = 0.1 * m_H * 6e18  #erg\n",
    "Lsun = 3.8e33  #erg/s\n",
    "\n",
    "t = E_H / Lsun / (1e9*yr)\n",
    "print(\"太阳在主序阶段的时间约为%.2f Gyr\" %t)\n",
    "\n",
    "R = 1.5e13  #日地距离 cm\n",
    "f_earth = (1-0.3) * Lsun / (4 * np.pi * R**2)\n",
    "\n",
    "sigma = 5.67e-5\n",
    "T = (f_earth / sigma) ** (1/4)\n",
    "print(\"假设地球为黑体时表面温度为%d K\" %T)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "一个星系的恒星形成历史可以用exp(-t/τ)来描述（其中τ=3Gyr），该星系的年龄为10Gy，请计算该星系的V波段的恒星质光比（随着时间的演化）。（选做）<p>\n",
    "• SSP的质量随着时间的变化 <p>\n",
    "• http://202.127.29.3/~shen/Lecture/bc2003_hr_m162_salp_ssp.4color"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "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.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}