PK EJBH mimetypetext/x-wxmathmlPK EJ@'v v content.xml
kill(all);definisco >0 proprietà geometriche e rigidezzeassume(a>0,b>0,k>0,k1>0,k2>0);definisco "non enormi" gli spostamentiassume(u<a,v<a,u<b,v<b,-u<a,-v<a,-u<b,-v<b);assume(l0>0);funzioni ausiliarie utilifunzione lunghezzalen(x1,y1,x2,y2):=sqrt((x2-x1)^2+(y2-y1)^2);componenti seno e coseno, definizione nonlineare P (x2,y2) ------- + ^ + | sin + | (x1,y1) O ----------------- | | | cos | |------>|attenti ai segni!!s(x1,y1,x2,y2):=(y2-y1)/len(x1,y1,x2,y2)$c(x1,y1,x2,y2):=(x2-x1)/len(x1,y1,x2,y2)$funzione rotazione angolare rispetto a indeformato (x3,y3) Q / P (x2,y2) /t + / + /+ (x1,y1) Ot(x1,y1,x2,y2,x3,y3):= atan2((x2*y3+x3*y1+x1*y2-x2*y1-x3*y2-x1*y3),(x3-x1)*(x2-x1)+(y3-y1)*(y2-y1));
controllot(1,2,1+1,2+1,1+0,2+1);t(1,2,1+1,2+1,1+1,2+0);t(0,0,0,1,1,-1);t(0,0,0,1,-1,-1);energia potenziale elastica dell'elemento barra incernierata a telaio al punto (x0,y0),connesso al punto mobile (x1,y1) con posizione deformata (x1+u,y1+v);k: rigidezza in estensione della molla lineares: rigidezza torsionale molla angolareU : 1/2*k*(len(x0,y0,x1+u,y1+v)-len(x0,y0,x1,y1))^2 + 1/2*s*t(x0,y0,x1,y1,x1+u,y1+v)^2;funzione potenziale, nel caso serva.U(x0,y0,x1,y1,u,v,k,s) := 1/2*k*(len(x0,y0,x1+u,y1+v)-len(x0,y0,x1,y1))^2 + 1/2*s*t(x0,y0,x1,y1,x1+u,y1+v)^2;matrice hessiana dell'energia potenziale elasticahessian(U,[u,v])$concide con la matrice di rigidezza tangenteforma nell'intorno di [u,v]=[0,0]ev(%,u=0,v=0);definisco matrice di rigidezza elementodefine(K(x0,y0,x1,y1,k,s),%);verifica su esempioK(0,0,0,a,k,s);suppongo ora l'elemento precaricato nella configurazione indeformata da uno sforzo normale N, positivo se trattivovado a considerare il lavoro svolto da tale forza sugli spostamenti (u,v), ovvero sull'allungamento da essi prodottoL : N * ( len(x0,y0,x1+u,y1+v) -len(x0,y0,x1 ,y1 ) );tale lavoro è da considerarsi un lavoro di forze interne, e si addiziona alla variazione di energia potenziale, da cui una correzione alla matrice di rigidezza nella formahessian(L,[u,v]);valuto nell'intorno di spostamenti nulliev(%,u=0,v=0);definisco la matrice di rigidezza geometrica = matrice di rigidezza da precaricodefine(G(x0,y0,x1,y1,N),%);FORMA ALTERNATIVA:in alternativa (e per controllo) posso definire tale lavoro come la differenza traa) l'energia potenziale elastica CON precarico iniziale (la condizione u=0,v=0 è precompressa di una lunghezza N/k), eb) l'energia potenziale elastica SENZA precarico inizialeLbis : fullratsimp( ( 1/2*k*(len(x0,y0,x1+u,y1+v)-len(x0,y0,x1,y1) + N/k )^2 + 1/2*s*t(x0,y0,x1,y1,x1+u,y1+v)^2 ) - ( 1/2*k*(len(x0,y0,x1+u,y1+v)-len(x0,y0,x1,y1) + 0/k )^2 + 1/2*s*t(x0,y0,x1,y1,x1+u,y1+v)^2 ));hessian(Lbis,[u,v]);FORMA ALTERNATIVA:valuto nell'intorno di spostamenti nulliev(%,u=0,v=0);FORMA ALTERNATIVA:definisco la matrice di rigidezza geometrica = matrice di rigidezza da precaricodefine(Gbis(x0,y0,x1,y1,N),%);FORMA ALTERNATIVA:le due forme sono EQUIVALENTIfullratsimp(Gbis(x0,y0,x1,y1,N)-G(x0,y0,x1,y1,N));esempio: | P | V O --------- | | | k | | | a |- | | \ s |____O___|____|_/////////////// A y||+---> xGradi di libertà:u : spostamento dir. x del punto di applicazione del carico v : spostamento dir. y del punto di applicazione del caricomatrice di rigidezza geometrica associata a precarico compressivo unitarioG(0,0,0,a,-1);matrice di rigidezza elasticaK(0,0,0,a, k,s);Linearized Pre Buckling analysisK(0,0,0,a, k,s)+P*G(0,0,0,a,-1);determinant(%);solve(%,P);in forma di estrazione autovalori/autovettoririduco a problema agli autovalori standard utilizzando(K + P*G).v =0(inv(K).K + P*inv(K).G).v = inv(K).0con K invertibile se correttamente vincolato( I + P*inv(K).G).v = 0( inv(K).G + (1/P) I ).v =0formalmente analogo a ( A - lambda I ).v =0con lambda = -1/P autovalorev autovettore[[vals,mults],vecs] : eigenvectors(invert(K(0,0,0,a, k,s)).G(0,0,0,a,-1));l'inverso degli autovalori coincide con l'opposto del carico critico -1/Pvals[1];mults[1];carico critico-1/vals[1];autovettori associati, potenzialmente multiplivecs[1];modo instabile con puro spostamento u, mentre v=0il secondo autovalore vals[2];è nullo ed è quindi associato a carico critico infinito.geometria di riferimento: |----- a ----| n2 n1|/P ==> o-----k1----o|/ ---- | |/ | | | k2 | | | b | | n3 o ----- ------- ///////nodo 2 spostamenti "u" in direzione x+ e "v" in direzione y+A y||+---> xGradi di libertà:u : spostamento dir. x del punto di applicazione del carico v : spostamento dir. y del punto di applicazione del caricomatrice di rigidezza elastica data dall'assemblaggio dei due elementiKel : K(0,0,0,b,k2,0) + K(0,b,a,b,k1,0);matrice di rigidezza geometrica associata a precompressione del solo elemento orizzontaleKG : G(0,0,0,b,0) + G(0,b,a,b,-P);Linearized Pre Buckling analysisdeterminant(Kel + KG);solve(%,P);Caso | P | V O -------------- k + + k | + + | b____O____ ____O____ ____|_///////// ///////// | | | | a | a | |-----|-----| A y||+---> xGradi di libertà:u : spostamento dir. x del punto di applicazione del carico v = -delta : spostamento dir. y del punto di applicazione del caricoangolo di inclinazione degli elementi rispetto all'orizzontale alpha t.c.sin(alpha)=b/sqrt(a^2+b^2)cos(alpha)=a/sqrt(a^2+b^2)dimensionamento : [b=1,a=2,k=1];matrice di rigidezza elastica data dall'assemblaggio dei due elementiKel : K(-a,0,0,b,k,0) + K(a,0,0,b,k,0);Kel:fullratsimp(Kel);matrice di rigidezza geometrica associata a precarico compressivo(P/2)/sin(alpha) su singola barraKG : ( G(-a,0,0,b,-1) + G(a,0,0,b,-1) ) * (P/2)/(b/sqrt(a^2+b^2));KG :fullratsimp(KG);Linearized Pre Buckling analysisdeterminant(Kel + KG);fullratsimp(%);solve(%,P);es. dimensionamento specifico, vedi sotto.Pcrit_LPBA : ev(%,dimensionamento,numer);rimuovo un'eventuale precedente assegnazione di Pkill(P);analisi nonlinearesomma contributi- energia potenziale elastica elemento (-a,0)-(0,b)- energia potenziale elastica elemento (+a,0)-(0,b)- energia potenziale carico FUtot:U(-a,0,0,b,u,v,k,0)+U(a,0,0,b,u,v,k,0)+P*v;definisco condizione di equilibrio come punto stazionario energia potenzialeequilibrio:[diff(Utot,u)=0,diff(Utot,v)=0];definisco P come equilibrante del caso u=0 (assumo spostamento orizzontale nullo)ev(equilibrio,u=0);linsolve(%[2],P),globalsolve=True;plotto P(delta) = P(-v) per un dimensionamento specificoqui delta è l'abbassamento del punto di applicazione del caricoplotme:ev(P,dimensionamento,v=-delta);wxplot2d(plotme,[delta,0,2]);osservo un punto stazionario di P(-v)=P(delta)deltacrit : find_root(diff(plotme,delta),delta,0.3,0.5);Pcrit: ev(plotme, delta=deltacrit),numer;noto che questo valore è notevolmente inferiore a quelli predetti dalla procedura Linearized pre-buckling analysis.Pcrit_LPBA;calcolo la matrice Hessiana del potenziale nel dimensionamento specifico e in v=-deltacrit,u=0hessian(Utot,[u,v]);ev(%,u=0,v=-deltacrit,dimensionamento),numer;si nota che l'hessiana della funzione potenziale (=matrice di rigidezza tangente) risulta singolare nell'intorno della condizione critica P=Pcrit, v=-deltacrit,u=0.Potrei quindi aspettarmi che il metodo proposto predica correttamente tale condizione.provo il metodo secantedefinisco un Hessiano "ridotto" 1x1, con u=0 e unica variabile vHred:diff( U(-a,0,0,b,0,v,k,0)+U(a,0,0,b,0,v,k,0), v, 2);definisco le due condizioni di riferimento per la secantecondizione 0; delta = -v = 0.3v0:-0.3;P0:ev(P,dimensionamento,v=v0,numer);Hred0 : ev(Hred,v=v0,dimensionamento,numer);condizione 1; delta = -v = 0.4v1:-0.4;P1:ev(P,dimensionamento,v=v1,numer);Hred1 : ev(Hred,v=v1,dimensionamento,numer);ricerco condizione di Hessiano ridotto singolarelinsolve( Hred0 + zeta * (Hred1-Hred0), zeta),globalsolve=True,numer;carico a cui ottengo tale condizioneP0 + zeta*(P1-P0);non lontanissimo dal massimo ottenuto su grafico se i punti di campionamento sono prossimi.
PK EJgn=> =>
image1.pngPNG
IHDR X W sBITO IDATx`Ms#3X36JgVE2Bu}T+a([
ǏL/Q&vOk{y>r;yϙ=uHj@DD$!i!i!i!i!i!i!i!i!iD`4[WfO>DG6FDDvA
a^^óǏĉRmƍgϞ.]IHDDA
aBBSXXN
.fԨQkB}#:u{yy-Rm
Sn݂+#"5rBh6]\\fs5SPP :ԺKt:U Fz}nnnbNN/Qe8H!NMM-ZLNNd weɒ%%!!!...$$@LLLIDDT\?r"""
k>ĉ4h ??VZk3f̢EJXwqZIKf ]vݷo_9rkvCz-
wuqdu[mҝ:RzXHXHXHXHXHXHXHXHXHXHXHXHXHXHXHgzRIVt3.\+/\
NNh͛
h|}ѢDD#