Sunday Engineer: NFCアンテナ(6)というかOpenEMS

NFCアンテナのシミュレーションをしてみたくなったのでOpenEMSでやってみた。

NFCで主流のISO14443Aってproximity card(近接カード)とかcontactless card(非接触カード)であって、タッチすることが前提の規格だと思うんだけど、なぜか通信距離を求められることが多い。だったらvicinity card(近傍カード=ISO15693)にすればいいのにってのは許されない。いろんなしがらみがあるんだと思う。ほんと大人の世界はいやらしい。は、さておき、こういうのなかなか実験できない場合にシミュレーションで何とかしたくなることもある(実験の方が早いことも多々あるので決断は慎重に)。
では、OpenEMSでシミュレーションしてみよう。OpenEMSの環境構築はまぁあっちこっち検索すればわかるっしょ。Octaveをインストールして、、、云々。
例によって唐突にソースコードはこちら。

gen_small_loop_coil.m



function [CSX port]=gen_small_loop_coil(
  CSX,
  loop_name,
  coil_width,coil_height,coil_nturn,
  metal_width,metal_pitch,metal_thickness,
  sub_thickness,
  series_cap_value,parallel_cap_value,port_n,port_excite,
  x_origin,y_origin,z_origin
)
 
metal_name=[loop_name '_metal'];
CSX=AddMetal(CSX,metal_name);
parallel_cap_name=[loop_name '_pcap_z'];
CSX=AddLumpedElement(CSX,parallel_cap_name,'z','C',parallel_cap_value);
series_cap_name=[loop_name '_scap_z'];
CSX=AddLumpedElement(CSX,series_cap_name,'z','C',series_cap_value);

x_start=coil_width/2-metal_width/2;
x_end=x_start-(coil_nturn)*metal_pitch;
x=x_start:-metal_pitch:x_end;

y_start=coil_height/2-metal_width/2;
y_end=y_start-(coil_nturn)*metal_pitch;
y=y_start:-metal_pitch:y_end;
 
z0=-sub_thickness/2;
z1=sub_thickness/2;
 
start_stop_vector=[metal_width/2 metal_width/2 metal_thickness/2];
origins=[x_origin y_origin z_origin];
 
for i=1:coil_nturn
  if i>1
    point1=[-x(i-1) -y(i) z1]+origins-start_stop_vector;
    point2=[x(i) -y(i) z1]+origins+start_stop_vector;
  else
    point1=[-x(coil_nturn) -y(i) z1]+origins-start_stop_vector;
    point2=[x(i) -y(i) z1]+origins+start_stop_vector;
  endif
  CSX=AddBox(CSX,metal_name,1,point1,point2);
 
  point1=[x(i) -y(i) z1]+origins-start_stop_vector;
  point2=[x(i) y(i) z1]+origins+start_stop_vector;
  CSX=AddBox(CSX,metal_name,1,point1,point2);
 
  point1=[-x(i) y(i) z1]+origins-start_stop_vector;
  point2=[x(i) y(i) z1]+origins+start_stop_vector;
  CSX=AddBox(CSX,metal_name,1,point1,point2);
 
  point1=[-x(i) -y(i+1) z1]+origins-start_stop_vector;
  point2=[-x(i) y(i) z1]+origins+start_stop_vector;
  CSX=AddBox(CSX,metal_name,1,point1,point2);
endfor
 
point1=[-x(coil_nturn) -y(1) z0]+origins-start_stop_vector;
point2=[-x(coil_nturn) -y(coil_nturn+1) z0]+origins+start_stop_vector;
CSX=AddBox(CSX,metal_name,1,point1,point2);
 
point1=[-x(coil_nturn) -y(coil_nturn+1) z0]+origins-start_stop_vector;
point2=[-x(coil_nturn) -y(coil_nturn+1) z1]+origins+start_stop_vector;
CSX=AddBox(CSX,series_cap_name,10,point1,point2);
 
point1=[-x(coil_nturn) -y(1) z0]+origins+start_stop_vector.*[-1 1 1];
point2=[-x(coil_nturn) -y(1) z1]+origins+start_stop_vector.*[1 1 -1];
CSX=AddBox(CSX,parallel_cap_name,10,point1,point2);
 
point1=[-x(coil_nturn) -y(1) z0]+origins+start_stop_vector.*[-1 -1 1];
point2=[-x(coil_nturn) -y(1) z1]+origins+start_stop_vector.*[1 -1 -1];
[CSX port]=AddLumpedPort(CSX,100,port_n,50,point1,point2,[0 0 1],port_excite);
 
 
endfunction

nfc_trx_analysis_main.m



function [freq s11 s21 s12 s22]=nfc_trx_analysis_main(
  port_direction,
  tx_coil_width,tx_coil_height,tx_coil_nturn,
  rx_coil_width,rx_coil_height,rx_coil_nturn,
  trx_distance,
  disgeomplot
)
 
physical_constants;
unit=1e-3;
 
FDTD=InitFDTD('NrTS',4e6,'EndCriteria',1e-3);
f0=14e6;
fc=12e6;
FDTD=SetGaussExcite(FDTD,f0,fc);
BC={'MUR','MUR','MUR','MUR','MUR','MUR'};
FDTD=SetBoundaryCond(FDTD,BC);
 
CSX=InitCSX();
if port_direction>0
  txport_excite=false;
  rxport_excite=true;
else
  txport_excite=true;
  rxport_excite=false;
endif
[CSX port1]=gen_small_loop_coil(
  CSX,
  'tx',
  tx_coil_width,tx_coil_height,tx_coil_nturn,
  0.2,0.4,0,
  0.2,
  1,1e-18,1,txport_excite,
  0,0,trx_distance/2
);
[CSX port2]=gen_small_loop_coil(
  CSX,
  'rx',
  rx_coil_width,rx_coil_height,rx_coil_nturn,
  0.2,0.4,0,
  0.2,
  1,1e-18,2,rxport_excite,
  0,0,-trx_distance/2
);
 
## magnetic material
##CSX = AddMaterial(CSX, 'mag_material');
##CSX = SetMaterialProperty(CSX, 'mag_material', 'Mue',45);
##x_mag=rx_coil_width/2-(rx_coil_nturn-1)*0.4-0.2-1;
##y_mag=rx_coil_height/2-(rx_coil_nturn-1)*0.4-0.2-1;
##z_mag=0.2/2;
##mag_origins=[0 0 -trx_distance/2];
##CSX=AddBox(CSX,'mag_material',1,[-x_mag -y_mag z_mag]+mag_origins,[x_mag y_mag z_mag]+mag_origins);
##
 
air_box=16;
mesh_box_size=[20 20 16];
mesh_box_resolution=16;
mesh=DetectEdges(CSX);
mesh.x=[-mesh_box_size(1) mesh.x +mesh_box_size(1)];
mesh.y=[-mesh_box_size(2) mesh.y +mesh_box_size(2)];
mesh.z=[-mesh_box_size(3) mesh.z +mesh_box_size(3)];
 
mesh=SmoothMesh(mesh,mesh_box_resolution);
 
mesh.x=[-air_box+mesh.x(1) mesh.x mesh.x(end)+air_box];
mesh.y=[-air_box+mesh.y(1) mesh.y mesh.y(end)+air_box];
mesh.z=[-air_box+mesh.z(1) mesh.z mesh.z(end)+air_box];
 
mesh=SmoothMesh(mesh,c0/(f0+fc)/unit/10,1.5,'algorithm',1);
 
mesh=AddPML(mesh,4);
 
disp(['number of cells: ' num2str(1e-6*numel(mesh.x)*numel(mesh.y)*numel(mesh.z)) 'Mcells'])
 
CSX=DefineRectGrid(CSX,unit,mesh);
 
Sim_Path=['tmp_' mfilename];
Sim_CSX=[mfilename '.xml'];
confirm_recursive_rmdir(0);
[status,message,messageid]=rmdir(Sim_Path,'s');
[status,message,messageid]=mkdir(Sim_Path);
WriteOpenEMS([Sim_Path '/' Sim_CSX],FDTD,CSX);
if disgeomplot==0
  CSXGeomPlot( [Sim_Path '/' Sim_CSX],['--export-polydata-vtk=' Sim_Path ' --RenderDiscMaterial -v']);
endif
RunOpenEMS(Sim_Path,Sim_CSX);
 
port={port1 port2};
 
freq=linspace(f0-fc,f0+fc,10001);
port_calc=calcPort(port,Sim_Path,freq,'RefImpedance',50);
 
s11=port_calc{1}.uf.ref./port_calc{1}.uf.inc;
s21=port_calc{2}.uf.ref./port_calc{1}.uf.inc;
s12=port_calc{1}.uf.ref./port_calc{2}.uf.inc;
s22=port_calc{2}.uf.ref./port_calc{2}.uf.inc;
 
##[fid,msg]=fopen('out.s2p','w');
##wdata=[freq' real(s11)' imag(s11)' real(s21)' imag(s21)' real(s12)' imag(s12)' real(s22)' imag(s22)'];
##dlmwrite(fid,wdata,'delimiter','\t');
##fclose(fid);
 
endfunction

nfc_trx_analysis.m



close all
clear all
clc
 
distance=[96 20 16 12 8 4 2 1];
 
for i=1:length(distance)
  [freq_1 s11_1 s21_1 s12_1 s22_1]=nfc_trx_analysis_main(0,36,32,4,26,12,4,distance(i),1);
  [freq_2 s11_2 s21_2 s12_2 s22_2]=nfc_trx_analysis_main(1,36,32,4,26,12,4,distance(i),1);
 
  fn=['out_' num2str(distance(i)) '.s2p']
  #fn=['out_m45_' num2str(distance(i)) '.s2p']
  [fid,msg]=fopen(fn,'w');
  fprintf(fid,"# Hz S RI R 50\n");
  wdata=[freq_1' real(s11_1)' imag(s11_1)' real(s21_1)' imag(s21_1)' real(s12_2)' imag(s12_2)' real(s22_2)' imag(s22_2)'];
  dlmwrite(fid,wdata,'delimiter','\t');
  fclose(fid);
endfor

コードの構成がまぁまぁいびつなんだけど、しゃーない。OpenEMSで集中定数も含めてシミュレーションできるので、そのようなコードにしているけど、今回は無効になるような値を入れて無効化している。OpenEMSで2portのシミュレーションをするには、port1をexcite=trueにしたケースとport2をexcite=trueにしたケースの２回のシミュレーションが必要になる。んで、コイル間距離を8パタンやってしまおうという上記のコードを我が家のCorei5 Ubuntu24.02マシンで実行すると6.5時間くらいで完了。まぁ思ったより短い。てか短くするためにいろいろごまかしている。メタルはperfect electric conductor(PEC)で、厚さを0としている。誘電体は入れてないので真空にメタルのコイルが浮かんでるモデル。
で、結果を見てみる。
コイル間距離96mm(まぁ十分遠いのでほぼ結合していないと見ていい状態)のSパラファイルをQucs-sで見てみて、で、同調も取ってみる。

ってなるので、TXコイルは1548nH、RXコイルは686nHってことになる。あってるかどうか、、、例のサイトによると、TXコイルは1530nH、RXコイルは670nHってなるっぽい。まぁまぁあってるすごい。ちなみに、同調とってインダクタンスを推定するのはシミュレーションデータだからこそたぶんやっても大丈夫なやつで、実測データには寄生容量が含まれているので、このやり方ではなく、過去、日曜技術者がやったように低周波で \[ L=\frac{X}{\omega} \] をインダクタンス値とすべきである(はず)。

で、実際のものはどんな感じなのかって思う、、、完全に同じものは無理だけど、こんな感じで大体で作っちゃう。

んで、
TXアンテナは

ってことで1540nH。
RXアンテナは

ってことで690nH。
すごいね、よくあってる。むしろ工作が正確なのにびっくり。
では、２つのコイルを結合させてみる。とりあえず、20mm。

こんなかんじ。スペーサーボルトを駆使してちょうどコイル間が20mmになっている。
これを我が家のNanoVNA H4で測定するにはport1とport2を入れ替えて２回測定する必要がある。、、、やっぱまともなネットアナが欲しいなって思っちゃうけど、がまんも大事。
んで、見てみるわけだが、まずはOpenEMSのシミュレーション結果がこちら(注:Qucs-s回路図内の同調Cは無効化している)。