This page was generated from unit-1.1-poisson/poisson.ipynb.

1.1 First NGSolve example

Let us solve the Poisson problem of finding \(u\) satisfying

\[\begin{split}\begin{aligned} -\Delta u & = f && \text { in the unit square}, \\ u & = 0 && \text{ on the bottom and right parts of the boundary}, \\ \frac{\partial u }{\partial n } & = 0 && \text{ on the remaining boundary parts}. \end{aligned}\end{split}\]

Quick steps to solution:

1. Import NGSolve and Netgen Python modules:

[1]:
from ngsolve import *
from ngsolve.webgui import Draw

2. Generate an unstructured mesh

[2]:
mesh = Mesh(unit_square.GenerateMesh(maxh=0.2))
mesh.nv, mesh.ne   # number of vertices & elements
[2]:
(39, 56)

Here we prescribed a maximal mesh-size of 0.2 using the maxh flag.

[3]:
Draw(mesh);

3. Declare a finite element space:

[4]:
fes = H1(mesh, order=2, dirichlet="bottom|right")
fes.ndof  # number of unknowns in this space
[4]:
133

Python’s help system displays further documentation.

[5]:
help(fes)
Help on H1 in module ngsolve.comp object:

class H1(FESpace)
 |  An H1-conforming finite element space.
 |
 |  The H1 finite element space consists of continuous and
 |  element-wise polynomial functions. It uses a hierarchical (=modal)
 |  basis built from integrated Legendre polynomials on tensor-product elements,
 |  and Jaboci polynomials on simplicial elements.
 |
 |  Boundary values are well defined. The function can be used directly on the
 |  boundary, using the trace operator is optional.
 |
 |  The H1 space supports variable order, which can be set individually for edges,
 |  faces and cells.
 |
 |  Internal degrees of freedom are declared as local dofs and are eliminated
 |  if static condensation is on.
 |
 |  The wirebasket consists of all vertex dofs. Optionally, one can include the
 |  first (the quadratic bubble) edge basis function, or all edge basis functions
 |  into the wirebasket.
 |
 |  Keyword arguments can be:
 |
 |  order: int = 1
 |    order of finite element space
 |  complex: bool = False
 |    Set if FESpace should be complex
 |  dirichlet: regexpr
 |    Regular expression string defining the dirichlet boundary.
 |    More than one boundary can be combined by the | operator,
 |    i.e.: dirichlet = 'top|right'
 |  dirichlet_bbnd: regexpr
 |    Regular expression string defining the dirichlet bboundary,
 |    i.e. points in 2D and edges in 3D.
 |    More than one boundary can be combined by the | operator,
 |    i.e.: dirichlet_bbnd = 'top|right'
 |  dirichlet_bbbnd: regexpr
 |    Regular expression string defining the dirichlet bbboundary,
 |    i.e. points in 3D.
 |    More than one boundary can be combined by the | operator,
 |    i.e.: dirichlet_bbbnd = 'top|right'
 |  definedon: Region or regexpr
 |    FESpace is only defined on specific Region, created with mesh.Materials('regexpr')
 |    or mesh.Boundaries('regexpr'). If given a regexpr, the region is assumed to be
 |    mesh.Materials('regexpr').
 |  dim: int = 1
 |    Create multi dimensional FESpace (i.e. [H1]^3)
 |  dgjumps: bool = False
 |    Enable discontinuous space for DG methods, this flag is needed for DG methods,
 |    since the dofs have a different coupling then and this changes the sparsity
 |    pattern of matrices.
 |  autoupdate: bool = False
 |    Automatically update on a change to the mesh.
 |  low_order_space: bool = True
 |    Generate a lowest order space together with the high-order space,
 |    needed for some preconditioners.
 |  hoprolongation: bool = False
 |    Create high order prolongation operators,
 |    only available for H1 and L2 on simplicial meshes
 |  order_policy: ORDER_POLICY = ORDER_POLICY.OLDSTYLE
 |    CONSTANT .. use the same fixed order for all elements,
 |    NODAL ..... use the same order for nodes of same shape,
 |    VARIABLE ... use an individual order for each edge, face and cell,
 |    OLDSTYLE .. as it used to be for the last decade
 |  print: bool = False
 |    (historic) print some output into file set by 'SetTestoutFile'
 |  wb_withedges: bool = true(3D) / false(2D)
 |    use lowest-order edge dofs for BDDC wirebasket
 |  wb_fulledges: bool = false
 |    use all edge dofs for BDDC wirebasket
 |  hoprolongation: bool = false
 |    (experimental, only trigs) creates high order prolongation,
 |    and switches off low-order space
 |
 |  Method resolution order:
 |      H1
 |      FESpace
 |      NGS_Object
 |      pybind11_builtins.pybind11_object
 |      builtins.object
 |
 |  Methods defined here:
 |
 |  __getstate__(...)
 |      __getstate__(self: ngsolve.comp.FESpace) -> tuple
 |
 |  __init__(...)
 |      __init__(self: ngsolve.comp.H1, mesh: ngsolve.comp.Mesh, **kwargs) -> None
 |
 |  __setstate__(...)
 |      __setstate__(self: ngsolve.comp.H1, arg0: tuple) -> None
 |
 |  ----------------------------------------------------------------------
 |  Static methods defined here:
 |
 |  __flags_doc__(...) from builtins.PyCapsule
 |      __flags_doc__() -> dict
 |
 |  ----------------------------------------------------------------------
 |  Data descriptors defined here:
 |
 |  __dict__
 |
 |  ----------------------------------------------------------------------
 |  Methods inherited from FESpace:
 |
 |  ApplyM(...)
 |      ApplyM(self: ngsolve.comp.FESpace, vec: ngsolve.la.BaseVector, rho: ngsolve.fem.CoefficientFunction = None, definedon: ngsolve.comp.Region = None) -> None
 |
 |      Apply mass-matrix. Available only for L2-like spaces
 |
 |  ConvertL2Operator(...)
 |      ConvertL2Operator(self: ngsolve.comp.FESpace, l2space: ngsolve.comp.FESpace) -> BaseMatrix
 |
 |  CouplingType(...)
 |      CouplingType(self: ngsolve.comp.FESpace, dofnr: int) -> ngsolve.comp.COUPLING_TYPE
 |
 |
 |               Get coupling type of a degree of freedom.
 |
 |      Parameters:
 |
 |      dofnr : int
 |        input dof number
 |
 |  CreateDirectSolverCluster(...)
 |      CreateDirectSolverCluster(self: ngsolve.comp.FESpace, **kwargs) -> list
 |
 |  CreateSmoothingBlocks(...)
 |      CreateSmoothingBlocks(self: ngsolve.comp.FESpace, **kwargs) -> pyngcore.pyngcore.Table_I
 |
 |
 |      Create table of smoothing blocks for block-Jacobi/block-Gauss-Seidel preconditioners.
 |
 |      Every table entry describes the set of dofs belonging a Jacobi/Gauss-Seidel block.
 |
 |      Paramters:
 |
 |      blocktype: string | [ string ] | int
 |          describes blocktype.
 |          string form ["vertex", "edge", "face",  "facet", "vertexedge", ....]
 |          or list of strings for combining multiple blocktypes
 |          int is for backward compatibility with old style blocktypes
 |
 |      condense: bool = False
 |          boolexclude dofs eliminated by static condensation
 |
 |  Elements(...)
 |      Elements(*args, **kwargs)
 |      Overloaded function.
 |
 |      1. Elements(self: ngsolve.comp.FESpace, VOL_or_BND: ngsolve.comp.VorB = <VorB.VOL: 0>) -> ngsolve.comp.FESpaceElementRange
 |
 |
 |      Returns an iterable range of elements.
 |
 |      Parameters:
 |
 |      VOL_or_BND : ngsolve.comp.VorB
 |        input VOL, BND, BBND,...
 |
 |
 |
 |      2. Elements(self: ngsolve.comp.FESpace, arg0: ngsolve.comp.Region) -> Iterator[ngsolve.comp.FESpaceElement]
 |
 |  FinalizeUpdate(...)
 |      FinalizeUpdate(self: ngsolve.comp.FESpace) -> None
 |
 |      finalize update
 |
 |  FreeDofs(...)
 |      FreeDofs(self: ngsolve.comp.FESpace, coupling: bool = False) -> pyngcore.pyngcore.BitArray
 |
 |
 |
 |      Return BitArray of free (non-Dirichlet) dofs\n
 |      coupling=False ... all free dofs including local dofs\n
 |      coupling=True .... only element-boundary free dofs
 |
 |      Parameters:
 |
 |      coupling : bool
 |        input coupling
 |
 |  GetDofNrs(...)
 |      GetDofNrs(*args, **kwargs)
 |      Overloaded function.
 |
 |      1. GetDofNrs(self: ngsolve.comp.FESpace, ei: ngsolve.comp.ElementId) -> tuple
 |
 |
 |
 |      Parameters:
 |
 |      ei : ngsolve.comp.ElementId
 |        input element id
 |
 |
 |
 |      2. GetDofNrs(self: ngsolve.comp.FESpace, ni: ngsolve.comp.NodeId) -> tuple
 |
 |
 |
 |      Parameters:
 |
 |      ni : ngsolve.comp.NodeId
 |        input node id
 |
 |  GetDofs(...)
 |      GetDofs(self: ngsolve.comp.FESpace, region: ngsolve.comp.Region) -> pyngcore.pyngcore.BitArray
 |
 |
 |      Returns all degrees of freedom in given region.
 |
 |      Parameters:
 |
 |      region : ngsolve.comp.Region
 |        input region
 |
 |  GetFE(...)
 |      GetFE(self: ngsolve.comp.FESpace, ei: ngsolve.comp.ElementId) -> object
 |
 |
 |      Get the finite element to corresponding element id.
 |
 |      Parameters:
 |
 |      ei : ngsolve.comp.ElementId
 |         input element id
 |
 |  GetOrder(...)
 |      GetOrder(self: ngsolve.comp.FESpace, nodeid: ngsolve.comp.NodeId) -> int
 |
 |      return order of node.
 |      by now, only isotropic order is supported here
 |
 |  GetTrace(...)
 |      GetTrace(self: ngsolve.comp.FESpace, arg0: ngsolve.comp.FESpace, arg1: ngsolve.la.BaseVector, arg2: ngsolve.la.BaseVector, arg3: bool) -> None
 |
 |  GetTraceTrans(...)
 |      GetTraceTrans(self: ngsolve.comp.FESpace, arg0: ngsolve.comp.FESpace, arg1: ngsolve.la.BaseVector, arg2: ngsolve.la.BaseVector, arg3: bool) -> None
 |
 |  HideAllDofs(...)
 |      HideAllDofs(self: ngsolve.comp.FESpace, component: object = <ngsolve.ngstd.DummyArgument>) -> None
 |
 |      set all visible coupling types to HIDDEN_DOFs (will be overwritten by any Update())
 |
 |  InvM(...)
 |      InvM(self: ngsolve.comp.FESpace, rho: ngsolve.fem.CoefficientFunction = None) -> BaseMatrix
 |
 |  Mass(...)
 |      Mass(self: ngsolve.comp.FESpace, rho: ngsolve.fem.CoefficientFunction = None, definedon: Optional[ngsolve.comp.Region] = None) -> BaseMatrix
 |
 |  ParallelDofs(...)
 |      ParallelDofs(self: ngsolve.comp.FESpace) -> ngsolve.la.ParallelDofs
 |
 |      Return dof-identification for MPI-distributed meshes
 |
 |  Prolongation(...)
 |      Prolongation(self: ngsolve.comp.FESpace) -> ngmg::Prolongation
 |
 |      Return prolongation operator for use in multi-grid
 |
 |  Range(...)
 |      Range(self: ngsolve.comp.FESpace, arg0: int) -> ngsolve.la.DofRange
 |
 |      deprecated, will be only available for ProductSpace
 |
 |  SetCouplingType(...)
 |      SetCouplingType(*args, **kwargs)
 |      Overloaded function.
 |
 |      1. SetCouplingType(self: ngsolve.comp.FESpace, dofnr: int, coupling_type: ngsolve.comp.COUPLING_TYPE) -> None
 |
 |
 |               Set coupling type of a degree of freedom.
 |
 |      Parameters:
 |
 |      dofnr : int
 |        input dof number
 |
 |      coupling_type : ngsolve.comp.COUPLING_TYPE
 |        input coupling type
 |
 |
 |
 |      2. SetCouplingType(self: ngsolve.comp.FESpace, dofnrs: ngsolve.ngstd.IntRange, coupling_type: ngsolve.comp.COUPLING_TYPE) -> None
 |
 |
 |               Set coupling type for interval of dofs.
 |
 |      Parameters:
 |
 |      dofnrs : Range
 |        range of dofs
 |
 |      coupling_type : ngsolve.comp.COUPLING_TYPE
 |        input coupling type
 |
 |  SetDefinedOn(...)
 |      SetDefinedOn(self: ngsolve.comp.FESpace, region: ngsolve.comp.Region) -> None
 |
 |
 |      Set the regions on which the FESpace is defined.
 |
 |      Parameters:
 |
 |      region : ngsolve.comp.Region
 |        input region
 |
 |  SetOrder(...)
 |      SetOrder(*args, **kwargs)
 |      Overloaded function.
 |
 |      1. SetOrder(self: ngsolve.comp.FESpace, element_type: ngsolve.fem.ET, order: int) -> None
 |
 |
 |
 |      Parameters:
 |
 |      element_type : ngsolve.fem.ET
 |        input element type
 |
 |      order : object
 |        input polynomial order
 |
 |
 |      2. SetOrder(self: ngsolve.comp.FESpace, nodeid: ngsolve.comp.NodeId, order: int) -> None
 |
 |
 |
 |      Parameters:
 |
 |      nodeid : ngsolve.comp.NodeId
 |        input node id
 |
 |      order : int
 |        input polynomial order
 |
 |  SolveM(...)
 |      SolveM(self: ngsolve.comp.FESpace, vec: ngsolve.la.BaseVector, rho: ngsolve.fem.CoefficientFunction = None, definedon: ngsolve.comp.Region = None) -> None
 |
 |
 |               Solve with the mass-matrix. Available only for L2-like spaces.
 |
 |      Parameters:
 |
 |      vec : ngsolve.la.BaseVector
 |        input right hand side vector
 |
 |      rho : ngsolve.fem.CoefficientFunction
 |        input CF
 |
 |  TestFunction(...)
 |      TestFunction(self: ngsolve.comp.FESpace) -> object
 |
 |      Return a proxy to be used as a testfunction for :any:`Symbolic Integrators<symbolic-integrators>`
 |
 |  TnT(...)
 |      TnT(self: ngsolve.comp.FESpace) -> tuple[object, object]
 |
 |      Return a tuple of trial and testfunction
 |
 |  TraceOperator(...)
 |      TraceOperator(self: ngsolve.comp.FESpace, tracespace: ngsolve.comp.FESpace, average: bool) -> BaseMatrix
 |
 |  TrialFunction(...)
 |      TrialFunction(self: ngsolve.comp.FESpace) -> object
 |
 |      Return a proxy to be used as a trialfunction in :any:`Symbolic Integrators<symbolic-integrators>`
 |
 |  Update(...)
 |      Update(self: ngsolve.comp.FESpace) -> None
 |
 |      update space after mesh-refinement
 |
 |  UpdateDofTables(...)
 |      UpdateDofTables(self: ngsolve.comp.FESpace) -> None
 |
 |      update dof-tables after changing polynomial order distribution
 |
 |  __eq__(...)
 |      __eq__(self: ngsolve.comp.FESpace, space: ngsolve.comp.FESpace) -> bool
 |
 |  __mul__(...)
 |      __mul__(self: ngsolve.comp.FESpace, arg0: ngsolve.comp.FESpace) -> ngcomp::CompoundFESpace
 |
 |  __pow__(...)
 |      __pow__(self: ngsolve.comp.FESpace, arg0: int) -> ngcomp::CompoundFESpaceAllSame
 |
 |  __str__(...)
 |      __str__(self: ngsolve.comp.FESpace) -> str
 |
 |  __timing__(...)
 |      __timing__(self: ngsolve.comp.FESpace) -> object
 |
 |  ----------------------------------------------------------------------
 |  Static methods inherited from FESpace:
 |
 |  __special_treated_flags__(...) from builtins.PyCapsule
 |      __special_treated_flags__() -> dict
 |
 |  ----------------------------------------------------------------------
 |  Readonly properties inherited from FESpace:
 |
 |  autoupdate
 |
 |  components
 |      deprecated, will be only available for ProductSpace
 |
 |  couplingtype
 |
 |  dim
 |      multi-dim of FESpace
 |
 |  globalorder
 |      query global order of space
 |
 |  is_complex
 |
 |  loembedding
 |
 |  lospace
 |
 |  mesh
 |      mesh on which the FESpace is created
 |
 |  ndof
 |      number of degrees of freedom
 |
 |  ndofglobal
 |      global number of dofs on MPI-distributed mesh
 |
 |  type
 |      type of finite element space
 |
 |  ----------------------------------------------------------------------
 |  Data and other attributes inherited from FESpace:
 |
 |  __hash__ = None
 |
 |  ----------------------------------------------------------------------
 |  Readonly properties inherited from NGS_Object:
 |
 |  __memory__
 |
 |  flags
 |
 |  ----------------------------------------------------------------------
 |  Data descriptors inherited from NGS_Object:
 |
 |  name
 |
 |  ----------------------------------------------------------------------
 |  Static methods inherited from pybind11_builtins.pybind11_object:
 |
 |  __new__(*args, **kwargs) from pybind11_builtins.pybind11_type
 |      Create and return a new object.  See help(type) for accurate signature.

4. Declare test function, trial function, and grid function

  • Test and trial function are symbolic objects - called ProxyFunctions - that help you construct bilinear forms (and have no space to hold solutions).

  • GridFunctions, on the other hand, represent functions in the finite element space and contains memory to hold coefficient vectors.

[6]:
u = fes.TrialFunction()  # symbolic object
v = fes.TestFunction()   # symbolic object
gfu = GridFunction(fes)  # solution

Alternately, you can get both the trial and test variables at once:

[7]:
u, v = fes.TnT()

5. Define and assemble linear and bilinear forms:

[8]:
a = BilinearForm(fes)
a += grad(u)*grad(v)*dx
a.Assemble()

f = LinearForm(fes)
f += x*v*dx
f.Assemble();

Alternately, we can do one-liners:

[9]:
a = BilinearForm(grad(u)*grad(v)*dx).Assemble()
f = LinearForm(x*v*dx).Assemble()

You can examine the linear system in more detail:

[10]:
print(f.vec)
 0.000333333
 0.00873506
 0.00633333
 0.000510926
 0.00430724
 0.00890537
 0.0106804
 0.0112774
 0.0148637
 0.0166495
 0.0150628
 0.0185037
 0.0216167
 0.0109997
 0.00797031
 0.00292185
 0.000628106
 0.0010144
 0.00172207
 0.00185741
 0.0142294
 0.01997
 0.0225059
 0.020015
 0.0311335
 0.0259928
 0.028898
 0.0191106
 0.0227498
 0.00992941
 0.00267856
 0.00497134
 0.0080662
 0.0203373
 0.0164825
 0.0204557
 0.0134827
 0.0158788
 0.018219
 -6.66667e-05
 -3.33333e-05
 -0.000516153
 -0.000566493
 -0.0010733
 -0.0008
 -0.000766667
 -5.90186e-05
 -1.43029e-05
 -7.99565e-05
 -0.000293781
 -0.000218551
 -0.000469969
 -0.000398925
 -0.000689288
 -0.000822998
 -0.000483647
 -0.000847511
 -0.00093399
 -0.000895046
 -0.000915843
 -0.00074941
 -0.00109658
 -0.00126105
 -0.000643908
 -0.00140539
 -0.00131192
 -0.000648404
 -0.00121351
 -0.00121755
 -0.00162374
 -0.00149151
 -0.00048888
 -0.00145143
 -0.00110173
 -0.000452064
 -0.000914624
 -0.00088675
 -0.000208119
 -0.0007941
 -0.000543629
 -0.000302341
 -0.000167379
 -2.0832e-05
 -6.88223e-05
 -8.44745e-05
 -2.94073e-05
 -0.000116264
 -0.000137817
 -4.79891e-05
 -0.000243969
 -0.000195256
 -0.000257349
 -0.000787648
 -0.00044492
 -0.000660183
 -0.000883158
 -0.000835026
 -0.000794567
 -0.000892677
 -0.0010382
 -0.000974429
 -0.00102675
 -0.00114594
 -0.000985686
 -0.000923147
 -0.00108379
 -0.0008854
 -0.00086299
 -0.00104694
 -0.000933801
 -0.000848623
 -0.000866758
 -0.000666847
 -0.000884373
 -0.000729732
 -0.000868298
 -0.000206371
 -0.000281267
 -0.000500567
 -0.000145299
 -0.000258405
 -0.000362166
 -0.000550055
 -0.000437637
 -0.00077891
 -0.000777094
 -0.000731175
 -0.000607743
 -0.000717949
 -0.000729496
 -0.000805781
 -0.000739587
 -0.000678956


[11]:
print(a.mat)
Row 0:   0: 1   4: -0.5   19: -0.5   39: -0.0833333   40: -0.0833333   50: 0.166667
Row 1:   1: 0.828732   7: -0.195774   8: -0.197948   23: -0.43501   41: -0.0360819   42: -0.0364197   43: -0.0656203   59: 0.0687109   61: 0.069411
Row 2:   2: 1   11: -0.5   12: -0.5   44: -0.0833333   45: -0.0833333   69: 0.166667
Row 3:   3: 0.834885   15: -0.164012   16: -0.157951   30: -0.512921   46: -0.0431063   47: -0.0423806   48: -0.0536605   81: 0.0704417   83: 0.0687058
Row 4:   0: -0.5   4: 1.89671   5: -0.546735   19: -0.180532   20: -0.66944   39: 5.83717e-17   40: 0.0833333   49: -0.0544706   50: -0.140436   51: -0.121211   53: 0.145593   91: 0.0871914
Row 5:   4: -0.546735   5: 1.88547   6: -0.352441   20: -0.196284   21: -0.790007   49: -0.0136061   51: 0.104729   52: -0.0481159   53: -0.174674   54: -0.0778481   56: 0.106856   92: 0.10266
Row 6:   5: -0.352441   6: 1.75375   7: -0.265716   21: -0.432048   22: -0.703548   52: -0.0385388   54: 0.097279   55: -0.0537377   56: -0.12226   57: -0.0777552   58: 0.0980238   95: 0.0969894
Row 7:   1: -0.195774   6: -0.265716   7: 1.822   22: -0.405418   23: -0.955095   41: -0.0839333   43: 0.116562   55: -0.0465664   57: 0.0908524   58: -0.119535   59: -0.0536322   98: 0.0962524
Row 8:   1: -0.197948   8: 1.87437   9: -0.341759   23: -1.05108   24: -0.283578   42: -0.0827377   43: 0.115729   60: -0.0333848   61: -0.0468695   62: -0.149403   64: 0.0903445   101: 0.106321
Row 9:   8: -0.341759   9: 1.74315   10: -0.259596   24: -0.651696   25: -0.490101   60: -0.055818   62: 0.112778   63: -0.0349443   64: -0.103699   65: -0.0960639   67: 0.0782102   102: 0.0995372
Row 10:   9: -0.259596   10: 1.76903   11: -0.259654   25: -0.823099   26: -0.426681   63: -0.0694609   65: 0.112727   66: -0.0367441   67: -0.0776353   68: -0.110998   70: 0.0800198   105: 0.102092
Row 11:   2: -0.5   10: -0.259654   11: 1.96958   12: -0.139204   26: -1.07072   44: 0   45: 0.0833333   66: -0.0669124   68: 0.110188   69: -0.194875   70: -0.0664764   72: 0.134742
Row 12:   2: -0.5   11: -0.139204   12: 1.81482   13: -0.275945   26: -0.20179   27: -0.697885   44: 0.0833333   45: 0   69: -0.115666   70: 0.0555335   71: -0.0418546   72: -0.0976603   73: -0.0472897   75: 0.0878455   108: 0.0757584
Row 13:   12: -0.275945   13: 1.79186   14: -0.426433   27: -0.746323   28: -0.343162   71: -0.0571403   73: 0.103131   74: -0.0399256   75: -0.063259   76: -0.138319   78: 0.110998   110: 0.084515
Row 14:   13: -0.426433   14: 1.77214   15: -0.271016   28: -0.405887   29: -0.668803   74: -0.0369995   76: 0.108072   77: -0.0450239   78: -0.137515   79: -0.0758177   80: 0.0901932   112: 0.0970916
Row 15:   3: -0.164012   14: -0.271016   15: 1.81591   29: -0.463066   30: -0.917816   46: -0.0818566   48: 0.109192   77: -0.054447   79: 0.0996163   80: -0.116282   81: -0.0500661   116: 0.0938434
Row 16:   3: -0.157951   16: 1.82092   17: -0.237555   30: -0.936175   31: -0.489235   47: -0.0848367   48: 0.111162   82: -0.0553671   83: -0.0524972   84: -0.110785   86: 0.0949597   119: 0.0973645
Row 17:   16: -0.237555   17: 1.74436   18: -0.336551   31: -0.666886   32: -0.50337   82: -0.0500457   84: 0.0896382   85: -0.0444322   86: -0.0790554   87: -0.117194   90: 0.100524   120: 0.100565
Row 18:   17: -0.336551   18: 1.8616   19: -0.527472   20: -0.241774   32: -0.755808   85: -0.0442896   87: 0.100381   88: -0.015565   89: -0.16959   90: -0.0808224   91: 0.103477   93: 0.106409
Row 19:   0: -0.5   4: -0.180532   18: -0.527472   19: 1.89097   20: -0.682967   39: 0.0833333   40: 6.245e-17   50: -0.143274   51: 0.0900292   88: -0.0538873   89: 0.141799   91: -0.118001
Row 20:   4: -0.66944   5: -0.196284   18: -0.241774   19: -0.682967   20: 3.47846   21: -0.529072   32: -0.463883   34: -0.695042   49: 0.0680767   50: 0.117043   51: -0.0735467   53: -0.106564   54: 0.0712017   88: 0.0694523   89: -0.100166   90: 0.0710095   91: -0.0726677   92: -0.0767684   93: -0.0829105   94: -0.0671198   97: 0.0937453   122: 0.0892148
Row 21:   5: -0.790007   6: -0.432048   20: -0.529072   21: 3.54163   22: -0.566892   33: -0.74118   34: -0.482434   52: 0.0866548   53: 0.135646   54: -0.0906326   56: -0.0977964   57: 0.0831496   92: -0.135435   94: 0.0879677   95: -0.0896637   96: -0.0733244   97: -0.10342   100: 0.100996   124: 0.0958583
Row 22:   6: -0.703548   7: -0.405418   21: -0.566892   22: 3.51389   23: -0.735409   24: -0.547842   33: -0.554777   55: 0.100304   56: 0.113201   57: -0.0962468   58: -0.109466   59: 0.0767314   95: -0.107352   96: 0.0886334   98: -0.0570951   99: -0.115471   100: -0.100017   101: 0.102932   103: 0.103846
Row 23:   1: -0.43501   7: -0.955095   8: -1.05108   22: -0.735409   23: 3.76027   24: -0.583672   41: 0.120015   42: 0.119157   43: -0.166671   58: 0.130977   59: -0.0918101   61: -0.0896687   62: 0.145692   98: -0.119279   99: 0.11087   101: -0.159283
Row 24:   8: -0.283578   9: -0.651696   22: -0.547842   23: -0.583672   24: 3.42538   25: -0.432896   33: -0.399388   37: -0.526313   60: 0.0892028   61: 0.0671272   62: -0.109067   64: -0.07836   65: 0.0977732   98: 0.0801216   99: -0.0850449   100: 0.0962304   101: -0.0499701   102: -0.0948866   103: -0.100845   104: -0.0527235   107: 0.0692628   125: 0.0711795
Row 25:   9: -0.490101   10: -0.823099   24: -0.432896   25: 3.58075   26: -0.51662   35: -0.398376   37: -0.919662   63: 0.104405   64: 0.0917145   65: -0.114436   67: -0.080167   68: 0.112945   102: -0.130489   104: 0.110924   105: -0.110081   106: -0.110408   107: -0.0512114   109: 0.0832393   129: 0.0935649
Row 26:   10: -0.426681   11: -1.07072   12: -0.20179   25: -0.51662   26: 3.72303   27: -0.801924   35: -0.705295   66: 0.103656   67: 0.0795922   68: -0.112135   69: 0.143875   70: -0.0690769   72: -0.201931   73: 0.0916879   105: -0.0969936   106: 0.103505   108: -0.0562238   109: -0.0841455   111: 0.09819
Row 27:   12: -0.697885   13: -0.746323   26: -0.801924   27: 3.70186   28: -0.78084   35: -0.674884   71: 0.0989949   72: 0.164849   73: -0.147529   75: -0.110283   76: 0.135675   108: -0.120983   109: 0.0897882   110: -0.133204   111: -0.104977   113: 0.127669
Row 28:   13: -0.343162   14: -0.405887   27: -0.78084   28: 3.4438   29: -0.579412   35: -0.387095   36: -0.569649   38: -0.377753   74: 0.076925   75: 0.0856965   76: -0.105428   78: -0.0904765   79: 0.0811993   110: -0.0397802   111: 0.0842238   112: -0.0751098   113: -0.0901223   114: -0.0805661   115: -0.0924837   118: 0.0904793   130: 0.0704143   131: 0.0850283
Row 29:   14: -0.668803   15: -0.463066   28: -0.579412   29: 3.51847   30: -0.682657   31: -0.574344   36: -0.550188   77: 0.0994709   78: 0.116994   79: -0.104998   80: -0.101799   81: 0.0795062   112: -0.120582   114: 0.100156   116: -0.0557497   117: -0.0945602   118: -0.108723   119: 0.0900197   121: 0.100265
Row 30:   3: -0.512921   15: -0.917816   16: -0.936175   29: -0.682657   30: 3.70011   31: -0.65054   46: 0.124963   47: 0.127217   48: -0.166693   80: 0.127888   81: -0.0998818   83: -0.0945679   84: 0.12338   116: -0.127349   117: 0.113237   119: -0.128193
Row 31:   16: -0.489235   17: -0.666886   29: -0.574344   30: -0.65054   31: 3.50264   32: -0.54324   36: -0.578395   82: 0.105413   83: 0.0783594   84: -0.102233   86: -0.100444   87: 0.106179   116: 0.0892549   117: -0.102008   118: 0.108477   119: -0.059191   120: -0.111729   121: -0.108168   123: 0.0960903
Row 32:   17: -0.50337   18: -0.755808   20: -0.463883   31: -0.54324   32: 3.52106   34: -0.590839   36: -0.663921   85: 0.0887218   86: 0.0845401   87: -0.0893669   89: 0.127957   90: -0.090711   93: -0.139279   94: 0.0886357   120: -0.0877797   121: 0.0937796   122: -0.0933717   123: -0.0863352   127: 0.103209
Row 33:   21: -0.74118   22: -0.554777   24: -0.399388   33: 3.5711   34: -0.567872   37: -0.917504   38: -0.390377   95: 0.100026   96: -0.096241   97: 0.119745   99: 0.0896463   100: -0.0972098   103: -0.129815   104: 0.106733   124: -0.113287   125: -0.0446607   126: -0.11397   128: 0.0881878   132: 0.0908447
Row 34:   20: -0.695042   21: -0.482434   32: -0.590839   33: -0.567872   34: 3.50498   36: -0.453329   38: -0.715467   92: 0.109543   93: 0.11578   94: -0.109484   96: 0.080932   97: -0.11007   122: -0.100779   123: 0.0834713   124: -0.0852611   126: 0.0989745   127: -0.103379   128: -0.0751922   131: 0.0954622
Row 35:   25: -0.398376   26: -0.705295   27: -0.674884   28: -0.387095   35: 3.55956   37: -0.661211   38: -0.732703   105: 0.104983   106: -0.118983   107: 0.080396   108: 0.101448   109: -0.088882   110: 0.0884695   111: -0.077437   113: -0.147633   115: 0.123679   129: -0.0644794   130: -0.0958469   132: 0.0942853
Row 36:   28: -0.569649   29: -0.550188   31: -0.578395   32: -0.663921   34: -0.453329   36: 3.48945   38: -0.673969   112: 0.0986001   114: -0.1153   115: 0.111641   117: 0.0833311   118: -0.0902331   120: 0.0989443   121: -0.0858763   122: 0.104936   123: -0.0932263   127: -0.113503   128: 0.0841228   131: -0.083436
Row 37:   24: -0.526313   25: -0.919662   33: -0.917504   35: -0.661211   37: 3.69615   38: -0.671464   102: 0.125838   103: 0.126814   104: -0.164933   106: 0.125886   107: -0.0984473   125: -0.106982   126: 0.133086   129: -0.120086   130: 0.104402   132: -0.125577
Row 38:   28: -0.377753   33: -0.390377   34: -0.715467   35: -0.732703   36: -0.673969   37: -0.671464   38: 3.56173   113: 0.110086   114: 0.0957097   115: -0.142836   124: 0.10269   125: 0.0804635   126: -0.118091   127: 0.113673   128: -0.0971183   129: 0.0910005   130: -0.0789689   131: -0.0970544   132: -0.0595534
Row 39:   0: -0.0833333   4: 5.83717e-17   19: 0.0833333   39: 0.0416667   40: 1.73472e-17   50: -0.0208333
Row 40:   0: -0.0833333   4: 0.0833333   19: 6.245e-17   39: 1.73472e-17   40: 0.0416667   50: -0.0208333
Row 41:   1: -0.0360819   7: -0.0839333   23: 0.120015   41: 0.038161   43: -0.0209833   59: -0.00902048
Row 42:   1: -0.0364197   8: -0.0827377   23: 0.119157   42: 0.0380372   43: -0.0206844   61: -0.00910492
Row 43:   1: -0.0656203   7: 0.116562   8: 0.115729   23: -0.166671   41: -0.0209833   42: -0.0206844   43: 0.0761982   59: -0.00815724   61: -0.00824784
Row 44:   2: -0.0833333   11: 0   12: 0.0833333   44: 0.0416667   45: 0   69: -0.0208333
Row 45:   2: -0.0833333   11: 0.0833333   12: 0   44: 0   45: 0.0416667   69: -0.0208333
Row 46:   3: -0.0431063   15: -0.0818566   30: 0.124963   46: 0.0380746   48: -0.0204642   81: -0.0107766
Row 47:   3: -0.0423806   16: -0.0848367   30: 0.127217   47: 0.0383856   48: -0.0212092   83: -0.0105951
Row 48:   3: -0.0536605   15: 0.109192   16: 0.111162   30: -0.166693   46: -0.0204642   47: -0.0212092   48: 0.0764602   81: -0.00683385   83: -0.00658129
Row 49:   4: -0.0544706   5: -0.0136061   20: 0.0680767   49: 0.0397998   51: -0.00340152   53: -0.0136176
Row 50:   0: 0.166667   4: -0.140436   19: -0.143274   20: 0.117043   39: -0.0208333   40: -0.0208333   50: 0.0784497   51: -0.0149851   91: -0.0142757
Row 51:   4: -0.121211   5: 0.104729   19: 0.0900292   20: -0.0735467   49: -0.00340152   50: -0.0149851   51: 0.0765828   53: -0.0227806   91: -0.00752216
Row 52:   5: -0.0481159   6: -0.0385388   21: 0.0866548   52: 0.0363487   54: -0.00963471   56: -0.012029
Row 53:   4: 0.145593   5: -0.174674   20: -0.106564   21: 0.135646   49: -0.0136176   51: -0.0227806   53: 0.0784882   54: -0.0130234   92: -0.020888
Row 54:   5: -0.0778481   6: 0.097279   20: 0.0712017   21: -0.0906326   52: -0.00963471   53: -0.0130234   54: 0.0750371   56: -0.014685   92: -0.00477698
Row 55:   6: -0.0537377   7: -0.0465664   22: 0.100304   55: 0.0361475   57: -0.0116416   58: -0.0134344
Row 56:   5: 0.106856   6: -0.12226   21: -0.0977964   22: 0.113201   52: -0.012029   54: -0.014685   56: 0.0730162   57: -0.0124201   95: -0.0158801
Row 57:   6: -0.0777552   7: 0.0908524   21: 0.0831496   22: -0.0962468   55: -0.0116416   56: -0.0124201   57: 0.072815   58: -0.0110715   95: -0.00836729
Row 58:   6: 0.0980238   7: -0.119535   22: -0.109466   23: 0.130977   55: -0.0134344   57: -0.0110715   58: 0.0741427   59: -0.013932   98: -0.0188123
Row 59:   1: 0.0687109   7: -0.0536322   22: 0.0767314   23: -0.0918101   41: -0.00902048   43: -0.00815724   58: -0.013932   59: 0.0761562   98: -0.0052508
Row 60:   8: -0.0333848   9: -0.055818   24: 0.0892028   60: 0.0365406   62: -0.0139545   64: -0.00834619
Row 61:   1: 0.069411   8: -0.0468695   23: -0.0896687   24: 0.0671272   42: -0.00910492   43: -0.00824784   61: 0.0779298   62: -0.0133123   101: -0.00346954
Row 62:   8: -0.149403   9: 0.112778   23: 0.145692   24: -0.109067   60: -0.0139545   61: -0.0133123   62: 0.0764332   64: -0.0142399   101: -0.0231108
Row 63:   9: -0.0349443   10: -0.0694609   25: 0.104405   63: 0.0369178   65: -0.0173652   67: -0.00873607
Row 64:   8: 0.0903445   9: -0.103699   24: -0.07836   25: 0.0917145   60: -0.00834619   62: -0.0142399   64: 0.0726688   65: -0.0112438   102: -0.0116848
Row 65:   9: -0.0960639   10: 0.112727   24: 0.0977732   25: -0.114436   63: -0.0173652   64: -0.0112438   65: 0.0730459   67: -0.0108165   102: -0.0131995
Row 66:   10: -0.0367441   11: -0.0669124   26: 0.103656   66: 0.0367331   68: -0.0167281   70: -0.00918602
Row 67:   9: 0.0782102   10: -0.0776353   25: -0.080167   26: 0.0795922   63: -0.00873607   65: -0.0108165   67: 0.0737464   68: -0.0113057   105: -0.00859235
Row 68:   10: -0.110998   11: 0.110188   25: 0.112945   26: -0.112135   66: -0.0167281   67: -0.0113057   68: 0.0735617   70: -0.0108189   105: -0.0169306
Row 69:   2: 0.166667   11: -0.194875   12: -0.115666   26: 0.143875   44: -0.0208333   45: -0.0208333   69: 0.0834355   70: -0.0080832   72: -0.0278854
Row 70:   10: 0.0800198   11: -0.0664764   12: 0.0555335   26: -0.0690769   66: -0.00918602   68: -0.0108189   69: -0.0080832   70: 0.0785019   72: -0.00580017
Row 71:   12: -0.0418546   13: -0.0571403   27: 0.0989949   71: 0.0362464   73: -0.0142851   75: -0.0104636
Row 72:   11: 0.134742   12: -0.0976603   26: -0.201931   27: 0.164849   69: -0.0278854   70: -0.00580017   72: 0.0833057   73: -0.0225973   108: -0.0186149
Row 73:   12: -0.0472897   13: 0.103131   26: 0.0916879   27: -0.147529   71: -0.0142851   72: -0.0225973   73: 0.0777833   75: -0.0114977   108: -0.00032471
Row 74:   13: -0.0399256   14: -0.0369995   28: 0.076925   74: 0.0369993   76: -0.00924986   78: -0.00998139
Row 75:   12: 0.0878455   13: -0.063259   27: -0.110283   28: 0.0856965   71: -0.0104636   73: -0.0114977   75: 0.0744823   76: -0.0171071   110: -0.00431704
Row 76:   13: -0.138319   14: 0.108072   27: 0.135675   28: -0.105428   74: -0.00924986   75: -0.0171071   76: 0.0752351   78: -0.017768   110: -0.0168117
Row 77:   14: -0.0450239   15: -0.054447   29: 0.0994709   77: 0.03616   79: -0.0136117   80: -0.011256
Row 78:   13: 0.110998   14: -0.137515   28: -0.0904765   29: 0.116994   74: -0.00998139   76: -0.017768   78: 0.0739099   79: -0.0126377   112: -0.0166108
Row 79:   14: -0.0758177   15: 0.0996163   28: 0.0811993   29: -0.104998   77: -0.0136117   78: -0.0126377   79: 0.0730707   80: -0.0112923   112: -0.00766209
Row 80:   14: 0.0901932   15: -0.116282   29: -0.101799   30: 0.127888   77: -0.011256   79: -0.0112923   80: 0.0738148   81: -0.0141939   116: -0.0177782
Row 81:   3: 0.0704417   15: -0.0500661   29: 0.0795062   30: -0.0998818   46: -0.0107766   48: -0.00683385   80: -0.0141939   81: 0.0757293   116: -0.00568268
Row 82:   16: -0.0553671   17: -0.0500457   31: 0.105413   82: 0.0362513   84: -0.0125114   86: -0.0138418
Row 83:   3: 0.0687058   16: -0.0524972   30: -0.0945679   31: 0.0783594   47: -0.0105951   48: -0.00658129   83: 0.0757736   84: -0.0130468   119: -0.00654301
Row 84:   16: -0.110785   17: 0.0896382   30: 0.12338   31: -0.102233   82: -0.0125114   83: -0.0130468   84: 0.0736393   86: -0.00989814   119: -0.0177981
Row 85:   17: -0.0444322   18: -0.0442896   32: 0.0887218   85: 0.0362034   87: -0.0110724   90: -0.011108
Row 86:   16: 0.0949597   17: -0.0790554   31: -0.100444   32: 0.0845401   82: -0.0138418   84: -0.00989814   86: 0.0726619   87: -0.0112693   120: -0.0098657
Row 87:   17: -0.117194   18: 0.100381   31: 0.106179   32: -0.0893669   85: -0.0110724   86: -0.0112693   87: 0.0726139   90: -0.0140229   120: -0.0152755
Row 88:   18: -0.015565   19: -0.0538873   20: 0.0694523   88: 0.0393411   89: -0.0134718   91: -0.00389124
Row 89:   18: -0.16959   19: 0.141799   20: -0.100166   32: 0.127957   88: -0.0134718   89: 0.077513   90: -0.0115697   91: -0.021978   93: -0.0204196
Row 90:   17: 0.100524   18: -0.0808224   20: 0.0710095   32: -0.090711   85: -0.011108   87: -0.0140229   89: -0.0115697   90: 0.0743753   93: -0.00618266
Row 91:   4: 0.0871914   18: 0.103477   19: -0.118001   20: -0.0726677   50: -0.0142757   51: -0.00752216   88: -0.00389124   89: -0.021978   91: 0.0761241
Row 92:   5: 0.10266   20: -0.0767684   21: -0.135435   34: 0.109543   53: -0.020888   54: -0.00477698   92: 0.0750955   94: -0.0129707   97: -0.0144151
Row 93:   18: 0.106409   20: -0.0829105   32: -0.139279   34: 0.11578   89: -0.0204196   90: -0.00618266   93: 0.0748758   94: -0.0144002   122: -0.014545
Row 94:   20: -0.0671198   21: 0.0879677   32: 0.0886357   34: -0.109484   92: -0.0129707   93: -0.0144002   94: 0.0731109   97: -0.0090212   122: -0.00775875
Row 95:   6: 0.0969894   21: -0.0896637   22: -0.107352   33: 0.100026   56: -0.0158801   57: -0.00836729   95: 0.0728744   96: -0.010958   100: -0.0140486
Row 96:   21: -0.0733244   22: 0.0886334   33: -0.096241   34: 0.080932   95: -0.010958   96: 0.0732738   97: -0.0131023   100: -0.0112004   124: -0.00713071
Row 97:   20: 0.0937453   21: -0.10342   33: 0.119745   34: -0.11007   92: -0.0144151   94: -0.0090212   96: -0.0131023   97: 0.0734739   124: -0.0168339
Row 98:   7: 0.0962524   22: -0.0570951   23: -0.119279   24: 0.0801216   58: -0.0188123   59: -0.0052508   98: 0.0747355   99: -0.0110074   101: -0.00902298
Row 99:   22: -0.115471   23: 0.11087   24: -0.0850449   33: 0.0896463   98: -0.0110074   99: 0.0729558   100: -0.0102538   101: -0.01671   103: -0.0121577
Row 100:   21: 0.100996   22: -0.100017   24: 0.0962304   33: -0.0972098   95: -0.0140486   96: -0.0112004   99: -0.0102538   100: 0.0724223   103: -0.0138038
Row 101:   8: 0.106321   22: 0.102932   23: -0.159283   24: -0.0499701   61: -0.00346954   62: -0.0231108   98: -0.00902298   99: -0.01671   101: 0.076633
Row 102:   9: 0.0995372   24: -0.0948866   25: -0.130489   37: 0.125838   64: -0.0116848   65: -0.0131995   102: 0.0743812   104: -0.0209374   107: -0.0105222
Row 103:   22: 0.103846   24: -0.100845   33: -0.129815   37: 0.126814   99: -0.0121577   100: -0.0138038   103: 0.0743062   104: -0.020296   125: -0.0114075
Row 104:   24: -0.0527235   25: 0.110924   33: 0.106733   37: -0.164933   102: -0.0209374   103: -0.020296   104: 0.0763439   107: -0.00679354   125: -0.00638735
Row 105:   10: 0.102092   25: -0.110081   26: -0.0969936   35: 0.104983   67: -0.00859235   68: -0.0169306   105: 0.0732945   106: -0.0156561   109: -0.0105897
Row 106:   25: -0.110408   26: 0.103505   35: -0.118983   37: 0.125886   105: -0.0156561   106: 0.0739468   107: -0.0140897   109: -0.0102201   129: -0.0173819
Row 107:   24: 0.0692628   25: -0.0512114   35: 0.080396   37: -0.0984473   102: -0.0105222   104: -0.00679354   106: -0.0140897   107: 0.0757339   129: -0.00600932
Row 108:   12: 0.0757584   26: -0.0562238   27: -0.120983   35: 0.101448   72: -0.0186149   73: -0.00032471   108: 0.0777152   109: -0.0116308   111: -0.0137312
Row 109:   25: 0.0832393   26: -0.0841455   27: 0.0897882   35: -0.088882   105: -0.0105897   106: -0.0102201   108: -0.0116308   109: 0.0726442   111: -0.0108162
Row 110:   13: 0.084515   27: -0.133204   28: -0.0397802   35: 0.0884695   75: -0.00431704   76: -0.0168117   110: 0.0757811   111: -0.00562801   113: -0.0164894
Row 111:   26: 0.09819   27: -0.104977   28: 0.0842238   35: -0.077437   108: -0.0137312   109: -0.0108162   110: -0.00562801   111: 0.0737236   113: -0.0154279
Row 112:   14: 0.0970916   28: -0.0751098   29: -0.120582   36: 0.0986001   78: -0.0166108   79: -0.00766209   112: 0.0730651   114: -0.0135347   118: -0.0111154
Row 113:   27: 0.127669   28: -0.0901223   35: -0.147633   38: 0.110086   110: -0.0164894   111: -0.0154279   113: 0.0755677   115: -0.0204188   130: -0.00710264
Row 114:   28: -0.0805661   29: 0.100156   36: -0.1153   38: 0.0957097   112: -0.0135347   114: 0.0727019   115: -0.0152903   118: -0.0115045   131: -0.00863707
Row 115:   28: -0.0924837   35: 0.123679   36: 0.111641   38: -0.142836   113: -0.0204188   114: -0.0152903   115: 0.0745698   130: -0.0105009   131: -0.01262
Row 116:   15: 0.0938434   29: -0.0557497   30: -0.127349   31: 0.0892549   80: -0.0177782   81: -0.00568268   116: 0.0742186   117: -0.014059   119: -0.00825474
Row 117:   29: -0.0945602   30: 0.113237   31: -0.102008   36: 0.0833311   116: -0.014059   117: 0.073073   118: -0.0114429   119: -0.0142502   121: -0.00938986
Row 118:   28: 0.0904793   29: -0.108723   31: 0.108477   36: -0.0902331   112: -0.0111154   114: -0.0115045   117: -0.0114429   118: 0.0726635   121: -0.0156763
Row 119:   16: 0.0973645   29: 0.0900197   30: -0.128193   31: -0.059191   83: -0.00654301   84: -0.0177981   116: -0.00825474   117: -0.0142502   119: 0.0739519
Row 120:   17: 0.100565   31: -0.111729   32: -0.0877797   36: 0.0989443   86: -0.0098657   87: -0.0152755   120: 0.0725123   121: -0.0120792   123: -0.0126569
Row 121:   29: 0.100265   31: -0.108168   32: 0.0937796   36: -0.0858763   117: -0.00938986   118: -0.0156763   120: -0.0120792   121: 0.0726108   123: -0.0113657
Row 122:   20: 0.0892148   32: -0.0933717   34: -0.100779   36: 0.104936   93: -0.014545   94: -0.00775875   122: 0.0731559   123: -0.0106497   127: -0.0155842
Row 123:   31: 0.0960903   32: -0.0863352   34: 0.0834713   36: -0.0932263   120: -0.0126569   121: -0.0113657   122: -0.0106497   123: 0.0725538   127: -0.0102181
Row 124:   21: 0.0958583   33: -0.113287   34: -0.0852611   38: 0.10269   96: -0.00713071   97: -0.0168339   124: 0.0732984   126: -0.0141846   128: -0.0114879
Row 125:   24: 0.0711795   33: -0.0446607   37: -0.106982   38: 0.0804635   103: -0.0114075   104: -0.00638735   125: 0.0761401   126: -0.0153381   132: -0.00477781
Row 126:   33: -0.11397   34: 0.0989745   37: 0.133086   38: -0.118091   124: -0.0141846   125: -0.0153381   126: 0.0742808   128: -0.0105591   132: -0.0179334
Row 127:   32: 0.103209   34: -0.103379   36: -0.113503   38: 0.113673   122: -0.0155842   123: -0.0102181   127: 0.0731092   128: -0.0127917   131: -0.0156265
Row 128:   33: 0.0881878   34: -0.0751922   36: 0.0841228   38: -0.0971183   124: -0.0114879   126: -0.0105591   127: -0.0127917   128: 0.0728888   131: -0.00823901
Row 129:   25: 0.0935649   35: -0.0644794   37: -0.120086   38: 0.0910005   106: -0.0173819   107: -0.00600932   129: 0.0736918   130: -0.0126396   132: -0.0101105
Row 130:   28: 0.0704143   35: -0.0958469   37: 0.104402   38: -0.0789689   113: -0.00710264   115: -0.0105009   129: -0.0126396   130: 0.0742333   132: -0.0134608
Row 131:   28: 0.0850283   34: 0.0954622   36: -0.083436   38: -0.0970544   114: -0.00863707   115: -0.01262   127: -0.0156265   128: -0.00823901   131: 0.0732047
Row 132:   33: 0.0908447   35: 0.0942853   37: -0.125577   38: -0.0595534   125: -0.00477781   126: -0.0179334   129: -0.0101105   130: -0.0134608   132: 0.0742602

6. Solve the system:

[12]:
gfu.vec.data = \
    a.mat.Inverse(freedofs=fes.FreeDofs()) * f.vec
Draw(gfu);

The Dirichlet boundary condition constrains some degrees of freedom. The argument fes.FreeDofs() indicates that only the remaining "free" degrees of freedom should participate in the linear solve.

You can examine the coefficient vector of solution if needed:

[13]:
print(gfu.vec)
       0
       0
       0
 0.0923179
       0
       0
       0
       0
       0
       0
       0
       0
 0.0578988
 0.0863391
 0.0954151
 0.0945036
 0.0888226
 0.0780489
 0.0595559
 0.0331334
 0.043177
 0.0384731
 0.0312766
 0.0180508
 0.0379191
 0.0426386
 0.0472136
 0.0758445
 0.0892927
 0.0932644
 0.0919339
  0.0863
 0.0722526
 0.0576401
 0.0666805
 0.0710587
 0.0846294
 0.0618883
 0.0766818
       0
 -0.00563337
       0
       0
 0.0203536
       0
 -0.0350728
 0.00280938
 -0.00728756
 -0.00137235
       0
 -0.00966674
 -0.0121644
       0
 -0.0238901
 -0.0131488
       0
 -0.0161892
 -0.0200118
 -0.00657505
 -0.0212905
       0
 -0.0255408
 -0.00925383
       0
 -0.0375247
 -0.0131487
       0
 -0.0261275
 -0.0183435
 -0.0333457
 -0.0247246
 -0.0242762
 -0.00473494
 -0.0120835
 -0.0145349
 -0.00640155
 -0.0108854
 -0.00549118
 -0.00937973
 -0.00705122
 -0.0052206
 -0.00304656
 -0.00758425
 0.00107111
 0.00110145
 -0.00824163
 0.000388213
 0.00186718
 -0.0082305
 -0.000562953
 0.00136394
 0.00464277
 -0.00573891
 -0.0141064
 -0.0092274
 -0.00723695
 -0.0140782
 -0.012669
 -0.0071093
 -0.0248443
 -0.00762662
 -0.00336514
 -0.00878309
 -0.018877
 -0.0045665
 -0.00280081
 -0.0126613
 -0.0165067
 -0.0151556
 -0.0190153
 -0.0164166
 -0.00385688
 -0.0106348
 -0.0142227
 -0.0122429
 -0.00847428
 -0.00104741
 -0.00422534
 -0.00710203
 -0.00456872
 -0.00792875
 -0.00279066
 -0.00418665
 -0.00609652
 -0.00917816
 -0.0111537
 -0.00714295
 -0.00950903
 -0.0114547
 -0.00372816
 -0.0141414
 -0.00887479
 -0.00923917


You can see the zeros coming from the zero boundary conditions.

Ways to interact with NGSolve

  • A jupyter notebook (like this one) gives you one way to interact with NGSolve. When you have a complex sequence of tasks to perform, the notebook may not be adequate.

  • You can write an entire python module in a text editor and call python on the command line. (A script of the above is provided in poisson.py.) python3 poisson.py

  • If you want the Netgen GUI, then use netgen on the command line: netgen poisson.py You can then ask for a python shell from the GUI’s menu options (Solve -> Python shell).