my''(t)+2ay'(t)+by(t)=g(t), 0<t<=T
y(0)=0
y'(0)=0

y(t) - смещение в момент времени, 0<t<=T

m- масса

b - коэффициент упругости пружины

2a - коэффициент вязкого трения

g(t) - сила, действующая на систему в момент времени t,

причем в начальный момент t=0 тело покоится.

Лучше всего это в Simulink'e сделать

Преобразуем уравнение в систему дуфуров 1-го порядка:
y2' = (1/m) * ( g(t) - 2*a*y2(t) - b*y1(t) )
y1' = y2
y2(0) = 0
y1(0) = 0.

Пишем 2 М-файла.
1-й М-файл reshenie.m:

function fun = reshenie
clc
clear
T = 60; % vremya
[X, Y] = ode45( 'MyDifEq1', [0, T], [0, 0]);
% ode45( 'название ф-ии', интервал_на_котором_решаем_ДУ, нач.условия)
plot( X, Y(:,1) );
grid on



2-й М-файл MyDifEq1.m:

function F = MyDifEq1(t, y)
a = 0.1;
b = 3;
m = 2;
g = my_g(t);

F = [y(2); (1/m) * (g - 2*a*y(2) - b*y(1)) ];

function g = my_g(t)
% здесь задаем ф-ю g(t)
g = sin(t); % например

Ну, типа, все!

Guest,Спасибо большое!

Граммотный код, хотя не понял, так мало, точно наверное.
Решу как могу сам, скину в форум m-файл. Сверь, если захочешь.

код решения ДУ на MatLab.

Другая версия - думаю сносно.

Правильная или нет, решать Вам.

[/CODE]close all;
global D;
% Zadacha 1
% x``(t)+ax`(t)+kx(t)=0;
% x(o)=xo; x`(o)=xo_1;
% v^2+av+v=0; x=exp.^v*dt;

% Zadacha 2
% z``(t)+az`(t)+kz(t)=f(t);
% f(t)=A*exp(w*dt);
% z = p*exp(i*w*dt);

m=2;
k=3; % koeff. uprug. pruzhini;
a=10;%koeff. vyzkogo trenia;
A=0.025;% amplituda vinuzhdennyh kolebaniy;
w= 1;% chastota vinuzhdaychei silu;
%a1=A/(-w*w+i*a+k);
a1=a;
dt = 1:60;
d=sqrt(a^2-4*k*m);
b=sqrt(-a^2+4*k*m);
v_1=(-a-d)/2*a;
v_2=(-a+d)/2*a;

xo=0.5;
xo_1=xo./dt(1);

y = exp(v_1*dt);
z = exp(v_2*dt);
D.x=xo;
D.tr_x=[D.x];
if( D.x > 0.05 & D.x <1)
if(d>0)
D.x = ((v_2*xo-xo_1)/v_2-v_1)*y...
+((v_1*xo-xo_1)/v_1-v_2)*z;

elseif(d<0)
D.x = xo+(exp(length(dt)).^(-a*dt*0.5))*cos(0.5*b*dt) + ...
(2/b)*(xo_1+0.5*a*xo)*(exp(length(dt)).^(-a*dt*0.5))*sin(0.5*b*dt);

elseif(d==0)
D.x = (exp(length(dt))).^(-a*dt*0.5)

end
end
end


D.x_neodn = A*((k-w^2)/(k-w^2)^2+a1*a1)*cos(w*dt)-A*(a1/(k-w^2)^2+a1*a1)*sin(w*dt);
D.x_neodn1 = D.x_neodn + D.x;


figure;
plot(D.x);% reshenie zadachu2: x``(t)+ax`(t)+kx(t)=0; x(o)=xo, x`(o)=xo_1;
grid on;
zoom;

figure;
hold on;
%plot(D.x_neodn, 'r');
plot(D.x_neodn1);grid on;% reshenie zadachu1: % z``(t)+az`(t)+kz(t)=f(t);
hold off;
zoom on;

Здравствуйте!!!

В продолжение темы,

Вопрос: как построить картинку опыта в цикле( координаты не в явном виде) , так, чтобы задать движение и наблюдать визуально (изменение координат в цикле)?

конкретика не нужна, если есть,то круто конечно
Объясните, пожалуйста, на качественном уровне или примере.

Спасибо

P.S. Лобовой вариант без цикла:

close all;
% q-podves
% p-porshen
% p1-porshen_1
% pk-porshen_korpus
% pk_p-porshen_korpus_1
% t-telo
% pr-pruzhina

q_x=[0,0];
q_y=[8,16];

p_x=[0,3];
p_y=[10,10];

p1_x=[6,10];
p1_y=[10,10];

pk_x=[7,3,3,7];
pk_y=[11,11,9,9];

pk_p_x=[6,6,5,5,6];
pk_p_y=[11,9,9,11,11];

t_x=[10,13,13,10,10];
t_y=[9,9,15,15,9];

pr_x=[10,9,8,7,6,5,4,3,2,1,0];
pr_y=[14,14,15,13,15,13,15,13,15,14,14];

figure;
hold on;
p_q=plot(q_x,q_y);
p_p=plot(p_x,p_y);
p_p1=plot(p1_x,p1_y);
pk=plot(pk_x,pk_y);
pk_p=plot(pk_p_x,pk_p_y);
t=plot(t_x,t_y);
pr=plot(pr_x,pr_y);
hold off;
%- the end of the file

см.также выше file"sistema.doc"

Немного не понял, ты хочешь. чтобы масса на пружинке двигалась на картинке в соответствии с законом Y(t)? Что значит "координаты не в явном виде"?

Здравствуй, Guest

Y(t) - это идеальный вариант.
Стандарт - перерисовка, таrим образом, что изменяется массштаб.
Например, точка передвигается:

function H=DunDrawOnLine2
% close all;
%---Initial acts:
axis([0 10 0 1]);
dx=0.001;x=0;
y=sinc(cos(x));
X=x;
Y=y;
h=line(X,Y,'EraseMode','none');
%---Cycle of dynamic graphing:---
for i =1:10000
xn=x+dx;
yn=sinc(cos(xn));
X=[x,xn];
Y=[y,yn];
set(h,'XData',X,'YData',Y);
x=xn;
y=yn;
drawnow;
end


В задаче не точка, а куча массивов, причем сгрупированных hold on/hold off
(перепишу в виде: x и y массивы координат ) .

Вопрос: Как устроенна перерисовка т.е. стерание предыдущего положения?
Ведь наслаиваются изображения.

Не явное задание координат - это я намудрил, мягко говоря.

Движение прямоугольника с увеличивающимися сторонами по траектории у=sinc(cos(х)).

clc
clear
d_x_ = 0.1;
XX = [0:d_x_:10];
YY = sinc(cos(XX));
plot(XX, YY);
axis([-0.4 10.4 -0.1 1.1]);
hold on

dx=0.01; x=0;
y=sinc(cos(x));
X=x;
Y=y;
h=line(X,Y,'EraseMode','none');

%---Cycle of dynamic graphing:---
for i =1:1000

xn=x+dx;
yn=sinc(cos(xn));
d_x = d_x_/2 + (i*d_x_)/300;
X = [xn-d_x/2, xn-d_x/2, xn+d_x/2, xn+d_x/2, xn-d_x/2];
Y = [yn-d_x/8, yn+d_x/8, yn+d_x/8, yn-d_x/8, yn-d_x/8];
set(h,'XData',X,'YData',Y, 'Color','r');
axis([-0.4 10.4 -0.1 1.1]);

set(h,'XData',XX,'YData',YY, 'Color','b'); % закрашиваем белым цветом

pause(0.00001);
set(h,'XData',X,'YData',Y, 'Color','w');

x=xn;
y=yn;

end


Т.е. численно решаешь диф. уравнение каким-нибудь методом. На i итерации находишь y(i). Закрашиваешь предыдущее изображение белым цветом. Рисуешь новое, в качестве расстояния массы на пружинке принимая y(i). Нужно наиболее внимательно подойти к перерисовке пружины. Как ее масштабировать при изменении y(t) понятно?

Guest, Большое Тебе Спасибо

Буду дальше разбираться, груз на пружинке и тд.

Масштабрование - будут Мысли напишу. Ботать буду, пока не в теме.