Ukázka: sériový RLC obvod

Zadání

Hodnoty součástek a zdroje napětí. Zdroj napětí dodává konstantní stejnosměrné napětí 1 V:

dos={R->6,c->(10^(-4)),L->0.1,u[t]->1};

Popis obvodu

Popis obvodu rovnicemi. Soustavy rovnic uložíme do seznamů (angl. Lists) rceSoucastek, rceNapeti, rceProudy a rcePocatecnichPodminek:

rceSoucastek={uR[t]==iR[t]*R,uL[t]==L*iL'[t],ic[t]==c*uc'[t]};
rceNapeti={0==u[t]-uR[t]-uL[t]-uc[t]};
rceProudy={iR[t]==iL[t],iL[t]==ic[t]};
rcePocatecnichPodminek={uc[0]==0,iL[0]==0};

Tvrdíme, že na počátku bylo napětí na kondenzátoru nulové a cívkou tekl nulový proud. Možná, že zdroj napětí byl až do času 0 odpojen a pak v čase 0 jsme jej k obvodu připojili...

Všechny rovnice sjednotíme do jednoho seznamu, dosadíme hodnoty součástek a zdroje napětí a s dosazenými hodnotami uložíme do seznamu vse:

vse=Union[rceSoucastek,rceProudy,rceNapeti,rcePocatecnichPodminek]/.dos;

Stanovíme interval řešení. Obvod bude kmitat, jednu periodu kmitu označme tau. Řešení chceme pro tři kmity:

Řešení

tau=2*Pi*(Sqrt[L*c])/.dos;
Print["tau= ",tau," s"];
tmax=3tau;

NDSolve vyžaduje 3 povinné parametry:
1) soustavu rovnic
2) soustavu neznámých
3) interval řešení.
Zbývá tedy sestavit seznam neznámých a obvod může být vyřešen. Řešení uložíme do proměnné res:

nezname={uR[t],iR[t],uL[t],iL[t],uc[t],ic[t]};
res=NDSolve[vse,nezname,{t,0,tmax}][[1]];

Zobrazení výsledků

Chceme-li zobrazit proud, musíme do výrazu iL[t] dosadit výsledek, který jsme obdrželi vyřešením diferenciálních rovnic. Podobně s napětím.

proud=iL[t]/.res;
napeti=uc[t]/.res;

Zobrazíme výsledky:

Plot[proud,{t,0,tmax},AxesOrigin->{0,0},PlotRange->All,PlotStyle->{Hue[0.63],Thickness[0.008]},AxesLabel->{"t [s]","i [A]"}]
Plot[napeti,{t,0,tmax},AxesOrigin->{0,0},PlotRange->All,PlotStyle->{Hue[0.63],Thickness[0.008]},AxesLabel->{"t [s]","uc [V]"}]

Totéž zobrazení, ale vychytaněji:

jak={AxesOrigin->{0,0},PlotRange->All,PlotStyle->{Hue[0.63],Thickness[0.008]}};
jakproud={AxesLabel->{"t [s]","i [A]"}};
jaknapeti={AxesLabel->{"t [s]","uc [V]"}};
Plot[proud,{t,0,tmax},Evaluate[Sequence@@Union[jak,jakproud]]]
Plot[napeti,{t,0,tmax},Evaluate[Sequence@@Union[jak,jaknapeti]]]