%22%22%22Variational%20Derivatives%20to%20Numerical%20Stencils%3A%20Inflation%20Simulations%0A%0ABased%20on%20arXiv%3A1608.04408%20-%20%22Robustness%20of%20Inflation%20to%20Inhomogeneous%20Initial%20Conditions%22%0Aby%20Clough%2C%20Lim%2C%20DiNunno%2C%20Fischler%2C%20Flauger%2C%20and%20Paban.%0A%0AThis%20notebook%20demonstrates%20how%20to%3A%0A1.%20Derive%20equations%20of%20motion%20from%20Lagrangians%20using%20variational%20calculus%0A2.%20Convert%20symbolic%20PDEs%20to%20finite%20difference%20stencils%0A3.%20Generate%20code%20for%20numerical%20simulations%20in%20multiple%20languages%0A%22%22%22%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.19.4%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%20%20%20%20from%20derive%20import%20(%0A%20%20%20%20%20%20%20%20Symbol%2C%20symbols%2C%20Function%2C%20Rational%2C%20R%2C%0A%20%20%20%20%20%20%20%20D%2C%20Simplify%2C%20Expand%2C%0A%20%20%20%20%20%20%20%20Exp%2C%20Sin%2C%20Cos%2C%20Sqrt%2C%20Pi%2C%0A%20%20%20%20%20%20%20%20TeXForm%2C%20Pipe%2C%0A%20%20%20%20)%0A%20%20%20%20from%20derive.calculus%20import%20VariationalDerivative%2C%20EulerLagrangeEquation%0A%20%20%20%20from%20derive.discretization%20import%20Discretize%2C%20ToStencil%2C%20StencilCodeGen%0A%20%20%20%20from%20derive.diffgeo%20import%20Metric%2C%20minkowski_metric%0A%20%20%20%20return%20(%0A%20%20%20%20%20%20%20%20D%2C%0A%20%20%20%20%20%20%20%20Discretize%2C%0A%20%20%20%20%20%20%20%20Function%2C%0A%20%20%20%20%20%20%20%20Metric%2C%0A%20%20%20%20%20%20%20%20Pipe%2C%0A%20%20%20%20%20%20%20%20R%2C%0A%20%20%20%20%20%20%20%20Simplify%2C%0A%20%20%20%20%20%20%20%20StencilCodeGen%2C%0A%20%20%20%20%20%20%20%20Symbol%2C%0A%20%20%20%20%20%20%20%20ToStencil%2C%0A%20%20%20%20%20%20%20%20VariationalDerivative%2C%0A%20%20%20%20%20%20%20%20mo%2C%0A%20%20%20%20%20%20%20%20symbols%2C%0A%20%20%20%20)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20From%20Lagrangians%20to%20Numerical%20Stencils%0A%0A%20%20%20%20This%20notebook%20demonstrates%20the%20complete%20pipeline%20from%20theoretical%20physics%20to%0A%20%20%20%20numerical%20simulation%20code%2C%20based%20on%20the%20scalar%20field%20dynamics%20in%0A%20%20%20%20%5BarXiv%3A1608.04408%5D(https%3A%2F%2Farxiv.org%2Fabs%2F1608.04408).%0A%0A%20%20%20%20%23%23%20The%20Physics%0A%0A%20%20%20%20Single-field%20inflation%20uses%20a%20scalar%20field%20with%20canonical%20kinetic%20term%3A%0A%0A%20%20%20%20%24%24%5Cmathcal%7BL%7D_%5Cphi%20%3D%20-%5Cfrac%7B1%7D%7B2%7Dg%5E%7B%5Cmu%5Cnu%7D%5Cpartial_%5Cmu%20%5Cphi%20%5Cpartial_%5Cnu%5Cphi%20-%20V(%5Cphi)%24%24%0A%0A%20%20%20%20The%20Klein-Gordon%20equation%20governing%20the%20field%20dynamics%20is%3A%0A%0A%20%20%20%20%24%24%5Cpartial_t%5E2%20%5Cphi%20-%20%5Cgamma%5E%7Bij%7D%5Cpartial_i%5Cpartial_j%5Cphi%20%2B%20%5Cfrac%7BdV%7D%7Bd%5Cphi%7D%20%3D%200%24%24%0A%0A%20%20%20%20We%20will%20derive%20this%20equation%20using%20variational%20calculus%2C%20then%20convert%20it%20to%0A%20%20%20%20finite%20difference%20form%20suitable%20for%20numerical%20simulation.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%201.%20Deriving%20the%20Klein-Gordon%20Equation%0A%0A%20%20%20%20Starting%20from%20the%20scalar%20field%20Lagrangian%20density%20in%20flat%20spacetime%3A%0A%0A%20%20%20%20%24%24%5Cmathcal%7BL%7D%20%3D%20%5Cfrac%7B1%7D%7B2%7D(%5Cpartial_t%5Cphi)%5E2%20-%20%5Cfrac%7B1%7D%7B2%7D(%5Cpartial_x%5Cphi)%5E2%20-%20V(%5Cphi)%24%24%0A%0A%20%20%20%20The%20Euler-Lagrange%20equation%20gives%20us%20the%20equation%20of%20motion.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Function%2C%20R%2C%20Simplify%2C%20Symbol%2C%20VariationalDerivative%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Define%20coordinates%20and%20field%0A%20%20%20%20x%2C%20t%20%3D%20symbols('x%20t')%0A%20%20%20%20phi%20%3D%20Function('phi')(x%2C%20t)%0A%20%20%20%20m%20%3D%20Symbol('m'%2C%20positive%3DTrue)%20%20%23%20mass%20parameter%0A%0A%20%20%20%20%23%20Klein-Gordon%20Lagrangian%3A%20L%20%3D%20(1%2F2)(d_t%20phi)%5E2%20-%20(1%2F2)(d_x%20phi)%5E2%20-%20(1%2F2)m%5E2%20phi%5E2%0A%20%20%20%20L_KG%20%3D%20R(1%2C%202)%20*%20D(phi%2C%20t)**2%20-%20R(1%2C%202)%20*%20D(phi%2C%20x)**2%20-%20R(1%2C%202)%20*%20m**2%20*%20phi**2%0A%0A%20%20%20%20%23%20Derive%20the%20equation%20of%20motion%0A%20%20%20%20eom_KG%20%3D%20VariationalDerivative(L_KG%2C%20phi%2C%20%5Bx%2C%20t%5D)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Klein-Gordon%20Lagrangian%3A**%0A%0A%20%20%20%20%24%5C%5Cmathcal%7B%7BL%7D%7D%20%3D%20%7BL_KG%7D%24%0A%0A%20%20%20%20**Equation%20of%20Motion**%20(%24%5C%5Cdelta%5C%5Cmathcal%7B%7BL%7D%7D%2F%5C%5Cdelta%5C%5Cphi%20%3D%200%24)%3A%0A%0A%20%20%20%20%24%7BSimplify(eom_KG)%7D%20%3D%200%24%0A%0A%20%20%20%20This%20is%20the%20Klein-Gordon%20equation%3A%20%24%5C%5Cpartial_t%5E2%5C%5Cphi%20-%20%5C%5Cpartial_x%5E2%5C%5Cphi%20%2B%20m%5E2%5C%5Cphi%20%3D%200%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%202.%20The%20Wave%20Equation%20(Massless%20Limit)%0A%0A%20%20%20%20Setting%20%24m%20%3D%200%24%20gives%20the%20wave%20equation%2C%20which%20appears%20in%20the%20gradient%20energy%0A%20%20%20%20dominated%20regime%20of%20inflation%20simulations.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Function%2C%20R%2C%20Simplify%2C%20VariationalDerivative%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Wave%20equation%20(massless%20Klein-Gordon)%0A%20%20%20%20x_w%2C%20t_w%20%3D%20symbols('x%20t')%0A%20%20%20%20phi_w%20%3D%20Function('phi')(x_w%2C%20t_w)%0A%0A%20%20%20%20L_wave%20%3D%20R(1%2C%202)%20*%20D(phi_w%2C%20t_w)**2%20-%20R(1%2C%202)%20*%20D(phi_w%2C%20x_w)**2%0A%0A%20%20%20%20eom_wave%20%3D%20VariationalDerivative(L_wave%2C%20phi_w%2C%20%5Bx_w%2C%20t_w%5D)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Wave%20Equation%20Lagrangian%3A**%0A%0A%20%20%20%20%24%5C%5Cmathcal%7B%7BL%7D%7D%20%3D%20%7BL_wave%7D%24%0A%0A%20%20%20%20**Equation%20of%20Motion%3A**%0A%0A%20%20%20%20%24%7BSimplify(eom_wave)%7D%20%3D%200%24%0A%0A%20%20%20%20This%20gives%3A%20%24%5C%5Cpartial_t%5E2%5C%5Cphi%20%3D%20%5C%5Cpartial_x%5E2%5C%5Cphi%24%20(the%20wave%20equation)%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%203.%20Converting%20to%20Finite%20Differences%0A%0A%20%20%20%20For%20numerical%20simulation%2C%20we%20need%20to%20discretize%20the%20spatial%20derivatives.%0A%20%20%20%20The%20%60Discretize%60%20function%20converts%20symbolic%20derivatives%20to%20finite%20difference%0A%20%20%20%20approximations%20using%20Taylor%20series%20matching.%0A%0A%20%20%20%20%23%23%23%20Central%20Difference%20Stencils%0A%0A%20%20%20%20For%20a%203-point%20central%20difference%3A%0A%20%20%20%20-%20First%20derivative%3A%20%24%5Cfrac%7B%5Cpartial%20f%7D%7B%5Cpartial%20x%7D%20%5Capprox%20%5Cfrac%7Bf(x%2Bh)%20-%20f(x-h)%7D%7B2h%7D%24%0A%20%20%20%20-%20Second%20derivative%3A%20%24%5Cfrac%7B%5Cpartial%5E2%20f%7D%7B%5Cpartial%20x%5E2%7D%20%5Capprox%20%5Cfrac%7Bf(x%2Bh)%20-%202f(x)%20%2B%20f(x-h)%7D%7Bh%5E2%7D%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20Symbol%2C%20mo)%3A%0A%20%20%20%20%23%20Discretization%20example%0A%20%20%20%20x_d%20%3D%20Symbol('x')%0A%20%20%20%20h%20%3D%20Symbol('h')%20%20%23%20grid%20spacing%0A%20%20%20%20f%20%3D%20Function('f')(x_d)%0A%0A%20%20%20%20%23%20Second%20derivative%0A%20%20%20%20d2f_dx2%20%3D%20D(f%2C%20(x_d%2C%202))%0A%0A%20%20%20%20%23%20Central%20difference%20discretization%0A%20%20%20%20stencil_3pt%20%3D%20Discretize(d2f_dx2%2C%20%7Bx_d%3A%20(%5Bx_d%20-%20h%2C%20x_d%2C%20x_d%20%2B%20h%5D%2C%20h)%7D)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Second%20Derivative%20Discretization%20(3-point%20central)%3A**%0A%0A%20%20%20%20Symbolic%3A%20%24%5C%5Cfrac%7B%7B%5C%5Cpartial%5E2%20f%7D%7D%7B%7B%5C%5Cpartial%20x%5E2%7D%7D%24%0A%0A%20%20%20%20Discretized%3A%20%24%7BSimplify(stencil_3pt)%7D%24%0A%0A%20%20%20%20This%20is%20the%20standard%203-point%20stencil%3A%20%24%5C%5Cfrac%7B%7Bf_%7B%7Bi%2B1%7D%7D%20-%202f_i%20%2B%20f_%7B%7Bi-1%7D%7D%7D%7D%7B%7Bh%5E2%7D%7D%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%204.%20Discretizing%20the%20Wave%20Equation%0A%0A%20%20%20%20Now%20let's%20discretize%20the%20full%20wave%20equation%20for%20numerical%20simulation.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Full%20wave%20equation%20discretization%0A%20%20%20%20x_full%2C%20t_full%20%3D%20symbols('x%20t')%0A%20%20%20%20hx%2C%20ht%20%3D%20symbols('h_x%20h_t')%20%20%23%20spatial%20and%20temporal%20grid%20spacing%0A%20%20%20%20u%20%3D%20Function('u')(x_full%2C%20t_full)%0A%0A%20%20%20%20%23%20Wave%20equation%3A%20d%5E2u%2Fdt%5E2%20-%20d%5E2u%2Fdx%5E2%20%3D%200%0A%20%20%20%20wave_eq%20%3D%20D(u%2C%20(t_full%2C%202))%20-%20D(u%2C%20(x_full%2C%202))%0A%0A%20%20%20%20%23%20Discretize%20both%20derivatives%0A%20%20%20%20step_map%20%3D%20%7B%0A%20%20%20%20%20%20%20%20x_full%3A%20(%5Bx_full%20-%20hx%2C%20x_full%2C%20x_full%20%2B%20hx%5D%2C%20hx)%2C%0A%20%20%20%20%20%20%20%20t_full%3A%20(%5Bt_full%20-%20ht%2C%20t_full%2C%20t_full%20%2B%20ht%5D%2C%20ht)%2C%0A%20%20%20%20%7D%0A%20%20%20%20wave_discrete%20%3D%20Discretize(wave_eq%2C%20step_map)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Wave%20Equation%3A**%0A%0A%20%20%20%20%24%5C%5Cfrac%7B%7B%5C%5Cpartial%5E2%20u%7D%7D%7B%7B%5C%5Cpartial%20t%5E2%7D%7D%20-%20%5C%5Cfrac%7B%7B%5C%5Cpartial%5E2%20u%7D%7D%7B%7B%5C%5Cpartial%20x%5E2%7D%7D%20%3D%200%24%0A%0A%20%20%20%20**Discretized%20Form%3A**%0A%0A%20%20%20%20%24%7BSimplify(wave_discrete)%7D%20%3D%200%24%0A%0A%20%20%20%20Rearranging%20for%20time-stepping%3A%20%24u(x%2C%20t%2Bh_t)%20%3D%202u(x%2Ct)%20-%20u(x%2C%20t-h_t)%20%2B%20%5C%5Cfrac%7B%7Bh_t%5E2%7D%7D%7B%7Bh_x%5E2%7D%7D%5Bu(x%2Bh_x%2Ct)%20-%202u(x%2Ct)%20%2B%20u(x-h_x%2Ct)%5D%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%205.%20Higher-Order%20Stencils%0A%0A%20%20%20%20For%20better%20accuracy%2C%20we%20can%20use%20wider%20stencils.%20The%20%60ToStencil%60%20function%0A%20%20%20%20automatically%20generates%20symmetric%20stencil%20points.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Function%2C%20Simplify%2C%20Symbol%2C%20ToStencil%2C%20mo)%3A%0A%20%20%20%20%23%20Higher-order%20stencils%0A%20%20%20%20x_ho%20%3D%20Symbol('x')%0A%20%20%20%20h_ho%20%3D%20Symbol('h')%0A%20%20%20%20f_ho%20%3D%20Function('f')(x_ho)%0A%0A%20%20%20%20%23%20Compare%203-point%20and%205-point%20stencils%20for%20first%20derivative%0A%20%20%20%20df_dx%20%3D%20D(f_ho%2C%20x_ho)%0A%0A%20%20%20%20stencil_3%20%3D%20ToStencil(df_dx%2C%20%7Bx_ho%3A%20h_ho%7D%2C%20width%3D3)%0A%20%20%20%20stencil_5%20%3D%20ToStencil(df_dx%2C%20%7Bx_ho%3A%20h_ho%7D%2C%20width%3D5)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**First%20Derivative%20Stencils%3A**%0A%0A%20%20%20%203-point%3A%20%24%7BSimplify(stencil_3)%7D%24%0A%0A%20%20%20%205-point%3A%20%24%7BSimplify(stencil_5)%7D%24%0A%0A%20%20%20%20The%205-point%20stencil%20provides%204th-order%20accuracy%20vs%202nd-order%20for%203-point.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%206.%20Code%20Generation%0A%0A%20%20%20%20The%20%60StencilCodeGen%60%20function%20converts%20discretized%20expressions%20to%20code%0A%20%20%20%20in%20multiple%20programming%20languages%20-%20useful%20for%20generating%20numerical%0A%20%20%20%20relativity%20simulation%20code.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20StencilCodeGen%2C%20Symbol%2C%20mo)%3A%0A%20%20%20%20%23%20Code%20generation%20example%0A%20%20%20%20x_cg%20%3D%20Symbol('x')%0A%20%20%20%20h_cg%20%3D%20Symbol('h')%0A%20%20%20%20phi_cg%20%3D%20Function('phi')(x_cg)%0A%0A%20%20%20%20%23%20Second%20derivative%20stencil%0A%20%20%20%20d2phi%20%3D%20D(phi_cg%2C%20(x_cg%2C%202))%0A%20%20%20%20stencil_cg%20%3D%20Discretize(d2phi%2C%20%7Bx_cg%3A%20(%5Bx_cg%20-%20h_cg%2C%20x_cg%2C%20x_cg%20%2B%20h_cg%5D%2C%20h_cg)%7D)%0A%0A%20%20%20%20%23%20Generate%20code%20in%20different%20languages%0A%20%20%20%20python_code%20%3D%20StencilCodeGen(stencil_cg%2C%20language%3D'python'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'phi'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%20%20%20%20c_code%20%3D%20StencilCodeGen(stencil_cg%2C%20language%3D'c'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'phi'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%20%20%20%20fortran_code%20%3D%20StencilCodeGen(stencil_cg%2C%20language%3D'fortran'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'phi'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%20%20%20%20latex_code%20%3D%20StencilCodeGen(stencil_cg%2C%20language%3D'latex')%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Code%20Generation%20for**%20%24%5C%5Cpartial%5E2%5C%5Cphi%2F%5C%5Cpartial%20x%5E2%24%3A%0A%0A%20%20%20%20**Python%3A**%0A%20%20%20%20%60%60%60python%0A%20%20%20%20d2phi_dx2%20%3D%20%7Bpython_code%7D%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20**C%3A**%0A%20%20%20%20%60%60%60c%0A%20%20%20%20double%20d2phi_dx2%20%3D%20%7Bc_code%7D%3B%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20**Fortran%3A**%0A%20%20%20%20%60%60%60fortran%0A%20%20%20%20d2phi_dx2%20%3D%20%7Bfortran_code%7D%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20**LaTeX%3A**%0A%20%20%20%20%24%7Blatex_code%7D%24%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%207.%20Complete%20Pipeline%3A%20Lagrangian%20to%20Simulation%20Code%0A%0A%20%20%20%20Using%20symderive's%20%60Pipe%60%20API%2C%20we%20can%20chain%20the%20entire%20workflow.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(%0A%20%20%20%20D%2C%0A%20%20%20%20Discretize%2C%0A%20%20%20%20Function%2C%0A%20%20%20%20Pipe%2C%0A%20%20%20%20R%2C%0A%20%20%20%20StencilCodeGen%2C%0A%20%20%20%20Symbol%2C%0A%20%20%20%20VariationalDerivative%2C%0A%20%20%20%20mo%2C%0A%20%20%20%20symbols%2C%0A)%3A%0A%20%20%20%20%23%20Complete%20pipeline%20example%0A%20%20%20%20x_pipe%2C%20t_pipe%20%3D%20symbols('x%20t')%0A%20%20%20%20h_pipe%20%3D%20Symbol('h')%0A%20%20%20%20psi%20%3D%20Function('psi')(x_pipe%2C%20t_pipe)%0A%0A%20%20%20%20%23%20Lagrangian%20for%20massless%20scalar%20field%20(1D)%0A%20%20%20%20L_pipe%20%3D%20R(1%2C%202)%20*%20D(psi%2C%20t_pipe)**2%20-%20R(1%2C%202)%20*%20D(psi%2C%20x_pipe)**2%0A%0A%20%20%20%20%23%20Pipeline%3A%20Lagrangian%20-%3E%20EoM%20-%3E%20Discretize%20spatial%20part%0A%20%20%20%20eom_spatial%20%3D%20(%0A%20%20%20%20%20%20%20%20Pipe(L_pipe)%0A%20%20%20%20%20%20%20%20.then(VariationalDerivative%2C%20psi%2C%20%5Bx_pipe%2C%20t_pipe%5D)%0A%20%20%20%20%20%20%20%20.then(Discretize%2C%20%7Bx_pipe%3A%20(%5Bx_pipe%20-%20h_pipe%2C%20x_pipe%2C%20x_pipe%20%2B%20h_pipe%5D%2C%20h_pipe)%7D)%0A%20%20%20%20%20%20%20%20.value%0A%20%20%20%20)%0A%0A%20%20%20%20%23%20Generate%20simulation%20code%0A%20%20%20%20sim_code%20%3D%20StencilCodeGen(eom_spatial%2C%20language%3D'python'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'psi'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Complete%20Pipeline%3A**%0A%0A%20%20%20%201.%20Start%20with%20Lagrangian%3A%20%24%5C%5Cmathcal%7B%7BL%7D%7D%20%3D%20%7BL_pipe%7D%24%0A%0A%20%20%20%202.%20Derive%20equation%20of%20motion%20via%20variational%20derivative%0A%0A%20%20%20%203.%20Discretize%20spatial%20derivatives%0A%0A%20%20%20%204.%20Generate%20Python%20code%3A%0A%0A%20%20%20%20%60%60%60python%0A%20%20%20%20%23%20Equation%20of%20motion%20(set%20to%20zero%20and%20solve%20for%20time%20evolution)%0A%20%20%20%20eom%20%3D%20%7Bsim_code%7D%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20This%20can%20be%20directly%20used%20in%20a%20time-stepping%20numerical%20scheme!%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%208.%20Application%3A%20Initial%20Conditions%20from%20arXiv%3A1608.04408%0A%0A%20%20%20%20The%20paper%20uses%20inhomogeneous%20initial%20conditions%20for%20the%20scalar%20field%3A%0A%0A%20%20%20%20%24%24%5Cphi(t%3D0%2C%20%5Cmathbf%7Bx%7D)%20%3D%20%5Cphi_0%20%2B%20%5Cfrac%7B%5CDelta%5Cphi%7D%7BN%7D%5Csum_%7Bn%3D1%7D%5E%7BN%7D%5Cleft(%5Ccos%5Cfrac%7B2%5Cpi%20nx%7D%7BL%7D%20%2B%20%5Ccos%5Cfrac%7B2%5Cpi%20ny%7D%7BL%7D%20%2B%20%5Ccos%5Cfrac%7B2%5Cpi%20nz%7D%7BL%7D%5Cright)%24%24%0A%0A%20%20%20%20With%20initial%20velocity%3A%0A%20%20%20%20%24%24%5Cfrac%7B%5Cpartial%5Cphi(t%3D0%2C%20%5Cmathbf%7Bx%7D)%7D%7B%5Cpartial%20t%7D%20%3D%200%24%24%0A%0A%20%20%20%20The%20gradient%20energy%20density%20is%3A%0A%20%20%20%20%24%24%5Crho_%7B%5Cmathrm%7Bgrad%7D%7D%20%3D%20%5Cfrac%7B1%7D%7B2%7D%5Cgamma%5E%7Bij%7D%5Cpartial_i%5Cphi%5Cpartial_j%5Cphi%24%24%0A%0A%20%20%20%20Let's%20compute%20the%20discretized%20gradient%20energy.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20R%2C%20Simplify%2C%20Symbol%2C%20mo)%3A%0A%20%20%20%20%23%20Gradient%20energy%20density%20discretization%0A%20%20%20%20x_grad%20%3D%20Symbol('x')%0A%20%20%20%20h_grad%20%3D%20Symbol('h')%0A%20%20%20%20phi_grad%20%3D%20Function('phi')(x_grad)%0A%0A%20%20%20%20%23%20Gradient%20energy%3A%20(1%2F2)(d%20phi%2Fdx)%5E2%0A%20%20%20%20rho_grad%20%3D%20R(1%2C%202)%20*%20D(phi_grad%2C%20x_grad)**2%0A%0A%20%20%20%20%23%20Discretize%0A%20%20%20%20rho_grad_discrete%20%3D%20Discretize(rho_grad%2C%20%7Bx_grad%3A%20(%5Bx_grad%20-%20h_grad%2C%20x_grad%2C%20x_grad%20%2B%20h_grad%5D%2C%20h_grad)%7D)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Gradient%20Energy%20Density%20(1D)%3A**%0A%0A%20%20%20%20Continuous%3A%20%24%5C%5Crho_%7B%7B%5C%5Cmathrm%7B%7Bgrad%7D%7D%7D%7D%20%3D%20%5C%5Cfrac%7B%7B1%7D%7D%7B%7B2%7D%7D%5C%5Cleft(%5C%5Cfrac%7B%7B%5C%5Cpartial%5C%5Cphi%7D%7D%7B%7B%5C%5Cpartial%20x%7D%7D%5C%5Cright)%5E2%24%0A%0A%20%20%20%20Discretized%3A%20%24%5C%5Crho_%7B%7B%5C%5Cmathrm%7B%7Bgrad%7D%7D%7D%7D%20%3D%20%7BSimplify(rho_grad_discrete)%7D%24%0A%0A%20%20%20%20This%20uses%20the%20central%20difference%20approximation%20for%20the%20first%20derivative.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%209.%20Multi-Dimensional%20Stencils%3A%20The%203D%20Laplacian%0A%0A%20%20%20%20Numerical%20relativity%20simulations%20operate%20in%203D.%20The%20spatial%20Laplacian%0A%20%20%20%20appears%20in%20the%20Hamiltonian%20constraint%20and%20evolution%20equations%3A%0A%0A%20%20%20%20%24%24%5Cnabla%5E2%20%5Cphi%20%3D%20%5Cfrac%7B%5Cpartial%5E2%5Cphi%7D%7B%5Cpartial%20x%5E2%7D%20%2B%20%5Cfrac%7B%5Cpartial%5E2%5Cphi%7D%7B%5Cpartial%20y%5E2%7D%20%2B%20%5Cfrac%7B%5Cpartial%5E2%5Cphi%7D%7B%5Cpartial%20z%5E2%7D%24%24%0A%0A%20%20%20%20Each%20term%20uses%20the%20standard%20second%20derivative%20stencil.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20StencilCodeGen%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%203D%20Laplacian%0A%20%20%20%20x_3d%2C%20y_3d%2C%20z_3d%20%3D%20symbols('x%20y%20z')%0A%20%20%20%20hx_3d%2C%20hy_3d%2C%20hz_3d%20%3D%20symbols('h_x%20h_y%20h_z')%0A%20%20%20%20phi_3d%20%3D%20Function('phi')(x_3d%2C%20y_3d%2C%20z_3d)%0A%0A%20%20%20%20%23%20Laplacian%20in%203D%0A%20%20%20%20laplacian_3d%20%3D%20D(phi_3d%2C%20(x_3d%2C%202))%20%2B%20D(phi_3d%2C%20(y_3d%2C%202))%20%2B%20D(phi_3d%2C%20(z_3d%2C%202))%0A%0A%20%20%20%20%23%20Discretize%20with%20uniform%20spacing%20h%0A%20%20%20%20h_uniform%20%3D%20symbols('h')%0A%20%20%20%20step_map_3d%20%3D%20%7B%0A%20%20%20%20%20%20%20%20x_3d%3A%20(%5Bx_3d%20-%20h_uniform%2C%20x_3d%2C%20x_3d%20%2B%20h_uniform%5D%2C%20h_uniform)%2C%0A%20%20%20%20%20%20%20%20y_3d%3A%20(%5By_3d%20-%20h_uniform%2C%20y_3d%2C%20y_3d%20%2B%20h_uniform%5D%2C%20h_uniform)%2C%0A%20%20%20%20%20%20%20%20z_3d%3A%20(%5Bz_3d%20-%20h_uniform%2C%20z_3d%2C%20z_3d%20%2B%20h_uniform%5D%2C%20h_uniform)%2C%0A%20%20%20%20%7D%0A%20%20%20%20laplacian_discrete%20%3D%20Discretize(laplacian_3d%2C%20step_map_3d)%0A%0A%20%20%20%20%23%20Generate%20C%20code%20for%20the%20stencil%0A%20%20%20%20c_laplacian%20%3D%20StencilCodeGen(%0A%20%20%20%20%20%20%20%20Discretize(D(phi_3d%2C%20(x_3d%2C%202))%2C%20%7Bx_3d%3A%20(%5Bx_3d%20-%20h_uniform%2C%20x_3d%2C%20x_3d%20%2B%20h_uniform%5D%2C%20h_uniform)%7D)%2C%0A%20%20%20%20%20%20%20%20language%3D'c'%2C%20array_name%3D'phi'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx'%0A%20%20%20%20)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**3D%20Laplacian%20Stencil%3A**%0A%0A%20%20%20%20Continuous%3A%20%24%5C%5Cnabla%5E2%5C%5Cphi%20%3D%20%7Blaplacian_3d%7D%24%0A%0A%20%20%20%20Discretized%20(uniform%20grid%20%24h%24)%3A%20%24%7BSimplify(laplacian_discrete)%7D%24%0A%0A%20%20%20%20This%20is%20the%20classic%207-point%20stencil%20used%20in%20numerical%20relativity%20codes.%0A%0A%20%20%20%20**C%20code%20for%20one%20dimension**%20(repeat%20for%20j%2C%20k%20indices)%3A%0A%20%20%20%20%60%60%60c%0A%20%20%20%20d2phi_dx2%20%3D%20%7Bc_laplacian%7D%3B%0A%20%20%20%20%60%60%60%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%2010.%20Stencil%20Accuracy%3A%202nd%20vs%204th%20vs%206th%20Order%0A%0A%20%20%20%20Higher-order%20stencils%20use%20more%20points%20but%20converge%20faster.%0A%20%20%20%20The%20truncation%20error%20scales%20as%20%24O(h%5En)%24%20where%20%24n%24%20is%20the%20accuracy%20order.%0A%0A%20%20%20%20%7C%20Order%20%7C%20Points%20%7C%20Error%20Scaling%20%7C%20Use%20Case%20%7C%0A%20%20%20%20%7C-------%7C--------%7C---------------%7C----------%7C%0A%20%20%20%20%7C%202nd%20%20%20%7C%203%20%20%20%20%20%20%7C%20%24O(h%5E2)%24%20%20%20%20%20%20%7C%20Simple%20problems%2C%20quick%20tests%20%7C%0A%20%20%20%20%7C%204th%20%20%20%7C%205%20%20%20%20%20%20%7C%20%24O(h%5E4)%24%20%20%20%20%20%20%7C%20Production%20simulations%20%7C%0A%20%20%20%20%7C%206th%20%20%20%7C%207%20%20%20%20%20%20%7C%20%24O(h%5E6)%24%20%20%20%20%20%20%7C%20High-precision%20work%20%7C%0A%20%20%20%20%7C%208th%20%20%20%7C%209%20%20%20%20%20%20%7C%20%24O(h%5E8)%24%20%20%20%20%20%20%7C%20Spectral-like%20accuracy%20%7C%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Function%2C%20Simplify%2C%20Symbol%2C%20ToStencil%2C%20mo)%3A%0A%20%20%20%20%23%20Compare%20stencil%20orders%20for%20second%20derivative%0A%20%20%20%20x_ord%20%3D%20Symbol('x')%0A%20%20%20%20h_ord%20%3D%20Symbol('h')%0A%20%20%20%20f_ord%20%3D%20Function('f')(x_ord)%0A%0A%20%20%20%20d2f%20%3D%20D(f_ord%2C%20(x_ord%2C%202))%0A%0A%20%20%20%20stencil_2nd%20%3D%20ToStencil(d2f%2C%20%7Bx_ord%3A%20h_ord%7D%2C%20width%3D3)%0A%20%20%20%20stencil_4th%20%3D%20ToStencil(d2f%2C%20%7Bx_ord%3A%20h_ord%7D%2C%20width%3D5)%0A%20%20%20%20stencil_6th%20%3D%20ToStencil(d2f%2C%20%7Bx_ord%3A%20h_ord%7D%2C%20width%3D7)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Second%20Derivative%20Stencils%20by%20Order%3A**%0A%0A%20%20%20%20**2nd%20order%20(3-point)%3A**%0A%20%20%20%20%24%7BSimplify(stencil_2nd)%7D%24%0A%0A%20%20%20%20**4th%20order%20(5-point)%3A**%0A%20%20%20%20%24%7BSimplify(stencil_4th)%7D%24%0A%0A%20%20%20%20**6th%20order%20(7-point)%3A**%0A%20%20%20%20%24%7BSimplify(stencil_6th)%7D%24%0A%0A%20%20%20%20Note%20how%20higher-order%20stencils%20have%20smaller%20coefficients%20on%20the%20outer%20points%2C%0A%20%20%20%20reducing%20numerical%20dispersion.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%2011.%20Mixed%20Derivatives%20and%20Cross%20Terms%0A%0A%20%20%20%20The%20BSSN%20formulation%20of%20general%20relativity%20includes%20mixed%20partial%20derivatives%0A%20%20%20%20like%20%24%5Cpartial_x%5Cpartial_y%5Cphi%24.%20These%20require%202D%20stencils.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Mixed%20derivative%0A%20%20%20%20x_mix%2C%20y_mix%20%3D%20symbols('x%20y')%0A%20%20%20%20h_mix%20%3D%20symbols('h')%0A%20%20%20%20phi_mix%20%3D%20Function('phi')(x_mix%2C%20y_mix)%0A%0A%20%20%20%20%23%20Mixed%20partial%20derivative%0A%20%20%20%20d2_dxdy%20%3D%20D(D(phi_mix%2C%20x_mix)%2C%20y_mix)%0A%0A%20%20%20%20%23%209-point%20stencil%20for%20mixed%20derivative%0A%20%20%20%20step_map_mix%20%3D%20%7B%0A%20%20%20%20%20%20%20%20x_mix%3A%20(%5Bx_mix%20-%20h_mix%2C%20x_mix%2C%20x_mix%20%2B%20h_mix%5D%2C%20h_mix)%2C%0A%20%20%20%20%20%20%20%20y_mix%3A%20(%5By_mix%20-%20h_mix%2C%20y_mix%2C%20y_mix%20%2B%20h_mix%5D%2C%20h_mix)%2C%0A%20%20%20%20%7D%0A%20%20%20%20mixed_stencil%20%3D%20Discretize(d2_dxdy%2C%20step_map_mix)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Mixed%20Partial%20Derivative%3A**%0A%0A%20%20%20%20Continuous%3A%20%24%5C%5Cfrac%7B%7B%5C%5Cpartial%5E2%5C%5Cphi%7D%7D%7B%7B%5C%5Cpartial%20x%5C%5Cpartial%20y%7D%7D%24%0A%0A%20%20%20%20Discretized%3A%20%24%7BSimplify(mixed_stencil)%7D%24%0A%0A%20%20%20%20This%20is%20the%20standard%204-corner%20stencil%3A%0A%20%20%20%20%60%60%60%0A%20%20%20%20(-1)-----(0)-----(%2B1)%0A%20%20%20%20%20%20%7C%20%20%20%20%20%20%20%7C%20%20%20%20%20%20%20%7C%0A%20%20%20%20(%2B1)%20%20%20%20(0%2C0)%20%20%20(%2B1)%0A%20%20%20%20%20%20%7C%20%20%20%20%20%20%20%7C%20%20%20%20%20%20%20%7C%0A%20%20%20%20(-1)-----(0)-----(%2B1)%0A%20%20%20%20%60%60%60%0A%20%20%20%20Only%20the%20four%20corners%20contribute%20(with%20alternating%20signs).%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%2012.%20The%20ADM%20Evolution%20Equations%0A%0A%20%20%20%20In%20numerical%20relativity%2C%20the%20ADM%20(Arnowitt-Deser-Misner)%20formalism%20evolves%0A%20%20%20%20the%203-metric%20%24%5Cgamma_%7Bij%7D%24%20and%20extrinsic%20curvature%20%24K_%7Bij%7D%24%3A%0A%0A%20%20%20%20%24%24%5Cpartial_t%20%5Cgamma_%7Bij%7D%20%3D%20-2%5Calpha%20K_%7Bij%7D%20%2B%20%5Cmathcal%7BL%7D_%5Cbeta%20%5Cgamma_%7Bij%7D%24%24%0A%0A%20%20%20%20%24%24%5Cpartial_t%20K_%7Bij%7D%20%3D%20-D_i%20D_j%20%5Calpha%20%2B%20%5Calpha(R_%7Bij%7D%20%2B%20K%20K_%7Bij%7D%20-%202K_%7Bik%7DK%5Ek%7B%7D_j)%20%2B%20%5Cmathcal%7BL%7D_%5Cbeta%20K_%7Bij%7D%24%24%0A%0A%20%20%20%20The%20spatial%20Ricci%20tensor%20%24R_%7Bij%7D%24%20contains%20second%20derivatives%20of%20the%20metric%2C%0A%20%20%20%20which%20we%20discretize%20using%20the%20stencils%20developed%20above.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20StencilCodeGen%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Simplified%20ADM-like%20term%3A%20second%20derivative%20of%20metric%20component%0A%20%20%20%20x_adm%2C%20y_adm%20%3D%20symbols('x%20y')%0A%20%20%20%20h_adm%20%3D%20symbols('h')%0A%20%20%20%20gamma_xx%20%3D%20Function('gamma_xx')(x_adm%2C%20y_adm)%0A%0A%20%20%20%20%23%20Part%20of%20the%20Ricci%20tensor%3A%20d%5E2(gamma_xx)%2Fdx%5E2%0A%20%20%20%20ricci_term%20%3D%20D(gamma_xx%2C%20(x_adm%2C%202))%0A%0A%20%20%20%20%23%20Discretize%0A%20%20%20%20step_adm%20%3D%20%7Bx_adm%3A%20(%5Bx_adm%20-%20h_adm%2C%20x_adm%2C%20x_adm%20%2B%20h_adm%5D%2C%20h_adm)%7D%0A%20%20%20%20ricci_discrete%20%3D%20Discretize(ricci_term%2C%20step_adm)%0A%0A%20%20%20%20%23%20Generate%20code%0A%20%20%20%20ricci_code%20%3D%20StencilCodeGen(ricci_discrete%2C%20language%3D'c'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'gamma_xx'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Ricci%20Tensor%20Component%20(simplified)%3A**%0A%0A%20%20%20%20One%20term%20in%20%24R_%7B%7Bxx%7D%7D%24%20involves%3A%20%24%5C%5Cfrac%7B%7B%5C%5Cpartial%5E2%20%5C%5Cgamma_%7B%7Bxx%7D%7D%7D%7D%7B%7B%5C%5Cpartial%20x%5E2%7D%7D%24%0A%0A%20%20%20%20Discretized%3A%20%24%7BSimplify(ricci_discrete)%7D%24%0A%0A%20%20%20%20**C%20code%3A**%0A%20%20%20%20%60%60%60c%0A%20%20%20%20double%20d2gamma_dx2%20%3D%20%7Bricci_code%7D%3B%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20The%20full%20Ricci%20tensor%20has%20many%20such%20terms%20for%20each%20metric%20component.%0A%20%20%20%20Derive%20automates%20the%20error-prone%20process%20of%20deriving%20and%20coding%20each%20stencil.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%2013.%20Using%20the%20Metric%20Class%20for%20Tensor%20Calculus%0A%0A%20%20%20%20symderive's%20%60diffgeo%60%20module%20provides%20a%20%60Metric%60%20class%20for%20full%20tensor%20calculus.%0A%20%20%20%20This%20computes%20Christoffel%20symbols%2C%20Riemann%20tensor%2C%20and%20Ricci%20tensor%20automatically.%0A%0A%20%20%20%20For%20numerical%20relativity%2C%20these%20tensors%20contain%20the%20derivatives%20we%20need%20to%20discretize.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(Metric%2C%20Simplify%2C%20Symbol%2C%20symbols)%3A%0A%20%20%20%20%23%20Define%20a%20general%202D%20metric%20(simplified%20example)%0A%20%20%20%20x_m%2C%20y_m%20%3D%20symbols('x%20y'%2C%20real%3DTrue)%0A%0A%20%20%20%20%23%20Metric%20components%20as%20symbols%20(would%20be%20functions%20in%20full%20NR)%0A%20%20%20%20g_xx%20%3D%20Symbol('g_xx'%2C%20positive%3DTrue)%0A%20%20%20%20g_yy%20%3D%20Symbol('g_yy'%2C%20positive%3DTrue)%0A%0A%20%20%20%20%23%20Create%20a%20diagonal%20metric%3A%20ds%5E2%20%3D%20g_xx%20dx%5E2%20%2B%20g_yy%20dy%5E2%0A%20%20%20%20%23%20Metric(coords%2C%20components)%0A%20%20%20%20metric_2d%20%3D%20Metric(%0A%20%20%20%20%20%20%20%20%5Bx_m%2C%20y_m%5D%2C%0A%20%20%20%20%20%20%20%20%5B%5Bg_xx%2C%200%5D%2C%0A%20%20%20%20%20%20%20%20%20%5B0%2C%20g_yy%5D%5D%0A%20%20%20%20)%0A%0A%20%20%20%20%23%20Christoffel%20symbols%20contain%20first%20derivatives%20of%20the%20metric%0A%20%20%20%20christoffels%20%3D%20metric_2d.christoffel_second_kind()%0A%0A%20%20%20%20%23%20Example%3A%20Gamma%5Ex_xx%20%3D%20(1%2F2)%20g%5Exx%20*%20d(g_xx)%2Fdx%0A%20%20%20%20Gamma_x_xx%20%3D%20Simplify(christoffels%5B0%2C%200%2C%200%5D)%0A%0A%20%20%20%20(metric_2d.g%2C%20Gamma_x_xx)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Discretize%2C%20Function%2C%20Simplify%2C%20StencilCodeGen%2C%20mo%2C%20symbols)%3A%0A%20%20%20%20%23%20Now%20discretize%20a%20Christoffel-like%20term%0A%20%20%20%20%23%20Gamma%5Ei_jk%20involves%20d(g_jk)%2Fdx%5Ei%2C%20so%20we%20discretize%20metric%20derivatives%0A%20%20%20%20x_chr%2C%20y_chr%20%3D%20symbols('x%20y')%0A%20%20%20%20h_chr%20%3D%20symbols('h')%0A%0A%20%20%20%20%23%20Metric%20component%20as%20a%20function%20of%20position%0A%20%20%20%20g_xx_func%20%3D%20Function('g_xx')(x_chr%2C%20y_chr)%0A%0A%20%20%20%20%23%20The%20derivative%20that%20appears%20in%20Christoffel%20symbols%0A%20%20%20%20dg_dx%20%3D%20D(g_xx_func%2C%20x_chr)%0A%0A%20%20%20%20%23%20Discretize%20using%20central%20difference%0A%20%20%20%20step_chr%20%3D%20%7Bx_chr%3A%20(%5Bx_chr%20-%20h_chr%2C%20x_chr%2C%20x_chr%20%2B%20h_chr%5D%2C%20h_chr)%7D%0A%20%20%20%20dg_discrete%20%3D%20Discretize(dg_dx%2C%20step_chr)%0A%0A%20%20%20%20%23%20Code%20for%20Christoffel%20computation%0A%20%20%20%20chr_code%20%3D%20StencilCodeGen(dg_discrete%2C%20language%3D'c'%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20array_name%3D'g_xx'%2C%20index_var%3D'i'%2C%20spacing_name%3D'dx')%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Discretizing%20Christoffel%20Symbol%20Terms%3A**%0A%0A%20%20%20%20Christoffel%20symbols%20involve%20metric%20derivatives%20like%20%24%5C%5Cpartial_x%20g_%7B%7Bxx%7D%7D%24%0A%0A%20%20%20%20Discretized%3A%20%24%7BSimplify(dg_discrete)%7D%24%0A%0A%20%20%20%20**C%20code%3A**%0A%20%20%20%20%60%60%60c%0A%20%20%20%20%2F%2F%20dg_xx%2Fdx%20for%20Christoffel%20computation%0A%20%20%20%20double%20dg_dx%20%3D%20%7Bchr_code%7D%3B%0A%0A%20%20%20%20%2F%2F%20Then%20Gamma%5Ex_xx%20%3D%200.5%20*%20g_xx_inv%20*%20dg_dx%0A%20%20%20%20double%20Gamma_x_xx%20%3D%200.5%20*%20g_xx_inv%20*%20dg_dx%3B%0A%20%20%20%20%60%60%60%0A%0A%20%20%20%20The%20full%20BSSN%20evolution%20uses%20these%20discretized%20Christoffels%20throughout.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%2014.%20Kreiss-Oliger%20Dissipation%0A%0A%20%20%20%20Numerical%20simulations%20often%20add%20artificial%20dissipation%20to%20control%0A%20%20%20%20high-frequency%20noise.%20The%20Kreiss-Oliger%20operator%20for%204th%20order%20schemes%3A%0A%0A%20%20%20%20%24%24%5Cepsilon%20(-1)%5E%7Br%2B1%7D%20h%5E%7B2r-1%7D%20D_%2B%5Er%20D_-%5Er%20%5Cphi%24%24%0A%0A%20%20%20%20where%20%24D_%2B%24%20and%20%24D_-%24%20are%20forward%2Fbackward%20difference%20operators.%0A%20%20%20%20For%20%24r%3D2%24%20this%20gives%20a%205-point%20dissipation%20stencil.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(D%2C%20Function%2C%20Simplify%2C%20Symbol%2C%20ToStencil%2C%20mo)%3A%0A%20%20%20%20%23%20Kreiss-Oliger%20dissipation%20(4th%20derivative%20for%204th%20order%20scheme)%0A%20%20%20%20x_ko%20%3D%20Symbol('x')%0A%20%20%20%20h_ko%20%3D%20Symbol('h')%0A%20%20%20%20eps%20%3D%20Symbol('epsilon')%0A%20%20%20%20phi_ko%20%3D%20Function('phi')(x_ko)%0A%0A%20%20%20%20%23%20Fourth%20derivative%20(appears%20in%20Kreiss-Oliger%20for%204th%20order%20methods)%0A%20%20%20%20d4_phi%20%3D%20D(phi_ko%2C%20(x_ko%2C%204))%0A%0A%20%20%20%20%23%20Discretize%20with%205-point%20stencil%0A%20%20%20%20ko_stencil%20%3D%20ToStencil(d4_phi%2C%20%7Bx_ko%3A%20h_ko%7D%2C%20width%3D5)%0A%0A%20%20%20%20mo.md(f%22%22%22%0A%20%20%20%20**Kreiss-Oliger%20Dissipation%20Term%3A**%0A%0A%20%20%20%20Fourth%20derivative%3A%20%24%5C%5Cfrac%7B%7B%5C%5Cpartial%5E4%5C%5Cphi%7D%7D%7B%7B%5C%5Cpartial%20x%5E4%7D%7D%24%0A%0A%20%20%20%20Discretized%3A%20%24%7BSimplify(ko_stencil)%7D%24%0A%0A%20%20%20%20Applied%20as%3A%20%24%5C%5Cphi%20%5C%5Cleftarrow%20%5C%5Cphi%20-%20%5C%5Cepsilon%20%5C%5Ccdot%20h%5E3%20%5C%5Ccdot%20(%5C%5Ctext%7B%7Bstencil%7D%7D)%24%0A%0A%20%20%20%20The%20coefficient%20%24%5C%5Cepsilon%20%5C%5Csim%200.1%24%20controls%20dissipation%20strength.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Summary%0A%0A%20%20%20%20This%20notebook%20demonstrated%20the%20complete%20pipeline%20from%20theoretical%20physics%0A%20%20%20%20to%20production-ready%20numerical%20simulation%20code%3A%0A%0A%20%20%20%20**Physics%20to%20Code%20Pipeline%3A**%0A%20%20%20%201.%20**Variational%20Calculus**%3A%20Derive%20equations%20of%20motion%20from%20Lagrangians%0A%20%20%20%202.%20**Discretization**%3A%20Convert%20PDEs%20to%20finite%20difference%20stencils%0A%20%20%20%203.%20**Multi-dimensional**%3A%20Handle%203D%20Laplacians%2C%20mixed%20derivatives%0A%20%20%20%204.%20**Code%20Generation**%3A%20Output%20C%2C%20Fortran%2C%20or%20Python%20for%20direct%20use%0A%0A%20%20%20%20**Numerical%20Relativity%20Applications%3A**%0A%20%20%20%20-%203D%20Laplacian%20(7-point%20stencil)%20for%20constraint%20equations%0A%20%20%20%20-%20Mixed%20derivatives%20for%20BSSN%20formulation%0A%20%20%20%20-%20ADM%20evolution%20equations%20with%20Ricci%20tensor%20terms%0A%20%20%20%20-%20Metric%20class%20for%20Christoffel%20symbols%20and%20curvature%20tensors%0A%20%20%20%20-%20Kreiss-Oliger%20dissipation%20for%20numerical%20stability%0A%0A%20%20%20%20**Accuracy%20Control%3A**%0A%20%20%20%20-%202nd%20order%20(3-point)%20to%208th%20order%20(9-point)%20stencils%0A%20%20%20%20-%20Error%20scaling%20from%20%24O(h%5E2)%24%20to%20%24O(h%5E8)%24%0A%0A%20%20%20%20This%20workflow%20is%20directly%20applicable%20to%20numerical%20relativity%20codes%20like%0A%20%20%20%20GRChombo%20(used%20in%20arXiv%3A1608.04408)%20for%20simulating%20inflation%20with%0A%20%20%20%20inhomogeneous%20initial%20conditions.%0A%0A%20%20%20%20%23%23%23%20Key%20Functions%3A%0A%20%20%20%20-%20%60VariationalDerivative(L%2C%20field%2C%20coords)%60%20-%20Euler-Lagrange%20equations%0A%20%20%20%20-%20%60Discretize(expr%2C%20step_map)%60%20-%20Finite%20difference%20conversion%0A%20%20%20%20-%20%60ToStencil(expr%2C%20spacing%2C%20width)%60%20-%20Auto-generate%20stencil%20points%0A%20%20%20%20-%20%60StencilCodeGen(expr%2C%20language)%60%20-%20C%2FFortran%2FPython%20code%20generation%0A%20%20%20%20-%20%60Metric(components%2C%20coords)%60%20-%20Tensor%20calculus%20with%20auto%20Christoffels%0A%20%20%20%20-%20%60Pipe(expr).then(...)%60%20-%20Composable%20API%20for%20chaining%20operations%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A
e603859ade21b58956ff4d13ffbe30d8