Matlab - Programa de metodo elemntos finitos libro Young W Kwon

 
Vista:
sin imagen de perfil
Val: 8
Ha aumentado 1 puesto en Matlab (en relación al último mes)
Gráfica de Matlab

Programa de metodo elemntos finitos libro Young W Kwon

Publicado por wil (3 intervenciones) el 21/02/2020 13:22:47
Hola, digite el programa del libro Young W Kwon y nao compila , donde estara el error???


Aqui el programa principal:


1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
% Exemplo 351.m
% para resolver a EDO dad por:
% au''+bu'+cu=1, 0<x<1
% u(0)= e u(1)=0
% usando 5 elementos
%
% descripcao de variaveis
% k= matriz elemento
% f= vetor elemento
% kk= matriz sistema
% ff= vetor sistema
% index= um vetor contendo um sistema associado dofs de cada elemento
% bcdof= um vetor contendo dofs associado com condicoes de fronteira
% bcval= um vetor contendo condicao de valores na fronteira associado com o dofs em bcdof
%
%--------------------------------------------------------------------------%
% entrada de datos
%-------------------------------------------------------------------------
nel=5;    % numero de elementos
nnel=2;    % numero de nodos por elemento
ndof=1;    % numero de dofs por nodo
nnode=6;     %numero total de nodos del sistema
sdof=nnode*ndof;    %dofs do sistema total
%
%--------------------------------------------
% datos de entrada de valores coordenados de cada nodo
%-------------------------------------------
gcoord(1)=0.0 ; gcoord(2)=0.2; gcoord(3)=0.4; gcoord(4)=0.6;
gcoord(5)=0.8; gcoord(6)=1.0;
%
%-------------------------------------------
% datos de entrada de conectividade nodal para cada elemento
%-------------------------------------------
nodes(1,1)=1; nodes(1,2)=2; nodes(2,1)=2; nodes(2,2)=3;
nodes(3,1)=3; nodes(3,2)=4; nodes(4,1)=4; nodes(4,2)=5;
nodes(5,1)=5; nodes(5,2)=6;
%
%------------------------------------------------------
% coeficientes de entrada da EDO
%-----------------------------------------------------
acoef=1;
bcoef=-3;
ccoef= 2;
%
%-----------------------------------------
%datos de entrada para a condicao de fronteira
%-------------------------------------------
bcdof(1)=1;
bcval(1)=0;
bcdof(2)=6;
bcval(2)=0;
%
%------------------------------------------------
% inicializacao de vetores e matrizes
%-----------------------------------------------
ff= zeros(sdof,1);     % inicializacao do sistema do vetor forca
kk=zeros(sdof,sdof);           % inicializacao do sistema matricial
index=zeros(nnel*ndof,1);           % inicializacao do vetor index
%
%------------------------------------------------
% calculo de elemntos das matrizes e vetores e seu emsamble
%---------------------------------------------------
for iel=1:nel                       % loop para o numero total de elementos
%
nl=nodes(iel,1); nr=nodes(iel,2);
xl=gcoord(nl); xr=gcoord(nr);
eleng=xr-xl;
index=feeldof1(iel,nnel,ndof);
%
k= feode2l(acoef,bcoef,ccoef,eleng);
f= fef1l(xl,xr);
[kk,ff]=feasmbl2(kk,ff,k,f,index);
%
end
%
%----------------------------------------
%aplicando condicoes de fronteira
%-----------------------------------------
[kk,ff]=feaplyc2(kk,ff,bcdof,bcval);
%kk, ff, k, f, index
%--------------------------------------
% resolve a equacao matricial
%------------------------------------
fsol=kk\ff;
%
%------------------------------------------------
c1= 0.5/exp(1);
c2=-0.5*(1+1/exp(1));
for i=1:nnode
    x=gcoord(i);
    esol(i)=c1*exp(2*x)+c2*exp(x)+1/2;
end
%
%------------------------------------------
% imprime ambas solucoes exata e por mef
%----------------------------------------
num=1:1:sdof;
results=[num' fsol  esol']
%-----------------------------------------------


Ahora las subrutinas:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
function [kk,ff]=feaplyc2(kk,ff,bcdof,bcval)
%--------------------------------------------------------
% aqui ha dof=2 e 10
%bcdof(1)=2 e bcdof(2)=10;
%bcval(1)=1.0 e bcval(2)=2.5
%---------------------------------------------
n=length(bcdof);
sdof=size(kk);
%
for i=1:n
    c=bcdof(i);
    for j=1:sdof
        kk(c,j)=0;
    end
    %
    kk(c,c)=1;
    ff(c)=bcval(i);
end
end


SEGUNDA SUB-RUTINA

1
2
3
4
5
6
7
8
function [index]=feeldof1(iel,nnel,ndof)
edof=nnel*ndof;
start=(iel-1)*(nnel-1)*ndof;
%
for i=1:edof
    index(i)=start+i;
end
end


TERCERA SUBRUTINA

1
2
3
4
5
6
function [f]=fef1l(xl,xr)
%----------------------------------
% elemnto vetor para f(x)=1
% usando elemento linear
 eleng=xr-xl;
f=[eleng/2;eleng/2];


QUARTA SUB RUTINA

1
2
3
4
5
6
function [k]=feode2l(acoef,bcoef,ccoef,eleng)
%---------------------------------------
a1=-(acoef/eleng); a2=bcoef/2; a3=ccoef*eleng/6;
k=[a1-a2+2*a3,  -a1+a2+a3;...
    -a1-a2+a3 , a1+a2+2*a3];
%----------------------------


QUINTA SUBRUTINA

1
2
3
4
5
6
7
8
9
10
11
function [kk,ff]=feasmbl2(kk,ff,k,f,index)
%--------------------------
edof=length(index);
for i=1:edof
    ii=index(i);
    ff(ii)=f(ii)+f(i);
    for j=1:edof
        jj=index(j);
        kk(ii,jj)=kk(ii,jj)+k(i,j);
    end
end


ME SALE ERROR EN LA QUINTA SUBRUTINA Y NOSE COMO ARREGRARLO PARA QUE COMPILE, COPIE EL CODIGO TAL Y COMO ESTA EN EL LIBRIO, ALGUNA AYUDA, GRACIAS .
Valora esta pregunta
Me gusta: Está pregunta es útil y esta claraNo me gusta: Está pregunta no esta clara o no es útil
0
Responder
Imágen de perfil de JOSE JEREMIAS CABALLERO
Val: 7.752
Oro
Ha mantenido su posición en Matlab (en relación al último mes)
Gráfica de Matlab

Programa de metodo elemntos finitos libro Young W Kwon

Publicado por JOSE JEREMIAS CABALLERO (5506 intervenciones) el 21/02/2020 13:38:26
1
2
3
4
5
6
7
8
9
10
11
12
13
function [kk,ff]=feasmbl2(kk,ff,k,f,index)
%--------------------------
edof=length(index);
for i=1:edof
    ii=index(i);
    ff(ii)=ff(ii)+f(i)   ;
    for j=1:edof
        jj=index(j);
        kk(ii,jj)=kk(ii,jj)+k(i,j);
    end
    end
end
% -----------------------------------------------

1
2
3
4
5
6
7
8
9
10
>> metodos_elementos_finitos_young
 
results =
 
    1.0000   -0.0000         0
    2.0000   -0.0621   -0.0610
    3.0000   -0.1133   -0.1110
    4.0000   -0.1388   -0.1355
    5.0000   -0.1142   -0.1111
    6.0000         0   -0.0000



Saludos
JOSE JEREMIAS CABALLERO
Asesor de Proyectos con Matlab
Servicios de programación matlab
Servicio de Asesoría Online en Matlab
jjcc94@hotmail.com


http://matlabcaballero.blogspot.com
https://www.facebook.com/matlabcaballero
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar
sin imagen de perfil
Val: 8
Ha aumentado 1 puesto en Matlab (en relación al último mes)
Gráfica de Matlab

Programa de metodo elemntos finitos libro Young W Kwon

Publicado por wil (3 intervenciones) el 21/02/2020 13:59:31
Hola, y donde estaba el erro??
pues lo coloco asi como lo digistaste y me sale :

Attempted to access f(3); index out of bounds because numel(f)=2.

Error in feasmbl2 (line 6)
ff(ii)=f(ii)+f(i);

Error in exemplo351 (line 77)
[kk,ff]=feasmbl2(kk,ff,k,f,index);
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar
Imágen de perfil de JOSE JEREMIAS CABALLERO
Val: 7.752
Oro
Ha mantenido su posición en Matlab (en relación al último mes)
Gráfica de Matlab

Programa de metodo elemntos finitos libro Young W Kwon

Publicado por JOSE JEREMIAS CABALLERO (5506 intervenciones) el 21/02/2020 14:02:07
Estimado no está colocando como lo he digitado. Esta expresión ff(ii)=f(ii)+f(i); no se encuentra en la rutina que he puesto al foro.
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar
sin imagen de perfil
Val: 8
Ha aumentado 1 puesto en Matlab (en relación al último mes)
Gráfica de Matlab

Programa de metodo elemntos finitos libro Young W Kwon

Publicado por wil (3 intervenciones) el 21/02/2020 14:16:15
muy bien gracias , necesito unos lentes.
Valora esta respuesta
Me gusta: Está respuesta es útil y esta claraNo me gusta: Está respuesta no esta clara o no es útil
0
Comentar