Hi Robert,

I compared the displacement calculated by sfepy and by Comsol (fine mesh, both for free boundary condition and fixed bottom boundary condition). No matter what boundary condition it is, sfepy result and Comsol result agree with each other perfectly. I mean sfepy free boundary result == Comsol free boundary result, sfepy fixed bottom result == Comsol fixed bottom result. But free boundary result =! fixed bottom result. So, I still have three questions, please forgive my ignorance on FEM:
1. you mentioned, with free boundary condition, there will be UmfpackWarning: (almost) singular matrix!. However, I do not have, why ?
2. In Comsol, the general solver's relative tolerance is 0.0001, in this case,  fixed bottom boundary condition can be solved successfully, but free boundary condition cannot. Only if I set relative tolerance to be over 1.5, then free boundary condition can be solved. So, the question is, if it is the case that, the free boundary condition actually cannot be solved theoretically, but if we do some tricks to the solver, it might be solved. So, can I say the free boundary condition result is kind of 'artificial', therefore not very physically meaningful?
3. Assuming I have a 3*3 quad mesh with following vertex coordinates
(0,3)  (1,3)  (2,3)  (3,3)
(0,2)  (1,2)  (2,2)  (3,2)
(0,0)  (1,1)  (2,1)  (3,1)
(0,0)  (1,0)  (2,0)  (3,0)
if the displacement result is the displacement at these vertices?

Assuming I use strain = pb.evaluate('ev_cauchy_strain.2.Omega(u)', mode='el_avg'), if the strain result is the strain at cell centers?

Regards
Ronghai


在 2014年11月17日星期一UTC+1下午3时47分20秒，Robert Cimrman写道：

On 11/17/2014 11:04 AM, Ronghai Wu wrote:
> Hi Robert,
>
> Thank you very much. But I do not understand what you mean by "Note that
> you still need to specify some boundary conditions so that the matrix is
> not singular". Because it seems no boundary condition is constrained in the
> modified codes, but it still runs, and the result looks correct. Would you
> please explain more about your saying?
>
> Regards
> Ronghai


sfepy: nls: iter: 0, residual: 7.497687e+01 (rel: 1.000000e+00)
/usr/local/lib/python/dist-packages/scikits/umfpack/umfpack.py:720:
UmfpackWarning: (almost) singular matrix! (estimated cond. number: 5.91e+15)
   warnings.warn(msg, UmfpackWarning)
sfepy:   rezidual:    0.00 [s]
sfepy:      solve:    0.01 [s]
sfepy:     matrix:    0.00 [s]
sfepy: nls: iter: 1, residual: 4.099472e-14 (rel: 5.467649e-16)

Yes, it is solved, because the direct solver somehow deals with the singular
matrix (see the warning above). But without any boundary conditions the matrix
coming from the linear elasticity is singular - in 2D, it has three zero
eigenvalues for the two rigid body displacements and one rotation. These need
to be fixed to have a unique solution. So, for example, fix the bottom left
corner vertex completely ({'u.all' : 0.0}), and fix the 'y' displacement
component of the bottom right vertex ({'u.1' : 0.0}). See [1], where two ways
(a constant and by a function) of specifying the Dirichlet boundary conditions
(EssentialBC) are shown for the interactive use. Let us know if you need more help.

r.
[1]
http://sfepy.org/doc-devel/examples/standalone/interactive/linear_elasticity.html