from sympy import *
init_printing()
t = 20
t
C = 200000/(1-0.26**2) * Matrix([[1, 0.26, 0],
[0.26, 1, 0],
[0, 0, (1-0.26)/2]])
C
Polinomio de Lagrange
\begin{equation*} L_{i}(r)=\prod_{\substack{j=1 \\ j \neq i}}^{n} \frac{r-r_{j}}{r_{i}-r_{j}} \end{equation*}# https://www.longqi.cf/python/2014/03/24/implement-of-lagrange-polynomial-in-sympy/
def LagrangPoly(x,order,i,xi=None):
if xi==None:
xi=symbols('x:%d'%(order+1))
index = range(order+1)
index.pop(i)
return prod([(x-xi[j])/(xi[i]-xi[j]) for j in index])
N1r = simplify(LagrangPoly(symbols('r'),1,0,[0,1]))
N2r = simplify(LagrangPoly(symbols('r'),1,1,[0,1]))
N1s = simplify(LagrangPoly(symbols('s'),1,0,[0,1]))
N2s = simplify(LagrangPoly(symbols('s'),1,1,[0,1]))
N1 = simplify(N1r*N1s)
N1
N2 = simplify(N2r*N1s)
N2
N3 = simplify(N2r*N2s)
N3
N4 = simplify(N1r*N2s)
N4
N = Matrix([[N1, N2, N3, N4]])
N
x1, y1 = 0, 0
x2, y2 = 250, 0
x3, y3 = 250, 250
x4, y4 = 0, 250
x = simplify(N1*x1 + N2*x2 + N3*x3 + N4*x4)
x
y = simplify(N1*y1 + N2*y2 + N3*y3 + N4*y4)
y
r = symbols('x')/250
r
s = symbols('y')/250
s
J = Matrix([[diff(x, symbols('r')), diff(y, symbols('r'))],
[diff(x, symbols('s')), diff(y, symbols('s'))]])
J
detJ = det(J)
detJ
M1 = Matrix([[1, 0, 0, 0],
[0, 0, 0, 1],
[0, 1, 1, 0]])
M1
M2 = Matrix([[diff(r, symbols('x')), diff(s, symbols('x')), 0, 0],
[diff(r, symbols('y')), diff(s, symbols('y')), 0, 0],
[0, 0, diff(r, symbols('x')), diff(s, symbols('x'))],
[0, 0, diff(r, symbols('y')), diff(s, symbols('y'))]])
M2
M3 = Matrix([[diff(N1, symbols('r')), 0, diff(N2, symbols('r')), 0, diff(N3, symbols('r')), 0, diff(N4, symbols('r')), 0],
[diff(N1, symbols('s')), 0, diff(N2, symbols('s')), 0, diff(N3, symbols('s')), 0, diff(N4, symbols('s')), 0],
[0, diff(N1, symbols('r')), 0, diff(N2, symbols('r')), 0, diff(N3, symbols('r')), 0, diff(N4, symbols('r'))],
[0, diff(N1, symbols('s')), 0, diff(N2, symbols('s')), 0, diff(N3, symbols('s')), 0, diff(N4, symbols('s'))]])
M3
B = M1*M2*M3
B
K = integrate(integrate(transpose(B)*C*B*detJ*t, (symbols('s'), 0, 1)), (symbols('r'), 0, 1))
K.evalf(9)