This page was generated from wta/adaptivity.ipynb.

Adaptivity

[1]:
from ngsolve import *
from ngsolve.webgui import Draw
from netgen.geom2d import SplineGeometry

Define a geometry by curves:

[2]:
#   point numbers 0, 1, ... 11
#   sub-domain numbers (1), (2), (3)
#
#
#             7-------------6
#             |             |
#             |     (2)     |
#             |             |
#      3------4-------------5------2
#      |                           |
#      |             11            |
#      |           /   \           |
#      |         10 (3) 9          |
#      |           \   /     (1)   |
#      |             8             |
#      |                           |
#      0---------------------------1
#

def MakeGeometry():
    geometry = SplineGeometry()

    # point coordinates ...
    pnts = [ (0,0), (1,0), (1,0.6), (0,0.6), \
             (0.2,0.6), (0.8,0.6), (0.8,0.8), (0.2,0.8), \
             (0.5,0.15), (0.65,0.3), (0.5,0.45), (0.35,0.3) ]
    pnums = [geometry.AppendPoint(*p) for p in pnts]

    # start-point, end-point, boundary-condition, left-domain, right-domain:
    lines = [ (0,1,"bot",1,0), (1,2,"outer",1,0), (2,5,"outer",1,0), (5,4,"inner",1,2), (4,3,"outer",1,0), (3,0,"outer",1,0), \
              (5,6,"outer",2,0), (6,7,"outer",2,0), (7,4,"outer",2,0), \
              (8,9,"inner",3,1), (9,10,"inner",3,1), (10,11,"inner",3,1), (11,8,"inner",3,1) ]

    for p1,p2,bc,left,right in lines:
        geometry.Append(["line", pnums[p1], pnums[p2]], bc=bc, leftdomain=left, rightdomain=right)

    geometry.SetMaterial(1,"base")
    geometry.SetMaterial(2,"top")
    geometry.SetMaterial(3,"chip")

    return geometry

geo = MakeGeometry()
# Draw(geo)

Piece-wise constant coefficients in sub-domains:

[3]:
mesh = Mesh(geo.GenerateMesh(maxh=0.2))

fes = H1(mesh, order=3, dirichlet="bot", autoupdate=True)
u, v = fes.TnT()

lam = CoefficientFunction([1, 1000, 10])
a = BilinearForm(fes)
a += lam*grad(u)*grad(v)*dx

# heat-source in inner subdomain
f = LinearForm(fes)
f += 1*v*dx(definedon="chip")

c = Preconditioner(a, type="multigrid", inverse="sparsecholesky")

gfu = GridFunction(fes, autoupdate=True)
 Generate Mesh from spline geometry
 Boundary mesh done, np = 24
 CalcLocalH: 24 Points 0 Elements 0 Surface Elements
 Meshing domain 1 / 3
 load internal triangle rules
 Surface meshing done
 Meshing domain 2 / 3
 Surface meshing done
 Meshing domain 3 / 3
 Surface meshing done
 Edgeswapping, topological
 Smoothing
 Split improve
 Combine improve
 Smoothing
 Edgeswapping, metric
 Smoothing
 Split improve
 Combine improve
 Smoothing
 Edgeswapping, metric
 Smoothing
 Split improve
 Combine improve
 Smoothing
 Update mesh topology
 Update clusters
creating low order biform on demand

Assemble and solve problem:

[4]:
def SolveBVP():
    a.Assemble()
    f.Assemble()
    inv = CGSolver(a.mat, c.mat)
    gfu.vec.data = inv * f.vec

SolveBVP()
Draw (gfu, mesh)
assemble VOL element 40/40
assemble VOL element 40/40
BlockJacobi Preconditioner, constructor called, #blocks = 70
BlockJacobi Preconditioner built
assemble VOL element 40/40
0 0.0243924
1 0.00169156
2 0.000314731
3 0.000109161
4 2.02152e-05
5 4.19968e-06
6 8.39502e-07
7 1.9635e-07
8 3.41343e-08
9 5.4817e-09
10 1.07336e-09
11 1.98063e-10
[4]:
BaseWebGuiScene

Gradient recovery error estimator: Interpolate finite element flux

\[q_h := I_h (\lambda \nabla u_h)\]

and take difference as element error indicator:

\[\eta_T := \tfrac{1}{\lambda} \| q_h - \lambda \nabla u_h \|_{L_2(T)}^2\]
[5]:
l = []    # l = list of estimated total error
space_flux = HDiv(mesh, order=2, autoupdate=True)
gf_flux = GridFunction(space_flux, "flux", autoupdate=True)

def CalcError():

    # FEM-flux
    flux = lam * grad(gfu)
    # interpolate into H(div)
    gf_flux.Set(flux)

    # compute estimator:
    err = 1/lam*(flux-gf_flux)*(flux-gf_flux)
    eta2 = Integrate(err, mesh, VOL, element_wise=True)
    l.append ((fes.ndof, sqrt(sum(eta2))))
    print("ndof =", fes.ndof, " toterr =", sqrt(sum(eta2)))

    # mark for refinement:
    maxerr = max(eta2)
    for el in mesh.Elements():
        mesh.SetRefinementFlag(el, eta2[el.nr] > 0.25*maxerr)

CalcError()
setvalues element 40/40
ndof = 208  toterr = 0.004786183577239372

Adaptive loop:

[6]:
level = 0
while fes.ndof < 50000:
    mesh.Refine()
    SolveBVP()
    CalcError()
    level = level+1
    if level%5 == 0:
        Draw (gfu)
 Mesh bisection
 resetting marked-element information
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 94/94
assemble VOL element 94/94
BlockJacobi Preconditioner, constructor called, #blocks = 152
BlockJacobi Preconditioner built
assemble VOL element 94/94
0 0.0243626
1 0.00199073
2 0.00035584
3 7.39726e-05
4 2.02669e-05
5 4.64731e-06
6 8.62037e-07
7 1.9379e-07
8 4.63641e-08
9 1.1328e-08
10 2.4138e-09
11 4.85511e-10
12 9.04996e-11
setvalues element 94/94
ndof = 454  toterr = 0.003608150849132442
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 114/114
assemble VOL element 114/114
BlockJacobi Preconditioner, constructor called, #blocks = 183
BlockJacobi Preconditioner built
assemble VOL element 114/114
0 0.0244433
1 0.00184925
2 0.000293038
3 7.68055e-05
4 1.52396e-05
5 3.23129e-06
6 5.84066e-07
7 1.09074e-07
8 2.24352e-08
9 4.39687e-09
10 8.72167e-10
11 1.82363e-10
setvalues element 114/114
ndof = 547  toterr = 0.00298640856132947
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 252/252
assemble VOL element 252/252
BlockJacobi Preconditioner, constructor called, #blocks = 390
BlockJacobi Preconditioner built
assemble VOL element 252/252
0 0.0243985
1 0.00192216
2 0.000285485
3 5.19631e-05
4 1.25468e-05
5 2.70349e-06
6 5.62839e-07
7 1.20442e-07
8 2.60699e-08
9 5.70759e-09
10 1.22883e-09
11 2.4496e-10
12 4.99458e-11
setvalues element 252/252
ndof = 1168  toterr = 0.0020124591768725672
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 432/432
assemble VOL element 432/432
BlockJacobi Preconditioner, constructor called, #blocks = 662
BlockJacobi Preconditioner built
assemble VOL element 432/432
0 0.0244232
1 0.00190806
2 0.000320351
3 4.54737e-05
4 1.1129e-05
5 2.29043e-06
6 4.85515e-07
7 9.83544e-08
8 2.1052e-08
9 4.31087e-09
10 1.00064e-09
11 2.02053e-10
setvalues element 432/432
ndof = 1984  toterr = 0.0012753055579931864
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 632/632
assemble VOL element 632/632
BlockJacobi Preconditioner, constructor called, #blocks = 964
BlockJacobi Preconditioner built
assemble VOL element 632/632
0 0.0244416
1 0.00185464
2 0.000349386
3 5.29323e-05
4 1.39633e-05
5 2.72756e-06
6 6.37959e-07
7 1.40696e-07
8 3.10923e-08
9 6.83372e-09
10 1.37585e-09
11 2.97984e-10
12 6.98171e-11
setvalues element 632/632
ndof = 2890  toterr = 0.0007999943134523451
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 824/824
assemble VOL element 824/824
BlockJacobi Preconditioner, constructor called, #blocks = 1254
BlockJacobi Preconditioner built
assemble VOL element 824/824
0 0.0244612
1 0.00175702
2 0.000350428
3 5.36635e-05
4 1.07209e-05
5 2.1686e-06
6 4.89429e-07
7 1.23308e-07
8 2.28523e-08
9 4.5923e-09
10 9.89557e-10
11 2.39916e-10
setvalues element 824/824
ndof = 3760  toterr = 0.0005088041990863475
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1002/1002
assemble VOL element 1002/1002
BlockJacobi Preconditioner, constructor called, #blocks = 1523
BlockJacobi Preconditioner built
assemble VOL element 1002/1002
0 0.0244803
1 0.00165973
2 0.00032907
3 5.9812e-05
4 1.38846e-05
5 3.04176e-06
6 6.73464e-07
7 1.64736e-07
8 3.38292e-08
9 7.69336e-09
10 1.70502e-09
11 2.95996e-10
12 6.11942e-11
setvalues element 1002/1002
ndof = 4567  toterr = 0.0003319535330953844
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1136/1136
assemble VOL element 1136/1136
BlockJacobi Preconditioner, constructor called, #blocks = 1726
BlockJacobi Preconditioner built
assemble VOL element 1136/1136
0 0.0244973
1 0.00156936
2 0.000297071
3 5.94095e-05
4 1.27087e-05
5 2.37632e-06
6 3.85673e-07
7 9.47142e-08
8 2.31579e-08
9 4.34395e-09
10 8.25615e-10
11 1.64937e-10
setvalues element 1136/1136
ndof = 5176  toterr = 0.0002370558087353895
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1278/1278
assemble VOL element 1278/1278
BlockJacobi Preconditioner, constructor called, #blocks = 1941
BlockJacobi Preconditioner built
assemble VOL element 1278/1278
0 0.0245123
1 0.00148848
2 0.000265622
3 5.74906e-05
4 1.33482e-05
5 2.70499e-06
6 6.053e-07
7 1.54969e-07
8 2.89651e-08
9 7.39195e-09
10 1.63036e-09
11 2.9092e-10
12 5.51084e-11
setvalues element 1278/1278
ndof = 5821  toterr = 0.0001733174243150628
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1426/1426
assemble VOL element 1426/1426
BlockJacobi Preconditioner, constructor called, #blocks = 2165
BlockJacobi Preconditioner built
assemble VOL element 1426/1426
0 0.0245277
1 0.00140869
2 0.000238603
3 5.46793e-05
4 1.20785e-05
5 2.09437e-06
6 3.92755e-07
7 1.05542e-07
8 2.43948e-08
9 4.54597e-09
10 9.38379e-10
11 1.84273e-10
setvalues element 1426/1426
ndof = 6493  toterr = 0.00012667678504619132
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1564/1564
assemble VOL element 1564/1564
BlockJacobi Preconditioner, constructor called, #blocks = 2374
BlockJacobi Preconditioner built
assemble VOL element 1564/1564
0 0.0245375
1 0.00134271
2 0.000215239
3 4.97489e-05
4 1.15651e-05
5 2.14975e-06
6 5.28613e-07
7 1.19138e-07
8 2.37197e-08
9 6.07934e-09
10 1.35272e-09
11 2.75511e-10
12 5.29853e-11
setvalues element 1564/1564
ndof = 7120  toterr = 0.00010246756794703842
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1720/1720
assemble VOL element 1720/1720
BlockJacobi Preconditioner, constructor called, #blocks = 2610
BlockJacobi Preconditioner built
assemble VOL element 1720/1720
0 0.0245441
1 0.00128233
2 0.000195272
3 4.47285e-05
4 1.03398e-05
5 1.64122e-06
6 2.56625e-07
7 5.54726e-08
8 1.16608e-08
9 2.56537e-09
10 5.2822e-10
11 1.05347e-10
setvalues element 1720/1720
ndof = 7828  toterr = 7.662553021200613e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 1895/1895
assemble VOL element 1895/1895
BlockJacobi Preconditioner, constructor called, #blocks = 2875
BlockJacobi Preconditioner built
assemble VOL element 1895/1895
0 0.0245534
1 0.00121685
2 0.000177665
3 4.07126e-05
4 9.60478e-06
5 1.93223e-06
6 4.96499e-07
7 1.21346e-07
8 2.16962e-08
9 4.44325e-09
10 8.67857e-10
11 1.99143e-10
setvalues element 1895/1895
ndof = 8623  toterr = 6.10522895701154e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 2083/2083
assemble VOL element 2083/2083
BlockJacobi Preconditioner, constructor called, #blocks = 3160
BlockJacobi Preconditioner built
assemble VOL element 2083/2083
0 0.0245614
1 0.00116158
2 0.000160018
3 3.58836e-05
4 8.42786e-06
5 1.58151e-06
6 2.49225e-07
7 5.2854e-08
8 1.31034e-08
9 2.87776e-09
10 6.68019e-10
11 1.41589e-10
setvalues element 2083/2083
ndof = 9478  toterr = 4.6717580837311625e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 2291/2291
assemble VOL element 2291/2291
BlockJacobi Preconditioner, constructor called, #blocks = 3476
BlockJacobi Preconditioner built
assemble VOL element 2291/2291
0 0.0245668
1 0.00111295
2 0.00014845
3 3.3192e-05
4 7.81181e-06
5 1.5827e-06
6 2.61757e-07
7 5.81601e-08
8 1.31129e-08
9 2.89267e-09
10 6.93718e-10
11 1.54129e-10
setvalues element 2291/2291
ndof = 10426  toterr = 3.678705750643948e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 2535/2535
assemble VOL element 2535/2535
BlockJacobi Preconditioner, constructor called, #blocks = 3845
BlockJacobi Preconditioner built
assemble VOL element 2535/2535
0 0.0245735
1 0.00106631
2 0.000136208
3 2.92387e-05
4 6.74726e-06
5 1.48215e-06
6 2.7995e-07
7 6.67623e-08
8 1.33643e-08
9 2.79308e-09
10 6.0621e-10
11 1.47059e-10
setvalues element 2535/2535
ndof = 11533  toterr = 2.9159799090287877e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 2814/2814
assemble VOL element 2814/2814
BlockJacobi Preconditioner, constructor called, #blocks = 4266
BlockJacobi Preconditioner built
assemble VOL element 2814/2814
0 0.0245791
1 0.00102457
2 0.000125873
3 2.62789e-05
4 6.07004e-06
5 1.35175e-06
6 2.30576e-07
7 5.11522e-08
8 1.21246e-08
9 2.58303e-09
10 6.47232e-10
11 1.48125e-10
setvalues element 2814/2814
ndof = 12796  toterr = 2.4079046532700772e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 3159/3159
assemble VOL element 3159/3159
BlockJacobi Preconditioner, constructor called, #blocks = 4786
BlockJacobi Preconditioner built
assemble VOL element 3159/3159
0 0.0245807
1 0.000988102
2 0.000115148
3 2.36182e-05
4 5.89632e-06
5 1.43815e-06
6 3.11396e-07
7 8.61875e-08
8 1.95442e-08
9 4.55646e-09
10 1.09217e-09
11 2.48427e-10
12 5.46581e-11
setvalues element 3159/3159
ndof = 14356  toterr = 1.908780843744537e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 3594/3594
assemble VOL element 3594/3594
BlockJacobi Preconditioner, constructor called, #blocks = 5441
BlockJacobi Preconditioner built
assemble VOL element 3594/3594
0 0.0245835
1 0.00095948
2 0.000105574
3 2.01372e-05
4 4.71568e-06
5 1.12272e-06
6 2.06175e-07
7 4.7191e-08
8 1.38603e-08
9 2.72121e-09
10 6.1004e-10
11 1.2681e-10
setvalues element 3594/3594
ndof = 16321  toterr = 1.5036749761423592e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 4182/4182
assemble VOL element 4182/4182
BlockJacobi Preconditioner, constructor called, #blocks = 6326
BlockJacobi Preconditioner built
assemble VOL element 4182/4182
0 0.0245845
1 0.000935179
2 9.8201e-05
3 1.85591e-05
4 4.53092e-06
5 1.15051e-06
6 2.53877e-07
7 7.11748e-08
8 1.83314e-08
9 3.60416e-09
10 9.16097e-10
11 2.04771e-10
setvalues element 4182/4182
ndof = 18976  toterr = 1.1511866941302552e-05
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 4779/4779
assemble VOL element 4779/4779
BlockJacobi Preconditioner, constructor called, #blocks = 7226
BlockJacobi Preconditioner built
assemble VOL element 4779/4779
0 0.0245872
1 0.000917381
2 9.22341e-05
3 1.66616e-05
4 3.93915e-06
5 9.48473e-07
6 1.88918e-07
7 3.92544e-08
8 1.06925e-08
9 2.22243e-09
10 5.63505e-10
11 1.28875e-10
setvalues element 4779/4779
ndof = 21676  toterr = 9.043238592090386e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 5566/5566
assemble VOL element 5566/5566
BlockJacobi Preconditioner, constructor called, #blocks = 8411
BlockJacobi Preconditioner built
assemble VOL element 5566/5566
0 0.0245883
1 0.000902709
2 8.80385e-05
3 1.58904e-05
4 4.06811e-06
5 9.88866e-07
6 2.02727e-07
7 4.43748e-08
8 1.16099e-08
9 2.17764e-09
10 4.98842e-10
11 1.16186e-10
setvalues element 5566/5566
ndof = 25231  toterr = 6.843042535436214e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 6353/6353
assemble VOL element 6353/6353
BlockJacobi Preconditioner, constructor called, #blocks = 9598
BlockJacobi Preconditioner built
assemble VOL element 6353/6353
0 0.0245898
1 0.000890542
2 8.45294e-05
3 1.51754e-05
4 4.27446e-06
5 1.06971e-06
6 2.24378e-07
7 5.05773e-08
8 1.29039e-08
9 2.37305e-09
10 5.18125e-10
11 1.16242e-10
setvalues element 6353/6353
ndof = 28792  toterr = 5.400326682173985e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 7300/7300
assemble VOL element 7300/7300
BlockJacobi Preconditioner, constructor called, #blocks = 11024
BlockJacobi Preconditioner built
assemble VOL element 7300/7300
0 0.0245913
1 0.000878475
2 8.15029e-05
3 1.40316e-05
4 3.99785e-06
5 9.87042e-07
6 2.06822e-07
7 4.52749e-08
8 1.16512e-08
9 2.20106e-09
10 4.65288e-10
11 1.09218e-10
setvalues element 7300/7300
ndof = 33070  toterr = 4.250477494356433e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 8632/8632
assemble VOL element 8632/8632
BlockJacobi Preconditioner, constructor called, #blocks = 13028
BlockJacobi Preconditioner built
assemble VOL element 8632/8632
0 0.0245922
1 0.000870296
2 7.87564e-05
3 1.28009e-05
4 3.54784e-06
5 8.57007e-07
6 1.79659e-07
7 3.79862e-08
8 9.65713e-09
9 1.95162e-09
10 4.22539e-10
11 1.12628e-10
setvalues element 8632/8632
ndof = 39082  toterr = 3.2304065705096825e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 10105/10105
assemble VOL element 10105/10105
BlockJacobi Preconditioner, constructor called, #blocks = 15241
BlockJacobi Preconditioner built
assemble VOL element 10105/10105
0 0.0245928
1 0.000861719
2 7.63198e-05
3 1.14752e-05
4 2.96794e-06
5 7.08207e-07
6 1.47492e-07
7 2.98205e-08
8 7.68038e-09
9 1.78268e-09
10 4.12258e-10
11 1.12947e-10
setvalues element 10105/10105
ndof = 45721  toterr = 2.461487672200506e-06
 Mesh bisection
 Update mesh topology
 Update clusters
 Bisection done
 Update clusters
assemble VOL element 11653/11653
assemble VOL element 11653/11653
BlockJacobi Preconditioner, constructor called, #blocks = 17569
BlockJacobi Preconditioner built
assemble VOL element 11653/11653
0 0.0245936
1 0.000854761
2 7.39737e-05
3 1.00461e-05
4 2.22604e-06
5 5.28993e-07
6 1.18174e-07
7 2.53788e-08
8 6.49663e-09
9 1.56909e-09
10 3.45042e-10
11 8.9469e-11
setvalues element 11653/11653
ndof = 52705  toterr = 1.9366962586682258e-06
[7]:
Draw (gfu)
[7]:
BaseWebGuiScene
[8]:
%matplotlib inline
import matplotlib.pyplot as plt
plt.yscale('log')
plt.xscale('log')
plt.xlabel("ndof")
plt.ylabel("H1 error-estimate")
ndof,err = zip(*l)
plt.plot(ndof,err, "-*")

plt.ion()
plt.show();
../../_images/i-tutorials_wta_adaptivity_13_0.png
[ ]:

[ ]: