まずは、どのような実行プログラムになるかを想定します。
以下は、抵抗とコンデンサからなる簡単な分圧回路です。
circuit=Circuit() # 各パーツを入れる回路オブジェクトの生成
junction1=Junction(name='junc1') # ジャンクションの生成
junction2=Junction(name='junc2') # ジャンクションの生成
circuit.add_oj(junction1) # Gndとなるジャンクションの回路への導入
circuit.add_j(junction2) # ジャンクションの回路への導入
C1=Capacitor(1.0e-6,name='C1') # キャパシタの生成
circuit.add_d(C1,junction1,junction2) # キャパシタの回路への導入。接続するジャンクションを指定する。
R1=Resister(1.0e5,name='R1') # 抵抗器の生成
junction3=Junction() # ジャンクションの生成
circuit.add_j(junction3) # ジャンクションの回路への導入
circuit.add_d(R1,junction3,junction1) # 抵抗器の回路への導入。接続するジャンクションを指定する。
V_Source=VoltageSource(amp=1.0,phase=0.0,name='V_Source') # 電圧源の生成。振幅と初期位相を設定します。
circuit.add_d(V_Source,junction2,junction3) # 電圧源の回路への導入。接続するジャンクションを指定する。
circuit.set_freq(1e3) # 周波数の設定
cuicuit.calc() # 計算
まずは、回路クラスを作ります。
In [1]:
class Circuit(object):
def __init__(self):
self.junction=[]
self.o_index=None
self.device=[]
self.freq=0.0
def add_j(self,junction):
self.junction.append(junction)
junction.circuit=self
junction.index=len(self.junction)-1
def add_oj(self,junction):
self.junction.append(junction)
junction.circuit=self
junction.index=len(self.junction)-1
self.o_index=len(self.junction)-1
def add_d(self,device,jun1,jun2):
if jun1.circuit!=self or jun2.circuit!=self:
print('add the junctions before adding this divice')
return
device.p_terminal=jun1
device.n_terminal=jun2
self.device.append(device)
device.circuit=self
device.index=len(self.device)-1
def set_matrix(self):# 上から順に、各ジャンクションの電圧,各デバイスの電流を並べたベクトルXに対し、成立する連立方程式を行列で表す。
N=len(self.junction)+len(self.device)# ベクトルのサイズはN、N
self.K=np.zeros((N,N),dtype='complex')
self.V=np.zeros(N,dtype='complex')
for j,junc in enumerate(self.junction):
if j==self.o_index:# 電圧基準点は0Volt
self.K[j,j]=1
else: # 各ジャンクションに流れる電流の和はゼロ
for i, dev in enumerate(self.device):
if dev.p_terminal==junc:
self.K[j,len(self.junction)+i]=1
if dev.n_terminal==junc:
self.K[j,len(self.junction)+i]=-1
for j,dev in enumerate(self.device):
if isinstance(dev,CurrentSource): # 電流源に流れる電流を指定する。
self.V[j+len(self.junction)]=dev.get_current()
self.K[j+len(self.junction),j+len(self.junction)]=1
else: # 電流源以外のデバイスは、端子間電圧を含む関係式
_ip=(dev.p_terminal).index
_in=(dev.n_terminal).index
self.K[j+len(self.junction),_ip]=1
self.K[j+len(self.junction),_in]=-1
if isinstance(dev,VoltageSource): # 電圧源なら、端子間電圧を指定する。
self.V[j+len(self.junction)]=dev.get_voltage()
else: # 普通の受動デバイスはインピーダンスを指定する。
self.K[j+len(self.junction),j+len(self.junction)]=dev.get_impedance()
def calc(self):
import numpy.linalg
self.set_matrix()
return numpy.linalg.solve(self.K,self.V)#K.X=VをXに対して解く。
次にジャンクションのクラス
In [2]:
class Junction(object):
def __init__(self,name=''):
self.circuit=None
self.index=None
self.name=name
次は様々なデバイスの元になるDeviceクラス。
In [3]:
class Device(object):
def __init__(self,name=''):
self.circuit=None
self.index=None
self.p_terminal=None
self.n_terminal=None
self.name=name
Deviceクラスを元に、抵抗器のクラスを作ります。
In [4]:
class Resistor(Device):
def __init__(self,R=0,name=''):
Device.__init__(self,name=name)
self.R=R
def get_impedance(self):
return self.R
Deviceクラスを元に、キャパシタのクラスを作ります。
In [5]:
class Capacitor(Device):
def __init__(self,C=1e-12,name=''):
Device.__init__(self,name=name)
self.C=C
def get_impedance(self):
omega=np.pi*2.0*self.circuit.freq
return 1.0/(1.0j*omega*self.C)
Deviceクラスを元に、インダクタのクラスを作ります。
In [6]:
class Inductor(Device):
def __init__(self,L=1e-9,name=''):
Device.__init__(self,name=name)
self.L=L
def get_impedance(self):
omega=np.pi*2.0*self.circuit.freq
return 1.0j*omega*self.L
Deviceクラスを元に、電圧源のクラスを作ります。
In [7]:
class VoltageSource(Device):
def __init__(self,amp=1,phase=0,name=''):
Device.__init__(self,name=name)
self.amp=amp
self.phase=phase
def get_voltage(self):
return self.amp*np.exp(1.0j*self.phase)
Deviceクラスを元に、電流源のクラスを作ります。
In [8]:
class CurrentSource(Device):
def __init__(self,amp=1,phase=0,name=''):
Device.__init__(self,name=name)
self.amp=amp
self.phase=phase
def get_current(self):
return self.amp*np.exp(1.0j*self.phase)
いよいよ、回路を作ります。
In [9]:
circuit=Circuit()
junction1=Junction(name='junc1')
junction2=Junction(name='junc2')
circuit.add_oj(junction1)# 電圧の基準点を示す必要があります。基準となるジャンクションはこの関数で指定します。
circuit.add_j(junction2)
C1=Capacitor(1.0e-9,name='C1')
circuit.add_d(C1,junction1,junction2)
R1=Resistor(1.0e5,name='R1')
junction3=Junction(name='junc3')
circuit.add_j(junction3)
circuit.add_d(R1,junction3,junction2)
V_Souce=VoltageSource(amp=1.0,phase=0.0,name='V_Source')
circuit.add_d(V_Souce,junction1,junction3)
様々な周波数で、junction2の電圧を見てみましょう。
In [10]:
freq_array=np.logspace(1,6)
vout=np.zeros(len(freq_array),dtype='complex')
for i,freq in enumerate(freq_array):
circuit.freq=freq
x=circuit.calc()
vout[i]=x[junction2.index]
In [11]:
plt.plot(freq_array,np.abs(vout),'b-',label='Amplitude')
plt.xlabel('frequency')
plt.ylabel('Voltage of junction2')
plt.xscale('log')
plt.yscale('log')
plt.grid(which='both')
plt.legend(loc='lower left')
plt.twinx()
plt.plot(freq_array,np.angle(vout)*180.0/np.pi,'r-',label='phase')
plt.legend(loc='upper right')
plt.ylabel('Phase(Degree)')
Out[11]:
<matplotlib.text.Text at 0x474e210>
こんどは、LCRの並列回路に電流源を付けた例です。
In [12]:
circuit2=Circuit()
junction1=Junction(name='junc1')
junction2=Junction(name='junc2')
circuit2.add_oj(junction1)# 電圧の基準点を示す必要があります。基準となるジャンクションはこの関数で指定します。
circuit2.add_j(junction2)
C1=Capacitor(1.0e-6,name='C1')
circuit2.add_d(C1,junction1,junction2)
R1=Resistor(5.0e3,name='R1')
circuit2.add_d(R1,junction1,junction2)
L1=Inductor(1.0e-2,name='L1')
circuit2.add_d(L1,junction1,junction2)
I_Source=CurrentSource(amp=1.0,phase=0.0,name='I_Source')
circuit2.add_d(I_Source,junction1,junction2)
circuit2.freq=1e3
circuit2.set_matrix()
In [13]:
freq_array2=np.logspace(2,4,num=500)
Iout=np.zeros(len(freq_array2),dtype='complex')
for i,freq in enumerate(freq_array2):
circuit2.freq=freq
x2=circuit2.calc()
Iout[i]=x2[C1.index+len(circuit2.junction)]# キャパシタに流れる電流
In [14]:
plt.plot(freq_array2,np.abs(Iout),'b-',label='Amplitude')
plt.xlabel('frequency')
plt.ylabel('Current of C1')
plt.xscale('log')
plt.yscale('log')
plt.grid(which='both')
plt.legend(loc='lower left')
plt.twinx()
plt.plot(freq_array2,np.angle(Iout)*180.0/np.pi,'r-',label='phase')
plt.legend(loc='upper right')
plt.ylabel('Phase(Degree)')
Out[14]:
<matplotlib.text.Text at 0x54bf1d0>
まとめ¶
現在、グラフィカルに回路を構成できる高機能の回路シミュレータソフトを無償で利用することが可能です。 それに比べると、今回のような方法だと複雑な回路をシミュレートするのは容易ではないでしょう。
しかしながら、自分で作れるということは、好きなようにカスタマイズできるという良さがあるでしょう。 例えば、任意のインピーダンス特性を持つデバイスを作ることも簡単です。 ただし、デバイスのインピーダンスを与える際、因果律を守るようにしないと、物理的にあり得ないシミュレーションになってしまう可能性がありますので、注意しましょう。