python电路仿真_蔡氏电路仿真实验
蔡氏電路仿真實(shí)驗(yàn)
模擬電路運(yùn)行
采用MATLAB進(jìn)行模擬。本電路為一常微分方程的初值問題。
取定參數(shù)
matlab代碼實(shí)現(xiàn)
首次模擬的圖像對比
參數(shù)調(diào)整后,模擬的圖像對比
分析頻譜
模擬非線性電路同步
非線性電路同步瞬態(tài)過程
取定參數(shù)
模擬實(shí)驗(yàn)中的電感電容都采用測量值,為了方便,取為定值。
非線性負(fù)阻的I-V特性取成奇對稱,數(shù)值和測量值稍作對稱性修正。
代碼實(shí)現(xiàn)
這里使用Matlab中ode45,四階龍格庫塔法四階龍格庫塔法求解常微分方程初值問題。C和C++的代碼
I-V函數(shù) myfun.m
function fun=myfun(x)
if x<0
fun=-myfun(-x);
elseif x<1.664849891
fun=x*(-7.44316E-04)-5.065E-06;
elseif x<11.29059205
fun=x*(-4.003E-04)-5.778E-04;
else
fun=x*3.574E-03-4.545E-02;
end
微分方程組 myode.m
function dy=myode(t,y)
global gg;
dy = zeros(3,1);
G=1/gg;
C1=9.91E-9;
C2=98.2E-9;
L=23E-3;
dy=[(G*(y(2)-y(1))-myfun(y(1)))/C1;
(G*(y(1)-y(2))+y(3))/C2;
-y(2)/L];
end
調(diào)節(jié)1/G輸出圖像nonlinear.m
clear;
global gg;
for gg=1800:1:2200
fs=100000;
[T,Y]=ode45('myode',0:1/fs:0.2,[0;0;0]);
X=Y(10002:end,1);
plot(Y(10002:end,1),Y(10002:end,2));
grid on
axis([-13 13 -3 3]);
xlabel('U1/V');
ylabel('U2/V');
str=[num2str(gg) '.jpg'];
saveas(gcf,str);
end
首次模擬的結(jié)果
改變
的由小到大數(shù)值,作出
相圖,可以觀察到一下幾種典型圖樣
之中的幾個穩(wěn)態(tài)
圖像在屏幕上收縮成一點(diǎn)
圖形變化的規(guī)律和真實(shí)實(shí)驗(yàn)觀察到的基本一致,但是圖像對應(yīng)的1/G和形狀與真實(shí)實(shí)驗(yàn)不太一致,究其原因,是我的L,C1,C2測量不準(zhǔn)確。
參數(shù)調(diào)整后,模擬的結(jié)果
將L、C1、C2換為參考值,再次計算。1/G在1530到2040之間變化。
C1=10E-9;
C2=100E-9;
L=18E-3;
改變參數(shù)以后圖形和實(shí)驗(yàn)更相似,但是變化的規(guī)律不變
分岔圖分析
分岔圖可以直觀得看到分岔,這里用U1的最小值來畫分岔圖,代碼如下:
clear;
fs=100000;
for gg=1992:-0.1:1940;
[T,Y]=ode23(@myode,0:1/fs:0.2,[0;0;0]);
data=Y(:,1);
n=length(data);
N=n-round(n/4);
for i=N:n-2
if data(i)
plot(gg,data(i));
hold on;
else
end
end
end
beep
分岔圖:
G從1940到1990G從1960到1990(放大)
G=1991為單倍周期G=1980為2倍周期G=1971.6為4倍周期
介于matlab提供的龍格庫塔法函數(shù)精度很難控制,且計算緩慢,建議使用C語言進(jìn)行數(shù)值計算,并提高精度。
從左邊大范圍的分岔圖上已經(jīng)清晰顯示了1→2→4→8→6→5→3的分岔過程。我們將分岔圖區(qū)域放大后(如右圖),我們發(fā)現(xiàn)分岔圖具有自相似性,即分形的最大特點(diǎn)。
處理與分析
利用快速傅立葉變換分析頻譜
X=fft(X);
N=length(X);
X(1) = [];
power = abs(X(1:N/2)).^2;
[mp,index] = max(power);
nyquist = 1/2;
freq = (1:N/2)/(N/2)*nyquist*fs;
figure(1);
subplot(211);
plot(freq,power), grid on
title (freq(index))
xlabel('Freqency');
ylabel('Power');
period = 1./freq;
figure(1);
subplot(212);
plot(period,power), grid on
title (period(index))
ylabel('Power');
xlabel('Period');
G=1500、1827、1946、1999采樣
G相圖功率譜最大周期20倍時長時域圖說明
1500
1827功率譜有5個峰,而吸引子有5層圈。最大峰位的頻率與其它相比顯著小
1946
1971.7
1980
1999
模擬非線性電路同步
隨機(jī)產(chǎn)生兩個混沌電路的初始參數(shù),只需要改變?nèi)缦乱恍小?/p>
[T,M]=ode45(@myode,0:1/fs:0.19,[rand();rand();rand()/1000]);
未同步,兩個電路的U1未同步,兩個電路的U2
同步U1,得到兩個U2的相圖同步U1,C2變大5%,得到兩個U2的相圖
非線性電路同步的瞬態(tài)過程
在極短時間內(nèi)兩個電路達(dá)到了同步
同步過程的相圖同步過程的時域圖
2016 nonlinear circuit with python implementation
Claim: This is a brief documentation of using python to understand the simulation of the Nonlinear circuit problem.
Background: Compared with Matlab , python is weigh more powerful in handling object-oriented simulation despite of simulink.
Also,python is really a hot languange which will be beneficial for our scientific computing.
這里面要求環(huán)境裝上scipy和numpy,一般windows用pip,而OS用戶就可以用brew,linux的話大家都懂
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
def Nlfunction(V_1):
Ga=-0.00076
Gb=-0.00049
E=15.0
retg=Gb*V_1+(Gb-Ga)/2*(abs(V_1-E)-abs(V_1+E))
return retg
def delta(y,t):
global R
C1=9.91e-9
C2=98.2e-9
l=23e-3
G=1/R
ret_delta =np.array([(G*(y[1]-y[0])-Nlfunction(y[0]))/C1,
(G*(y[0]-y[1])+y[2])/C2,
-y[1]/l])
return ret_delta
time=np.linspace(0,0.2,10000)
r=range(1930,1960,1)
yinit=np.array([0.0,0.0,0.1])
global R
for R in r:
y = odeint(delta,yinit,time)
plt.figure()
plt.plot(y[1000:,0],y[1000:,1])
plt.xlabel('U_1')
plt.ylabel('U_2')
plt.legend()
plt.savefig('%domega'%R)
plt.clf()
張貼一些相圖
R=1800
R=1841
R=1859
to be continued —Yufeng Ma
討論區(qū)
羅頁,向你致敬!
接下來我苦了,要學(xué)習(xí)的內(nèi)容太多!
— 樂永康 2008/12/15 23:14就感覺東西好多好雜啊!
— 羅頁 2008/12/16 00:05
很贊馬同學(xué)的嘗試。 — 樂永康 2016/05/17 22:30之前用了ubuntu所以沒有用中文寫,終于可以切換到windows了,說正題,老師我想問一下這個非線性元件里的E是怎么計算的,我感覺跟學(xué)長模擬出來的區(qū)間不一樣
我懷疑是這個E的值不太對,另外我查文獻(xiàn)貌似15V輸入到里面會有飽和電壓,不知道是什么問題。
樂老師,我不知道怎么把圖貼上來啊,為什么點(diǎn)插入圖片毫無用處? — 馬雨楓 2016/5/18 21:13編輯欄的倒數(shù)第五個按鈕是上傳附件。圖片上傳到服務(wù)器了,就可以插入到頁面了。有可能你沒有足夠的權(quán)限。
電壓15V時飽和,那是因?yàn)槟莻€負(fù)阻的電源電壓是15V。余下問題來面談吧。 — 樂永康 2016/05/19 00:16
非線性負(fù)阻的I-V曲線必須包含Gb外側(cè)的另外兩個區(qū)域,不然模擬的過程中,雙吸引子出不來。另外,odeint使用的不是RK4,而是‘Isoda’。應(yīng)該使用scipy.integrate.RK45這個類才是4階RK,上面的matlab代碼使用的ode23可能只算到2階,所以精度上不去也是正常的。 — 戴植銳 2019/04/05 20:31
總結(jié)
以上是生活随笔為你收集整理的python电路仿真_蔡氏电路仿真实验的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: eclipse设置代码格式化(详解)
- 下一篇: python中lt是什么意思_pytho