21.5. Iterative Solvers#
So far we have used direct solvers to solve the linear system of equations. Although a direct solver can profit from the sparse matrix, it’s arithmetic complexity is sub-optimal. For large-scale problems iterative solvers are a must.
The conjugate gradient (cg) method is the standard method for symmetric and positive definite matrices. It’s convergence rate depends on a preconditioner, what is a cheap approximative inverse to the matrix.
from ngsolve import *
from ngsolve.webgui import Draw
We generate a 3D geometry and mesh using the OCC constructive solid geometry (CSG) modeler:
from netgen.occ import *
cube = Box((0,0,0),(1,1,1))
cyl = Cylinder((0,0.5,0.5),X, r=0.2, h=1)
cube.faces.name = "outer"
cyl.faces.name = "cyl"
shape = cube-cyl
ea = { "euler_angles" : [-180, -80, 165] }
Draw(shape, **ea);
ngmesh = OCCGeometry(shape).GenerateMesh(maxh=0.1)
mesh = Mesh(ngmesh).Curve(3)
for l in range(0):
mesh.Refine()
Draw (mesh, **ea);
fes = H1(mesh, order=3, dirichlet="outer", wb_withedges=False)
print ("we have", fes.ndof, "unknowns")
u,v = fes.TnT()
a = BilinearForm(grad(u)*grad(v)*dx)
f = LinearForm(1*v*dx)
# pre = preconditioners.MultiGrid(a)
pre = preconditioners.Local(a)
# pre = preconditioners.BDDC(a)
gfu = GridFunction(fes)
we have 23787 unknowns
assemble system and setup preconditioner in parallel:
with TaskManager():
a.Assemble()
f.Assemble()
solve the system using the preconditioned conjugate gradient method:
from ngsolve.krylovspace import CGSolver
from time import time
ts = time()
with TaskManager():
# inv = a.mat.Inverse(inverse="sparsecholesky", freedofs=fes.FreeDofs())
inv = CGSolver(mat=a.mat, pre=pre, printrates='\r', maxiter=1000)
gfu.vec.data = inv * f.vec
te = time()
print ("needed", te-ts, "seconds for", fes.ndof, "unknowns")
CG iteration 1, residual = 0.047739949075105856
CG iteration 2, residual = 0.05294727938543792
CG iteration 3, residual = 0.05500960698375211
CG iteration 4, residual = 0.04606877270694542
CG iteration 5, residual = 0.029876107300354116
CG iteration 6, residual = 0.023382706739268
CG iteration 7, residual = 0.018047146455209783
CG iteration 8, residual = 0.012650995831141893
CG iteration 9, residual = 0.008947539939651404
CG iteration 10, residual = 0.007134209747036224
CG iteration 11, residual = 0.006486603723518199
CG iteration 12, residual = 0.006445337337709603
CG iteration 13, residual = 0.006679235208108392
CG iteration 14, residual = 0.006938067131387422
CG iteration 15, residual = 0.00670502457467798
CG iteration 16, residual = 0.005583967413619987
CG iteration 17, residual = 0.004207964726985903
CG iteration 18, residual = 0.0030777087432416264
CG iteration 19, residual = 0.002292553710133819
CG iteration 20, residual = 0.0017210286974316954
CG iteration 21, residual = 0.0013918891862596495
CG iteration 22, residual = 0.0012867256293030037
CG iteration 23, residual = 0.0012207824425805743
CG iteration 24, residual = 0.0011508422442719894
CG iteration 25, residual = 0.0010326259360892076
CG iteration 26, residual = 0.0008650872416098794
CG iteration 27, residual = 0.0006834122553817224
CG iteration 28, residual = 0.0005203227711120672
CG iteration 29, residual = 0.00038568935816493373
CG iteration 30, residual = 0.0003028901691673216
CG iteration 31, residual = 0.0002468748185900818
CG iteration 32, residual = 0.0002005777049766653
CG iteration 33, residual = 0.00016399701908886814
CG iteration 34, residual = 0.00013686401056622853
CG iteration 35, residual = 0.0001141336130356841
CG iteration 36, residual = 9.45136943161346e-05
CG iteration 37, residual = 8.065214654944001e-05
CG iteration 38, residual = 7.15081833603998e-05
CG iteration 39, residual = 6.345362194690671e-05
CG iteration 40, residual = 5.5982876737655796e-05
CG iteration 41, residual = 4.817075802670549e-05
CG iteration 42, residual = 3.927699677793063e-05
CG iteration 43, residual = 3.084502899259499e-05
CG iteration 44, residual = 2.3441691859581664e-05
CG iteration 45, residual = 1.8068091668968058e-05
CG iteration 46, residual = 1.434983360644625e-05
CG iteration 47, residual = 1.1766739852763844e-05
CG iteration 48, residual = 9.964987268310976e-06
CG iteration 49, residual = 9.056719106388729e-06
CG iteration 50, residual = 8.561745619842757e-06
CG iteration 51, residual = 8.026554471089204e-06
CG iteration 52, residual = 7.09753603319757e-06
CG iteration 53, residual = 5.944747573670302e-06
CG iteration 54, residual = 4.55591009987537e-06
CG iteration 55, residual = 3.381907062916226e-06
CG iteration 56, residual = 2.4711107258406186e-06
CG iteration 57, residual = 1.8477611726611023e-06
CG iteration 58, residual = 1.4857145554661993e-06
CG iteration 59, residual = 1.2825964453132414e-06
CG iteration 60, residual = 1.1972381968942671e-06
CG iteration 61, residual = 1.1728671909923243e-06
CG iteration 62, residual = 1.126179049197469e-06
CG iteration 63, residual = 9.894511688700751e-07
CG iteration 64, residual = 7.958772966268723e-07
CG iteration 65, residual = 5.90160246772615e-07
CG iteration 66, residual = 4.353425978780637e-07
CG iteration 67, residual = 3.323706447740356e-07
CG iteration 68, residual = 2.6987925058119707e-07
CG iteration 69, residual = 2.302861222710536e-07
CG iteration 70, residual = 2.083853711174417e-07
CG iteration 71, residual = 1.9228490301356692e-07
CG iteration 72, residual = 1.7336077681121774e-07
CG iteration 73, residual = 1.4925977570475138e-07
CG iteration 74, residual = 1.2342270369710318e-07
CG iteration 75, residual = 1.0181258757675431e-07
CG iteration 76, residual = 8.45073856458963e-08
CG iteration 77, residual = 7.132426090697321e-08
CG iteration 78, residual = 6.050788344741731e-08
CG iteration 79, residual = 4.996349738876612e-08
CG iteration 80, residual = 4.0097491801673234e-08
CG iteration 81, residual = 3.144267233626898e-08
CG iteration 82, residual = 2.5066837296132902e-08
CG iteration 83, residual = 2.032771136151076e-08
CG iteration 84, residual = 1.750910849868475e-08
CG iteration 85, residual = 1.6082330771704167e-08
CG iteration 86, residual = 1.507852595383362e-08
CG iteration 87, residual = 1.3724650482122064e-08
CG iteration 88, residual = 1.2121199971090606e-08
CG iteration 89, residual = 1.008733464423273e-08
CG iteration 90, residual = 7.887163683063633e-09
CG iteration 91, residual = 5.86826803264321e-09
CG iteration 92, residual = 4.382454388703578e-09
CG iteration 93, residual = 3.387794835746971e-09
CG iteration 94, residual = 2.8340591227470708e-09
CG iteration 95, residual = 2.525311709921613e-09
CG iteration 96, residual = 2.409119551094e-09
CG iteration 97, residual = 2.339112405901695e-09
CG iteration 98, residual = 2.175314801856461e-09
CG iteration 99, residual = 1.903411323767574e-09
CG iteration 100, residual = 1.5439785580918631e-09
CG iteration 101, residual = 1.1763410478301127e-09
CG iteration 102, residual = 8.702406596584805e-10
CG iteration 103, residual = 6.43310398914805e-10
CG iteration 104, residual = 4.873018740096167e-10
CG iteration 105, residual = 3.9075295014215546e-10
CG iteration 106, residual = 3.3945265176337165e-10
CG iteration 107, residual = 3.151626969794772e-10
CG iteration 108, residual = 3.020145149625444e-10
CG iteration 109, residual = 2.8052314362462527e-10
CG iteration 110, residual = 2.442468495530773e-10
CG iteration 111, residual = 1.9591035782560552e-10
CG iteration 112, residual = 1.5149546389735437e-10
CG iteration 113, residual = 1.154754124727927e-10
CG iteration 114, residual = 8.664783149793581e-11
CG iteration 115, residual = 6.77217996187556e-11
CG iteration 116, residual = 5.466584010371964e-11
CG iteration 117, residual = 4.73614644500883e-11
CG iteration 118, residual = 4.2935538816882465e-11
CG iteration 119, residual = 3.9497856353158464e-11
CG iteration 120, residual = 3.619754117340719e-11
CG iteration 121, residual = 3.207356014002957e-11
CG iteration 122, residual = 2.702375129240235e-11
CG iteration 123, residual = 2.2041581493077788e-11
CG iteration 124, residual = 1.764558901491903e-11
CG iteration 125, residual = 1.3711072509145261e-11
CG iteration 126, residual = 1.0698819315453693e-11
CG iteration 127, residual = 8.390185438907755e-12
CG iteration 128, residual = 6.6420888354799266e-12
CG iteration 129, residual = 5.670529486054028e-12
CG iteration 130, residual = 5.195182434110249e-12
CG iteration 131, residual = 4.987197782596553e-12
CG iteration 132, residual = 4.683204131847481e-12
CG iteration 133, residual = 4.1873935225983635e-12
CG iteration 134, residual = 3.444689911333205e-12
CG iteration 135, residual = 2.61789558212488e-12
CG iteration 136, residual = 1.900146994721986e-12
CG iteration 137, residual = 1.3673930147381892e-12
CG iteration 138, residual = 9.948828724576634e-13
CG iteration 139, residual = 7.754671050264385e-13
CG iteration 140, residual = 6.747101210732056e-13
CG iteration 141, residual = 6.420861432933858e-13
CG iteration 142, residual = 6.333470586398766e-13
CG iteration 143, residual = 5.859148464051114e-13
CG iteration 144, residual = 5.100313624886648e-13
CG iteration 145, residual = 4.076752248328444e-13
CG iteration 146, residual = 3.043123211880837e-13
CG iteration 147, residual = 2.1879197346135391e-13
CG iteration 148, residual = 1.6124693349052223e-13
CG iteration 149, residual = 1.2394436474754466e-13
CG iteration 150, residual = 1.0401909929593768e-13
CG iteration 151, residual = 9.537471937105421e-14
CG iteration 152, residual = 9.073660846892042e-14
CG iteration 153, residual = 8.756141320012737e-14
CG iteration 154, residual = 8.170659984649713e-14
CG iteration 155, residual = 7.020705678530655e-14
CG iteration 156, residual = 5.5938970902465504e-14
CG iteration 157, residual = 4.200594152666575e-14
CG converged in 157 iterations to residual 4.200594152666575e-14
needed 0.01515507698059082 seconds for 23787 unknowns
Draw (gfu, **ea);