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();
[ ]:
[ ]: