2b:null 2c:["$","title",null,{"children":"Elastic-Pendulum-Solver with Runge-Kutta | Chuncheng | Observable"}] 31:T1073,// 使用 Runge-Kutta 4阶方法求解 class ElasticPendulumSolver { constructor(params) { // 物理参数 this.m = params.m || 1.0; // 质量 (kg) this.g = params.g || 9.8; // 重力加速度 (m/s²) this.k = params.k || 20.0; // 弹簧系数 (N/m) this.r0 = params.r0 || 1.0; // 弹簧原长 (m) // 初始条件 this.initialConditions = params.initialConditions || { r: 1.2, // 初始长度 (m) theta: 0.5, // 初始角度 (rad) v: 0.0, // 初始径向速度 (m/s) omega: 0.0 // 初始角速度 (rad/s) }; // 数值求解参数 this.dt = params.dt || 0.01; // 时间步长 (s) this.maxTime = params.maxTime || 10; // 总模拟时间 (s) } // 定义微分方程组 derivatives(t, state) { const [r, theta, v, omega] = state; // 方程1: m*d²r/dt² - m*r*(dθ/dt)² - m*g*cosθ + k*(r - r0) = 0 const r_ddot = r * omega * omega + this.g * Math.cos(theta) - (this.k / this.m) * (r - this.r0); // 方程2: r*d²θ/dt² + 2*(dr/dt)*(dθ/dt) + g*sinθ = 0 const theta_ddot = (-2 * v * omega - this.g * Math.sin(theta)) / r; return [v, omega, r_ddot, theta_ddot]; } // 4阶Runge-Kutta方法 rungeKuttaStep(t, state, dt) { const k1 = this.derivatives(t, state); const k2 = this.derivatives(t + dt/2, state.map((s, i) => s + k1[i] * dt/2)); const k3 = this.derivatives(t + dt/2, state.map((s, i) => s + k2[i] * dt/2)); const k4 = this.derivatives(t + dt, state.map((s, i) => s + k3[i] * dt)); return state.map((s, i) => s + (k1[i] + 2*k2[i] + 2*k3[i] + k4[i]) * dt / 6); } // 求解整个时间序列 solve() { const results = { time: [], r: [], theta: [], v: [], omega: [], x: [], y: [] }; // 初始状态向量: [r, theta, v, omega] let state = [ this.initialConditions.r, this.initialConditions.theta, this.initialConditions.v, this.initialConditions.omega ]; let t = 0; const numSteps = Math.floor(this.maxTime / this.dt); for (let i = 0; i <= numSteps; i++) { // 保存当前状态 results.time.push(t); results.r.push(state[0]); results.theta.push(state[1]); results.v.push(state[2]); results.omega.push(state[3]); // 计算笛卡尔坐标（用于可视化） results.x.push(state[0] * Math.sin(state[1])); results.y.push(-state[0] * Math.cos(state[1])); // y向上为负 // 如果到达最大时间，停止计算 if (t >= this.maxTime) break; // RK4步进 state = this.rungeKuttaStep(t, state, this.dt); t += this.dt; } return results; } // 能量计算（用于验证结果） calculateEnergy(results) { const kineticEnergy = []; const potentialEnergy = []; const totalEnergy = []; for (let i = 0; i < results.time.length; i++) { const r = results.r[i]; const theta = results.theta[i]; const v = results.v[i]; const omega = results.omega[i]; // 动能: 0.5*m*(v_r² + v_theta²) const KE = 0.5 * this.m * (v*v + (r * omega)*(r * omega)); // 势能: 重力势能 + 弹性势能 const PE_gravity = -this.m * this.g * r * Math.cos(theta); // 参考点设在支点 const PE_spring = 0.5 * this.k * Math.pow(r - this.r0, 2); const PE = PE_gravity + PE_spring; kineticEnergy.push(KE); potentialEnergy.push(PE); totalEnergy.push(KE + PE); } return { kineticEnergy, potentialEnergy, totalEnergy }; } } 32:[] 2f:{"id":"3c1d0e036806f3b7","slug":"elastic-pendulum-solver-with-runge-kutta","title":"Elastic-Pendulum-Solver with Runge-Kutta","thumbnail":"39213f55c14b8e929146ec9f3d9405799f7da85994024ecfef74b006b24d39de","access_level":"public","create_time":"$D2025-09-10T07:50:16.660Z","update_time":"$D2025-09-10T09:03:53.596Z","view_time":null,"first_public_time":"$D2025-09-10T08:56:37.554Z","fork":false,"stars":2,"starred":false,"owner":{"type":"team","id":"bcd65cef1cb4d2ed","login":"listenzcc","name":"Chuncheng","avatar_url":"https://avatars-next.observableusercontent.com/bb3c68e77a708a587ffc2d663bb83421f20d67991dbe6c1623c6e5d2fb73d15e","tier":"starter_2024"},"version":213,"publish_version":null,"role":null,"authored":false,"authors":[{"type":"individual","id":"16328058352682a8","login":"listenzcc","name":"Chuncheng","avatar_url":"https://avatars-next.observableusercontent.com/bb3c68e77a708a587ffc2d663bb83421f20d67991dbe6c1623c6e5d2fb73d15e","tier":null}],"members":[{"type":"individual","id":"16328058352682a8","login":"listenzcc","name":"Chuncheng","avatar_url":"https://avatars-next.observableusercontent.com/bb3c68e77a708a587ffc2d663bb83421f20d67991dbe6c1623c6e5d2fb73d15e","tier":null,"role":"editor"},{"type":"team","id":"bcd65cef1cb4d2ed","login":"listenzcc","name":"Chuncheng","avatar_url":"https://avatars-next.observableusercontent.com/bb3c68e77a708a587ffc2d663bb83421f20d67991dbe6c1623c6e5d2fb73d15e","tier":"starter_2024","role":"editor"}],"comments":[],"custom_thumbnail":false,"trashed":false,"subscribed":null,"forks":0,"collection":null,"type":"notebook","body":{"title":"Elastic-Pendulum-Solver with Runge-Kutta","cells":[{"id":0,"value":"# Elastic-Pendulum-Solver with Runge-Kutta\n\n\n---\n\nThe input parameters are:\n\n- experimentSetup(Const):\n - m: 1.0, // 质量 1kg\n - g: 9.8, // 重力加速度\n \n- experimentSetup(Variable):\n - k: 25.0, // 弹簧刚度\n - r0: 1.0, // 弹簧原长\n - dt: 0.005, // 时间步长\n - maxTime: 15, // 模拟15秒\n\n- initialConditions:\n - r: 1.2, // 初始拉伸长度\n - theta: 0.8, // 初始角度\n - v: 0.0, // 初始径向速度\n - omega: 0.0 // 初始角速度\n","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":173,"value":"

Experiment setup

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":131,"value":"viewof setup = Inputs.form({\n k: Inputs.range([20, 50], { label: \"k\", step: 0.1, value: 31.2 }),\n r0: Inputs.range([0.5, 2.5], { label: \"r0\", step: 0.01, value: 1.5 }),\n dt: Inputs.range([0.001, 0.01], { label: \"dt\", step: 0.001, value: 0.005 }),\n maxTime: Inputs.range([10, 20], { label: \"maxTime\", step: 1, value: 15 })\n})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":175,"value":"

Initial conditions

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":122,"value":"viewof initialConditions = Inputs.form({\n r: Inputs.range([1, 3], { label: \"r\", step: 0.01, value: 1.54 }),\n theta: Inputs.range([0.5, 1.2], { label: \"theta\", step: 0.01, value: 0.8 }),\n v: Inputs.range([-2, 2], { label: \"v\", step: 0.1, value: -1.1 }),\n omega: Inputs.range([-1.5, 1.5], { label: \"omega\", step: 0.1, value: -0.3 })\n})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":192,"value":"

Plot for how long?

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":44,"value":"viewof t = Inputs.range([0, d3.max(results, (d) => d.time)], {\n label: \"T(seconds)\",\n step: 0.1,\n value: d3.max(results, (d) => d.time)\n})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":104,"value":"

Trace of x, y

\n ${plt1}\n

Trace of v & omega

\n ${plt2}\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":196,"value":"---","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":14,"value":"plt1 = {\n let xDomain = d3.extent(results, (d) => d.x),\n yDomain = d3.extent(results, (d) => d.y);\n\n return Plot.plot({\n x: { nice: true, domain: xDomain },\n y: { nice: true, domain: yDomain },\n r: { range: [1, 5] },\n grid: true,\n color: {\n legend: true,\n scheme: \"Inferno\",\n reverse: true,\n range: [0.1, 0.9]\n },\n marks: [\n Plot.dot(\n results.filter((d) => d.time < t),\n {\n x: \"x\",\n y: \"y\",\n fill: \"time\",\n r: (d) => Math.abs(d.v * d.omega),\n tip: true\n }\n )\n ]\n });\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":93,"value":"plt2 = {\n let xDomain = d3.extent(results, (d) => d.v),\n yDomain = d3.extent(results, (d) => d.omega);\n\n return Plot.plot({\n x: { nice: true, domain: xDomain },\n y: { nice: true, domain: yDomain },\n grid: true,\n color: {\n legend: true,\n scheme: \"Inferno\",\n reverse: true,\n range: [0.1, 0.9]\n },\n marks: [\n Plot.dot(\n results.filter((d) => d.time < t),\n { x: \"v\", y: \"omega\", fill: \"time\", tip: true }\n )\n ]\n });\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":16,"value":"results = {\n const data = [],\n { time, r, theta, v, omega, x, y } = simulation.results;\n\n time.map((_, i) => {\n data.push({\n time: time[i],\n r: r[i],\n theta: theta[i],\n v: v[i],\n omega: omega[i],\n x: x[i],\n y: y[i]\n });\n });\n\n return data;\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":36,"value":"results","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":9,"value":"simulation = simulateElasticPendulum()","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":138,"value":"defaultParams = {\n return {\n m: 1.0, // 质量 1kg\n g: 9.8, // 重力加速度\n k: 25.0, // 弹簧刚度\n r0: 1.0, // 弹簧原长\n dt: 0.005, // 时间步长\n maxTime: 15, // 模拟15秒\n initialConditions: {\n r: 1.2, // 初始拉伸长度\n theta: 0.8, // 初始角度\n v: 0.0, // 初始径向速度\n omega: 0.0 // 初始角速度\n }\n };\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":83,"value":"function simulateElasticPendulum() {\n // 创建求解器实例\n const params = Object.assign({}, setup);\n\n Object.assign(params, { initialConditions });\n\n const solver = new ElasticPendulumSolver(params);\n\n // 求解微分方程\n const results = solver.solve();\n\n // 计算能量（用于验证）\n const energy = solver.calculateEnergy(results);\n\n return { results, energy, solver };\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":4,"value":"$31","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":53,"value":"d3 = require(\"d3\")","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"}],"files":"$W32","stdlib":"1","titling":"auto"}} 30:"$2f:body"