How to evaluate a symbolic matrix integration using three-point gauss rule?

Hello everyone I'm working on a routine to solve a stiffness matrix about hexahedral finite element. In order to improve my routines, I'm using three-point Gauss Quadrature Rule. The matrix to want solve has syms variable. So, when I want to use "eval" comand, matlab says something like this..

??? Error using ==> evalin
Undefined function or method 'conjugate' for input arguments of type 'double'.

Error in ==> sym.eval at 15
    s = evalin('caller',vectorize(map2mat(char(x))));

Error in ==> Khexaedro8nodos at 63
    int2=eval(int1)*w(i)+int2;

The code is as follow

% STIFFNESS MATRIX
clc
clear all
syms e n c
 coor_global=[0            0            0
            3            0            0
            3        0.125            0
            0        0.125            0
            0            0          0.2
            3            0          0.2
            3        0.125          0.2
            0        0.125          0.2];
        E=1;
        nu=0.2;
coor_local=[-1 -1 -1;1 -1 -1;1 1 -1;-1 1 -1;-1 -1 1;1 -1 1;1 1 1;-1 1 1];
aux=[1+coor_local(:,1)*n 1+coor_local(:,2)*c 1+coor_local(:,3)*e];
N=prod(aux')';
dN_local=[diff(N,n) diff(N,c) diff(N,e)];
J=double(expand(dN_local'*coor_global));
detJ=det(J);
dN=J\dN_local';
lambda=E*nu/((1+nu)*(1-2*nu));
mu=E/(2*(1+nu));
D=vpa([lambda+2*mu lambda lambda 0 0 0;...
    lambda lambda+2*mu lambda 0 0 0;...
    lambda lambda lambda+2*mu 0 0 0;...
    0 0 0 mu 0 0;0 0 0 0 mu 0;0 0 0 0 0 mu]);
B=vpa(zeros(24,6));
for j=1:8
    Bi=expand([dN(1,j) 0 0;...
        0 dN(2,j) 0;...
        0 0 dN(3,j);...
        dN(2,j) dN(1,j) 0;...
        0 dN(3,j) dN(2,j);...
        dN(3,j) 0 dN(1,j)]);
    B(3*j-2:3*j,:)=Bi';
end
% Very slow. This calculating consumes high memory. Instead, to use Gauss-Legendre.
% Klocal=detJ*double(int(int(int(B*D*B',n,-1,1),c,-1,1),e,-1,1));

% Five-point Gauss Quadrature Rule
% x=[-0.90617984 -0.53846931 0 0.53846931 0.90617984];
% w=[0.23692688509 0.4786286705 0.56888888 0.4786286705 0.23692688509];

% Three-point Gauss Quadrature Rule
x=[-0.7745966 0 0.7745966];
w=[0.55555 0.88888 0.55555];
npi=size(x,2);

aux1=0;
aux2=0;
aux3=0;

tic
int1=0;
for i=1:npi
    n=x(i);
    int1=eval(B*D*B')*w(i)+int1;
end
toc

int2=0;
for i=1:npi
    c=x(i);
    int2=eval(int1)*w(i)+int2;
end

int3=0;
for i=1:npi
    e=x(i);
    int3=eval(int2)*w(i)+int3;
end

Klocal=detJ*double(aux3);

pd: Sorry about my english. I try.

Answer 1

[Walter Roberson]

Do not eval() symbolic expressions. You can use subs() and vpa() and double(), and matlabFunction()

Symbolic expressions are written in terms of the MuPAD language, which has some differences in function names and syntax compared to MATLAB itself.