Helmholtz线圈磁场的仿真研究.docx
- 文档编号:7513422
- 上传时间:2023-01-24
- 格式:DOCX
- 页数:13
- 大小:217.58KB
Helmholtz线圈磁场的仿真研究.docx
《Helmholtz线圈磁场的仿真研究.docx》由会员分享,可在线阅读,更多相关《Helmholtz线圈磁场的仿真研究.docx(13页珍藏版)》请在冰豆网上搜索。
Helmholtz线圈磁场的仿真研究
一、仿真题目要求
为了获得一定区域上的匀强磁场,可采用多组Helmholtz线圈结构。
一种两对线圈的结构如下图所示。
线圈半径a1,a2,线圈间距h1,h2,以及线圈中通过电流i1,i2可变化量,如图1(a)所示。
为了定量衡量关注区域的磁场均压程度,过轴线做截面ABo1o2,取CD=0.8×AB和Eo3=0.8×o1o2,在CD和Eo3线段上每边均匀取20采样点,从而形成如图1(b)所示的采样节点,定义z方向B的不均压系数为:
其中,
为所有采样点的z方向磁感应强度平均值;
为第n个采样点的z方向磁感应强度值,N为采样点总数。
如给定a1=1m,i1=1A,且h1>0.2m,
情况下,如果要使得δ达到最小,线圈半径a2,线圈间距h1,h2,以及线圈中通过电流i2应如何取值?
(a)线圈结构示意图
(b)磁场采样节点示意图
二、公式推导
1.对于单线圈情况推导磁感应强度
由于仿真要求在xoz平面上取点,对磁感应强度z轴分量
进行不均压度分析,则只需关注在xoz平面上点的
,推导过程如下:
如上图所示,位于
处的电流元
,对场点P(r,z)处产生的磁感应强度z轴分量
则整个线圈电流元
中
从0到
,进行积分得到整个线圈对场点P产生的磁感应强度z轴分量
2.不均匀度公式的简化
题目中给定的不均匀度为所有采样点的磁感应强度标准差除以平均值的绝对值,而由数学推导可得如下化简形式,即采样点磁感应强度平方的均值除以均值的平方再减一,即得不均匀度的平方。
此处使用不均匀度的平方是因为可以消去绝对值运算和开方运算,减少运算次数、同时也使目标函数变得可导。
三、线圈分布规律猜测
根据之前的推导易于得出,对一个中心位于柱坐标系中z=h处、半径为R、通电流为I的线圈,对于轴线上、高度为z处的场点,磁感应强度z向分量为
由于该函数是关于z的偶函数,奇数阶导数一定都恒等于0,对其求偶数阶导数如下。
而对于4个线圈,分别按照题目中定义的参数如下表:
线圈
1
2
3
4
中心位置
-h2
-h1
h1
h2
半径
R2
R1
R1
R2
电流
I2
I1
I1
I2
可以求得4线圈同时作用时,原点处磁感应强度的z向分量各阶导数为
而对同一个上述线圈,仍根据之前得出的式子,可以推导出对于与轴线相距x、高度为0处的场点,磁感应强度z向分量为
由于问题本身的物理背景,仍然可以看出此函数是关于x的偶函数,奇数阶导数一定都恒等于0。
并且根据高等数学的知识,可以证明此函数对于x的各阶导数相当于积分内的函数先对x求偏导,再做对θ的积分。
求导工作量很大,借助MATLAB的符号工具箱对其求原点(0,0)处偶数阶导数的值,如下所示。
仍与刚才一样代入4个线圈的参数,得到4个线圈共同作用时原点处的各阶导数如下
至此,得到了
(1)~(4)式。
令它们均等于0,可以联立求解出4个未知变量R2、h1、h2、I2。
联立的方程组如下。
得到一组解:
R2=2.153********0254
h1=0.293332954474377
h2=1.967170747875403
I2=3.000000000000000
可以认为按照此规律分布的磁场比较均匀,故可以用此作为优化初值。
四、仿真程序编写
1.圆周上小段直电流产生的磁场函数dBz(θ,R,r)
在柱坐标系中,按照上述推导,编写了在(R,θ)处长为Rdθ的一段电流,方向为
,对于位于(r,z)处的场点,产生的磁场大小在z轴上分量的计算函数如下(没有乘系数μ/4π)。
dBz=@(theta)(R-r.*cos(theta))./sqrt((r.^2-2.*r.*R.*cos(theta)+R.^2+z.^2).^3)
2.一个圆产生的磁感应强度函数Bz_single(r,z,R,I)
在柱坐标系中对位于(r,z)处的场点,在原点处的半径为R的线圈产生的磁场z轴分量计算函数如下。
由于每次计算磁场都会乘系数μ/4π,这部分可以留到函数外面计算,故函数的返回值实际上是没有乘系数μ/4π的磁场强度。
函数采用自己编写的梯形公式进行积分,其中step为积分步长,last为上一积分步的被积函数值,now为当前积分步的被积函数值。
functionBz=Bz_single(r,z,R,I)
%tocalculateBz*4*pi/muofpoint(x,z)inmagneticfield,
%withintegralstepsvectortheta
dBz=@(theta)(R-r.*cos(theta))./sqrt((r.^2-2.*r.*R.*cos(theta)+R.^2+z.^2).^3);
step=pi/20;
Bz=0;
last=dBz(0);%tostoredBzoflaststep
fortheta=step:
step:
2*pi;
now=dBz(theta);%tostoredBzofthisstep
Bz=Bz+(now+last);
last=now;
Bz=I*R*0.5*step*Bz;
end
对此函数进行正确性测试。
在《大学物理学(A版)·电磁学》中有一道习题17.6,求一个半圆形电流I在半圆轴线上距离圆心z处的磁场强度。
习题答案中给出的磁场强度轴线方向分量为
,为了适应此函数,也即
。
由于题目中是半圆,相当于积分只从0积分到π,在不同的step下得到的积分结果如下表
step
π/1000
π/100
π/10
准确值
积分结果
z=0
3.141592653589793
3.141592653589794
3.141592653589793
3.141592653589793
z=1
1.110720734539600
1.110720734539591
1.110720734539592
1.110720734539592
z=2
0.280992589241630
0.280992589241630
0.280992589241629
0.280992589241629
z=3
0.099345882657959
0.099345882657961
0.099345882657961
0.099345882657961
可以看出首先函数编写的是准确的,但在这种情况下,积分步长对函数值变化的影响并不是很大,甚至出现了积分步长越粗大,结果反而更精确的反常现象。
这有可能是因为π除以一个大数之后变得不精确的原因。
为了找到一个经济合适的积分步长,之后再取一个更有一般性的例子,如线圈半径R=1,电流I=2,场点距轴线距离r=0.3,场点距圆心高度z=0.5,这个例子中的各个参数和最终要解决的问题是一个数量级的。
测试不同步长情况下的原函数(从0积分到2π)计算一个线圈磁感应强度的积分结果如下表
step
π/10000
π/1000
π/100
π/20
π/10
π/5
积分结果
8.957661590821102
8.957661590821102
8.957661590821099
8.957661590821099
8.957661590831195
8.957675492861105
误差
<1e-15
<1e-15
1e-14
1e-14
1e-11
1e-4
此处对θ的步长选取相当于是将圆等效成了多边形,可以看出π/20、π/10的步长已经很能满足一般的计算精度要求了,也很显然可以从直观上感觉,40边形和20边形对于圆的近似程度是很高的。
因此最后的程序中采用步长π/20进行计算。
再对函数进行效率优化。
由于此积分函数在最终的程序中会被反复调用,此函数耗时一定要短。
最初想用MATLAB自带的高斯-拉布拉托积分quad()或quadl(),但由于其中涉及高阶差分和自适应步长等过程,精度虽高,却不如简单的梯形积分快。
且由之前过程可看出,只要步长合适,梯形积分也可以满足精度要求。
此函数中在自己编写梯形积分时,每个θ对应的被积函数仅计算一遍,相同的乘法部分均提出到for循环外,没有任何冗余运算,同时也在保证效率前提下将占用变量(即系统内存)尽可能的进行了缩减。
3.两对圆电流总磁感应强度函数Bz_total(R2,h1,h2,I2,r,z)
其中R2为原图中外侧线圈半径,h1、h2分别为内、外侧线圈到原点的距离,I2为外侧线圈电流,以和I1同向为正,r、z分别为采样点的半径坐标向量、纵坐标向量。
本函数4次调用单个圆电流产生磁场的函数,在柱坐标系中将4个线圈的总效应相加,调用一次即可计算所有采样点的z向磁感应强度,结果返回给矩阵变量Bz。
编写的程序如下。
functionBz=Bz_total(R2,h1,h2,I2,r,z)
R1=1;I1=1;
%ifuseparallelcomputing
[X,Z]=meshgrid(r,z);
Z1=Z+h2;%coil1
Z2=Z+h1;%coil2
Z3=Z-h1;%coil3
Z4=Z-h2;%coil4
Bz1=Bz_single(X,Z1,R2,I2);
Bz2=Bz_single(X,Z2,R1,I1);
Bz3=Bz_single(X,Z3,R1,I1);
Bz4=Bz_single(X,Z4,R2,I2);
Bz=Bz1+Bz2+Bz3+Bz4;
end
值得一提的是,本函数没有采用大家普遍采用的for循环方式来计算400个采样点,而是采用了meshgrid函数生成400个点的坐标矩阵网格,然后对坐标矩阵进行操作,并行的同时计算400个采样点。
在不利用其他并行工具的前提下,MATLAB中for循环是一种低效率的串行方式,而通过上网查找资料发现较高版本的MATLAB中矩阵三角函数运算、矩阵开方运算、矩阵点乘除乘方运算都是自动将运算量分配到多个计算单元上并行执行的,能够显著提速。
由于我们用到的运算也就是以上几种的组合,因此本函数中先生成了关于r、z坐标向量的两个二维网格矩阵,再直接对这两个网格矩阵进行运算,以达到并行计算的目的。
经过实测,其他条件相同的情况下,采用for循环方式做一次优化至3分钟后仍未出结果(其他组同学的程序平均时间一般也在10分钟左右);采用坐标网格矩阵的并行方式做一次优化只需要2秒多。
4.不均匀度函数differ(x)
由于题中所给采样区的磁场分布一定是对称的,故此处不均匀度函数仅选取了题目中一半仿真区域11×21个点,半径方向20等分取21点,纵轴方向10等分取11点。
函数返回值d为所有采样点不均匀度的平方(这对于找其最小值并无影响),输入变量x为一个4维列向量,分别是4个待定参数R2、h1、h2、I2。
按照之前推导的简化计算量的公式,编写程序如下。
functiond=differ(x)
%x=[R2;h1;h2;I2]
%d为不均匀系数delta的平方
R2=x
(1);h1=x
(2);h2=x(3);I2=x(4);
%每边采样点数设置
r=0:
0.04:
0.8;
z=0:
0.08*h1:
0.8*h1;%z向仅对上半部分进行采样
Bz_half=Bz_total(R2,h1,h2,I2,r,z);
Bz=[Bz_half;Bz_half(2:
end,:
)];
d=mean2(Bz.^2)/mean2(Bz)^2-1;
end
因为此函数中没有绝对值,之前的其他函数部分也没有超越函数,所以至此,待优化的目标函数是可导的,有利于后续优化算法的选择。
5.优化主函数
选用优化工具箱中的函数fmincon进行优化,找到目标函数的最小值。
其中约束条件的设置在矩阵A1中和向量v1、v2中,限定了h2>h1、h2>0.2、-3 向量v1、v2的其他部分为人为设定的另外约束,为了方便找到最优解。 Clearall;clc; tic; mu=4*pi*1e-7; R1=1;I1=1; %variablex=[R2;h1;h2;I2]; x0=[1;0.5;1;3]; A1=[0,1,-1,0]; b1=0; v1=[0.1;0.2;0.3;-3]; v2=[4;0.7;4;3]; opt=optimset(‘largescale’,’off’,’MaxIter’,200,’MaxFun’,1000); [x,fv,ef,out,grad,hess]=fmincon(@differ,x0,A1,b1,[],[],v1,v2,[],opt) R2=x (1),h1=x (2),h2=x(3),I2=x(4) delta=sqrt(fv) toc; 五、数值试验 1.按之前导数为零的推导选取优化初值 首先选取了之前按照原点处各阶导数为0推导出的一组参数作为初值,变量范围设在其附近,得到了一组优化结果,画出其磁感应强度分布如下图(a)。 但随后发现当初值取为其他值时会出现不均匀度更小的解,如下图(b)。 这说明之前的推导是不正确的,它仅仅保证了原点附近磁场足够均匀,却不能保证整个采样区域的不均匀度最小。 (a)按照导数为零得到的结果 (b)按照另一组初值得到的结果 2.按数值试验逐渐缩小范围的方法选取初值 我们希望能通过数值实验先计算几组参数值下的不均匀系数平方,再进一步缩小范围寻找其最小值,确定所优化目标的大致位置。 由于要优化四个参数,若参数值一起改变,不太方便观察,可先固定两个参数不变,改变其他两个,再反过来改变另两个即可。 根据之前几次优化发现,一般最后结果都会有h1=0.2,I2=3,而另外两个量会有一些变化,因此先固定h1=0.2,I2=3。 由于只是要得到大致初值,故不必非常精确。 (1)首先固定h1=0.2,I2=3,改变h2,R2: 当h2由0.2变化到2,取步长为0.1,R2由0.2变化到2,取步长为0.1;得到h2=1,R2=1时differ(min)= ;进一步缩小范围,当h2由0.8变化到1.1,取步长为0.01,R2由0.9变化到1.1,取步长为0.01;得到h2=0.97,R2=1.04时differ(min)= 。 至此,h2和R2的位置已经大致确定。 (2)将h2,R2固定为第一步求得的值,即h2=0.97,R2=1.04,改变h1和I2的值: 当h1由0变化到4,取步长为0.01,I2由1变化到3,取步长为0.1;得到h1=0.19,I2=3.0时differ(min)= 。 最终选取的初值为h1=0.20,h2=0.97,R2=1.04,I2=3.0。 各个量的范围为h1=[0.20,0.50],h2=[0.80,1.20],R2=[0.80,1.20],I2=[2.50,3.0]。 六、结论 按照之前得出的初值和范围进行优化,最终得到的最优解如下: R2=1.0027,h1=0.2000,h2=0.9412,I2=2.9961,不均匀度delta=0.008539。 绘制出的图像如下: 整个优化程序耗时共: Elapsedtimeis2.150968seconds;迭代次数5次。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- Helmholtz 线圈 磁场 仿真 研究
![提示](https://static.bdocx.com/images/bang_tan.gif)