Progetto Serbatoio

Definizione dati

Es=200000000;(*[(kN/(m^2))]*)
s1=0.010;(*[m]*)
s2=0.007;(*[m]*)
s3=0.005;(*[m]*)
R=6.5;(*[m]*)
ν=0.3;(*[-]*)
h1=7;(*[m]*)
h2=7;(*[m]*)
h3=7;(*[m]*)(* parte di serbatoio vuota*)
γ1=9.3157*1;(*[(kN/(m^3))]*)
γ2=9.3157*1;(*[(kN/(m^3))]*)
γ3=9.3157*1;(*[(kN/(m^3))]*)

Calcolo Delle Grandezze Caratteristiche

H=h1+h2+h3

D1=(((Es (s1^3)))/(12 (1-(ν^2))))(*[kNm]*)

D2=(((Es (s2^3)))/(12 (1-(ν^2))))(*[kNm]*)

D3=(((Es (s3^3)))/(12 (1-(ν^2))))(*[kNm]*)

c1=(((Es s1))/(R^2))(*[(kN/(m^3))]*)(*costante di sottofondo equivalente*)

c2=(((Es s2))/(R^2))(*[(kN/(m^3))]*)(*costante di sottofondo equivalente*)

c3=(((Es s3))/(R^2))(*[(kN/(m^3))]*)(*costante di sottofondo equivalente*)

β1=(((((Es*s1)/(4 D1 (R^2)))))^(1/4))(*[(1/m)]*)

β2=(((((Es*s2)/(4 D2 (R^2)))))^(1/4))(*[(1/m)]*)

β3=(((((Es*s3)/(4 D3 (R^2)))))^(1/4))(*[(1/m)]*)

λ1=2Pi/ β1(*[m]*)

λ2=2Pi/ β2(*[m]*)

λ3=2Pi/ β3(*[m]*)

Lc1=1/ β1(*[m]*)(*Lunghezza caratteristica del tubo*)

Lc2=1/ β2(*[m]*)(*Lunghezza caratteristica del tubo*)

Lc3=1/ β3(*[m]*)(*Lunghezza caratteristica del tubo*)

tipotubo1= (h1/(λ1)) (* se l>>λ tubo lungo, altrimenti tubo corto *)

tipotubo2= (h2/(λ2)) (* se l>>λ tubo lungo, altrimenti tubo corto *)

tipotubo3= (h3/(λ3)) (* se l>>λ tubo lungo, altrimenti tubo corto *)

Definizione Carichi

P1[x_]=γ1 (h1-x)+γ2 (h2)+γ3 (h3)
P2[x_]=γ2 (h2+h1-x)+γ3 (h3)
P3[x_]=γ3 (h3+h2+h1-x)

Equazioni di Campo

EQ1=D1 w1''''[x]+c1 w1[x]==P1[x];
EQ2=D2 w2''''[x]+c2 w2[x]==P2[x];
EQ3=D3 w3''''[x]+c3 w3[x]==P3[x];

sol=DSolve[{EQ1,EQ2,EQ3},{w1[x],w2[x],w3[x]},x];

Equazioni delle restanti grandezze

φ1[x_]=w1'[x];
M1[x_]=D1 w1''[x];
T1[x_]=-D1 w1'''[x];
N1[x_]=((Es s1)/R)w1[x];
φ2[x_]=w2'[x];
M2[x_]=D2 w2''[x];
T2[x_]=-D2 w2'''[x];
N2[x_]=((Es s2)/R) w2[x];
φ3[x_]=w3'[x];
M3[x_]=D3 w3''[x];
T3[x_]=-D3 w3'''[x];
N3[x_]=((Es s3)/R) w3[x];

w1[x_]=w1[x]/.sol[[1]];
w2[x_]=w2[x]/.sol[[1]];
w3[x_]=w3[x]/.sol[[1]];

Condizioni al Contorno

cc1=w1[0]==0;
cc2=φ1[0]==0;
cc3=w1[h1]-w2[h1]==0;
cc4=φ1[h1]-φ2[h1]==0;
cc5=T1[h1]-T2[h1]==0;
cc6=M1[h1]-M2[h1]==0;
cc7=w2[h2+h1]-w3[h2+h1]==0;
cc8=φ2[h2+h1]-φ3[h2+h1]==0;
cc9=T2[h2+h1]-T3[h2+h1]==0;
cc10=M2[h2+h1]-M3[h2+h1]==0;
cc11=w3[h3+h2+h1]==0;
cc12=φ3[h3+h2+h1]==0;

sol1=NSolve[{cc1,cc2,cc3,cc4,cc5,cc6,cc7,cc8,cc9,cc10,cc11,cc12},{C[1],C[2],C[3],C[4],C[5],C[6],C[7],C[8],C[9],C[10],C[11],C[12]},WorkingPrecision->1000];

w1[x_]=w1[x]/.sol1[[1]];
w2[x_]=w2[x]/.sol1[[1]];
w3[x_]=w3[x]/.sol1[[1]];

Grafici

(*Spostamento della parete*)

s1=Plot[w1[x],{x,0,h1},PlotRange->All,Filling->Axis];

s2=Plot[w2[x],{x,h1,h2+h1},PlotRange->All,Filling->Axis];

s3=Plot[w3[x],{x,h1+h2,h2+h1+h3},PlotRange->All,Filling->Axis];

Show[s1,s2,s3,PlotRange->All]

(*Rotazione del generico meridiano*)

r1=Plot[φ1[x],{x,0,h1},PlotRange->All,Filling->Axis];
r2=Plot[φ2[x],{x,h1,h1+h2},PlotRange->All,Filling->Axis];
r3=Plot[φ3[x],{x,h1+h2,h1+h2+h3},PlotRange->All,Filling->Axis];

Show[r1,r2,r3,PlotRange->All]

(*taglio del generico meridiano*)

t1=Plot[T1[x],{x,0,h1},PlotRange->All,Filling->Axis];
t2=Plot[T2[x],{x,h1,h1+h2},PlotRange->All,Filling->Axis];
t3=Plot[T3[x],{x,h1+h2,h1+h2+h3},PlotRange->All,Filling->Axis];

Show[t1,t2,t3,PlotRange->All]

T1[x_]=-D1  w1'''[0]

T2[x_]=-D2  w2'''[h1]

T3[x_]=-D3  w3'''[h1+h2]

T3[x_]=-D3  w3'''[h1+h2+h3]

(*momento del generico meridiano*)

m1=Plot[M1[x],{x,0,h1},PlotRange->All,Filling->Axis]
m2=Plot[M2[x],{x,h1,h1+h2},PlotRange->All,Filling->Axis]
m3=Plot[M3[x],{x,h1+h2,h1+h2+h3},PlotRange->All,Filling->Axis]

Show[m1,m2,m3,PlotRange->All]

M1[x_]=D1  w1''[0]

M2[x_]=D2  w2''[h1]

M3[x_]=D3  w3''[h1+h2]

M3[x_]=D3  w3''[h1+h2+h3]

(*sforzo normale di parallelo*)

n1=Plot[N1[x],{x,0,h1},PlotRange->All,Filling->Axis];
n2=Plot[N2[x],{x,h1,h1+h2},PlotRange->All,Filling->Axis];
n3=Plot[N3[x],{x,h1+h2,h1+h2+h3},PlotRange->All,Filling->Axis];

Show[n1,n2,n3,PlotRange->All]

(*si ha una discontinuità quando c'è una variazione di spessore*)