Continue to Site

Welcome to EDAboard.com

Welcome to our site! EDAboard.com is an international Electronics Discussion Forum focused on EDA software, circuits, schematics, books, theory, papers, asic, pld, 8051, DSP, Network, RF, Analog Design, PCB, Service Manuals... and a whole lot more! To participate you need to register. Registration is free. Click here to register now.

Need urgent help regarding Finite element method

Status
Not open for further replies.

kakar133

Member level 2
Joined
Sep 16, 2010
Messages
49
Helped
0
Reputation
0
Reaction score
0
Trophy points
1,286
Location
USA
Activity points
1,579
hi every one I am trying calculate stiffness matrix for 2D fem and using quadrilateral element(2 elements). but i am facing problem with the integration.
it is giving me warnings with out calculating stiffness matrix. the warning is
Warning: Cannot attach the property of being close to the limit point to limit
variable. [limit]
please help me
syms z n;
epsilon=4*8.86e-12;
u=4*pi*1e-7;
w=2*pi*300e6;
si(1)=1/4*(1-z).*(1-n);
si(2)=1/4*(1+z).*(1-n);
si(3)=1/4*(1+z).*(1+n);
si(4)=1/4*(1-z).*(1+n);
for e=1:1:2;
J11(e) =1/4*(-1*(1-n)*x(e)+(1-n)*x(e+1)+(1+n)*x(e+2)-(1+n)*x(e+3));
J12(e) =1/4*(-1*(1-n)*y(e)+(1-n)*y(e+1)+(1+n)*y(e+2)-(1+n)*y(e+3));
J21(e) =1/4*(-1*(1-z)*x(e)-(1-z)*x(e+1)+(1+z)*x(e+2)+(1+z)*x(e+3));
J22(e) =1/4*(-1*(1-z)*y(e)-(1-z)*y(e+1)+(1+z)*y(e+2)+(1+z)*y(e+3));
j_d(e)=J11(e)*J22(e)-J12(e)*J21(e);
si_x(e,1)=(1/(4.*j_d(e))).*(-1*J22(e).*(1-n)+J12(e).*(1-z));
si_y(e,1)=(1/(4.*j_d(e))).*(J21(e).*(1-n)-J11(e).*(1-z));
si_x(e,2)=(1/(4.*j_d(e))).*(J22(e).*(1-n)+J12(e).*(1+z));
si_y(e,2)=(1/(4.*j_d(e))).*(-1*J21(e).*(1-n)-J11(e).*(1+z));
si_x(e,3)=(1/(4.*j_d(e))).*(J22(e).*(1+n)-J12(e).*(1+z));
si_y(e,3)=(1/(4.*j_d(e))).*(-1*J21(e).*(1+n)+J11(e).*(1+z));
si_x(e,4)=(1/(4.*j_d(e))).*(-1*J22(e).*(1+n)-J12(e).*(1-z));
si_y(e,4)=(1/(4.*j_d(e))).*(J21(e).*(1+n)+J11(e).*(1-z));
for j=1:1:4;
for i=1:1:4
K(j,i)= double(int(int((((1/epsilon)*(si_x(e,j).*si_x(e,i)+si_y(e,j).*si_y(e,i))-((w^2)*u*si(j).*si(i))).*j_d(e)),z,-1,1),n,-1,1));
%K(j,i)=dblquad(f,-1,1,-1,1);
end
end
end
 

Status
Not open for further replies.

Part and Inventory Search

Welcome to EDABoard.com

Sponsor

Back
Top