I'm getting some very strange behaviour in my "interactive" solution.  For the time being I'm just solving a linear Poisson equation on a square region with a constant source term everywhere and Dirichelet BCs on the left and right boundaries.  I'm set up and solved the same problem using a problem description file with the simple.py routine (see myPoisson_Soln.pdf) yet the interactive form has massive oscillations (see myPoissonInteractive_solution.
pdf) in the solution.  The error seems to be in that the matrix for the interactive case is too large (1599x1599 elements for a domain of 21x21 cells, myPoisson.py gives 399x399 elements for the same number of cells in the domain) and so the inversion is singular, yet I can't find where the error in my code could be.  Would someone please mind pointing it out to me?

Thanks.
D

Output of ./simple.py myPoisson.py:
sfepy: left over: ['verbose', '__builtins__', 'n_step', 'dims', 'shape', '__file__', '__name__', 't1', 'center', 'UserMeshIO', 'gen_block_mesh', 't0', '__package__', 'output_dir', '_filename', 'np', 'output', '__doc__', 'mesh_hook']
sfepy
: reading mesh [user] (function:mesh_hook)...
sfepy
: ...done in 0.00 s
sfepy
: creating regions...
sfepy
:     Right
sfepy
:     Top
sfepy
:     Bottom
sfepy
:     Omega
sfepy
:     Left
sfepy
: ...done in 0.00 s
sfepy
: equation "Temperature":
sfepy
:
           dw_laplace
.i.Omega(cond.val, s, T)
         
- dw_surface_integrate.2.Top(insulated.val, s)
         
- dw_volume_lvf.2.Omega(G.val, s)

sfepy
: using solvers:
                ts
: ts
               nls
: newton
                ls
: ls
sfepy
: updating variables...
sfepy
: ...done
sfepy
: setting up dof connectivities...
sfepy
: ...done in 0.00 s
sfepy
: matrix shape: (399, 399)
sfepy
: assembling matrix graph...
sfepy
: ...done in 0.00 s
sfepy
: matrix structural nonzeros: 3355 (2.11e-02% fill)
sfepy
: ====== time 0.000000e+00 (step 1 of 2) =====
sfepy
: updating materials...
sfepy
:     G
sfepy
:     cond
sfepy
:     insulated
sfepy
: ...done in 0.00 s
sfepy
: nls: iter: 0, residual: 2.530862e+01 (rel: 1.000000e+00)
sfepy
:   rezidual:    0.00 [s]
sfepy
:      solve:    0.00 [s]
sfepy
:     matrix:    0.00 [s]
sfepy
: nls: iter: 1, residual: 6.515061e-14 (rel: 2.574246e-15)
sfepy
: ====== time 1.000000e-01 (step 2 of 2) =====
sfepy
: updating variables...
sfepy
: ...done
sfepy
: updating materials...
sfepy
:     G
sfepy
:     cond
sfepy
:     insulated
sfepy
: ...done in 0.00 s
sfepy
: nls: iter: 0, residual: 6.515061e-14 (rel: 1.000000e+00)

Output of python myPoissonInteractive.py:
Enter code here...sfepy: saving regions as groups...
sfepy
:   Omega
sfepy
:   Left
sfepy
:   Right
sfepy
:   Bottom
sfepy
:   Top
sfepy
: ...done
sfepy
: updating variables...
sfepy
: ...done
sfepy
: setting up dof connectivities...
sfepy
: ...done in 0.00 s
sfepy
: matrix shape: (1599, 1599)
sfepy
: assembling matrix graph...
sfepy
: ...done in 0.00 s
sfepy
: matrix structural nonzeros: 24311 (9.51e-03% fill)
sfepy
: updating materials...
sfepy
:     cond
sfepy
:     insulated
sfepy
:     G
sfepy
: ...done in 0.00 s
sfepy
: nls: iter: 0, residual: 6.169742e+01 (rel: 1.000000e+00)
sfepy
:   rezidual:    0.00 [s]
sfepy
:      solve:    0.01 [s]
sfepy
:     matrix:    0.00 [s]
sfepy
: warning: linear system solution precision is lower
sfepy
: then the value set in solver options! (err = 2.856021e+01 < 1.000000e-10)
sfepy
: nls: iter: 1, residual: 2.946994e+01 (rel: 4.776527e-01)
IndexedStruct
  condition
:
   
1
  err
:
   
29.4699421316
  err0
:
   
61.6974210411
  n_iter
:
   
1
  time_stats
:
    dict
with keys: ['rezidual', 'solve', 'matrix']




On Friday, 16 December 2016 10:54:39 UTC+1, Robert Cimrman wrote:
Hi David,

The main problem is having 'vector' in the field creation (line 119) - this
means, that your variables are vectors, and not scalars, as required by the
Laplace term. Try replacing that with 'scalar' or 1.

As for materials, you can use simply:

coef = Material('coef', val=2.0)

etc.

r.

On 12/16/2016 10:46 AM, David Jessop wrote:
> Here's the complete problem description file.
>
> On Tuesday, 6 December 2016 15:26:44 UTC+1, David Jessop wrote:
>>
>>
>> <https://lh3.googleusercontent.com/-nqPKBodss14/WEbGCb-rgXI/AAAAAAAAA70/sQ8dedCI87g1P1S984VDxaC8wnJ26NsLwCLcB/s1600/non-linear_Poisson.png>
>> Hi,
>>
>> I have a non-linear form of the Poisson equation where the "diffusion"
>> coefficient depends on the derivatives of the state variable (see image).
>> There are no Dirichlet-type boundary conditions (i.e. no ebcs) so I really
>> need to specify a good initial guess at the solution.
>>
>> My questions are:
>> 1) How are initial guesses passed to the solver?  From the documentation,
>> it doesn't seem like initial guesses can be directly passed to the
>> nls.newton solver as an argument (though they are an argument of the
>> subroutine __call__).
>> 2) In what form should the initial guess vector/array be written?  Can a
>> function be passed?  Do numerical values have be be specified at the mesh
>> nodes?
>>
>> Thanks for any help on this matter.
>>
>> D
>>
>