Bistable
Z Skrypty dla studentów Ekonofizyki UPGOW
(Różnice między wersjami)
m (→Układ) |
(→Układ) |
||
| Linia 1: | Linia 1: | ||
==Oscylator bistabliny== | ==Oscylator bistabliny== | ||
| + | |||
| + | ===TEST=== | ||
| + | |||
| + | :{| class="toccolours collapsible collapsed" width="60%" style="text-align:left" | ||
| + | !Derivation of equilibrium equations | ||
| + | |- | ||
| + | |Consider a continuum body (see Figure 4) occupying a volume <math>V\,\!</math>, having a surface area <math>S\,\!</math>, with defined traction or surface forces <math>T_i^{(n)}\,\!</math> per unit area acting on every point of the body surface, and body forces <math>F_i\,\!</math> per unit of volume on every point within the volume <math>V\,\!</math>. Thus, if the body is in [[Momentum#Relating_to_force_-_General_equations_of_motion|equilibrium]] the resultant force acting on the volume is zero, thus: | ||
| + | :<math>\int_S T_i^{(n)}dS + \int_V F_i dV = 0\,\!</math> | ||
| + | By definition the stress vector is <math>T_i^{(n)} =\sigma_{ji}n_j\,\!</math>, then | ||
| + | :<math>\int_S \sigma_{ji}n_j\, dS + \int_V F_i\, dV = 0\,\!</math> | ||
| + | Using the [[Divergence theorem|Gauss's divergence theorem]] to convert a surface integral to a volume integral gives | ||
| + | :<math>\int_V \sigma_{ji,j}\, dV + \int_V F_i\, dV = 0\,\!</math> | ||
| + | :<math>\int_V (\sigma_{ji,j} + F_i\,) dV = 0\,\!</math> | ||
| + | For an arbitrary volume the integral vanishes, and we have the ''equilibrium equations'' | ||
| + | :<math>\sigma_{ji,j} + F_i = 0\,\!</math> | ||
| + | |} | ||
===Układ === | ===Układ === | ||
Wersja z 08:56, 8 lis 2010
Spis treści |
Oscylator bistabliny
TEST
-
Derivation of equilibrium equations Consider a continuum body (see Figure 4) occupying a volume \(V\,\!\), having a surface area \(S\,\!\), with defined traction or surface forces \(T_i^{(n)}\,\!\) per unit area acting on every point of the body surface, and body forces \(F_i\,\!\) per unit of volume on every point within the volume \(V\,\!\). Thus, if the body is in equilibrium the resultant force acting on the volume is zero, thus: \[\int_S T_i^{(n)}dS + \int_V F_i dV = 0\,\!\] By definition the stress vector is \(T_i^{(n)} =\sigma_{ji}n_j\,\!\), then \[\int_S \sigma_{ji}n_j\, dS + \int_V F_i\, dV = 0\,\!\] Using the Gauss's divergence theorem to convert a surface integral to a volume integral gives \[\int_V \sigma_{ji,j}\, dV + \int_V F_i\, dV = 0\,\!\] \[\int_V (\sigma_{ji,j} + F_i\,) dV = 0\,\!\] For an arbitrary volume the integral vanishes, and we have the equilibrium equations \[\sigma_{ji,j} + F_i = 0\,\!\]
Układ
Rozważmy oscylator nieliniowy (usaa):
\(m \ddot x = -\frac{dU(x)}{dx}-\gamma x ,\)
równoważnie:
\(\dot x = v\)
\(m \dot v = -\frac{dU(x)}{dx}-\gamma x ,\)
jako potencjal weżmy funkcję z dwoma minimami np.:
\(U(x) = \displaystyle x^4-4 x^2,\)
Wykres U(x) można otrzymać poleceniem:
fplot(@(x) x.^4-4*x.^2,[-2.1,2.1],200)
\(f(x) =-\frac{dU(x)}{dx} = -4 x^3+8 x^2,\)
close all clear all f =@(x,E) sqrt(8*x.^2-2*x.^4-E); xa=linspace(-2,2,200); xb=linspace(2,-2,200); hold on for E=-8:2:8 Rts=roots([-2,0,8,0,-E]); AllRealRoots=sort(Rts(find (imag(Rts)==0))); if (length(AllRealRoots)==2) xminmax=sort([AllRealRoots',-0.00001,0.00001]); else xminmax=AllRealRoots'; endif xa=linspace(xminmax(1),xminmax(2),200); xb=linspace(xminmax(3),xminmax(4),200); if (E==0) plot([rotdim(xa,2),xa],[-f(rotdim(xa,2),E),f(xa,E)],"r-;E=0;"); plot([xb,rotdim(xb,2)],[f(xb,E),-f(rotdim(xb,2),E)],"r-"); else plot([rotdim(xa,2),xa],[-f(rotdim(xa,2),E),f(xa,E)],"g-"); plot([xb,rotdim(xb,2)],[f(xb,E),-f(rotdim(xb,2),E)],"g-"); endif endfor hold off
Analiza
lsode_options("absolute tolerance",1e-6) lsode_options("relative tolerance",1e-6) n=500; tmax=10; t = linspace(0,tmax,n); global gama; gama=0.9 ny=10; y=linspace(0,5,ny); x0=-sqrt(2); trajsxP=zeros(ny,n); trajsyP=zeros(ny,n); trajsxN=zeros(ny,n); trajsyN=zeros(ny,n); for iy=1:ny X0 = [x0,y(iy)]; [X,T,MSG]=lsode(@ODEbistable,X0,t); if(X(length(X),1)>0) trajsxP(iy,:)=X(:,1); trajsyP(iy,:)=X(:,2); else trajsxN(iy,:)=X(:,1); trajsyN(iy,:)=X(:,2); endif endfor disp("Integration done"); hold on every=20 for i=(every+1):every:size(trajsxP)(2) plot(trajsxP(:,(i-every):i)',trajsyP(:,(i-every):i)',"r-",\ trajsxN(:,(i-every):i)',trajsyN(:,(i-every):i)',"b-") drawnow() disp(t(i)); fflush(stdout); print ( sprintf("Traj_anim_%05d.png",i), "-S300,220"); endfor hold off
który możemy zaimplementować jako funkcję w matlabie:
function dx = ODEbistable(X,T) global gama; dx = zeros(2,1); dx(1) = X(2); dx(2) = -gama*X(2)-4*X(1).^3+8*X(1); return end
