You need to enable JavaScript to run this app.
优惠活动
大模型
产品
解决方案
定价
更多

如何修正Runge-Kutta实现，正确模拟带电粒子间的电场线？

avatar
安全阿岳
2026-6-13
问题描述

我用Runge-Kutta算法模拟两个带电粒子的电场线，在正电荷周围采样点，通过小步长近似电场方向，但遇到两个问题：

  • 电场线无法延伸至左侧负电荷的后方
  • 从负电荷向正电荷生成电场线时，线条碰撞、视觉杂乱

现有效果：现有电场线效果
期望效果：期望电场线效果

以下是我的代码：

<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="UTF-8">
</head>
<body>
<canvas id="canvas" width="700" height="500"></canvas>
<script>
const canvas = document.getElementById("canvas");
const ctx = canvas.getContext("2d");

function vec3_add(a, b) {
   return [a[0] + b[0], a[1] + b[1], a[2] + b[2]];
}

function vec3_sub(a, b) {
   return [a[0] - b[0], a[1] - b[1], a[2] - b[2]];
}

function vec3_scale(v, s) {
   return [v[0] * s, v[1] * s, v[2] * s];
}

function vec3_magnitude(v) {
   return Math.sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
}

function vec3_normalize(v) {
   const mag = vec3_magnitude(v);
   if (mag === 0) return [0, 0, 0];
   return [v[0] / mag, v[1] / mag, v[2] / mag];
}
const charges = [{
       pos: [250, 250, 0],
       id: "negative"
   },
   {
       pos: [450, 250, 0],
       id: "positive"
   }
];
const lines = 128;
const step_size = 4;
const max_steps = 1000;

function calculate_field_at_point(point, charges) {
   let field = [0, 0, 0];
   for (let charge of charges) {
       let dir = vec3_sub(charge.pos, point);
       let dist = vec3_magnitude(dir);
       if (dist < 5) continue;
       let strength = charge.id.includes("negative") ? 1 : -1;
       strength /= (dist * dist);
       field = vec3_add(field, vec3_scale(dir, strength));
   }
   return vec3_normalize(field);
}

function rk4_step(point, charges) {
   let k1 = calculate_field_at_point(point, charges);
   let temp = vec3_add(point, vec3_scale(k1, step_size * 0.5));
   let k2 = calculate_field_at_point(temp, charges);
   temp = vec3_add(point, vec3_scale(k2, step_size * 0.5));
   let k3 = calculate_field_at_point(temp, charges);
   temp = vec3_add(point, vec3_scale(k3, step_size));
   let k4 = calculate_field_at_point(temp, charges);
   return vec3_scale(
       vec3_add(
           vec3_add(
               vec3_scale(k1, 1 / 6),
               vec3_scale(k2, 1 / 3)
           ),
           vec3_add(
               vec3_scale(k3, 1 / 3),
               vec3_scale(k4, 1 / 6)
           )
       ),
       step_size
   );
}

function generate_circle_points(center, radius, num_points) {
   let points = [];
   for (let i = 0; i < num_points; i++) {
       let angle = (i / num_points) * Math.PI * 2.0;
       let x = center[0] + radius * Math.cos(angle);
       let y = center[1] + radius * Math.sin(angle);
       points.push([x, y, 0]);
   }
   return points;
}

function integrate_field_line(start_point, charges) {
   let points = [start_point];
   let current_point = [...start_point];
   for (let i = 0; i < max_steps; i++) {
       let step = rk4_step(current_point, charges);
       if (vec3_magnitude(step) < 0.1) break;
       current_point = vec3_add(current_point, step);
       points.push([...current_point]);
       let too_close = charges.some(charge => vec3_magnitude(vec3_sub(current_point, charge.pos)) < 15);
       // if (too_close) break;
       // if (current_point[0] < 0 || current_point[0] > canvas.width ||
       //     current_point[1] < 0 || current_point[1] > canvas.height) break;
   }
   return points;
}
const positive_charges = charges.filter(c => c.id.includes("positive"));
for (let charge of positive_charges) {
   let start_points = generate_circle_points(charge.pos, 35, lines);
   for (let start_point of start_points) {
       let points = integrate_field_line(start_point, charges);
       if (points.length > 5) {
           ctx.beginPath();
           ctx.moveTo(points[0][0], points[0][1]);
           for (let i = 1; i < points.length; i++) {
               ctx.lineTo(points[i][0], points[i][1]);
           }
           ctx.strokeStyle = "rgba(0, 0, 0, 0.5)";
           ctx.lineWidth = 2;
           ctx.stroke();
       }
   }
}
ctx.beginPath();
ctx.arc(charges[0].pos[0], charges[0].pos[1], 30, 0, Math.PI * 2);
ctx.fillStyle = "blue";
ctx.fill();
ctx.lineWidth = 4;
ctx.strokeStyle = "black";
ctx.stroke();
ctx.beginPath();
ctx.moveTo(charges[0].pos[0] - 15, charges[0].pos[1]);
ctx.lineTo(charges[0].pos[0] + 15, charges[0].pos[1]);
ctx.lineWidth = 5;
ctx.strokeStyle = "black";
ctx.stroke();
ctx.beginPath();
ctx.arc(charges[1].pos[0], charges[1].pos[1], 30, 0, Math.PI * 2);
ctx.fillStyle = "red";
ctx.fill();
ctx.lineWidth = 4;
ctx.strokeStyle = "black";
ctx.stroke();
ctx.beginPath();
ctx.moveTo(charges[1].pos[0] - 15, charges[1].pos[1]);
ctx.lineTo(charges[1].pos[0] + 15, charges[1].pos[1]);
ctx.moveTo(charges[1].pos[0], charges[1].pos[1] - 15);
ctx.lineTo(charges[1].pos[0], charges[1].pos[1] + 15);
ctx.lineWidth = 5;
ctx.strokeStyle = "black";
ctx.stroke();
</script>
</body>
</html>
解决方案

核心问题分析

  1. 电场线无法延伸至负电荷后方：当前仅从正电荷正向生成电场线，负电荷后方的电场线需要从负电荷反向延伸（模拟电场线终止于负电荷的路径）。
  2. 负电荷生成线条杂乱：直接正向积分会导致线条全部指向负电荷，互相重叠，需改为反向积分电场方向。

具体修改步骤

  1. 修正电场强度计算的符号逻辑，确保正/负电荷的电场方向正确。
  2. 给integrate_field_line添加反向步进参数，支持反向积分电场线。
  3. 恢复边界与近电荷终止检测，避免线条混乱。
  4. 添加负电荷的反向电场线绘制逻辑。

修正后的完整代码 

<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="UTF-8">
</head>
<body>
<canvas id="canvas" width="700" height="500"></canvas>
<script>
const canvas = document.getElementById("canvas");
const ctx = canvas.getContext("2d");

function vec3_add(a, b) {
   return [a[0] + b[0], a[1] + b[1], a[2] + b[2]];
}

function vec3_sub(a, b) {
   return [a[0] - b[0], a[1] - b[1], a[2] - b[2]];
}

function vec3_scale(v, s) {
   return [v[0] * s, v[1] * s, v[2] * s];
}

function vec3_magnitude(v) {
   return Math.sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
}

function vec3_normalize(v) {
   const mag = vec3_magnitude(v);
   if (mag === 0) return [0, 0, 0];
   return [v[0] / mag, v[1] / mag, v[2] / mag];
}
const charges = [{
       pos: [250, 250, 0],
       id: "negative"
   },
   {
       pos: [450, 250, 0],
       id: "positive"
   }
];
const lines = 128;
const step_size = 4;
const max_steps = 1000;

function calculate_field_at_point(point, charges) {
   let field = [0, 0, 0];
   for (let charge of charges) {
       let dir = vec3_sub(charge.pos, point);
       let dist = vec3_magnitude(dir);
       if (dist < 5) continue;
       // 修正电场方向：正电荷向外，负电荷向内
       let strength = charge.id.includes("negative") ? -1 : 1;
       strength /= (dist * dist);
       field = vec3_add(field, vec3_scale(vec3_normalize(dir), strength));
   }
   return vec3_normalize(field);
}

function rk4_step(point, charges) {
   let k1 = calculate_field_at_point(point, charges);
   let temp = vec3_add(point, vec3_scale(k1, step_size * 0.5));
   let k2 = calculate_field_at_point(temp, charges);
   temp = vec3_add(point, vec3_scale(k2, step_size * 0.5));
   let k3 = calculate_field_at_point(temp, charges);
   temp = vec3_add(point, vec3_scale(k3, step_size));
   let k4 = calculate_field_at_point(temp, charges);
   return vec3_scale(
       vec3_add(
           vec3_add(
               vec3_scale(k1, 1 / 6),
               vec3_scale(k2, 1 / 3)
           ),
           vec3_add(
               vec3_scale(k3, 1 / 3),
               vec3_scale(k4, 1 / 6)
           )
       ),
       step_size
   );
}

function generate_circle_points(center, radius, num_points) {
   let points = [];
   for (let i = 0; i < num_points; i++) {
       let angle = (i / num_points) * Math.PI * 2.0;
       let x = center[0] + radius * Math.cos(angle);
       let y = center[1] + radius * Math.sin(angle);
       points.push([x, y, 0]);
   }
   return points;
}

function integrate_field_line(start_point, charges, reverse = false) {
   let points = [start_point];
   let current_point = [...start_point];
   for (let i = 0; i < max_steps; i++) {
       let step = rk4_step(current_point, charges);
       if (reverse) step = vec3_scale(step, -1);
       if (vec3_magnitude(step) < 0.1) break;
       current_point = vec3_add(current_point, step);
       points.push([...current_point]);
       
       // 终止条件：接近电荷或超出画布
       let too_close = charges.some(charge => vec3_magnitude(vec3_sub(current_point, charge.pos)) < 15);
       if (too_close) break;
       if (current_point[0] < 0 || current_point[0] > canvas.width ||
           current_point[1] < 0 || current_point[1] > canvas.height) break;
   }
   return points;
}

// 绘制正电荷出发的电场线（正向积分）
const positive_charges = charges.filter(c => c.id.includes("positive"));
for (let charge of positive_charges) {
   let start_points = generate_circle_points(charge.pos, 35, lines);
   for (let start_point of start_points) {
       let points = integrate_field_line(start_point, charges);
       if (points.length > 5) {
           ctx.beginPath();
           ctx.moveTo(points[0][0], points[0][1]);
           for (let i = 1; i < points.length; i++) {
               ctx.lineTo(points[i][0], points[i][1]);
           }
           ctx.strokeStyle = "rgba(0, 0, 0, 0.5)";
           ctx.lineWidth = 2;
           ctx.stroke();
       }
   }
}

// 绘制负电荷出发的电场线（反向积分，模拟终止于负电荷的路径）
const negative_charges = charges.filter(c => c.id.includes("negative"));
for (let charge of negative_charges) {
   let start_points = generate_circle_points(charge.pos, 35, lines);
   for (let start_point of start_points) {
       let points = integrate_field_line(start_point, charges, true);
       if (points.length > 5) {
           ctx.beginPath();
           ctx.moveTo(points[0][0], points[0][1]);
           for (let i = 1; i < points.length; i++) {
               ctx.lineTo(points[i][0], points[i][1]);
           }
           ctx.strokeStyle = "rgba(0, 0, 0, 0.5)";
           ctx.lineWidth = 2;
           ctx.stroke();
       }
   }
}

// 绘制电荷图形
ctx.beginPath();
ctx.arc(charges[0].pos[0], charges[0].pos[1], 30, 0, Math.PI * 2);
ctx.fillStyle = "blue";
ctx.fill();
ctx.lineWidth = 4;
ctx.strokeStyle = "black";
ctx.stroke();
ctx.beginPath();
ctx.moveTo(charges[0].pos[0] - 15, charges[0].pos[1]);
ctx.lineTo(charges[0].pos[0] + 15, charges[0].pos[1]);
ctx.lineWidth = 5;
ctx.strokeStyle = "black";
ctx.stroke();

ctx.beginPath();
ctx.arc(charges[1].pos[0], charges[1].pos[1], 30, 0, Math.PI * 2);
ctx.fillStyle = "red";
ctx.fill();
ctx.lineWidth = 4;
ctx.strokeStyle = "black";
ctx.stroke();
ctx.beginPath();
ctx.moveTo(charges[1].pos[0] - 15, charges[1].pos[1]);
ctx.lineTo(charges[1].pos[0] + 15, charges[1].pos[1]);
ctx.moveTo(charges[1].pos[0], charges[1].pos[1] - 15);
ctx.lineTo(charges[1].pos[0], charges[1].pos[1] + 15);
ctx.lineWidth = 5;
ctx.strokeStyle = "black";
ctx.stroke();
</script>
</body>
</html>

关键修改说明

  1. 电场方向修正：在calculate_field_at_point中，正电荷的电场强度设为1（向外），负电荷设为-1（向内），先归一化方向向量再缩放，确保计算准确。
  2. 反向积分支持：integrate_field_line新增reverse参数，反转步进方向，实现从负电荷向外绘制电场线（模拟电场线终止于负电荷的效果）。
  3. 终止条件恢复：重新启用接近电荷和画布边界的终止逻辑，避免线条重叠或超出画布。
  4. 负电荷电场线绘制：单独处理负电荷，调用反向积分生成电场线，覆盖负电荷后方区域，匹配期望效果。

内容的提问来源于stack exchange，提问作者imad

相关产品推荐
方舟 Agent Plan

超全模态模型 × Harness 升级，最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列、Doubao-Seedream-5.0-pro、Kimi-K3 (部分)， 限时 9.9 元起

产品介绍相关文档
相关产品推荐
方舟 Agent Plan

超全模态模型 × Harness 升级，最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列、Doubao-Seedream-5.0-pro、Kimi-K3 (部分)， 限时 9.9 元起

产品介绍相关文档
火山引擎 最新活动
全面开放
Doubao-Seedance-2.5全新模型升级

模型再升级,30秒超长叙事, 模态参考扩容

HOT
方舟Codingplan

模型自由，工具不限，最新支持 Deepseek-V4 系列、GLM-5.3 系列

HOT
方舟Agent plan

超全模态模型 × Harness 升级，最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列

HOT
新用户特惠专场

大模型19元起，Al应用9.9元畅享，新人首购爆款尽享优惠

最近更新时间：2026.06.13 23:27:29