[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [Getfem-users] Badly Scaled Bilaplacian
From: |
Yves Renard |
Subject: |
Re: [Getfem-users] Badly Scaled Bilaplacian |
Date: |
Tue, 17 Mar 2015 19:18:16 +0100 (CET) |
Dear Oladayo Eluyefa,
Your stiffness matrix is ill-conditionned because you use C0 finite elements
instead of hermite elements. For this problem you do have to use hermite
element such as argyris one (FEM_ARGYRIS) but for triangles or FVS element
(FEM_QUADC1_COMPOSITE) on quadrilaterals (but with a specific integration
method). Additionally, of course, the stiffness matrix will be singular unless
a Dirichlet condition is prescribed.
Note also that the Kirchhoff-love term is better conditionned that the
bi-laplacian term and equivalent but with different boundary conditions.
Yves.
----- Original Message -----
From: "OLADAYO ELUYEFA (RIT Student)" <address@hidden>
To: address@hidden
Sent: Tuesday, March 17, 2015 6:17:22 PM
Subject: [Getfem-users] Badly Scaled Bilaplacian
I am working on a 4th order 2D Inverse Problem using the matlab OO interface in
getfem++ . As a result, I calculate the bilaplacian, but I keep getting a badly
scaled matrix with a close to 0 condition number. Below is the code that
produces this. Is there a way to correct this.
Thank you
Code:
%% Note equation is -(au'')''=f
N = 20;
h = 1/(N+1);
m = gf_mesh('cartesian', 0:h:1, 0:h:1);
mf = gf_mesh_fem(m,1);
mfd = gf_mesh_fem(m,1);
[xs, ys] = meshgrid(0:1*h:1, 0:1*h:1);
gf_mesh_fem_set(mf, 'fem', gf_fem('FEM_QK(2,4)'));
gf_mesh_fem_set(mfd, 'fem', gf_fem('FEM_QK(2,4)'));
mim = gf_mesh_im(m, gf_integ('IM_GAUSS_PARALLELEPIPED(2,10)'));
% Boundary conditions:
border = gf_mesh_get(m, 'outer_faces');
gf_mesh_set(m, 'boundary', 1, border);
% Degrees of freedom
dof = gf_mesh_fem_get(mf, 'nbdof');
dofa = gf_mesh_fem_get(mfd, 'nbdof');
u = '256.*(x.^2).*(y.^2).*((1-x).^2).*(1-y).^2';
Us = gf_mesh_fem_get(mfd, 'eval', {u});
a = '1+0.25.*sin(pi.*x).*sin(pi.*y)';
avec = gf_mesh_fem_get(mfd, 'eval', {a});
a0 = gf_mesh_fem_get(mfd, 'eval', {1.5})';
K = @(a) gf_asm('bilaplacian', mim, mf, mfd, a);
K1 = gf_asm('laplacian', mim, mfd, mfd, ones(dofa,1));
L = @(u) gf_asm('volumic',...
['u=data(#2);', ...
't=comp(Hess(#2).Hess(#2).Base(#1))(:,i,i,:,i,i,:);',...
'M(#2,#1)+=t(j,:,:).u(j)'],...
mim, ...
mfd, ...
mf, ...
u);
M = gf_asm('mass matrix', mim, mf);
Ma = gf_asm('mass matrix', mim, mfd);
Oladayo Eluyefa,
Rochester Institute of Technology '15,
B.S., Computational Mathematics,
M.S., Applied and Computational Mathematics.
( address@hidden , address@hidden )
(585)-287-1955
CONFIDENTIALITY NOTE : The information transmitted, including attachments, is
intended only for the person(s) or entity to which it is addressed and may
contain confidential and/or privileged material. Any review, re-transmission,
dissemination or other use of, or taking of any action in reliance upon this
information by persons or entities other than the intended recipient is
prohibited. If you received this in error, please contact the sender and
destroy any copies of this information.
Please consider the environment before printing this email.
_______________________________________________
Getfem-users mailing list
address@hidden
https://mail.gna.org/listinfo/getfem-users