Plate defs 2023

Series solver defs

Clear[SOLVE];

δ[exp_]:=exp/.{a[any__] -> δa[any],b[any__] -> δb[any],c[any__] -> δc[any],
(Subscript[u, X])[any__]->(Subscript[δu, X])[any],(Subscript[u, Y])[any__]->(Subscript[δu, Y])[any],(Subscript[u, Z])[any__]->(Subscript[δu, Z])[any],
(Subscript[θ, X])[any__]->(Subscript[δθ, X])[any],(Subscript[θ, Y])[any__]->(Subscript[δθ, Y])[any],(Subscript[θ, Z])[any__]->(Subscript[δθ, Z])[any],
(Subscript[u, x])[any__]->(Subscript[δu, x])[any],(Subscript[u, y])[any__]->(Subscript[δu, y])[any],(Subscript[u, z])[any__]->(Subscript[δu, z])[any],
(Subscript[θ, x])[any__]->(Subscript[δθ, x])[any],(Subscript[θ, y])[any__]->(Subscript[δθ, y])[any],(Subscript[θ, z])[any__]->(Subscript[δθ, z])[any],
ux[any__]->δux[any],uy[any__]->δuy[any],uz[any__]->δuz[any],
θx[any__]->δθx[any],θy[any__]->δθy[any],θz[any__]->δθz[any],
uX[any__]->δuX[any],uY[any__]->δuY[any],uZ[any__]->δuZ[any],
θX[any__]->δθX[any],θY[any__]->δθY[any],θZ[any__]->δθZ[any],
Fx[any__]->δFx[any],Fy[any__]->δFy[any],Fz[any__]->δFz[any],
Mx[any__]->δMx[any],My[any__]->δMy[any],Mz[any__]->δMz[any],
FX[any__]->δFX[any],FY[any__]->δFY[any],FZ[any__]->δFZ[any],
MX[any__]->δMX[any],MY[any__]->δMY[any],MZ[any__]->δMZ[any],
ϑ[any__]->δϑ[any]};

DOFS[exp_] :=Union[Cases[exp,a[__]|b[__]|c[__] |(Subscript[u, X])[__]|(Subscript[u, Y])[__]|(Subscript[u, Z])[__]| (Subscript[θ, X])[__]| (Subscript[θ, Y])[__]|(Subscript[θ, Z])[__]|(Subscript[u, x])[__]|(Subscript[u, y])[__]|(Subscript[u, z])[__]| (Subscript[θ, x])[__]| (Subscript[θ, y])[__]|(Subscript[θ, z])[__]|uX[__]|uY[__]|uZ[__]| θX[__]| θY[__]|θZ[__]|ux[__]|uy[__]|uz[__]| θx[__]| θy[__]|θz[__]|FX[__]|FY[__]|FZ[__]| MX[__]| MY[__]|MZ[__]|Fx[__]|Fy[__]|Fz[__]| Mx[__]| My[__]|Mz[__]| ϑ[__],Infinity]];

δDOFS[exp_] :=Union[Cases[exp,δa[__]|δb[__]|δc[__] |(Subscript[δu, X])[__]|(Subscript[δu, Y])[__]|(Subscript[δu, Z])[__]| (Subscript[δθ, X])[__]| (Subscript[δθ, Y])[__]|(Subscript[δθ, Z])[__]|(Subscript[δu, x])[__]|(Subscript[δu, y])[__]|(Subscript[δu, z])[__]| (Subscript[δθ, x])[__]| (Subscript[δθ, y])[__]|(Subscript[δθ, z])[__]|δuX[__]|δuY[__]|δuZ[__]| δθX[__]| δθY[__]|δθZ[__]|δux[__]|δuy[__]|δuz[__]| δθx[__]| δθy[__]|δθz[__]|δFX[__]|δFY[__]|δFZ[__]| δMX[__]| δMY[__]|δMZ[__]|δFx[__]|δFy[__]|δFz[__]| δMx[__]| δMy[__]|δMz[__]| δϑ[__],Infinity]];

SOLVE[δWW_]:=Module[{equ,sol},
equ=D[δWW,{δDOFS[δWW]}];
sol=Solve[Thread[equ==0],DOFS[δWW]];
Return[First[sol]]];

Plate strip xyn defs

equilibrium[xyn]:={
Nxx'[x]+bx,
Nxy'[x]+by,
Qx'[x]+bn,
Mxy'[x]-Qy[x]+cy,
Mxx'[x]-Qx[x]+cx};

constitutive[xyn]:={
Nxx[x]-EE*t/(1-ν^2)*u'[x],
Nyy[x]-EE*t/(1-ν^2)*ν*u'[x],
Nxy[x]-EE*t/(1-ν^2)*(1-ν)*v'[x]/2,
Qx[x]-GG*t*(w'[x]+θ[x])*α,
Qy[x]-GG*t*(-ϕ[x])*α,
Mxx[x]-DD*θ'[x],
Myy[x]-DD*ν*θ'[x],
Mxy[x]-DD*(1-ν)*(-ϕ'[x])/2
};

(*α=1 or α=5/6 Reissner-Mindlin*)
(*α->Infinity Kirchhoff*)

clamped[xyn][x_]:={u[x],v[x],w[x],ϕ[x],θ[x]};
simple[xyn][x_]:={u[x],v[x],w[x],ϕ[x],Mxx[x]};
slide[xyn][x_]:={Nxx[x],v[x],w[x],ϕ[x],Mxx[x]};
free[xyn][x_]:={Nxx[x],v[x],Qx[x],ϕ[x],Mxx[x]};

functions[xyn]:={u[x],v[x],w[x],ϕ[x],θ[x],Nxx[x],Nyy[x],Nxy[x],Qx[x],Qy[x],Mxx[x],Myy[x],Mxy[x]};

Plate strip rϕn defs

wn[r_]:=((bn (r^4))/(64 DD))+a+b (r^2) +c (r^2) (1-Log[r])+d Log[r]

θϕ[r_]:=-(∂_r)wn[r];

Mrr[r_]:=-DD((∂_(r,r))wn[r]+ν(1/r)(∂_r)wn[r]);

Mϕϕ[r_]:=-DD(ν(∂_(r,r))wn[r]+(1/r)(∂_r)wn[r]);

Qr[r_]:=(1/r)((∂_r)(r Mrr[r])-Mϕϕ[r]);

Examples

Lecture example 5.1 (Double sine series solution)

rule={L->10,H->10, t->0.01,ν->1/3, EE->210*10^9, 
DD -> EE*t^3/(12*(1-ν^2)),b->90};

wij=((16 b (H^4) (L^4) ((Sin[((i π)/2)])^2) ((Sin[((j π)/2)])^2) Sin[((i π x)/L)] Sin[((j π y)/H)])/(DD (H^4) (i^5) j (π^6)+2 DD (H^2) (i^3) (j^3) (L^2) (π^6)+2 DD i (j^5) (L^4) (π^6))) //.rule;

(*derived terms*)
Mxxij=-(EE*t^3/(12*(1-ν^2)))*((∂_(x,x))wij+ν*(∂_(y,y))wij) //.rule;
Myyij=-(EE*t^3/(12*(1-ν^2)))*((∂_(y,y))wij+ν*(∂_(x,x))wij)//.rule;
Mxyij=-(EE*t^3/(12*(1-ν^2)))*(1-ν)*(∂_(x,y))wij//.rule;
Qxij=(∂_x)Mxxij+(∂_y)Mxyij;
Qyij=(∂_y)Myyij+(∂_x)Mxyij;

(*quantities*)
n=25;
ws=ParallelSum[wij,{i,1,n,2},{j,1,n,2}]; 
Mxxs=ParallelSum[Mxxij,{i,1,n,2},{j,1,n,2}]; 
Myys=ParallelSum[Myyij,{i,1,n,2},{j,1,n,2}]; 
Mxys=ParallelSum[Mxyij,{i,1,n,2},{j,1,n,2}]; 
Qxs=ParallelSum[Qxij,{i,1,n,2},{j,1,n,2}]; 
Qys=ParallelSum[Qyij,{i,1,n,2},{j,1,n,2}]; 

fsize=14;
Grid[Partition[MapIndexed[Plot3D[#1,Evaluate[{x,0,L} /.rule],Evaluate[{y,0,H} /.rule],Boxed->False,Axes->True,PlotStyle->{GrayLevel[.25],Specularity[White,10]},
AxesLabel->{"x","y"},PlotLabel->{"w","Mxx","Myy","Mxy","Qx","Qy"}[[#2[[1]]]],Lighting->"Neutral",SphericalRegion->True,BaseStyle->{FontFamily->"Times",FontSize->fsize},PlotRange->{Automatic,Automatic,All}]&,{ws,Mxxs,Myys,Mxys,Qxs,Qys}],3]]

Lecture example 5.2 (Plate strip xyn)

rule ={bx->0,by->0,bn->-ρ*g*t,cx->0,cy->0,cz->0};

equs=Join[equilibrium[xyn],constitutive[xyn],clamped[xyn][0],clamped[xyn][L]] /.rule;  

sol=DSolve[Thread[equs==0],functions[xyn],x][[1]];

Print["w(x) = ",Simplify[Limit[w[x]/. sol ,α->1]]]
Print["θ(x) = ",Simplify[Limit[θ[x]/. sol ,α->1]]]
Print["Mxx(x) = ",Simplify[Limit[Mxx[x]/. sol ,α->1]]]

Lecture example 5.3 (Plate strip xyn)

rule ={bx->0,by->0,bn->-ρ*g*t,cx->0,cy->0,cz->0};

equs=Join[equilibrium[xyn],constitutive[xyn],clamped[xyn][0],clamped[xyn][L]] /.rule;  

sol=DSolve[Thread[equs==0],functions[xyn],x][[1]];
Print["w(x) = ",Simplify[Limit[w[x]/. sol ,α->Infinity]]]
Print["θ(x) = ",Simplify[Limit[θ[x]/. sol ,α->Infinity]]]

Lecture example 5.4 (Plate strip rϕn)

conditions:={Mrr[R],wn[R],-Qr[ϵ],-Mrr[ϵ]};

sol=Simplify[Solve[Thread[conditions==0],{a,b,c,d}]][[1]];

Limit[Simplify[wn[ϵ] /.sol],ϵ->0]

Lecture example 5.5 (Plate strip rϕn)

rule={bn->0};

conditions:={Mrr[R],wn[R],-Mrr[ϵ],-Qr[ϵ]+P/(2 Pi ϵ)} /.rule;

sol=Simplify[Solve[Thread[conditions==0],{a,b,c,d}]][[1]];

Limit[Simplify[wn[ϵ] /.sol],ϵ->0]

Lecture example 5.6 (Double sine series solution)

(*approximation*)
w[ser] = a[i,j]*Sin[Pi*i(x/L)]*Sin[Pi*j(y/L)];

(*virtual work expression*)
δW=Integrate[-DD*(δ[{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]}].{{1,ν,0},{ν,1,0},{0,0,(1-ν)/2}}.
{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]})+bn*δ[w[ser]],{x,0,L},{y,0,H}] ;

sol=Simplify[SOLVE[δW],{i∈Integers,j∈Integers}]

Example problem 4 (Plate strip xyn)

rule ={bx->0,by->0,bn->-ρ*g*t,cx->0,cy->0,cz->0};

equs=Join[equilibrium[xyn],constitutive[xyn],simple[xyn][0],slide[xyn][L]] /.rule;  

sol=DSolve[Thread[equs==0],functions[xyn],x][[1]];

Print["w(x) = ",Simplify[Limit[w[x]/. sol ,α->Infinity]]]
Print["θ(x) = ",Simplify[Limit[θ[x]/. sol ,α->Infinity]]]

Example problem 5 (Plate strip series solution)

(*approximation*)
nt=1;
w[ser]= Sum[a[i]*(1-x/L)^i*(x/L),{i,1,nt}];

(*virtual work expression*)
Ω =Rectangle[{0,0},{L,H}];
δW=Integrate[-DD*(δ[{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]}].{{1,ν,0},{ν,1,0},{0,0,(1-ν)/2}}.
{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]})+bn*δ[w[ser]],{x,y}∈Ω ,Assumptions->L>0&&H>0] ;

sol=SOLVE[δW];

w[ser] /.sol

Example problem 6 (Rectangle plate series solution)

(*approximation*)
w[ser]= a[0]*(1-x/L)*(x/L)*(1-y/H)*(y/H);

(*virtual work expression*)
Ω =Rectangle[{0,0},{L,H}];
δW=Integrate[-DD*(δ[{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]}].{{1,ν,0},{ν,1,0},{0,0,(1-ν)/2}}.
{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]})+bn*δ[w[ser]],{x,y}∈Ω ,Assumptions->L>0&&H>0] ;

sol=SOLVE[δW];

w[ser] /.sol

Example problem 7 (Rectangle plate series solution)

(*approximation*)
w[ser]= a[0]*Sin[Pi(x/L)]*Sin[Pi(y/H)];

(*virtual work expression*)
δW=Integrate[-DD*(δ[{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]}].{{1,ν,0},{ν,1,0},{0,0,(1-ν)/2}}.
{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]})+bn*δ[w[ser]],{x,0,L},{y,0,H}] ;

sol=SOLVE[δW];

w[ser]/.sol

Example problem 10 (Plate strip rϕn)

rule={bn->-ρ g t,DD->Ε*t^3/(12*(1-ν^2))};

conditions:={Mrr[R],wn[R],Qr[ϵ],Mrr[ϵ]}  /.rule;

sol=Simplify[Solve[Thread[conditions==0],{C[1],C[2],C[3],C[4]}]][[1]];

Limit[Simplify[wn[ϵ] /.sol],ϵ->0]

Assignment 1 (Plate strip xyn)

rule ={bx->0,by->0,bn->ρ*g*t,cx->0,cy->0,cz->0};

equs=Join[equilibrium[xyn],constitutive[xyn],clamped[xyn][0],free[xyn][L]] /.rule;  

sol=DSolve[Thread[equs==0],functions[xyn],x][[1]];

Print["Nxx(x) = ",Simplify[Limit[Nxx[x]/. sol ,α->1]]]
Print["Qx(x) = ",Simplify[Limit[Qx[x]/. sol ,α->1]]]
Print["Mxx(x) = ",Simplify[Limit[Mxx[x]/. sol ,α->1]]]

Assignment 5 (Rectangle xyn)

rule={bn->p,DD-> EE*(t^3)/(12 (1-((νν)^2)))};

(*approximation*)
w[ser]= a[0]*(x/L)*(1-x/L);

(*virtual work expression*)
Ω =Rectangle[{0,0},{L,H}];
δW=Integrate[-DD*(δ[{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]}].{{1,ν,0},{ν,1,0},{0,0,(1-ν)/2}}.
{(∂_(x,x))w[ser],(∂_(y,y))w[ser],2(∂_(x,y))w[ser]})+bn*δ[w[ser]]/.rule ,{x,y}∈Ω ,Assumptions->L>0&&H>0] ;

sol=SOLVE[δW]