电磁场有限元方法EX2.2-里兹法求解泊松方程控制的边值问题
简单学习一下有限元法的基础理论，书本为电磁场有限元经典教材：
THE FINITE ELEMENT METHOD IN ELECTROMAGNETICS, JIAN-MING JIN
目录
- 电磁场有限元方法EX2.2-里兹法求解泊松方程控制的边值问题
- 1、问题定义
- 1.1、英文描述
- 1.2、中文描述
- 2、解题
- 2.1、解题思路
- 2.2、Matlab代码实现
- 2.2.1 定义试探函数和变分问题
- 2.2.2 变分F对c1、c2变量求偏导
- 2.2.3 F对c1、c2变量偏导为0，建立方程并求解
- 2.2.4 对比理论解和数值解，并绘图
- 2.3、结果对比
- 3、全部代码
- 4、变分公式的证明
1、问题定义
1.1、英文描述
1.2、中文描述
给定由方程定义的边值问题：
d 2 ϕ d x 2 + 2 ϕ = 0 ( 0 < x < 1 ) \frac{\mathrm{d}^{2}\phi}{\mathrm{d}x^{2}} + 2\phi = 0 \quad (0 < x < 1) dx2d2ϕ+2ϕ=0(0<x<1)
边界条件：
ϕ ∣ x = 0 = 0 和 ϕ ∣ x = 1 = 1 \left.\phi\right|_{x = 0} = 0 \quad \text{和} \quad \left.\phi\right|_{x = 1} = 1 ϕ∣x=0=0和ϕ∣x=1=1
证明该问题的泛函 F F F为：
F ( ϕ ~ ) = 1 2 ∫ 0 1 [ ( d ϕ ~ d x ) 2 − 2 ϕ ~ 2 ] d x F(\tilde{\phi}) = \frac{1}{2} \int_{0}^{1} \left[ \left( \frac{\mathrm{d} \tilde{\phi}}{\mathrm{d} x} \right)^{2} - 2 \tilde{\phi}^{2} \right] \mathrm{d}x F(ϕ~)=21∫01 (dxdϕ~)2−2ϕ~2 dx
应用里茨方法，使用试探函数：
ϕ ~ ( x ) = x + ∑ n = 1 N c n sin ( n π x ) ( 取 N = 2 ) \tilde{\phi}(x) = x + \sum_{n=1}^{N} c_{n} \sin(n\pi x) \quad (\text{取 } N = 2) ϕ~(x)=x+n=1∑Ncnsin(nπx)(取 N=2)
求近似解，并将结果与精确解：
ϕ = sin ( 2 x ) sin 2 \phi = \frac{\sin(\sqrt{2}\, x)}{\sin\sqrt{2}} ϕ=sin2sin(2x)
进行比较。
2、解题
2.1、解题思路
这题比较简单，使用经典里兹方法，将试探函数 ϕ ~ ( x ) \tilde{\phi}(x) ϕ~(x)带入泛函 F F F中，然后对泛函 F F F求c1和c2的偏导即可（ F F F对c1，c2的偏导为0，以此求得极小值）。变分公式的证明见最后。
下面使用matlab实现代码思路
2.2、Matlab代码实现
2.2.1 定义试探函数和变分问题
定义试探函数和变分问题，并将试探函数带入变分：
%% Solve
syms x c1 c2
% define the trail function and variational functional
phi = x+c1*sin(pi*x) +c2*sin(2*pi*x);
F=0.5*int((diff(phi, x))^2-2*(phi^2),x,0,1);
disp('F=');
pretty(simplify(F))
2.2.2 变分F对c1、c2变量求偏导
% Find the extreme value of the functional. calculate the partial derivative of c
diff_c1 = diff(F, c1); % ∂F/∂c1
diff_c2 = diff(F, c2); % ∂F/∂c2
2.2.3 F对c1、c2变量偏导为0，建立方程并求解
% Solve the system of equations for variables [c1, c2]
equations = [diff_c1 == 0, diff_c2 == 0];
solutions = solve(equations, [c1, c2]);
disp('Solution:');
% Extract all solutions for c1 and c2
c1_val = double(solutions.c1);
c2_val = double(solutions.c2);
% Display results (both symbolic and numeric forms)
disp(['c1 = ', num2str(c1_val)]);
disp(['c2 = ', num2str(c2_val)]);
2.2.4 对比理论解和数值解，并绘图
phi_curve=subs(phi, {c1, c2}, {c1_val, c2_val});
phi_curve_exact=sin(sqrt(2)*x)/sin(sqrt(2));%% Draw figure
figure(1)
fplot(phi_curve, [0, 1], ...'Color', [0, 0.4470, 0.7410], ... % MATLAB经典蓝'LineStyle', '--', ... % 虚线'LineWidth', 2.5, ...'Marker', 'o', ... % 圆形标记'MarkerSize', 7, ...'MarkerFaceColor', 'auto', ... % 自动填充'MarkerEdgeColor', 'k', ... % 黑色边框'DisplayName', sprintf('Numerical: $\\phi(x) = x + %.2f\\sin(\\pi x) + %.2f\\sin(2\\pi x)$', c1_val, c2_val));hold on
fplot(phi_curve_exact, [0, 1], ...'Color', [0.8500, 0.3250, 0.0980], ... % MATLAB经典红'LineStyle', '-', ... % 实线'LineWidth', 2.5, ...'Marker', 's', ... % 正方形标记'MarkerSize', 8, ...'MarkerFaceColor', 'auto', ... % 自动填充'DisplayName', 'Exact: $\frac{\sin(\sqrt{2} x)}{\sin(\sqrt{2})}$');xlabel('Position (x)', 'FontSize', 14);
ylabel('Function Value \phi(x)', 'FontSize', 14);
title('Comparison of Numerical and Exact Solutions', 'FontSize', 16);
xlim([0, 1]);legend('Interpreter', 'latex', 'Location', 'best', 'FontSize', 12);set(gca, 'FontSize', 12, 'GridAlpha', 0.3, 'Box', 'on');
grid minor;
2.3、结果对比
可以看到数值解和解析解是几乎一致的：
3、全部代码
% Exercise 2.2 THE FINITE ELEMENT METHOD IN ELECTROMAGNETICS JIN JIAN MINGclc
clear
%% Solve
syms x c1 c2
% define the trail function and variational functional
phi = x+c1*sin(pi*x) +c2*sin(2*pi*x);
F=0.5*int((diff(phi, x))^2-2*(phi^2),x,0,1);
disp('F=');
pretty(simplify(F))% Find the extreme value of the functional. calculate the partial derivative of c
diff_c1 = diff(F, c1); % ∂F/∂c1
diff_c2 = diff(F, c2); % ∂F/∂c2% Solve the system of equations for variables [c1, c2]
equations = [diff_c1 == 0, diff_c2 == 0];
solutions = solve(equations, [c1, c2]);
disp('Solution:');
% Extract all solutions for c1 and c2
c1_val = double(solutions.c1);
c2_val = double(solutions.c2);
% Display results (both symbolic and numeric forms)
disp(['c1 = ', num2str(c1_val)]);
disp(['c2 = ', num2str(c2_val)]);phi_curve=subs(phi, {c1, c2}, {c1_val, c2_val});
phi_curve_exact=sin(sqrt(2)*x)/sin(sqrt(2));%% Draw figure
figure(1)
fplot(phi_curve, [0, 1], ...'Color', [0, 0.4470, 0.7410], ... % MATLAB经典蓝'LineStyle', '--', ... % 虚线'LineWidth', 2.5, ...'Marker', 'o', ... % 圆形标记'MarkerSize', 7, ...'MarkerFaceColor', 'auto', ... % 自动填充'MarkerEdgeColor', 'k', ... % 黑色边框'DisplayName', sprintf('Numerical: $\\phi(x) = x + %.2f\\sin(\\pi x) + %.2f\\sin(2\\pi x)$', c1_val, c2_val));hold on
fplot(phi_curve_exact, [0, 1], ...'Color', [0.8500, 0.3250, 0.0980], ... % MATLAB经典红'LineStyle', '-', ... % 实线'LineWidth', 2.5, ...'Marker', 's', ... % 正方形标记'MarkerSize', 8, ...'MarkerFaceColor', 'auto', ... % 自动填充'DisplayName', 'Exact: $\frac{\sin(\sqrt{2} x)}{\sin(\sqrt{2})}$');xlabel('Position (x)', 'FontSize', 14);
ylabel('Function Value \phi(x)', 'FontSize', 14);
title('Comparison of Numerical and Exact Solutions', 'FontSize', 16);
xlim([0, 1]);legend('Interpreter', 'latex', 'Location', 'best', 'FontSize', 12);set(gca, 'FontSize', 12, 'GridAlpha', 0.3, 'Box', 'on');
grid minor;4、变分公式的证明
)
