I wrote my codes in sfepy-2014.4, but recently I upgraded to sfepy-2015.1. If every other things are the same except for solvers, I got some cases work and some do not. Do you have any idea about the reason? (1) If I keep the problem as nonlinear, it works(original code) ls = ScipyDirect({'presolve':True}) nls_status = IndexedStruct() nls = Newton({'problem' : 'nonlinear'}, lin_solver=ls, status=nls_status) pb = Problem('elasticity', equations=eqs, nls=nls, ls=ls) ############################################################# sfepy: nls: iter: 0, residual: 2.785678e+02 (rel: 1.000000e+00) sfepy: rezidual: 0.46 [s] sfepy: solve: 57.95 [s] sfepy: matrix: 0.67 [s] sfepy: nls: iter: 1, residual: 2.625559e-11 (rel: 9.425207e-14)
(2) if I change to linear, it does not work. ls = ScipyDirect({'presolve':True}) nls_status = IndexedStruct() nls = Newton({'is_linear' : True}, lin_solver=ls, status=nls_status) pb = Problem('elasticity', equations=eqs, nls=nls, ls=ls) ################################################# sfepy: nls: iter: 0, residual: 2.785678e+02 (rel: 1.000000e+00) virtual201503/local/lib/python2.7/site-packages/scipy/sparse/linalg/dsolve/linsolve.py:145: MatrixRankWarning: Matrix is exactly singular warn("Matrix is exactly singular", MatrixRankWarning) sfepy: rezidual: 0.47 [s] sfepy: solve: 19.58 [s] sfepy: matrix: 0.00 [s] sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-01) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-02) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-03) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-04) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-05) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-06) sfepy: linesearch: iter 1, (nan < 2.78565e+02) (new ls: 1.000000e-07) sfepy: linesearch failed, continuing anyway sfepy: nls: iter: 1, residual: nan (rel: nan)
(3)
ls = PETScKrylovSolver({'method' : 'cg', 'precond' : 'icc', 'eps_r' : 1e-8, 'i_max' : 100}) nls_status = IndexedStruct() nls = Newton({'is_linear' : True}, lin_solver=ls, status=nls_status) pb = Problem('elasticity', equations=eqs, nls=nls, ls=ls) ####################################################### sfepy: nls: iter: 0, residual: 2.785678e+02 (rel: 1.000000e+00) sfepy: cg(icc) convergence: -10 (DIVERGED_INDEFINITE_MAT) sfepy: rezidual: 0.47 [s] sfepy: solve: 0.26 [s] sfepy: matrix: 0.00 [s] sfepy: warning: linear system solution precision is lower sfepy: then the value set in solver options! (err = 2.785678e+02 < 1.000000e-10) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-01) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-02) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-03) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-04) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-05) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-06) sfepy: linesearch: iter 1, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-07) sfepy: linesearch failed, continuing anyway sfepy: nls: iter: 1, residual: 2.785678e+02 (rel: 1.000000e+00) sfepy: cg(icc) convergence: -10 (DIVERGED_INDEFINITE_MAT) sfepy: rezidual: 0.29 [s] sfepy: solve: 0.26 [s] sfepy: matrix: 0.00 [s] sfepy: warning: linear system solution precision is lower sfepy: then the value set in solver options! (err = 2.785678e+02 < 1.000000e-10) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-01) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-02) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-03) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-04) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-05) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-06) sfepy: linesearch: iter 2, (2.78568e+02 < 2.78565e+02) (new ls: 1.000000e-07) sfepy: linesearch failed, continuing anyway sfepy: nls: iter: 2, residual: 2.785678e+02 (rel: 1.000000e+00) sfepy: equation "tmp":
在 2015年3月12日星期四 UTC+1下午2:50:25，Robert Cimrman写道：
On 03/11/2015 08:02 PM, Ronghai Wu wrote:
Yes, it is linear elasticity problem. If using linear Newton solver, do you exactly mean ? ls = ScipyDirect({}) nls_status = IndexedStruct() nls = Newton({'problem' : 'linear'}, lin_solver=ls, status=nls_status) pb = Problem('elasticity', equations=eqs, nls=nls, ls=ls)
Yes, sorry, I forgot your script is interactive. In SfePy prior to 2015.1 the option is {'problem' : 'linear'}. In 2015.1 it reads {'is_linear' : True}, as 'problem' was ambiguous.
and if using PETScKrylovSolver, do you exactly mean ? ls = PETScKrylovSolver({})
Essentially, yes. The above uses the default solver, preconditioner, tolerances and other options. You may need to experiment with those. Try, for example:
ls = PETScKrylovSolver({ 'method' : 'cg', 'precond' : 'icc' 'eps_r' : 1e-8, 'i_max' : 100, })
and use this as in your snippet above instead of ScipyDirect.
Another question in the ''Interactive Example: Linear Elasticity'' : In [23]: bc_fun = Function(’shift_u_fun’, shift_u_fun, extra_args={’shift’ : 0.01}) In [24]: shift_u = EssentialBC(’shift_u’, gamma2, {’u.0’ : bc_fun}) is this shift a necessity to get correct result? I mean without this, if the result will be ridiculous or basically right? Thanks
Here, the shift is the amplitude of the imposed Dirichlet conditions on gamma2. It serves as an example of how to pass extra arguments to the ebc functions. If you set shift to zero, it would be equivalent to set {'u.0' : 0.0} on the gamma2 boundary. Did I answer your question?
r.