Integrazione di Newmark in Matlab...
Ciao ragazzi, mi trovo a dover implementare in matlab la soluzione di Newmark per ricavare spostamento e velocità di un sistema accoppiato a 2gdl.
Ho implementato il seguente codice, fino ad ottenere lo spostamento all'istante (tk+1) partento da condizioni iniziali (tk) nulle (stato di quiete). La forzante sollecitante il sistema primario è di tipo sinusoidale.
Le equazioni che ho utilizzato sono:

Uploaded with ImageShack.us
(scusatemi per le dimensioni dell'img)
Il procedimento che dovrei implementare:
1)risolvere la prima eq. e ricavare u2dot(tk+1), le quantità in funzione di (tk) sono all'istante iniziale=zero e quindi tutte note; la forzante p(tk+1) ha il valore della funzione sen(t) all'istante t(2).
2) sostituire il risultato trovato nell'eq. 2 e nell'eq. 3, trovando così rispettivamente spostamento e velocità allìistante (tk+1);
3)sostituire le quantita trovate nuovamente nell'eq. 1, sostituendole alle quantità in funzione di (tk) trovando stavolta la u2dot(tk+2)
4) Iterara il procedimento tramite ciclo for (che non ho saputo implementare) e ricavare così spostamento,velocità e accelerazione al variare del tempo
P.S: credo che i u0; udot0;u2dot0; debbano essere vettori 2x1 solo che sorgono problemi quando li definisco con tali dimensioni.
Trascrivo il codice:
clear all
close all
clc
%.......................DATI...........................
zita_p=0.05; %rapp smorzamento sistema Principale
zita_s=0.1; %rapp smorzamento sistema Secondario
gamma_=0.1; %rapp tra ms e mp; più è alto più il TMD è efficace
wp=1; %ATTENZIONE-?
ws=1; %ATTENZIONE-?
t=[0:0.1:100]; %vettore dei tempi **
u0=0; %è la u(tk)
udot0=0; %è la udot(tk)
u2dot0=0;%è la u2dot(tk)
dt=t(2)-t(1);
%M K C rispettivamente M. Masse, M. Rigidezza, M. Dissipazione (vedi
%quaderno tesi)
M=blkdiag(1,1);
K=[wp^2+ws^2*gamma_,-ws^2*gamma_;-ws^2,ws^2];
C=[(2*zita_p*wp+2*zita_s*ws*gamma_),-2*zita_s*ws*gamma_;-2*zita_s*ws,2*zita_s*ws];
%devo definire una forzante sinusoidale
p=[sin(t),0]';
%Procedimento di integrazione....
acc=(inv(M+(dt/2)*C+0.25*dt^2*K))*(p(2)-C*(udot0+dt/2*u2dot0)+K*(u0+dt*udot0+(0.5-0.25)*dt^2*u2dot0));
%Spostamento all'istante tk+1
u2=u0+dt*udot0+((0.5-0.25)*u2dot0+0.25*acc)*dt^2;
%Accelerazione all'istante tk+1
u2dot=udot0+dt/2*(u2dot0+acc);
.........................................................FINE Codice
GRAZIE A TUTTI COLORO CHE VORRANO DARMI UNA MANO!!!
Ho implementato il seguente codice, fino ad ottenere lo spostamento all'istante (tk+1) partento da condizioni iniziali (tk) nulle (stato di quiete). La forzante sollecitante il sistema primario è di tipo sinusoidale.
Le equazioni che ho utilizzato sono:

Uploaded with ImageShack.us
(scusatemi per le dimensioni dell'img)
Il procedimento che dovrei implementare:
1)risolvere la prima eq. e ricavare u2dot(tk+1), le quantità in funzione di (tk) sono all'istante iniziale=zero e quindi tutte note; la forzante p(tk+1) ha il valore della funzione sen(t) all'istante t(2).
2) sostituire il risultato trovato nell'eq. 2 e nell'eq. 3, trovando così rispettivamente spostamento e velocità allìistante (tk+1);
3)sostituire le quantita trovate nuovamente nell'eq. 1, sostituendole alle quantità in funzione di (tk) trovando stavolta la u2dot(tk+2)
4) Iterara il procedimento tramite ciclo for (che non ho saputo implementare) e ricavare così spostamento,velocità e accelerazione al variare del tempo
P.S: credo che i u0; udot0;u2dot0; debbano essere vettori 2x1 solo che sorgono problemi quando li definisco con tali dimensioni.
Trascrivo il codice:
clear all
close all
clc
%.......................DATI...........................
zita_p=0.05; %rapp smorzamento sistema Principale
zita_s=0.1; %rapp smorzamento sistema Secondario
gamma_=0.1; %rapp tra ms e mp; più è alto più il TMD è efficace
wp=1; %ATTENZIONE-?
ws=1; %ATTENZIONE-?
t=[0:0.1:100]; %vettore dei tempi **
u0=0; %è la u(tk)
udot0=0; %è la udot(tk)
u2dot0=0;%è la u2dot(tk)
dt=t(2)-t(1);
%M K C rispettivamente M. Masse, M. Rigidezza, M. Dissipazione (vedi
%quaderno tesi)
M=blkdiag(1,1);
K=[wp^2+ws^2*gamma_,-ws^2*gamma_;-ws^2,ws^2];
C=[(2*zita_p*wp+2*zita_s*ws*gamma_),-2*zita_s*ws*gamma_;-2*zita_s*ws,2*zita_s*ws];
%devo definire una forzante sinusoidale
p=[sin(t),0]';
%Procedimento di integrazione....
acc=(inv(M+(dt/2)*C+0.25*dt^2*K))*(p(2)-C*(udot0+dt/2*u2dot0)+K*(u0+dt*udot0+(0.5-0.25)*dt^2*u2dot0));
%Spostamento all'istante tk+1
u2=u0+dt*udot0+((0.5-0.25)*u2dot0+0.25*acc)*dt^2;
%Accelerazione all'istante tk+1
u2dot=udot0+dt/2*(u2dot0+acc);
.........................................................FINE Codice
GRAZIE A TUTTI COLORO CHE VORRANO DARMI UNA MANO!!!
Risposte
nessuno può darmi una mano? ne avrei davvero bisogno...
grazie
grazie